Linear Algebra and Programming Skills

Dr Jon Shiach, Department of Computing and Mathematics, Manchester Metropolitan University

---

Symbolic Mathematics

Learning Outcomes

On successful completion of this page readers will be able to:

• use commands from the Python library SymPy to perform algebraic manipulation of equations;
• perform symbolic differentiation, integration and definite integration using Python;
• solve algebraic and simple differential equations symbolically using Python.

---

SymPy

The Python library SymPy allows us to use a suite of functions for symbolically evaluating mathematical expressions. Before we can use SymPy functions we need to import the libary using the following

from sympy import *

This command will import all SymPy functions so that we can use them without the sympy. prefix. Execute the cell cell above to import all functions from the SymPy library.

---

Algebra

Symbolic variables

The first thing we need to do when working with symbolic mathematics is to define any variable and function that we want Python to treat symbolically. This is done using the symbols and function commands

variable1, variable2 = symbols('variable1, variable2')
function = Function('function')(variable1, variable2)

This defines the variables variable1 and variable2 as symbolic variables and function as a symbolic function expressed in terms of variable1 and variable2. Once symbolic variables and functions are defined they can be used to create symbolic expressions.

Example 1

The commands below declare symbolic variable x and the function $f(x) = e^x + \cos(x)$

x = symbols('x')
f = Function('f')(x)

f = exp(x) + cos(x)
f

Enter them into the code cell below and execute it.

x = symbols('x')
f = Function('f')(x)

f = exp(x) + cos(x)
f

$\displaystyle e^{x} + \cos{\left(x \right)}$

We can also define the function $g(x) = x^2 + 1$ and use the subs command to substitute functions into one another to form composite functions.

g = Function('g')(x)
g = x ** 2 + 1
g

g = Function('g')(x)
g = x ** 2 + 1
g

$\displaystyle x^{2} + 1$

To evaluate the composite function $f\circ g$ we substitute g in place of the x variable in the function f

f.subs(x, g)

f.subs(x, g)

$\displaystyle e^{x^{2} + 1} + \cos{\left(x^{2} + 1 \right)}$

Or to evaluate $g\circ f$ we substitute f in place of x in g

g.subs(x, f)

g.subs(x, f)

$\displaystyle \left(e^{x} + \cos{\left(x \right)}\right)^{2} + 1$

---

Expanding expressions

To expand an algebraic expression we can use the expand command.

expand(expression)

Example 2

The commands below simplify the expression $(x + y)^4$

y = symbols('y')
expand((x + y) ** 4)

Enter them into the code cell below and run the live script.

y = symbols('y')
expand((x + y) ** 4)

$\displaystyle x^{4} + 4 x^{3} y + 6 x^{2} y^{2} + 4 x y^{3} + y^{4}$

---

Factorising expressions

The factorise an expression we can use the factor command

factor(expression)

Example 3

The command below factorises the expression $x^3 - 3x^2 + 3x - 1$

factor(x ** 3 - 3 * x ** 2 + 3 * x - 1)

Enter it into the code cell below and run the live script.

factor(x ** 3 - 3 * x ** 2 + 3 * x - 1)

$\displaystyle \left(x - 1\right)^{3}$

---

Partial fractions

To split up a fraction of the form $\dfrac{f(x)}{g(x)}$ where $f$ and $g$ are polynomial functions we can use the apart command

apart(expression)

Example 4

The commands below determine the partial fraction decomposition of the following expression

\begin{align*} \frac{3x+1}{x^3 - 3x+ 2}. \end{align*}

apart((3 * x + 1) / (x ** 3 - 3 * x + 2))

apart((3 * x + 1) / (x ** 3 - 3 * x + 2))

$\displaystyle - \frac{5}{9 \left(x + 2\right)} + \frac{5}{9 \left(x - 1\right)} + \frac{4}{3 \left(x - 1\right)^{2}}$

---

Exercise 1 - Algebraic expressions

1. Use the expand command to expand the following expressions:

    (a)   $(x + 2)^2$;

expand((x + 2) ** 2)

$\displaystyle x^{2} + 4 x + 4$

    (b)   $(x + 1)(ax^2 + aby + cxy)$

a, b, c = symbols("a, b, c")
expand((x + 1) * (a * x ** 2 + a * b * y + c * x * y))

$\displaystyle a b x y + a b y + a x^{3} + a x^{2} + c x^{2} y + c x y$

    (c)   $(x + y + z)^3$.

