#!/usr/bin/perl

package Fraction;
use Exporter;

@ISA = qw(Exporter);
@EXPORT_OK = qw(decimal_to_continued_fraction
		continued_fraction_to_fraction
		decimal_to_fraction
		convergents
	       );

$DEFAULT_PRECISION = 2**(-24);

sub decimal_to_continued_fraction {
  my ($n, $e) = @_;
  $e = $DEFAULT_PRECISION unless defined $e;
  my $i = int $n;
  my $f = $n - $i;
#  print ">> $i + $f ($e)\n";
  if ($f <= $e) {
    return $i;
  } elsif (1-$f <= $e) {
    return $i + 1;
  }
  my @result = decimal_to_continued_fraction(1/$f, $e*$i);
  return ($i, @result);
}

sub continued_fraction_to_fraction {
  my @cf = reverse @_;
  my ($n, $d) = (shift @cf,1);
  while (@cf) {
    ($n, $d) = ($d, $n);
    $n += $d * shift(@cf);
  }
  ($n, $d);
}

sub convergents {
  my @cf = decimal_to_continued_fraction(@_);
  my @c;
  my @result;
  for (@cf) {
    push @c, $_;
    push @result, [continued_fraction_to_fraction(@c)];
  }
  @result;
}

unless (caller) {
  while ($r = <>) {
    @l = decimal_to_continued_fraction($r, 1e-14);
    ($n, $d) = continued_fraction_to_fraction(@l);
    print "That looks like $n/$d.\nLet me check.\n";
    $r2 = $n/$d;
    if ($r2 == $r) {
      print "Yup!\n";
    } elsif (abs($r2 - $r) < $DEFAULT_PRECISION) {
      print "Close enough.\n";
    } else {
      print "No, $n/$d is actually $r2.  Something is wrong.\n"
    }
  }
}

1;