1
program cosseno;
var fat1, fat2, n1, n2, cont1, cont2, s :real;
begin
  s := 0.5;
  fat1 := 24;
  fat2 := 720;
  n1 := 1/fat1;
  n2 := 1/fat2;
  cont1 := 8;
  cont2 := 10;

  while (n1 - n2) > 0.000001 do
  begin
    fat1 := (cont1) * (cont1-1) * (cont1-2) * (cont1-3) * fat1;
    fat2 := (cont2) * (cont2-1) * (cont2-2) * (cont2-3) * fat2;
    cont1 := cont1 + 4;
    cont2 := cont2 + 4;
    n1 := n1 + 1/fat1;
    n2 := n2 + 1/fat2;
  end;

  s := s + n1 - n2;
  writeln(s);
end.

This is a program to estimate e value of cos(1 rad), apparently, the error occurs within the while section.

lurker
  • 56,987
  • 9
  • 69
  • 103
  • Welcome to SO:SE. Please use `code` formatting to display the code and make it legible. See [ask] and [How do I format my posts using Markdown or HTML?](http://stackoverflow.com/help/formatting). – mins Apr 07 '15 at 00:09
  • 2
    By running this code in mind only, I think fat2 will have some overlflow after a few loops. Put an WriteLn Statement in the loop to see what your last turn is. – Thomas Krojer Apr 07 '15 at 05:39
  • When you step through the code using the debugger, what does it show you? What line in the `while` loop causes the error? What are the values of the variables used at that point? – Ken White Apr 07 '15 at 22:00
  • You should check the source of the algorithm because it doesn't work. if you [ignore the FPU errors](http://docs.embarcadero.com/products/rad_studio/delphiAndcpp2009/HelpUpdate2/EN/html/delphivclwin32/System_Set8087CW.html) (add `Set8087CW($133f);` before entering the loop) the program will run endlessly. That's mean that the condition `while (n1 - n2) > 0.000001` is never met. – Abstract type Apr 12 '15 at 12:27

1 Answers1

0

The Problem

The issue is that your variable fat1 is getting too big and overflowing. With free pascal, your code causes a 205 error (Floating point overflow). This is similar to a 207 Invalid floating point operation - free pascal run-time error codes. The real types that are relevant to your code have the following constraints (from the free pascal data types wiki page):

Type       Range                   Significant digits   Bytes
Real       platform dependent      ???                  4 or 8
Single     1.5E-45 .. 3.4E38       7-8                  4
Double     5.0E-324 .. 1.7E308     15-16                8
Extended   1.9E-4932 .. 1.1E4932   19-20                10

The Work Around

You simply need to account for the possibility of overflow with your largest variable. It is not possible to simply check if the variable is above it's max value so you should do the reverse of the operation that could cause the overflow. I also used Extended to have a better estimate.

var max_fat : Extended = 1.1 * power(10,4932);
var should_step : Extended= max_fat / (cont1) / (cont1-1) / (cont1-2) / (cont1-3);

then in your while loop include the check:

fat1 < should_step

The result of this reverse operation is the maximum that the variable fat1 can be before it overflows. I use the Math library to help calculate the maximum Extended. The following test code works for me:

program cosseno;

Uses Math;

var fat1, fat2, n1, n2, cont1, cont2, s, max_fat, should_step : Extended;
   begin
      s := 0.5;
      fat1 := 24;
      fat2 := 720;
      n1 := 1/fat1;
      n2 := 1/fat2;
      cont1 := 8;
      cont2 := 10;

      max_fat := 1.1 * power(10,4932);

      should_step := max_fat;

      while ((n1 - n2) > 0.000001) AND (fat1 < should_step) do
       begin
          should_step := max_fat / (cont1) / (cont1-1) / (cont1-2) / (cont1-3);

          fat1 := (cont1) * (cont1-1) * (cont1-2) * (cont1-3) * fat1;
          fat2 := (cont2) * (cont2-1) * (cont2-2) * (cont2-3) * fat2;
          cont1 := cont1 + 4;
          cont2 := cont2 + 4;
          n1 := n1 + 1/fat1;
          n2 := n2 + 1/fat2;
          writeln('fat1  > ', fat1);
          writeln('fat2  > ', fat2);
          writeln('cont1 > ', cont1);
          writeln('cont2 > ', cont2);
          writeln('n1    > ', n1);
          writeln('n2    > ', n2);
       end;

      s := s + n1 - n2;
      writeln('<==================>');
      writeln(s);
   end.

The code yields the following result:

5.40302305868139717414E-0001

which I checked to be correct for at least the first 10 digits.

Caleb Adams
  • 4,445
  • 1
  • 26
  • 26