z = symbols("z")
expand((x + y + z) ** 3)

$\displaystyle x^{3} + 3 x^{2} y + 3 x^{2} z + 3 x y^{2} + 6 x y z + 3 x z^{2} + y^{3} + 3 y^{2} z + 3 y z^{2} + z^{3}$

1. Use the factor command to simplify the following expressions:

    (a)   $x^3 + 4x^2 + 4x$;

factor(x ** 3 + 4 * x ** 2 + 4 * x)

$\displaystyle x \left(x + 2\right)^{2}$

    (b)  $x^4 - 12x^3 + 54x^2 - 108x + 81$;

factor(x ** 4 - 12 * x ** 3 + 54 * x ** 2 - 108 * x + 81)

$\displaystyle \left(x - 3\right)^{4}$

    (c)   $\dfrac{2x}{x^2-4x+4} + \dfrac{1}{x^2 - 4x + 4} + \dfrac{x^2}{x^2 - 4x + 4}$

factor(2 * x / (x ** 2 - 4 * x + 4) + 1 / (x ** 2 - 4 * x + 4) + x ** 2 / (x ** 2 - 4 * x + 4))

$\displaystyle \frac{\left(x + 1\right)^{2}}{\left(x - 2\right)^{2}}$

1. Consider the functions $f$, $g$, $h$ all $\mathbb{R} \to \mathbb{R}$ defined by
   
   \begin{align*} f(x) &= - 2x - 3, \\ g(x) &= x^2 + 3x + 1, \\ h(x) &= \frac{x}{2}. \end{align*}
   Use Python commands to find the formulas, simplified as far as possible using the simplify command, that represent the functions:
   
   (a) $f\circ g \circ h$;

f = Function("f")(x)
g = Function("g")(x)
h = Function("h")(x)

f = -2 * x - 3
g = x ** 2 + 3 * x + 1
h = x / 2

f.subs(x, g.subs(x, h))

$\displaystyle - \frac{x^{2}}{2} - 3 x - 5$

    (b)   $g \circ f \circ h$;

g.subs(x, f.subs(x, h))

$\displaystyle - 3 x + \left(- x - 3\right)^{2} - 8$

    (c)   $h \circ g \circ f$;

h.subs(x, g.subs(x, f))

$\displaystyle - 3 x + \frac{\left(- 2 x - 3\right)^{2}}{2} - 4$

    (d)   $f \circ f$;

f.subs(x, f)

$\displaystyle 4 x + 3$

    (e)  $g \circ g$;

g.subs(x, g)

$\displaystyle 3 x^{2} + 9 x + \left(x^{2} + 3 x + 1\right)^{2} + 4$

    (f)   $h \circ h$.

h.subs(x, h)

$\displaystyle \frac{x}{4}$

1. Obtain the partial fraction expansions of the following rational function expressions.

    (a)   $\dfrac{13t - 23}{2(3t - 1)(t - 3)}$;

t = symbols("t")

apart((13 * t - 23) / (2 * (3 * t - 1) * (t - 3)))

$\displaystyle \frac{7}{2 \left(3 t - 1\right)} + \frac{1}{t - 3}$

    (b)   $-\dfrac{9t^2 + 10t + 15}{3(t^2 + 1)(3t - 5)}$;

apart(-(9 * t ** 2 + 10 * t + 15) / (3 * (t ** 2 + 1) * (3 * t - 5)))

$\displaystyle \frac{2 t}{3 \left(t^{2} + 1\right)} - \frac{5}{3 t - 5}$

    (c)   $\dfrac{2t^2 - 11 t + 14}{(5t - 4) (t + 1)^2}$;

apart((2 * t ** 2 - 11 * t + 14) / ((5 * t - 4) * (t + 1) ** 2))

$\displaystyle \frac{2}{5 t - 4} - \frac{3}{\left(t + 1\right)^{2}}$

    (d)   $-\dfrac{25 t + 14}{ (2t^2 + 5t + 5)(t-4)}$.

apart(-(25 * t + 14) / ((2 * t ** 2 + 5 * t + 5) * (t - 4)))

$\displaystyle \frac{4 t + 1}{2 t^{2} + 5 t + 5} - \frac{2}{t - 4}$

---

Calculus

Differentiation

To evaluate the derivative of the function $f(x)$ with respect to $x$

\begin{align*} \dfrac{d}{dx}f(x) = \lim_{\Delta x \to 0} \frac{f(x + \Delta x) - f(x)}{\Delta x}, \end{align*}

