【Python】Chapter 4: Sympy - "Scientific Computing with Python"

Title: "【Python】Chapter 4: Sympy - 'Scientific Computing with Python'"
Date: 2017-06-08
Permalink: /posts/2017/06/【Python】Chapter 4: Sympy - 'Scientific Computing with Python'/
Tags:

Python

SymPy is a symbolic mathematics library for Python that allows for symbolic manipulation of mathematical formulas. For convenience, all the example programs below assume that all the contents of the sympy library have been imported beforehand:

>>> from sympy import *

The asterisk (*) represents "all".

1. Classic Formula#

This is called Euler's identity, where e is the base of the natural logarithm, i is the imaginary unit, and pi is the mathematical constant pi. This formula is considered to be the most remarkable formula in mathematics, as it connects five fundamental mathematical constants through addition, multiplication, and exponentiation. Let's use SymPy to verify this formula.

In the loaded symbols, E represents the base of the natural logarithm, I represents the imaginary unit, and pi represents the mathematical constant pi. Therefore, the above formula can be calculated as follows:

>>> E**(I*pi)+1
0

Euler's identity can also be calculated using the following formula:

To prove the above formula using SymPy, we need to introduce the variable x. In SymPy, mathematical symbols are objects of the Symbol class, so they must be created before they can be used:

>>> x = Symbol('x')

The expand function can expand formulas. Let's try expanding E**(I*x):

>>> expand( E**(I*x) )
exp(I*x)

It didn't work, it just changed the notation. Here, exp is not math.exp or numpy.exp, but sympy.exp. It is a class used to represent the exponential function.

math.exp(x): Returns the exponential of x, where x must be a number.
numpy.exp(x): Returns the exponential of x, where x can be an array or a number.

>>> import math
>>> import numpy as np

>>> x = 2
>>> y = np.array([1,2,3,4])
>>> z = np.array([2])

>>> math.exp(x)
7.38905609893065
>>> math.exp(z)
7.38905609893065
>>> math.exp(y)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: only length-1 arrays can be converted to Python scalars

>>> np.exp(x)
7.3890560989306504
>>> np.exp(z)
7.3890560989306504
>>> np.exp(y)
array([  2.71828183,   7.3890561 ,  20.08553692,  54.59815003])

numpy is a package for array operations in Python!

The expand function has a keyword argument complex (a parameter type with an imaginary part, essentially a complex number, as discussed in NumPy). When it is set to True, expand will split the formula into real and imaginary parts:

>>> expand(exp(I*x), complex=True)
I*exp(-im(x))*sin(re(x)) + cos(re(x))*exp(-im(x))

The result obtained this time is quite complex. Here, sin, cos, re, and im are classes defined in sympy. re represents the real part and im represents the imaginary part. Clearly, the calculation here treats the symbol x as a complex number. In order to specify that the symbol x must be a real number, we need to redefine the symbol x as follows:

>>> x = Symbol("x", real=True)
>>> expand(exp(I*x), complex=True)
I*sin(x) + cos(x)

Usage of single quotes '' and double quotes "" in Python:
In most cases, there is no clear distinction or special usage between them, and they can be used interchangeably. However, there is one case to note:

>>> str = "Hello,world!"
>>> str1 = 'Hello,world!'
>>> str2 = "Let's go!"
>>> str3 = 'Let\'s go!' # \ is an escape character, telling the system that the following ' is not a quotation mark.
>>> str4 = "Say:\"Hello!\""
>>> str5 = 'Say:"Hello"'
>>> print str
Hello,world!
>>> print str1
Hello,world!
>>> print str2
Let's go!
>>> print str3
Let's go!
>>> print str4
Say:"Hello!"
>>> print str5
Say:"Hello"

Without the escape character, an error message will be displayed:

>>> str6 ='Let's go!'
  File "<stdin>", line 1
    str6 ='Let's go!'
               ^
SyntaxError: invalid syntax

2. Integration#

SymPy's symbolic integration function integrate can help us perform symbolic integration. integrate can perform indefinite integration:

>>> integrate(x*sin(x), x)
-x*cos(x) + sin(x)

If the range of x is specified, integrate performs definite integration:

>>> integrate(x*sin(x), (x, 0, 2*pi))
-2*pi

To calculate the volume of a sphere, let's first see how to calculate the area of a circle. Assuming the radius of the circle is r, the Y-coordinate function of any point on the circle is:

>>> x, y, r = symbols('x,y,r')
>>> 2 * integrate(sqrt(r*r-x**2), (x, -r, r))
2*Integral((r**2 - x**2)**(1/2), (x, -r, r))

Unfortunately, the integrate function did not compute the result, but instead returned the input expression. This is because SymPy does not know that r is greater than 0. By redefining r as follows, we can get the correct answer:

>>> r = symbols('r', positive=True)
>>> circle_area = 2 * integrate(sqrt(r**2-x**2), (x, -r, r))
>>> circle_area
pi*r**2

Next, we need to substitute the variable r in circle_area. The radius of the cross-section changes as the X-coordinate changes. Now, assuming the X-coordinate is x and the radius of the sphere is r, the radius of the cross-section at x can be calculated using the previous formula y(x). Therefore, we need to substitute the variable r in circle_area as follows:

>>> circle_area = circle_area.subs(r, sqrt(r**2-x**2))
>>> circle_area
pi*(r**2 - x**2)

Then, integrate the variable x in circle_area over the interval -r to r to obtain the volume of the sphere:

>>> integrate(circle_area, (x, -r, r))
4*pi*r**3/3

Using subs for expression substitution
The subs function can be used to substitute symbols in an expression. It has three calling methods:
expression.subs(x, y): Replaces x with y in the expression.
expression.subs({x,u}): Uses a dictionary for multiple substitutions.
expression.subs([(x,y),(u,v)]): Uses a list for multiple substitutions.
Please note that multiple substitutions are executed in order, so:
expression.sub([(x,y),(y,x)])
will not swap the two symbols x and y.

Here are some examples to make it clear:

>>> from sympy import *
>>> a,b,c = symbols('a,b,c') # lowercase letters
>>> y = a + 10*b + 100*c
>>> y = y.subs(a,b)
>>> y
11*b + 100*c
>>> y = a + 10*b + 100*c
>>> y = y.subs([(a,b),(b,c)])
>>> y
111*c
>>> y = a + 10*b + 100*c
>>> y = y.subs([(a,b),(b,a)])
>>> y
11*a + 100*c
>>> x=Symbol('x') # note the uppercase letter
>>> y = a + 10*b + 100*c
>>> y = y.subs([(a,x),(b,a),(x,b)])
>>> y
10*a + b + 100*c

Finally, write it down with a pen in order.

To read the complete content, please visit the original article: http://old.sebug.net/paper/books/scipydoc/sympy_intro.html