HiveBrain v1.2.0
Get Started
← Back to all entries
patternMinor

Rational Approximation for e with Matlab

Submitted by: @import:stackexchange-codereview··
0
Viewed 0 times
matlabapproximationwithrationalfor

Problem

I've written a program that calculates the best rational approximation for e with a 1, 2, 3, 4, 5, and 6 digit denominator respectively (on Matlab). For a 5-digit denominator it takes about a minute. For a 6-digit denominator it takes 2 hours. I've tried to improve the code, but I am a beginner. Suggestions are kindly appreciated.

tic
clear
p=1;
q=1;
e=exp(1);
digits=1;
trial=1;
while numel(num2str(q))e
            q=q+1;
        end
        A(trial)=p/q;
        B(trial,1:2)=[p,q];
        trial=trial+1;
    end
    for i=1:trial-1
        error(i)=abs(A(i)-e);
    end
    [b,r]=min(error);
    p_q(digits,1:2)=B(r,1:2);
    digits=digits+1;
end
p_q
toc


This is the output:

p_q =
          19           7
         106          39
        1264         465
       25946        9545
      271801       99990
     1084483      398959

Elapsed time is 7232.153581 seconds.


Edit: Thanks to all who answered. Dennis and Jonas were right. It was the strings and growing matrices that slowed everything down. Besides I got two wrong values for p and q. The new code is this:

tic
p_best=1;
q_best=1;
e=exp(1);
for i=1:6
    for q=10^(i-1):10^i-1
        p=round(q*e);
        if abs(p/q-e)<abs(p_best/q_best-e)
            p_best=p;
            q_best=q;
        end
    end
    disp([p_best,q_best])
end
toc


The output is:

19           7
     193          71
    1457         536
   25946        9545
  271801       99990
 1084483      398959


Elapsed time is 0.054986 seconds.

It's faster and actually correct. I assumed I knew e like Dennis said. To Joni and hardmath: I can't use continued fractions because the teacher wanted us to use brute force, but thanks.

Solution

Continued fractions are an easy way to calculate good rational approximations to real numbers. You can find a disussion of an algorithm to convert a floating point number into a fraction in my blog, with code in JavaScript which should be easy to convert to matlab:

function float2rat(x) {
  var tolerance = 1.0E-6;
  var h1=1; var h2=0;
  var k1=0; var k2=1;
  var b = x;
  do {
     var a = Math.floor(b);
     var aux = h1; h1 = a*h1+h2; h2 = aux;
     aux = k1; k1 = a*k1+k2; k2 = aux;
     b = 1/(b-a);
  } while (Math.abs(x-h1/k1) > x*tolerance);
  return h1+"/"+k1;
}


If you're specifically interested in e you can take a shortcut and use the regularities in its continued fraction representation to produce rational approximations that are far more accurate than the machine epsilon. This method is described here.

Code Snippets

function float2rat(x) {
  var tolerance = 1.0E-6;
  var h1=1; var h2=0;
  var k1=0; var k2=1;
  var b = x;
  do {
     var a = Math.floor(b);
     var aux = h1; h1 = a*h1+h2; h2 = aux;
     aux = k1; k1 = a*k1+k2; k2 = aux;
     b = 1/(b-a);
  } while (Math.abs(x-h1/k1) > x*tolerance);
  return h1+"/"+k1;
}

Context

StackExchange Code Review Q#21379, answer score: 2

Revisions (0)

No revisions yet.