Python language has both good support of numerical calculation as well as symbolic calculation. SymPy module actually is the one for symbolic calculation in python [1][2]. Here we give some brief introduction to its usage in symbolic calculation.

Basic symbolic operations

Define symbolic expressions

1
2
import sympy
x, y, c1, c2 = sympy.symbols('x, y, c1, c2', real=True)

Substitute value to symbolic expressions

Substitute number

1
2
3
x_new = x.subs(3)
z = c1*x + c2*y
z_new = z.subs({c1:1, c2:3})

Notice the difference between the value substitution for expressions with only one variable and that with more than one variables.

Substitute other symbolic expressions

1
y_new = sympy.Subs(x**2, y-1)

What you will get is: (y1)2(y-1)**2.

Simplify the displayed results

Sometimes, after symbolic calculation, the results will be a bunch of mass. You could make results displayed more elegantly with simplify commands [3].

1
s_simple = sympy.simplify(s)

Complex numbers variables

Define complex symbolic variables

The complex symbolic variables are a combination of real and imaginary part. The imaginary unit: ii can be defined as:

1
2
import sympy
sympy.I

Calculate absolute value of a complex variable

1
2
3
4
5
# test complex symbolic calculation
x, y = sympy.symbols('x, y', real=True)
z = x + sympy.I*y
abs_z = sympy.Abs(z)
print('|z| = ', abs_z)

Plot with sympy

Basic plot

1
2
3
4
5
6
7
# test equation solve with multi-variables
c1, c2, c3 = sympy.symbols('c1, c2, c3', real=True)
# symbolic plot figures
f0 = sympy.cos(x)
f1 = sympy.cos(x) + c1
f2 = f1.subs({c1:3})
sympy.plot(f0)

Plot more curves in one figure

1
sympy.plot(f0, f2)

Notice that, sympy and pylab are separated modules, thus we do not need to create a new with pylab.figure().

Solve equations

Define an equation

1
2
3
# test equation solve with multi-variables
c1, c2, c3 = sympy.symbols('c1, c2, c3', real=True)
eq1 = sympy.Eq(c1*x + c2*y + c3, 0)

This defines the equation: c1x+c2y+c3=0c1 x + c2 y + c3 = 0.

Solve equation for a particular variable

1
2
sol1 = sympy.solve(eq1, x)
sol2 = sympy.solve(e11, y)

Notice that for an equation with more than one variable, eg: f(x,y)=0f(x, y) = 0, you should specify which variable to solve. Otherwise the solver will solve the first variable by default.

Symbolic calculus

Limitation

With sympy, you can even calculate limitations in calculus.

1
2
3
4
oo = sympy.oo # define infinity
fx = 1/(x**2+1)
Lim1 = sympy.limit(fx, x, 0)
Lim2 = sympy.limit(fx, x, oo)

Differential

Notice that for functions with more than one variable, you shall specify which variable to conduct differential.

1
2
3
f = c1*x + c2*y**2
prime_fx = sympy.diff(f, x)
prime_fy = sympy.diff(f, y)

Integration

You can calculate indefinite integration like this [3]:

1
Ii = sympy.integrate(f, x)

You can also calculate definite integration by setting upper and lower limits of integral range.

1
2
3
# declare the infinity large
oo = sympy.oo
Id = sympy.integrate(f, (x, -oo, 5))

Taylor’s expansion

You can even apply Taylor’s expansion on any given functions [3].

1
sympy.series(sympy.cos(2*x), x)

Practice example

Eg: Solve a complex circuit

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
# test the solution to a simple circuit
C, L, R,omega = sympy.symbols('C, L, R, omega', real=True)
Z1 = R
Z2 = 1/(sympy.I*omega*C)
Z3 = sympy.I*omega*L

Z = Z1 + Z2 + Z3
Z_simple = sympy.simplify(Z.expand(complex=True))
abs_Z = sympy.Abs(Z)

# calculate extrem values
prime_abs_Z = sympy.diff(abs_Z, omega)
# define the equation to be solved
eq2 = sympy.Eq(prime_abs_Z, 0)
# solve the equation about omega
sol2 = sympy.solve(eq2, omega)
# output of sol2
# [-1/(sqrt(Abs(C))*sqrt(Abs(L))),
# 1/(sqrt(Abs(C))*sqrt(Abs(L))),
# -I/(sqrt(Abs(C))*sqrt(Abs(L))),
# I/(sqrt(Abs(C))*sqrt(Abs(L)))]

root2 = 1/(sympy.sqrt(sympy.Abs(C))*sympy.sqrt(sympy.Abs(L)))
root2_simple = sympy.simplify(root2)

abs_Z1 = abs_Z.subs({C:1, L:1, R:1})
sympy.plot(abs_Z1, (omega, -15, 15))

This is just the analytic solution of alternating circuit with resonance frequency at: ω=1LC\omega = \frac{1}{\sqrt{LC}}.

Reference:

[1] https://zhuanlan.zhihu.com/p/96738286

[2] https://docs.sympy.org/latest/index.html

[3] https://scipy-lectures.org/packages/sympy.html#integration