Ask Your Question
2

Calculation with arbitrary precision

asked 2021-07-29 02:45:37 +0100

Apprentice gravatar image

updated 2023-01-09 23:59:55 +0100

tmonteil gravatar image

I know I can ask Sage to display the digits of Pi with arbitrary precision with n(pi, digits= 202) for example. Is it possible to ask it to perform an operation -- like a numerical integral -- to an arbitrary precision? That is, by specifying the number of significative digits.

I tried:

n(numerical_integral(exp(-1/x)/x, 0, 1 ),digits=200)

but that obviously produce the wrong answer.

I also know I can do:

integrate(e * exp(-1/x)/x,x, 0, 1 ).n(300)

but that integral seems to be computed symbolically.

edit retag flag offensive close merge delete

2 Answers

Sort by » oldest newest most voted
1

answered 2021-07-29 14:10:04 +0100

tmonteil gravatar image

updated 2021-07-29 15:17:44 +0100

You can first solve the integral symbolically :

sage: i = integral(exp(-1/x)/x,x,0,1)
sage: i
-Ei(-1)

You can get information about the Ei function with:

sage: Ei?

Then, you can convert the result into some real field approximation with high precision:

sage: F = RealIntervalField(200)
sage: F
Real Interval Field with 200 bits of precision
sage: F(i)
0.21938393439552027367716377546012164903104729340690820757797849?

Note that, compared to RealField(200), you have the guaranty that the actual value belongs to some interval:

sage: F(i).endpoints()
(0.21938393439552027367716377546012164903104729340690820757797,
 0.21938393439552027367716377546012164903104729340690820757798)

Note : RealBallField(200) which is supposed to be faster fails with a RecursionError: maximum recursion depth exceeded This issue is now tracked at trac ticket 32301

edit flag offensive delete link more
1

answered 2021-07-29 14:18:27 +0100

Emmanuel Charpentier gravatar image

but that integral seems to be computed symbolically.

Nope. It's just that :

  1. it can't be expressed in terms of "elementary" functions, but that

  2. it occurs frequently enough to warrant the creation of a "special function" (the curse of engineering maths...).

Try :

sage: (e*exp(-1/x)/x).integrate(x)
-Ei(-1/x)*e

And read Ei?...

This integral can be checked :

sage: bool((-Ei(-1/x)*e).diff(x)==(e*exp(-1/x)/x))
True

Therefore, one can trust :

sage: (e*exp(-1/x)/x).integrate(x, 0, 1).n(digits=200)

0.5963473623231940743410784993692793760741778601525487815734849104823272191148744174704304970936
127603442370347484286236898120782995290571966173692226658940243185135143682937632962547711879740
2524323021

as far as one trusts the mutliprécision arithmetic used by Sage.

Now, if you want high precision on a numerical integral, numerical_integralhas parameterseps_absandeps_rel` allowing you to set the desired precision, as long as this precision can be achieved with (machine) double precision arithmetic. Therefore, this routine cannot give you 200-digit precision.

mpmath being part of Sagemath, you may try to use mpmath's quad routine after setting the relevant precision parameters (mp.prec and mp.dps), but I am not familiar with this direct use of mpmath and can't possibly comment...

Final note of caution : there is no "magic bullet" in numerical analysis. Analyse your problem and try to find "reasonable" algorithms. For example, compute a small quantity as the difference of two large quantities is probably not a good idea...

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2021-07-29 02:21:07 +0100

Seen: 356 times

Last updated: Jul 29 '21