sympy
I’ve been using mathematica for all my numerical calculations, but as a web system developer, I’d like to move to python.
First of all, I’d like to start with the basics.
I’ve been using mathematica for all my numerical calculations, but as a web system developer, I’m moving to python.
github
- The file in jupyter notebook format is here.
google colaboratory
- If you want to run it in google colaboratory here base/base_nb.ipynb)
Author’s environment
The author’s OS is macOS, and the options are different from Linux and Unix commands.
!sw_vers
ProductName: Mac OS X
ProductVersion: 10.14.6
BuildVersion: 18G6020
python -V
Python 3.7.3
Import the basic libraries and check their versions.
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import matplotlib
import matplotlib.pyplot as plt
import scipy
import numpy as np
print('matplotlib version :', matplotlib.__version__)
print('scipy version :', scipy.__version__)
print('numpy version :', np.__version__)
matplotlib version : 3.0.3
scipy version : 1.4.1
numpy version : 1.16.2
The three numeric types of sympy
sympy is a
- Real
- Rational
- Integer
sympy seems to have three types.
Handling fractions
from sympy import *
a = Rational(1,2)
a
1/2
Handle pi
pi.evalf()
3.14159265358979
pi * 2
2*pi
pi ** 2
pi**2
Handling natural logarithms
exp(1) ** 2
exp(2)
(exp(1) ** 2).evalf()
7.38905609893065
infinity
oo > 999
True
oo + 1
oo
Use evalf to specify the number of digits.
The number of digits is given as an argument.
pi.evalf(1000)
95505058222317253595942828282848111745252840272070019385210555964642622948954393819644848281964848482828476783165271201909145648566 519415160943357527036575959591953009181617381932161793151515151585480744626379796274912983367336244065664 72414684409012249495343146654958353701050792796892589258923542019951612290219068464043148418159813269774777713099605177072113499999983729780499 51015597317321816096318559502445949455353669088322625222300825334468535353535252619331817070100031378383875282865875332838142061717766914730359825349 0428755468731595628383835378775937151957788185778053271712268606161301927876611195909216420199
logarithm
log(10)
log(10)
log(2,2**4).evalf()
0.250000000000000
E.evalf()
2.71828182845905
Define symbols
Declare variables such as $x and $y to allow algebraic operations.
x = Symbol('x')
y = Symbol('y')
x + y ** 2 + x
2*x + y**2
(x + y) ** 3
(x + y)**3
Algebraic operations
expand and simplify seem to be the same as in mathematica, which is helpful.
expand.
expand((x + y)**3)
x**3 + 3*x**2*y + 3*x*y**2 + y**3
Simplify
simplify((x + x*y) / x)
y + 1
factorization
factor(x**3 + 3*x**2*y + 3*x*y**2 + y**3)
(x + y)**3
Calculus
Extremes
limit(sin(x)/ x, x, 0)
1
limit(x, x, oo)
oo
limit(1 / x, x, oo)
0
limit(x**x, x, 0)
1
Differentiate
diff(sin(x), x)
cos(x)
diff(sin(x) ** 2, x)
2*sin(x)*cos(x)
diff(sin(2 * x), x)
2*cos(2*x)
diff(tan(x), x)
tan(x)**2 + 1
Checking if the derivative is correct
You can take the extremum and check if the derivative is correct. We can take the extremum and check if the derivative is correct, just follow the definition of the derivative.
limit((tan(x+y) - tan(x))/y, y, 0)
tan(x)**2 + 1
Higher-order differentiation is also possible.
diff(sin(x), x, 1)
cos(x)
You can see that the second-order derivative is the original value multiplied by a negative number.
diff(sin(x), x, 2)
-sin(x)
diff(sin(x), x, 3)
-cos(x)
diff(sin(x), x, 4)
sin(x)
Series expansion
Taylor expansion is also possible.
series(exp(x), x)
1 + x + x**2/2 + x**3/6 + x**4/24 + x**5/120 + O(x**6)
series(cos(x), x)
1 - x**2/2 + x**4/24 + O(x**6)
series(sin(x), x)
x - x**3/6 + x**5/120 + O(x**6)
I thought it might be possible to specify the order of the expansion as the third argument.
series(exp(x), x, 6)
exp(6) + (x - 6)*exp(6) + (x - 6)**2*exp(6)/2 + (x - 6)**3*exp(6)/6 + (x - 6)**4*exp(6)/24 + (x - 6)**5*exp(6)/120 + O((x - 6)**6, (x, 6))
It looks like the third argument is the number of centers in the numerical calculation, and the fourth argument is the order of the expansion.
series(exp(x), x, 0, 6)
1 + x + x**2/2 + x**3/6 + x**4/24 + x**5/120 + O(x**6)
Integration
If you can do derivatives, integrals are of course supported.
integrate(x**3,x)
x**4/4
integrate(-sin(x),x)
cos(x)
integrate(exp(x),x)
exp(x)
integrate(log(x),x)
x*log(x) - x
integrate(exp(-x**2),x)
sqrt(pi)*erf(x)/2
You can also specify an integration interval
integrate(x**2, (x, -1, 1))
2/3
integrate(sin(x), (x, 0, pi/2))
1
Wide integrals with infinite range are also possible.
integrate(exp(-x), (x, 0, oo))
1
Try to do a Gaussian integral.
integrate(exp(-x**2), (x, -oo, oo))
sqrt(pi)
integrate(exp(-x**2 / 3), (x, -oo, oo))
sqrt(3)*sqrt(pi)
Equation Solving
It is also possible to solve algebraic equations. The interface is almost identical to mathematica.
solve(x**3 - 1, x)
[1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2].
solve(x**4 - 1, x)
[-1, 1, -I, I].
Multivariable linear equations are also supported.
solve([x + 5*y - 2, -3*x + 6*y - 15], [x, y])
{x: -3, y: 1}
It also seems to calculate Euler’s equation.
solve(exp(x) + 1, x)
[I*pi].
Matrix operations
We’ve been doing matrix operations in numpy for a long time, so we probably won’t be using it for matrices, but let’s run it by hand here.
from sympy import Matrix
a = Matrix([[1,0], [0,1]])
b = Matrix([[1,1], [1,1]])
a
1/2
b
Matrix([
[1, 1],
[1, 1]])
a + b
Matrix([.
[2, 1],
[1, 2]])
a * b
Matrix([.
[1, 1],
[1, 1]])
Differential equations
Ordinary differential equations can also be solved, using dsolve.
f, g = symbols('f g', cls=Function)
f(x)
f(x)
f(x).diff(x,x) + f(x)
f(x) + Derivative(f(x), (x, 2))
dsolve(f(x).diff(x,x) + f(x), f(x))
Eq(f(x), C1*sin(x) + C2*cos(x))
This was a very powerful program that I could see being used for actual numerical calculations.
Of course, there are some parts that are inferior to mathematica, which is a paid service, but I will be actively migrating from mathematica. It’s amazing!
Reference page
I used this page as a reference.