17 SymPy

The module sympy in Python allows you to carry out symbolic manipulations. The short examples included below illustrate some of the more basic capabilities of this module (that are relevant to what we've been covering so far in CFM2103).

We are going to look at expanding brackets, differentiating and integrating functions, as well as solving linear ODEs. You can do much more -- if you are interested check out the link at the end of this document.

17.1 General stuff

In [67]:
from sympy import *
init_session()

# algebraic manipulations:
f1 = (x+1)*(x-2)*(x+6)        # this is a symbolic expression 
expand(f1)                 

IPython console for SymPy 1.0 (Python 3.5.2-64-bit) (ground types: python)

These commands were executed:
>>> from __future__ import division
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> k, m, n = symbols('k m n', integer=True)
>>> f, g, h = symbols('f g h', cls=Function)
>>> init_printing()

Documentation can be found at http://docs.sympy.org/1.0/
Out[67]:
$$x^{3} + 5 x^{2} - 8 x - 12$$
In [68]:
# another expansion:
f2 = (x+1)**4
expand(f2)
Out[68]:
$$x^{4} + 4 x^{3} + 6 x^{2} + 4 x + 1$$
In [84]:
# you can solve algebraic equations (e.g., quadratics)
eqn = x**2-5*x + 6
print('\nThe solutions are:')
solveset(eqn)                     # this solves eqn == 0

The solutions are:
Out[84]:
$$\left\{2, 3\right\}$$
In [83]:
# another quadratic:
eqn = x**2 + 3*x + 7
print('\nThe solutions are:')
solveset(eqn)                   # this solves eqn == 0

The solutions are:
Out[83]:
$$\left\{- \frac{3}{2} - \frac{\sqrt{19} i}{2}, - \frac{3}{2} + \frac{\sqrt{19} i}{2}\right\}$$
In [93]:
# and here's a cubic:
eqn = x**3 - 2*x + 1
print('\nThe cubic is:')
eqn

The cubic is:
Out[93]:
$$x^{3} - 2 x + 1$$
In [91]:
print('\nThe solutions are:')
solveset(eqn)                      # this solves eqn == 0 

The solutions are:
Out[91]:
$$\left\{1, - \frac{1}{2} + \frac{\sqrt{5}}{2}, - \frac{\sqrt{5}}{2} - \frac{1}{2}\right\}$$
In [78]:
# 1st order derivative for some function of x:
f3 = (x**2)*exp(-3*x)
print('\nThe function is:')
f3

The function is:
Out[78]:
$$x^{2} e^{- 3 x}$$
In [79]:
print('\n...and its derivative is:')
diff(f3, x)

...and its derivative is:
Out[79]:
$$- 3 x^{2} e^{- 3 x} + 2 x e^{- 3 x}$$
In [80]:
# 3rd order derivative ....
print('\n...its third-order derivative:')
diff(f3, x, 3)

...its third-order derivative:
Out[80]:
$$9 \left(- 3 x^{2} + 6 x - 2\right) e^{- 3 x}$$
In [10]:
# (indefinite) integration:
f4 = (x**3)*exp(-7*x)
integrate(f4, x)
Out[10]:
$$\frac{1}{2401} \left(- 343 x^{3} - 147 x^{2} - 42 x - 6\right) e^{- 7 x}$$
In [12]:
# another integral:
f5 = sqrt(x**2+1)
integrate(f5, x)
Out[12]:
$$\frac{x}{2} \sqrt{x^{2} + 1} + \frac{1}{2} \operatorname{asinh}{\left (x \right )}$$
In [14]:
# definite integral:
integrate(f4, (x, 1, 3))
Out[14]:
$$- \frac{10716}{2401 e^{21}} + \frac{538}{2401 e^{7}}$$
In [15]:
# Maclaurin series:
# f(x) = cosh(x)
cosh(x).series(x, 0, 8)
Out[15]:
$$1 + \frac{x^{2}}{2} + \frac{x^{4}}{24} + \frac{x^{6}}{720} + \mathcal{O}\left(x^{8}\right)$$
In [16]:
# another Maclaurin series
# inverse tan function:
atan(x).series(x, 0, 10)
Out[16]:
$$x - \frac{x^{3}}{3} + \frac{x^{5}}{5} - \frac{x^{7}}{7} + \frac{x^{9}}{9} + \mathcal{O}\left(x^{10}\right)$$

17.2 Ordinary Differential Equations (ODEs)

SymPy can provide closed-form solutions for any first or second-order ODE that admits such a solution. Let's start with some simple second-order ODEs with constant coefficients. Such equations are typically of the form:

$$ y'' + ay' + by = F(x)\,, $$

where $a$, $b\in\mathbb{R}$ are some given constants, $y=y(x)$ is the unknown function, and $f(x)$ is a given function (the so-called "forcing term"); here, the dash denotes differentiation with respect to $x$, i.e.

$$ y' = \dfrac{dy}{dx}\,,\qquad\quad y'' = \dfrac{d^2y}{dx^2}\,,\qquad\mbox{etc}. $$

In [45]:
# EX 1:
# example for solving a 2nd order ODE
f = symbols('y', cls=Function)
diffeq = Eq(y(x).diff(x,x)+2*y(x).diff(x)+y(x), sin(x))
print('\nThis is the equation that we are solving:')
diffeq

This is the equation that we are solving:
Out[45]:
$$y{\left (x \right )} + 2 \frac{d}{d x} y{\left (x \right )} + \frac{d^{2}}{d x^{2}} y{\left (x \right )} = \sin{\left (x \right )}$$
In [46]:
# solve the above ODE:
print('\nThe solution of this ODE is:')
dsolve(diffeq, y(x))

The solution of this ODE is:
Out[46]:
$$y{\left (x \right )} = \left(C_{1} + C_{2} x\right) e^{- x} - \frac{1}{2} \cos{\left (x \right )}$$
In [47]:
# EX 2:
# example for solving a 2nd order ODE
f = symbols('y', cls=Function)
diffeq = Eq(y(x).diff(x,x)+4*y(x).diff(x)+20*y(x), sin(x))
print('\nThis is the equation that we are solving:')
diffeq

This is the equation that we are solving:
Out[47]:
$$20 y{\left (x \right )} + 4 \frac{d}{d x} y{\left (x \right )} + \frac{d^{2}}{d x^{2}} y{\left (x \right )} = \sin{\left (x \right )}$$
In [53]:
# solve the above ODE:
print('\nThe solution of this ODE is:')
dsolve(diffeq, y(x))

The solution of this ODE is:
Out[53]:
$$y{\left (x \right )} = \left(C_{1} + 3 e^{x}\right) e^{x}$$

Next, we look at a couple of examples involving linear first-order ODEs. The solutions of these equations usually require the calculation of the so-called integrating factor. There is a brief handout posted on Brightspace (this topic is supposed to be covered in one of your other modules).




In [49]:
# linear ODEs (the ones that require integrating factors):

f = symbols('y', cls=Function) 
P = -1/(x+1)
Q = x/(x+1)
diffeq = Eq(y(x).diff(x) + P*y(x), Q)
print('\nThis is the equation that we are solving:')
diffeq

This is the equation that we are solving:
Out[49]:
$$\frac{d}{d x} y{\left (x \right )} - \frac{y{\left (x \right )}}{x + 1} = \frac{x}{x + 1}$$
In [50]:
print('\nThe solution of this ODE is:')
dsolve(diffeq, y(x))

The solution of this ODE is:
Out[50]:
$$y{\left (x \right )} = C_{1} x + C_{1} + x \log{\left (x + 1 \right )} + \log{\left (x + 1 \right )} + 1$$
In [51]:
# another example of linear ODE:
# 

f = symbols('y', cls=Function)
P = -1
Q = 3*exp(2*x)
diffeq = Eq(y(x).diff(x) + P*y(x), Q)
print('\nThis is the equation that we are solving:')
diffeq

This is the equation that we are solving:
Out[51]:
$$- y{\left (x \right )} + \frac{d}{d x} y{\left (x \right )} = 3 e^{2 x}$$
In [52]:
print('\nThe solution of this ODE is:')
dsolve(diffeq, y(x))

The solution of this ODE is:
Out[52]:
$$y{\left (x \right )} = \left(C_{1} + 3 e^{x}\right) e^{x}$$

You can discover more features of SymPy by visiting the official documentation page:

https://docs.sympy.org/latest/index.html