we can use the diff command.

diff(function, variable)

Where function is a symbolic function which is differentiated with respect to the variable and n is an optional argument which is the order of the derivative.

Example 5

The commands below evaluate the first, second and third order derivatives of the function $f(x) = 2x^3$. Enter them into the code cells below and execute each one.

f = Function('f')(x)
f = 2 * x ** 3
diff(f, x)

f = Function('f')(x)
f = 2 * x ** 3
diff(f, x)

$\displaystyle 6 x^{2}$

diff(f, x, 2)

diff(f, x, 2)

$\displaystyle 12 x$

diff(f, x, 3)

diff(f, x, 3)

$\displaystyle 12$

Example 6

The commands below produce the product and quotient rules for differentiation of single variable functions

\begin{align*} \frac{d}{dx}f(x)g(x) &= f(x)\frac{d}{dx}g(x) + g(x) \frac{d}{dx}f(x), \\ \frac{d}{dx}\left(\frac{f(x)}{g(x)}\right) &= \frac{g(x) \dfrac{d}{dx}f(x) - f(x) \dfrac{d}{dx}g(x)}{g(x)^2}, \\ \frac{d}{dx}f(g(x)) &= \frac{d}{dg(x)}f(x) \frac{d}{dx}g(x). \end{align*}

Enter each one into the code cell below and execute each one.

f = Function('f')(x)
g = Function('g')(x)

diff(f * g, x)

f = Function('f')(x)
g = Function('g')(x)

diff(f * g, x)

$\displaystyle f{\left(x \right)} \frac{d}{d x} g{\left(x \right)} + g{\left(x \right)} \frac{d}{d x} f{\left(x \right)}$

simplify(diff(f / g, x))

simplify(diff(f / g, x))

$\displaystyle \frac{- f{\left(x \right)} \frac{d}{d x} g{\left(x \right)} + g{\left(x \right)} \frac{d}{d x} f{\left(x \right)}}{g^{2}{\left(x \right)}}$

diff(f.subs(x, g))

diff(f.subs(x, g))

$\displaystyle \frac{d}{d g{\left(x \right)}} f{\left(g{\left(x \right)} \right)} \frac{d}{d x} g{\left(x \right)}$

---

Integration

To evaluate the integral of the function $f(x)$ with respect to $x$,

\begin{align*} F(x) = \int f(x) dx, \end{align*}

we can use the int command

integral(function, variable)

Example 7

The command below evaluates the integral

\begin{align*} \int 2x^3 dx = \frac{x^4}{2} + c. \end{align*}

integrate(2 * x ** 3, x)

integrate(2 * x ** 3, x)

$\displaystyle \frac{x^{4}}{2}$

Note that Python does not include the constant of integration in the output.

---

Definite integration

To calculate the value of the definite integral

\begin{align*} \int_a^b f(x) dx = [F(x)]_a^b = F(b) - F(a), \end{align*}

we specify the limits of integration $a$ and $b$ in the integrate command

integrate(function, (variable, a, b) )

Example 8

The commands below evaluates the definite integral

\begin{align*} \int_{\pi/4}^{\pi/2} \cos(x) dx = \sin\left(\frac{\pi}{2}\right) - \sin\left(\frac{\pi}{4}\right) = 1 - \frac{\sqrt{2}}{2}. \end{align*}

integrate(cos(x), (x, pi/4, pi/2) )

integrate(cos(x), (x, pi/4, pi/2) )

$\displaystyle 1 - \frac{\sqrt{2}}{2}$

---

Exercise 2 - Differentiation and integration

1. Use the diff command to obtain the derivatives of the following functions

    (a)   $f(x) = x^5 + 3x^4 - 5x^2 + 10$;

diff(x ** 5 + 3 * x ** 4 - 5 * x ** 2 + 10, x)

$\displaystyle 5 x^{4} + 12 x^{3} - 10 x$

    (b)   $g(x) = 3x^{1/2} - x^{3/2} + 2x^{-1/2}$; (hint: use Rational(a, b) to ensure fractional values in the output)

diff(3 * x ** Rational(1, 2) - x ** Rational(3, 2) + 2 * x ** -Rational(1, 2), x)

$\displaystyle - \frac{3 \sqrt{x}}{2} + \frac{3}{2 \sqrt{x}} - \frac{1}{x^{\frac{3}{2}}}$

    (c)   $\lambda(t) = \sqrt{3t} + 3\sqrt{t}$;

