ResearcherDaily
ResearcherDaily

Reputation: 183

Writeln output extended data type as floating-point not as Scientific Expression

Good night, I'm having a issue with writeln on Embarcadero 10.2, I'm trying to write the sin() function as a Taylor expansion infinite serie.

Everything it's running fine, but the output is being in Scientific Expression like: 3.60448486921676E-0158 when the correct was 0.912945250727627654376099983845.

I need a 30 digits precision. My code below, It's a console program.

program Project1;

{$APPTYPE CONSOLE}

uses
  System.SysUtils,
  Windows,
  Math;

function fact(n: LongInt): extended;
begin
  if (n = 0) then
     fact := 1
  else
     fact := n * fact(n - 1);
end;

var
  one:integer;
  fin:extended;
  number:LongInt;
  i:integer;
  cons:LongInt;
  cons1:extended;
begin
  one := -1;
  writeln('Digite o angulo em radianos para o Seno: ');
  readln(number);
  for i := 1 to 120 do
  begin
    cons := (2*i)+1;
    if(i mod 2) = 0 then 
      one := 1
    else 
      one := -1;
    cons1 := fact(cons);
    fin := (one/cons1)*(power(number,cons));
    cons := 0;
  end;
  writeln(fin);
  readln;

end.

Upvotes: 3

Views: 232

Answers (3)

Victor Melo
Victor Melo

Reputation: 310

program Project1;

{$APPTYPE CONSOLE}

uses
  System.SysUtils;


function sine(n: real): real;
asm
fld n
fsin
fstp Result
end;

var
n:real;
begin
write('Write the number: ');
readln(n);
writeln(sine(n):10:11);
readln;
end.

Upvotes: -1

MBo
MBo

Reputation: 80197

To perform calculation with very long precision, you need some library for arbitrary precision arithmetics - use Delphi wrapper/interface to GMP or some Delphi library. Example. Another one.

Note that for libraries with only long integers support (not floats) you still can calculate series like this Python code (it has internal support for big integers. ** is power, // is div).
The only division (slow operation) is performed once.

d30 = 10**30
d60 = d30 * d30
x = 3141592653589793238462643383279 // 4   #Pi/4 * d30
xx = - x * x
denom = d30
nom = x
sum = nom
i = 2
while abs(nom) * d30 > denom:
    nom =  nom * xx
    mul = i * (i + 1) * d60
    denom = denom * mul
    sum = sum * mul + nom
    i += 2

print('0.' + str(d30 * sum // denom))

>> 0.707106780551955811469384356540

Approach is based on this: x - x^3/6 + x^5/120 = (x*120 - x^3 * 20 + x^5) / 120

Upvotes: 3

gammatester
gammatester

Reputation: 1141

Computing sin(x) with Taylor expansion is stable only for small arguments, say less than Pi/2. For larger arguments you suffer from catastrophic cancelation, because intermediate terms are very large, e.g. for x=100 the largest terms is the 49'th with a value 0.107151028812546692318354675952E43. All these large terms are added together and give a sin() value with magnitude <= 1. You would need about 70 decimal digits, do compute the result with about 30 digits.

The solutions is: Use range reduction and multiprecision floating-point routines (see mpf_sin in my opensource Pascal MPArith package). Using the interactive calculator you get

sin(100) = -0.5063656411097587936565576104597854321 in 0.2 ms.

Upvotes: 3

Related Questions