AUTHORS:
This module provides easy access to many exponential integral special functions. It utilizes Maxima’s special functions package and the mpmath library.
REFERENCES:
[AS] | (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13) ‘Handbook of Mathematical Functions’, Milton Abramowitz and Irene A. Stegun, National Bureau of Standards Applied Mathematics Series, 55. See also http://www.math.sfu.ca/~cbm/aands/. |
AUTHORS:
Benjamin Jones
Implementations of the classes Function_exp_integral_*.
David Joyner and William Stein
Authors of the code which was moved from special.py and trans.py. Implementation of exp_int() (from sage/functions/special.py). Implementation of exponential_integral_1() (from sage/functions/transcendental.py).
Bases: sage.symbolic.function.BuiltinFunction
The trigonometric integral defined by
where is the Euler gamma constant (euler_gamma in Sage),
see [AS] 5.2.1.
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: cos_integral(3.0)
0.119629786008000
The alias can be used instead of
:
sage: Ci(3.0)
0.119629786008000
Compare cos_integral(3.0) to the definition of the value using numerical integration:
sage: N(euler_gamma + log(3.0) + integrate((cos(x)-1)/x, x, 0, 3.0) - cos_integral(3.0)) < 1e-14
True
Arbitrary precision and complex arguments are handled:
sage: N(cos_integral(3), digits=30)
0.119629786008000327626472281177
sage: cos_integral(ComplexField(100)(3+I))
0.078134230477495714401983633057 - 0.37814733904787920181190368789*I
The limit as
is zero:
sage: N(cos_integral(1e23))
-3.24053937643003e-24
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = cos_integral(x)
sage: f.diff(x)
cos(x)/x
sage: f.integrate(x)
x*cos_integral(x) - sin(x)
The Nielsen spiral is the parametric plot of (Si(t), Ci(t)):
sage: t=var('t')
sage: f(t) = sin_integral(t)
sage: g(t) = cos_integral(t)
sage: P = parametric_plot([f, g], (t, 0.5 ,20))
sage: show(P, frame=True, axes=False)
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Bases: sage.symbolic.function.BuiltinFunction
The trigonometric integral defined by
see [AS] 5.2.4.
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: cosh_integral(1.0)
0.837866940980208
The alias can be used instead of
:
sage: Chi(1.0)
0.837866940980208
Here is an example from the mpmath documentation:
sage: f(x) = cosh_integral(x)
sage: find_root(f, 0.1, 1.0)
0.5238225713894826
Compare cosh_integral(3.0) to the definition of the value using numerical integration:
sage: N(euler_gamma + log(3.0) + integrate((cosh(x)-1)/x, x, 0, 3.0) -
... cosh_integral(3.0)) < 1e-14
True
Arbitrary precision and complex arguments are handled:
sage: N(cosh_integral(3), digits=30)
4.96039209476560976029791763669
sage: cosh_integral(ComplexField(100)(3+I))
3.9096723099686417127843516794 + 3.0547519627014217273323873274*I
The limit of as
is
:
sage: N(cosh_integral(Infinity))
+infinity
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = cosh_integral(x)
sage: f.diff(x)
cosh(x)/x
sage: f.integrate(x)
x*cosh_integral(x) - sinh(x)
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Bases: sage.symbolic.function.BuiltinFunction
The generalized complex exponential integral Ei(z) defined by
for x > 0 and for complex arguments by analytic continuation, see [AS] 5.1.2.
EXAMPLES:
sage: Ei(10)
Ei(10)
sage: Ei(I)
Ei(I)
sage: Ei(3+I)
Ei(I + 3)
sage: Ei(1.3)
2.72139888023202
The branch cut for this function is along the negative real axis:
sage: Ei(-3 + 0.1*I)
-0.0129379427181693 + 3.13993830250942*I
sage: Ei(-3 - 0.1*I)
-0.0129379427181693 - 3.13993830250942*I
ALGORITHM: Uses mpmath.
TESTS:
Show that the evaluation and limit issue in trac ticket #13271 is fixed:
sage: var('Z')
Z
sage: (Ei(-Z)).limit(Z=oo)
0
sage: (Ei(-Z)).limit(Z=1000)
Ei(-1000)
sage: (Ei(-Z)).limit(Z=1000).n()
-5.07089306023517e-438
Bases: sage.symbolic.function.BuiltinFunction
The generalized complex exponential integral defined by
for complex numbers and
, see [AS] 5.1.4.
The special case where is denoted in Sage by
exp_integral_e1.
EXAMPLES:
Numerical evaluation is handled using mpmath:
sage: N(exp_integral_e(1,1))
0.219383934395520
sage: exp_integral_e(1, RealField(100)(1))
0.21938393439552027367716377546
We can compare this to PARI’s evaluation of exponential_integral_1():
sage: N(exponential_integral_1(1))
0.219383934395520
We can verify one case of [AS] 5.1.45, i.e.
:
sage: N(exp_integral_e(2, 3+I))
0.00354575823814662 - 0.00973200528288687*I
sage: N((3+I)*gamma(-1, 3+I))
0.00354575823814662 - 0.00973200528288687*I
Maxima returns the following improper integral as a multiple of exp_integral_e(1,1):
sage: uu = integral(e^(-x)*log(x+1),x,0,oo)
sage: uu
e*exp_integral_e(1, 1)
sage: uu.n(digits=30)
0.596347362323194074341078499369
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = exp_integral_e(2,x)
sage: f.diff(x)
-exp_integral_e(1, x)
sage: f.integrate(x)
-exp_integral_e(3, x)
sage: f = exp_integral_e(-1,x)
sage: f.integrate(x)
Ei(-x) - gamma(-1, x)
Some special values of exp_integral_e can be simplified. [AS] 5.1.23:
sage: exp_integral_e(0,x)
e^(-x)/x
[AS] 5.1.24:
sage: exp_integral_e(6,0)
1/5
sage: nn = var('nn')
sage: assume(nn > 1)
sage: f = exp_integral_e(nn,0)
sage: f.simplify()
1/(nn - 1)
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
Bases: sage.symbolic.function.BuiltinFunction
The generalized complex exponential integral defined by
see [AS] 5.1.4.
EXAMPLES:
Numerical evaluation is handled using mpmath:
sage: N(exp_integral_e1(1))
0.219383934395520
sage: exp_integral_e1(RealField(100)(1))
0.21938393439552027367716377546
We can compare this to PARI’s evaluation of exponential_integral_1():
sage: N(exp_integral_e1(2.0))
0.0489005107080611
sage: N(exponential_integral_1(2.0))
0.0489005107080611
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = exp_integral_e1(x)
sage: f.diff(x)
-e^(-x)/x
sage: f.integrate(x)
-exp_integral_e(2, x)
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
Bases: sage.symbolic.function.BuiltinFunction
The logarithmic integral defined by
for x > 1 and by analytic continuation for complex arguments z (see [AS] 5.1.3).
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: N(log_integral(3))
2.16358859466719
sage: N(log_integral(3), digits=30)
2.16358859466719197287692236735
sage: log_integral(ComplexField(100)(3+I))
2.2879892769816826157078450911 + 0.87232935488528370139883806779*I
sage: log_integral(0)
0
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = log_integral(x)
sage: f.diff(x)
1/log(x)
sage: f.integrate(x)
x*log_integral(x) - Ei(2*log(x))
Here is a test from the mpmath documentation. There are 1,925,320,391,606,803,968,923 many prime numbers less than 1e23. The value of log_integral(1e23) is very close to this:
sage: log_integral(1e23)
1.92532039161405e21
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Bases: sage.symbolic.function.BuiltinFunction
The offset logarithmic integral, or Eulerian logarithmic integral,
is defined by
for .
The offset logarithmic integral should also not be confused with the
polylogarithm (also denoted by ), which is
implemented as sage.functions.log.Function_polylog.
is identical to
except that
the lower limit of integration is
rather than
to avoid the
singularity at
of
See Function_log_integral for details of .
Thus
can also be represented by
So we have:
sage: li(4.5)-li(2.0)-Li(4.5)
0.000000000000000
is extended to complex arguments
by analytic continuation (see [AS] 5.1.3):
sage: Li(6.6+5.4*I)
3.97032201503632 + 2.62311237593572*I
The function is an approximation for the number of
primes up to
. In fact, the famous Riemann Hypothesis is
For “small” ,
is always slightly bigger
than
. However it is a theorem that there are very
large values of
(e.g., around
), such that
. See “A new bound for the
smallest x with
”,
Bays and Hudson, Mathematics of Computation, 69 (2000) 1285-1296.
Here is a test from the mpmath documentation. There are 1,925,320,391,606,803,968,923 prime numbers less than 1e23. The value of log_integral_offset(1e23) is very close to this:
sage: log_integral_offset(1e23)
1.92532039161405e21
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: N(log_integral_offset(3))
1.11842481454970
sage: N(log_integral_offset(3), digits=30)
1.11842481454969918803233347815
sage: log_integral_offset(ComplexField(100)(3+I))
1.2428254968641898308632562019 + 0.87232935488528370139883806779*I
sage: log_integral_offset(2)
0
sage: for n in range(1,7):
... print '%-10s%-10s%-20s'%(10^n, prime_pi(10^n), N(Li(10^n)))
10 4 5.12043572466980
100 25 29.0809778039621
1000 168 176.564494210035
10000 1229 1245.09205211927
100000 9592 9628.76383727068
1000000 78498 78626.5039956821
Symbolic derivatives are handled by Sage and integration by Maxima:
sage: x = var('x')
sage: f = log_integral_offset(x)
sage: f.diff(x)
1/log(x)
sage: f.integrate(x)
-x*log_integral(2) + x*log_integral(x) - Ei(2*log(x))
sage: Li(x).integrate(x,2.0,4.5)
-2.5*log_integral(2) + 5.79932114741
sage: N(f.integrate(x,2.0,3.0))
0.601621785860587
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Bases: sage.symbolic.function.BuiltinFunction
The trigonometric integral defined by
see [AS] 5.2.1.
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: sin_integral(0)
0
sage: sin_integral(0.0)
0.000000000000000
sage: sin_integral(3.0)
1.84865252799947
sage: N(sin_integral(3), digits=30)
1.84865252799946825639773025111
sage: sin_integral(ComplexField(100)(3+I))
2.0277151656451253616038525998 + 0.015210926166954211913653130271*I
The alias can be used instead of
:
sage: Si(3.0)
1.84865252799947
The limit of as
is
:
sage: N(sin_integral(1e23))
1.57079632679490
sage: N(pi/2)
1.57079632679490
At 200 bits of precision agrees with
up to
:
sage: sin_integral(RealField(200)(1e23))
1.5707963267948966192313288218697837425815368604836679189519
sage: N(pi/2, prec=200)
1.5707963267948966192313216916397514420985846996875529104875
The exponential sine integral is analytic everywhere:
sage: sin_integral(-1.0)
-0.946083070367183
sage: sin_integral(-2.0)
-1.60541297680269
sage: sin_integral(-1e23)
-1.57079632679490
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = sin_integral(x)
sage: f.diff(x)
sin(x)/x
sage: f.integrate(x)
x*sin_integral(x) + cos(x)
sage: integrate(sin(x)/x, x)
-1/2*I*Ei(I*x) + 1/2*I*Ei(-I*x)
Compare values of the functions and
, which are both anti-derivatives of
, at some random positive real numbers:
sage: f(x) = 1/2*I*Ei(-I*x) - 1/2*I*Ei(I*x) - pi/2
sage: g(x) = sin_integral(x)
sage: R = [ abs(RDF.random_element()) for i in range(100) ]
sage: all(abs(f(x) - g(x)) < 1e-10 for x in R)
True
The Nielsen spiral is the parametric plot of (Si(t), Ci(t)):
sage: x=var('x')
sage: f(x) = sin_integral(x)
sage: g(x) = cos_integral(x)
sage: P = parametric_plot([f, g], (x, 0.5 ,20))
sage: show(P, frame=True, axes=False)
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Bases: sage.symbolic.function.BuiltinFunction
The trigonometric integral defined by
see [AS] 5.2.3.
EXAMPLES:
Numerical evaluation for real and complex arguments is handled using mpmath:
sage: sinh_integral(3.0)
4.97344047585981
sage: sinh_integral(1.0)
1.05725087537573
sage: sinh_integral(-1.0)
-1.05725087537573
The alias can be used instead of
:
sage: Shi(3.0)
4.97344047585981
Compare sinh_integral(3.0) to the definition of the value using numerical integration:
sage: N(integrate((sinh(x))/x, x, 0, 3.0) - sinh_integral(3.0)) < 1e-14
True
Arbitrary precision and complex arguments are handled:
sage: N(sinh_integral(3), digits=30)
4.97344047585980679771041838252
sage: sinh_integral(ComplexField(100)(3+I))
3.9134623660329374406788354078 + 3.0427678212908839256360163759*I
The limit as
is
:
sage: N(sinh_integral(Infinity))
+infinity
Symbolic derivatives and integrals are handled by Sage and Maxima:
sage: x = var('x')
sage: f = sinh_integral(x)
sage: f.diff(x)
sinh(x)/x
sage: f.integrate(x)
x*sinh_integral(x) - cosh(x)
Note that due to some problems with the way Maxima handles these expressions, definite integrals can sometimes give unexpected results (typically when using inexact endpoints) due to inconsistent branching:
sage: integrate(sinh_integral(x), x, 0, 1/2)
-cosh(1/2) + 1/2*sinh_integral(1/2) + 1
sage: integrate(sinh_integral(x), x, 0, 1/2).n() # correct
0.125872409703453
sage: integrate(sinh_integral(x), x, 0, 0.5).n() # fixed in maxima 5.29.1
0.125872409703453
ALGORITHM:
Numerical evaluation is handled using mpmath, but symbolics are handled by Sage and Maxima.
REFERENCES:
Returns the exponential integral . If the optional
argument
is given, computes list of the first
values of the exponential integral
.
The exponential integral is
INPUT:
OUTPUT:
A real number if n is 0 (the default) or a list of reals if n > 0. The precision is the same as the input, with a default of 53 bits in case the input is exact.
EXAMPLES:
sage: exponential_integral_1(2)
0.0489005107080611
sage: exponential_integral_1(2,4) # abs tol 1e-18
[0.0489005107080611, 0.00377935240984891, 0.000360082452162659, 0.0000376656228439245]
sage: exponential_integral_1(40,5)
[1.03677326145166e-19, 2.22854325868847e-37, 6.33732515501151e-55, 2.02336191509997e-72, 6.88522610630764e-90]
sage: exponential_integral_1(0)
+Infinity
sage: r = exponential_integral_1(RealField(150)(1))
sage: r
0.21938393439552027367716377546012164903104729
sage: parent(r)
Real Field with 150 bits of precision
sage: exponential_integral_1(RealField(150)(100))
3.6835977616820321802351926205081189876552201e-46
TESTS:
The relative error for a single value should be less than 1 ulp:
sage: for prec in [20..1000]: # long time (22s on sage.math, 2013)
....: R = RealField(prec)
....: S = RealField(prec+64)
....: for t in range(8): # Try 8 values for each precision
....: a = R.random_element(-15,10).exp()
....: x = exponential_integral_1(a)
....: y = exponential_integral_1(S(a))
....: e = float(abs(S(x) - y)/x.ulp())
....: if e >= 1.0:
....: print "exponential_integral_1(%s) with precision %s has error of %s ulp"%(a, prec, e)
The absolute error for a vector should be less than , where
is the precision in bits of
and
:
sage: for prec in [20..128]: # long time (15s on sage.math, 2013)
....: R = RealField(prec)
....: S = RealField(prec+64)
....: a = R.random_element(-15,10).exp()
....: n = 2^ZZ.random_element(14)
....: x = exponential_integral_1(a, n)
....: y = exponential_integral_1(S(a), n)
....: c = RDF(2 * max(1.0, y[0]))
....: for i in range(n):
....: e = float(abs(S(x[i]) - y[i]) << prec)
....: if e >= c:
....: print "exponential_integral_1(%s, %s)[%s] with precision %s has error of %s >= %s"%(a, n, i, prec, e, c)
ALGORITHM: use the PARI C-library function eint1.
REFERENCE: