Mpmath provides a rich collection of mathematical functions.
The predefined objects j (imaginary unit), inf (positive infinity) and nan (not-a-number) are shortcuts to mpc and mpf instances with these fixed values.
Mpmath supports arbitrary-precision computation of various common (and less common) mathematical constants.
Symbol | Value | Description |
pi | 3.141592654 | Area of the unit circle |
degree | 0.017453292 | 1 degree of angle = pi/180 |
e | 2.718281828 | Base of the natural logarithm |
phi | 1.618033989 | The golden ratio, (1+sqrt(5))/2 |
ln2 | 0.693147180 | The natural logarithm of 2 |
ln10 | 2.302585093 | The natural logarithm of 10 |
euler | 0.577215664 | Euler’s constant (gamma) |
catalan | 0.915965594 | Catalan’s constant |
apery | 1.202056903 | Apery’s constant, zeta(3) |
khinchin | 2.685452001 | Khinchin’s constant |
glaisher | 1.282427129 | Glaisher’s constant |
These constants are implemented as lazy objects that can be computed with any precision. Whenever the objects are used as function arguments or as operands in arithmetic operations, they automagically evaluate to the current working precision. A lazy number can be converted to a regular mpf using the unary + operator:
>>> from mpmath import *
>>> mp.dps = 15
>>> pi
<pi: 3.14159~>
>>> 2*pi
mpf('6.2831853071795862')
>>> +pi
mpf('3.1415926535897931')
>>> mp.dps = 40
>>> pi
<pi: 3.14159~>
>>> 2*pi
mpf('6.283185307179586476925286766559005768394338')
>>> +pi
mpf('3.141592653589793238462643383279502884197169')
To show many digits of a number, you can just set the desired precision and print it:
>>> mp.dps = 1000
>>> print pi # doctest:+SKIP
3.141592653589793238462643383279502884197169399375105820974944592307816406286208
99862803482534211706798214808651328230664709384460955058223172535940812848111745
02841027019385211055596446229489549303819644288109756659334461284756482337867831
65271201909145648566923460348610454326648213393607260249141273724587006606315588
17488152092096282925409171536436789259036001133053054882046652138414695194151160
94330572703657595919530921861173819326117931051185480744623799627495673518857527
24891227938183011949129833673362440656643086021394946395224737190702179860943702
77053921717629317675238467481846766940513200056812714526356082778577134275778960
91736371787214684409012249534301465495853710507922796892589235420199561121290219
60864034418159813629774771309960518707211349999998372978049951059731732816096318
59502445945534690830264252230825334468503526193118817101000313783875288658753320
83814206171776691473035982534904287554687311595628638823537875937519577818577805
32171226806613001927876611195909216420198
The time needed to compute many digits of the various constants differs greatly. On the author’s computer, the golden ratio can be evaluated to 100,000 digits in one second, while obtaining just 500 digits of Khinchin’s constant takes several seconds. However, mpmath caches computed values of constants so that subsequent calculations at the same or lower precision are instant. This important for general efficiency reasons since for example the value of pi is required internally by mpmath to compute trigonometric functions.
Mpmath implements the standard functions available in Python’s math and cmath modules, for both real and complex numbers and with arbitrary precision. Many other named elementary functions are also available, including the reciprocal trigonometric and hyperbolic functions sec, csc, cot, sech, etc. Here are some simple evaluations using elementary functions:
>>> mp.dps = 25
>>> print cosh('1.234')
1.863033801698422589073644
>>> print asin(1)
1.570796326794896619231322
>>> print log(1+2j)
(0.8047189562170501873003797 + 1.107148717794090503017065j)
>>> print exp(2+3j)
(-7.315110094901102517486536 + 1.042743656235904414101504j)
Trigonometric functions can be evaluated at degree angles using the degree constant:
>>> mp.dps = 15
>>> print cos(45*degree)
0.707106781186548
The functions cospi and sinpi allow accurate evaluation at multiples of pi:
>>> print sin(3*pi)
3.67394039744206e-16
>>> print sinpi(3)
0.0
An example function evaluation with very high precision:
>>> mp.dps = 500
>>> print tanh(pi) # doctest:+SKIP
0.996272076220749944264690580012536711896899190804587614362612415978541298982821
71768528355120691683454615502559430681804165984973439112096298906391022044704505
48122748988701224416244977555668537821981195835666118518137460301886473587466464
21806234986385738703533092864596585365203757599676939623117433432914420827886553
76127620464757884951425425306744563300307924663159394900310518383521418687640452
46705057548209662896867568997817626620037806847842642167105610253517933039999796
6863567179835779017033
The list of available elementary functions (and related simple functions) is given in the following table.
Function | Description |
---|---|
sqrt(x) | Square root, x^(1/2) |
cbrt(x) | Cube root, x^(1/3) |
nthroot(x) | Nth root, x^(1/n) |
hypot(x,y) | Euclidean norm, sqrt(x^2+y^2) |
exp(x) | Exponential function |
log(x,b) | Natural logarithm (optionally base-b logarithm) |
ln(x) | Natural logarithm |
log10(x) | Base-10 logarithm |
power(x,y) | Power, x**y |
modf(x,y) | Modulo, x mod y |
ldexp(x,n) | Fast calculation of x * 2**n |
frexp(x) | Scales x to the range [0.5, 1), returning value and exponent |
floor(x) | Floor function (round to integer in the direction of -inf) |
ceil(x) | Ceiling function (round to integer in the direction of +inf) |
arg(z) | Complex argument of z |
degrees(x) | Convert radians to degrees |
radians(x) | Convert degrees to radians |
cos(x) | Cosine |
sin(x) | Sine |
tan(x) | Tangent |
sec(x) | Secant |
csc(x) | Cosecant |
cot(x) | Cotangent |
cosh(x) | Hyperbolic cosine |
sinh(x) | Hyperbolic sine |
tanh(x) | Hyperbolic tangent |
sech(x) | Hyperbolic secant |
csch(x) | Hyperbolic cosecant |
coth(x) | Hyperbolic cotangent |
acos(x) | Inverse cosine |
asin(x) | Inverse sine |
atan(x) | Inverse tangent |
atan2(y,x) | Inverse tangent atan(y/x) with attention to signs of both x and y |
asec(x) | Inverse secant |
acsc(x) | Inverse cosecant |
acot(x) | Inverse cotangent |
acosh(x) | Inverse hyperbolic cosine |
asinh(x) | Inverse hyperbolic sine |
atanh(x) | Inverse hyperbolic tangent |
asech(x) | Inverse hyperbolic secant |
acsch(x) | Inverse hyperbolic cosecant |
acoth(x) | Inverse hyperbolic cotangent |
Various useful nonelementary functions that do not exist in the standard Python math library are available, such as factorials (with support for noninteger arguments):
>>> mp.dps = 20
>>> print factorial(10)
3628800.0
>>> print factorial(0.25)
0.90640247705547707798
>>> print factorial(2+3j)
(-0.44011340763700171113 - 0.06363724312631702183j)
The nonelementary functions can all be evaluated to high precision:
>>> mp.dps = 200
>>> print j0(ellipe(1/pi)) + erf(1) # doctest:+SKIP
1.389312963715905348881873975556967070995660558533781734716705228600460043954957
61618134488095800072249215436991243117271403513444828929928306553194119494900584
6834728102001723155155631593301614083021
The Lambert W function is the inverse function of w*exp(w):
>>> mp.dps = 35
>>> w = lambertw(1)
>>> print w
0.56714329040978387299996866221035555
>>> print w*exp(w)
1.0
The Lambert W function has infinitely many complex branches. By default, the principal branch (k = 0) is evaluated. You can choose a custom branch:
>>> w = lambertw(1, k=3)
>>> print w
(-2.8535817554090378072068187234910812 + 17.113535539412145912607826671159289j)
>>> print w*exp(w)
(1.0 + 3.5075477124212226194278700785075126e-36j)
The function hyper sums a general hypergeometric series (this works only for argument where the series converges):
>>> mp.dps = 50
>>> print hyper([1,2,3],[4,5,6], 1/pi)
1.0162105742486967023020961610585170205537831496322
Function | Description |
---|---|
lambertw(x,k) | Lambert W function |
factorial(x) | Factorial |
gamma(x) | Gamma function |
binomial(x,y) | Binomial coefficient |
gammaprod(x,y) | Gamma function product G(x[0])*G(x[1])... / (G(y[0])*G(y[1])...) |
rf(x,n) | Rising factorial (Pochhammer symbol), x^(n) |
ff(x,n) | Falling factorial, x_(n) |
lower_gamma(a,x) | Lower incomplete gamma function |
upper_gamma(a,x) | Upper incomplete gamma function |
psi(n,x) | Polygamma function of order n |
stieltjes(n) | Nth Stieltjes constant |
erf(x) | Error function |
erfc(x) | Complementary error function |
erfi(x) | Imaginary error function |
npdf(x, mu, sigma) | Normal probability density function |
ncdf(x, mu, sigma) | Normal cumulative distribution function |
zeta(x) | Riemann zeta function |
bernoulli(n) | Bernoulli number, B_n |
bernfrac(n) | Bernoulli number as an exact fraction |
ei(x) | Exponential integral, Ei(x) |
li(x) | Logarithmic integral, li(x) |
ci(x) | Cosine integral, Ci(x) |
si(x) | Sine integral, Si(x) |
chi(x) | Hyperbolic cosine integral, Ci(x) |
shi(x) | Hyperbolic sine integral, Si(x) |
j0(x) | Bessel function J_0(x) |
j1(x) | Bessel function J_1(x) |
jn(n,x) | Bessel function J_n(x) |
ellipe(x) | Complete elliptic integral E |
ellipk(x) | Complete elliptic integral K |
airyai(x) | Airy function Ai(x) |
airybi(x) | Airy function Bi(x) |
fresnels(x) | Fresnel integral S(x) |
fresnelc(x) | Fresnel integral C(x) |
agm(x,y) | Arithmetic geometric mean |
hyper(as,bs,z) | Hypergeometric function PFQ |
hyp0f1(a,z) | Hypergeometric function 0F1 |
hyp1f1(a,b,z) | Hypergeometric function 1F1 |
hyp2f1(a,b,c,z) | Hypergeometric function 2F1 |
jacobi(n,a,b,x) | Jacobi polynomial / Jacobi function |
legendre(n,x) | Legendre polynomial / Legendre function |
chebyt(n,x) | Chebyshev polynomial of the first kind, T_n(x) |
chebyu(n,x) | Chebyshev polynomial of the second kind, U_n(x) |
jacobi_elliptic_xx(u,m) | Jacobi elliptic functions, xx = cn, sn, dn |
jacobi_theta_x(u,m) | Jacobi theta functions, x = 1, 2, 3, 4 |
arange, with the same argument order as Python’s built-in range function, generates a range of floating-point numbers, not including the endpoint:
>>> nprint(arange(0, 5, 0.5))
[0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5]
>>> nprint(arange(pi))
[0.0, 1.0, 2.0, 3.0]
>>> nprint(arange(1.5, pi))
[1.5, 2.5]
>>> nprint(arange(1.5, pi, 0.4))
[1.5, 1.9, 2.3, 2.7, 3.1]
The special number eps is defined as the difference between 1 and the smallest floating-point number after 1 that can be represented with the current working precision:
>>> mp.dps = 15
>>> eps
<epsilon of working precision: 2.22045e-16~>
>>> 1 + eps
mpf('1.0000000000000002')
>>> 1 + eps/2 # Too small to make a difference
mpf('1.0')
>>>
>>> mp.dps = 100
>>> eps
<epsilon of working precision: 1.42873e-101~>
An useful application of eps is to perform approximate comparisons that work at any precision level, for example to check for convergence of iterative algorithms:
>>> def a_series():
... s = 0
... n = 1
... while 1:
... term = mpf(5) ** (-n)
... s += term
... if term < eps:
... print "added", n, "terms"
... return s
... n += 1
...
>>> mp.dps = 15
>>> a_series()
added 23 terms
mpf('0.25000000000000011')
>>>
>>> mp.dps = 40
>>> a_series()
added 59 terms
mpf('0.2500000000000000000000000000000000000000057')
rand() returns a random number chosen from the uniform distribution on the interval [0, 1). All bits up to the working precision are random (this is not true if generating a float directly with Python’s random module and then converting it into an mpf).
The output might look like this:
>>> mp.dps = 30
>>> for i in range(10): print rand() # doctest:+SKIP
...
0.391978954448554897401435871089
0.0807882783964168543644048960474
0.703930940186206635051721420436
0.626158757167772268326252881759
0.662757493867073525281224309434
0.756414034497179767347019515565
0.927214900855958669681240990671
0.244655678804068448731650701685
0.34560622640274419514821879678
0.353366977589549657670073479468
Polynomial functions can be evaluated using polyval, which takes as input a list of coefficients and the desired evaluation point. The following example evaluates x^3 + 5x + 2 at x = 3.5:
>>> mp.dps = 20
>>> polyval([1, 0, 5, 2], mpf('3.5'))
mpf('62.375')
With derivative=True, both the polynomial and its derivative are evaluated at the same point:
>>> polyval([1, 0, 5, 2], mpf('3.5'), derivative=True)
(mpf('62.375'), mpf('41.75'))
The point and coefficients may be complex numbers.