diff(sqrt(3 * t) + 3 * sqrt(t), t)

$\displaystyle \frac{\sqrt{3}}{2 \sqrt{t}} + \frac{3}{2 \sqrt{t}}$

    (d)   $\mu(y) = \left(2 + 5y + \dfrac{1}{2}y^2\right)^{1/2}$;

diff(sqrt(2 + 5 * y + 1 / 2 * y ** 2), y)

$\displaystyle \frac{\sqrt{5} \left(0.1 y + \frac{1}{2}\right)}{\sqrt{0.1 y^{2} + y + \frac{2}{5}}}$

    (e)   $\phi(z) = (z + 1) \sqrt{z^4 - 2z + 2}$

diff((z + 1) * sqrt(z ** 4 - 2 * z + 2), z)

$\displaystyle \frac{\left(z + 1\right) \left(2 z^{3} - 1\right)}{\sqrt{z^{4} - 2 z + 2}} + \sqrt{z^{4} - 2 z + 2}$

    (f)   $h(x) = \dfrac{x}{\sqrt{1 - 4x^2}}$;

diff(x / sqrt(1 - 4 * x ** 2))

$\displaystyle \frac{4 x^{2}}{\left(1 - 4 x^{2}\right)^{\frac{3}{2}}} + \frac{1}{\sqrt{1 - 4 x^{2}}}$

1. Use the integrate command to determine the following integrals:

    (a)   $\int 4x^3 + 2x^2 - 3x + 4 \, dx$;

integrate(4 * x ** 3 + 2 * x ** 2 - 4 * x + 4, x)

$\displaystyle x^{4} + \frac{2 x^{3}}{3} - 2 x^{2} + 4 x$

    (b)   $\int \cos(x) \sin(x)\, dx$;

integrate(cos(x) * sin(x), x)

$\displaystyle \frac{\sin^{2}{\left(x \right)}}{2}$

    (c)   $\int x + e^{x^{1/3}} \, dx$;

integrate(x + exp(x ** Rational(1, 3)), x)

$\displaystyle 3 x^{\frac{2}{3}} e^{\sqrt[3]{x}} - 6 \sqrt[3]{x} e^{\sqrt[3]{x}} + \frac{x^{2}}{2} + 6 e^{\sqrt[3]{x}}$

    (d)   $\displaystyle \int_1^3 x^2 + \cos(x) \, dx$;

integrate(x ** 2 + cos(x), (x, 1, 3))

$\displaystyle - \sin{\left(1 \right)} + \sin{\left(3 \right)} + \frac{26}{3}$

    (e)   $\displaystyle \int_0^{\pi/2} x^3 \cos(x) \, dx$;

integrate(x ** 3 * cos(x), (x, 0, pi / 2))

$\displaystyle - 3 \pi + \frac{\pi^{3}}{8} + 6$

---

Solving equations

Solving algebraic equations

We can use Python to solve equations or systems of equations symbolically using the solve command

solve((equation1, equation2, ...))

Where equation1, equation2 etc. are equations written using symbolic variables and variable1, variable2 etc. are the variables that are being solved. The equations being solved are written such that they are equal to zero, e.g., for the equation $x + y = 1$

x + y - 1

Example 9

The commands below solve the equation $y = 2x - 5$ for $x$.

solve(y - 2 * x + 5)

solve(y - 2 * x + 5, x)

[y/2 + 5/2]

Example 10

The commands below solves the following linear system of equations

\begin{align*} x + y + z &= 0, \\ 2x + z &= -1, \\ x - 2y + 3z &= -12. \end{align*}

eqn1 = x + y + z
eqn2 = 2 * x + z + 1
eqn3 = x - 2 * y + 3 * z + 12
solve( (eqn1, eqn2, eqn3))

Enter them into the code cell below and execute it

eqn1 = x + y + z
eqn2 = 2 * x + z + 1
eqn3 = x - 2 * y + 3 * z + 12
solve( (eqn1, eqn2, eqn3), (x, y, z))

{x: 1, y: 2, z: -3}

Example 11

The solve command can also solve non-linear equations. The commands below derive the quadratic formula used to find the roots of the quadratic equation $f(x) = ax^2 + bx + c$ which is

\begin{align*} x = \frac{-b \pm \sqrt{b^2 - 4ac}}{2a}. \end{align*}

solve(a * x ** 2 + b * x + c, x)

solve(a * x ** 2 + b * x + c, x)

[(-b + sqrt(-4*a*c + b**2))/(2*a), -(b + sqrt(-4*a*c + b**2))/(2*a)]

---

Solving differential equations

To solve differential equations (where solutions exist) we can use the dsolve command.

dsolve(diffeq, function)

Where diffeq is a differential equation defined using the Eq command and function is the solution function. Eq is used specify an equivalance relationship between two objects for example

Eq(a, b)

specifies that a and b are equivalent. So to define a differential equation we can use

diffeq = Eq(lhs, rhs)

Example 12

The commands below determine the solution to the simplest differential equation

\begin{align*} \frac{d}{dt}y(t) &= ay & &\implies & y(t) &= C e^{at}. \end{align*}

y = Function('y')(t)
diffeq = Eq(diff(y, t), a * y)
dsolve(diffeq, y)

Enter them into the code cell below and execute it.

y = Function('y')(t)
diffeq = Eq(diff(y, t), a * y)
dsolve(diffeq, y)

$\displaystyle y{\left(t \right)} = C_{1} e^{a t}$

Example 13

The commands below use dsolve to solve the second-order ODE

\begin{align*} \frac{d^2}{dx^2}y(x) + \frac{d}{dx}y(x) + y(x) = 0. \end{align*}

y = Function('y')(x)
diffeq = Eq(diff(y, x, 2) + diff(y, x) + y, 0)
expand(dsolve(diffeq, y))

y = Function('y')(x)
diffeq = Eq(diff(y, x, 2) + diff(y, x) + y, 0)
expand(dsolve(diffeq, y))

$\displaystyle y{\left(x \right)} = C_{1} e^{- \frac{x}{2}} \sin{\left(\frac{\sqrt{3} x}{2} \right)} + C_{2} e^{- \frac{x}{2}} \cos{\left(\frac{\sqrt{3} x}{2} \right)}$

Example 14

The commands below use dsolve to solve the following system of ODEs

\begin{align*} \frac{d}{dt}y(t) &= z, \\ \frac{d}{dt}z(t) &= -y. \end{align*}

Enter them into the code cell below and run the live script.

y = Function('y')(t)
z = Function('z')(t)
diffeq1 = diff(y, t, 1) - z
diffeq2 = diff(z, t, 1) + y
dsolve((diffeq1, diffeq2))

y = Function('y')(t)
z = Function('z')(t)

diffeq1 = diff(y, t, 1) - z
diffeq2 = diff(z, t, 1) + y
dsolve((diffeq1, diffeq2))

[Eq(y(t), C1*sin(t) + C2*cos(t)), Eq(z(t), C1*cos(t) - C2*sin(t))]

---

Exercise 3 - Solving equations

1. Use the solve command to find the solution to the following linear system of equations

\begin{align*} x + y + z &= 0, \\ 2x + 3y + z &= 1, \\ -x -2y + 3z &= 2. \end{align*}

x, y, z = symbols('x, y, z')
eq1 = x + y + z
eq2 = 2 * x + 3 * y + z - 1
eq3 = -x - 2 * y + 3 * z - 2

solve((eq1, eq2, eq3))

{x: -3, y: 2, z: 1}

1. Use the solve command to find the solutions to the following system of nonlinear equations

\begin{align*} p^2 + q^2 &= 10,\\ 2 p + q &= 1. \end{align*}

p, q = symbols('p, q')
eq1 = p ** 2 + q ** 2 - 10
eq2 = 2 * p + q - 1

solve((eq1, eq2))

[{p: -1, q: 3}, {p: 9/5, q: -13/5}]

1. Use the dsolve command to find the solutions to the following system of differential equations

\begin{align*} \frac{dy}{dt} &= 4y + 7z, \\ \frac{dz}{dt} &= -2y - 5z. \end{align*}

y = Function('y')(t)
z = Function('z')(t)

diffeq1 = diff(y, t, 1) - 4 * y - 7 * z
diffeq2 = diff(z, t, 1) + 2 * y - 5 * z
dsolve((diffeq1, diffeq2))

[Eq(y(t), (C1/4 + sqrt(55)*C2/4)*exp(9*t/2)*cos(sqrt(55)*t/2) + (sqrt(55)*C1/4 - C2/4)*exp(9*t/2)*sin(sqrt(55)*t/2)),
 Eq(z(t), C1*exp(9*t/2)*cos(sqrt(55)*t/2) - C2*exp(9*t/2)*sin(sqrt(55)*t/2))]