Overview of the Equation class

from sympy import *
from sympy_equation import Equation, Eqn, solve
init_printing()

In this overview we are going to explore the capabilities of the Equation class. Note that:

1. Eqn is an alias for Equation.

2. the second line of code imports solve from sympy_equation: this function is aware of objects of type Equation. It is just a wrapper to SymPy’s solve.

Defining equations

var("a:d")
# first mode: instantiate the Equation class directly
e1 = Eqn(a, b + 1)
# second mode: short syntax
e2 =@ c + 1 = d
display(e1, e2)

$\displaystyle a = b + 1$

$\displaystyle c + 1 = d$

The above cell created two equations:

• e1 was created by explicitly instanting the Equation class. While it might be longer to type, it is more readable as it clearly indicates what the code is doing.

• e2 was created with the short syntax =@. This symbol combination was chosen to avoid conflicts with reserved python symbols while minimizing impacts on syntax highlighting and autoformatting. The short syntax is only available on interactive shells (like IPython, Jupyter Notebook).

Mathematical operators

Objects of type Equation supports the following mathematical operations:

e1 + e2

$\displaystyle a + c + 1 = b + d + 1$

e1 - e2

$\displaystyle a - c - 1 = b - d + 1$

e1 * e2

$\displaystyle a \left(c + 1\right) = d \left(b + 1\right)$

e1 / e2

$\displaystyle \frac{a}{c + 1} = \frac{b + 1}{d}$

e1**2

$\displaystyle a^{2} = \left(b + 1\right)^{2}$

2**e1

$\displaystyle 2^{a} = 2^{b + 1}$

e1**e2

$\displaystyle a^{c + 1} = \left(b + 1\right)^{d}$

Applying functions to Equation

We can apply mathematical functions (or user defined functions) to an Equation (or to a particular side) with the following methods:

• eq.apply(func, *args, **kwargs): apply to both side of the equation, where *args, **kwargs are additional positional and keyword arguments to the function.

• eq.applylhs(func, *args, **kwargs): apply to the LHS of the equation.

• eq.applyrhs(func, *args, **kwargs): apply to the RHS of the equation.

e1

$\displaystyle a = b + 1$

e1.apply(sin)

$\displaystyle \sin{\left(a \right)} = \sin{\left(b + 1 \right)}$

e1.applylhs(sin)

$\displaystyle \sin{\left(a \right)} = b + 1$

e1.applyrhs(sin)

$\displaystyle a = \sin{\left(b + 1 \right)}$

eq3 = Eqn(a, b / c)
eq3

$\displaystyle a = \frac{b}{c}$

Compute the first cubic root:

eq3.apply(root, 3)

$\displaystyle \sqrt[3]{a} = \sqrt[3]{\frac{b}{c}}$

Compute the second cubic root, only on the RHS:

eq3.applyrhs(root, 3, 1)

$\displaystyle a = \left(-1\right)^{\frac{2}{3}} \sqrt[3]{\frac{b}{c}}$

Compute the third cubic root, only on the LHS:

eq3.applylhs(root, 3, 2)

$\displaystyle - \sqrt[3]{-1} \sqrt[3]{a} = \frac{b}{c}$

Apply a custom function:

def add_square(expr):
    return expr + expr**2

eq3.apply(add_square)

$\displaystyle a^{2} + a = \frac{b^{2}}{c^{2}} + \frac{b}{c}$

add_square(eq3)

$\displaystyle a^{2} + a = \frac{b^{2}}{c^{2}} + \frac{b}{c}$

Expression manipulation

factor(eq), simplify(eq), collect(eq) works on both sides of the Equation. The same is true for eq.factor(), eq.simplify(), eq.collect().

var("x, y")
eq4 = Eqn(2**(x**2 + 2*x + 1), x * y**2 + 2*y*x)
eq4

$\displaystyle 2^{x^{2} + 2 x + 1} = x y^{2} + 2 x y$

eq4.factor()

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

eq4.factor(deep=True)

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

simplify(eq4)

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

eq4.collect(x)

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

In addition to the apply(), applylhs(), applyrhs() methods, it is also possible to work directly on a specified side of the equation and perform a complex sequence of manipulations by chaining different methods of the Expr class. This can be done with the following machinery:

• do: work on both sides of the Equation.

• dolhs: work on the LHS.

• dorhs: work on the RHS.

These are not attributes, in the sense that they don’t return anything useful to the user. They are just machinery that allows to work on a particular side. For example:

eq5 = Eqn(sin(x), series(sin(x), x, 0))
eq5

$\displaystyle \sin{\left(x \right)} = x - \frac{x^{3}}{6} + \frac{x^{5}}{120} + O\left(x^{6}\right)$

eq5 = eq5.dorhs.removeO().factor()
eq5

$\displaystyle \sin{\left(x \right)} = \frac{x \left(x^{4} - 20 x^{2} + 120\right)}{120}$

plot(*eq5.args, (x, -pi, pi), legend=True)

<sympy.plotting.backends.matplotlibbackend.matplotlib.MatplotlibBackend at 0x755696cda450>

Substitution and numerical evaluation

Equation.subs() is just a wrapper to SymPy’s Basic.subs() in order to make it aware of objects of type Equation. The following modes of operation are supported:

• eq.subs(old, new)

• eq.subs([(old1, new1), (old2, new2), ...])

• eq.subs({old1: new1, old2: new2, ...})

• eq.subs(other_eq): here, other_eq.lhs is the old value to replace, while other_eq.rhs is the new value to use.

• eq.subs(eq2, eq3): substitute multiple equations. Here, the order of the arguments is important.

var("p V n R T")
ideal_gas_law = Eqn(p * V, n * R * T)
pressure = ideal_gas_law / V
pressure

$\displaystyle p = \frac{R T n}{V}$

var("L atm mol K") # units
d = {R: 0.08206*L*atm/mol/K, T: 273*K, V: 24*L, n: 1*mol}
pressure.subs(d)

$\displaystyle p = 0.9334325 atm$

pressure.evalf(subs=d, n=3)

$\displaystyle p = 0.933 atm$

An example showing the substitution of an equation into another:

var("a:d")
eq_a = Eqn(a * (b + c), b / d)
eq_b = Eqn(b + c, 3 * d)
display(eq_a, eq_b)

$\displaystyle a \left(b + c\right) = \frac{b}{d}$

$\displaystyle b + c = 3 d$

eq_a.subs(eq_b)

$\displaystyle 3 a d = \frac{b}{d}$

Differentiation and Integration

diff(eq, *args, **kwargs) and integrate(eq, *args, **kwargs) works on both sides of the Equation. The same is true for eq.diff(*args, **kwargs) and eq.integrate(*args, **kwargs).

Use eq.applylhs or eq.applyrhs (or eq.dolhs, eq.dorhs) to differentiate only one side of the equation.

e = Eqn(a * x**2 * y + b * x * y**2, x*y - c)
e

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

e.diff(y)

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

e.diff(x, 2)

$\displaystyle 2 a y = 0$

e.integrate((y, 1, 2))

$\displaystyle \frac{3 a x^{2}}{2} + \frac{7 b x}{3} = - c + \frac{3 x}{2}$

e.applylhs(diff, x)

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

e.dorhs.diff(y)

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

e.dorhs.integrate(x)

$\displaystyle a x^{2} y + b x y^{2} = - c x + \frac{x^{2} y}{2}$

Attributes and Methods

• eq.swap: swap the LHS with RHS.

• eq.cross_multiply()

• eq.as_expr(): convert the Equation to a symbolic expression by writing LHS - RHS.

• eq.as_Boolean(): convert the Equation to a SymPy’s Equality (which does not support mathematical operations).

• eq.check(): verify if the LHS is mathematically equivalent to the RHS.

var("a:d")
e = Eqn(a/b, c/d)
e

$\displaystyle \frac{a}{b} = \frac{c}{d}$

e.swap

$\displaystyle \frac{c}{d} = \frac{a}{b}$

e.cross_multiply()

$\displaystyle a d = b c$

e.as_expr()

$\displaystyle \frac{a}{b} - \frac{c}{d}$

equality = e.as_Boolean()
print(type(equality))
equality

<class 'sympy.core.relational.Equality'>

$\displaystyle \frac{a}{b} = \frac{c}{d}$

Note that there is nothing stopping us from writing “wrong” equations. For example:

f = Eqn(0, 1)
f

$\displaystyle 0 = 1$

To verify if the two sides are mathematically equivalent let’s execute:

f.check()

$\displaystyle \text{False}$

If the equivalence can’t be determited, an object of type Equality will be returned:

res = e.check()
print(type(res))
res

<class 'sympy.core.relational.Equality'>

$\displaystyle \frac{a}{b} = \frac{c}{d}$

Solving equations

If one or more objects of type Equation are provided to solve(), then the output will contains object of the same type:

e = Eqn(a*x**2, -b*c - c)
e

$\displaystyle a x^{2} = - b c - c$

solve(e, x)

$\displaystyle \left[ x = - \sqrt{- \frac{c \left(b + 1\right)}{a}}, \ x = \sqrt{- \frac{c \left(b + 1\right)}{a}}\right]$

Differently, if equations are converted to expressions, than solve() will behave exactly like SymPy’s solve():

solve(e.as_expr(), x)

$\displaystyle \left[ - \sqrt{- \frac{c \left(b + 1\right)}{a}}, \ \sqrt{- \frac{c \left(b + 1\right)}{a}}\right]$

Solving multiple equations:

e1 = Eqn(x**2 * y, -2 + x)
e2 = Eqn(y**2 * x, -3 * y)
solve([e1, e2], [x, y])

$\displaystyle \left[ \left[ x = \frac{1}{2}, \ y = -6\right], \ \left[ x = 2, \ y = 0\right]\right]$

solve([e.as_expr() for e in [e1, e2]], [x, y])

$\displaystyle \left[ \left( \frac{1}{2}, \ -6\right), \ \left( 2, \ 0\right)\right]$

solve([e.as_expr() for e in [e1, e2]], [x, y], dict=True)

$\displaystyle \left[ \left\{ x : \frac{1}{2}, \ y : -6\right\}, \ \left\{ x : 2, \ y : 0\right\}\right]$

Configuration options

equation_config is an object containing a few properties to customize the behaviour of the module:

from sympy_equation import equation_config

Arguably the most useful options are :

• equation_config.show_label

• equation_config.integers_as_exact

Hide/Show the equation label

By default, the equation’s label will be hidden:

equation_config.show_label

False

var("a:d")
e = Eqn(a / b, c / d)
e

$\displaystyle \frac{a}{b} = \frac{c}{d}$

When it’s True, a label with the name of the equation in the python environment will be shown on the screen.

equation_config.show_label = True
e

$\displaystyle \frac{a}{b} = \frac{c}{d}\qquad (e)$

Integer as exact

When it’s True and we are running in an IPython/Jupyter environment, it preparses the content of a code line in order to convert integer numbers to sympy’s Integer. In doing so, we can write 2/3, which will be converted to Integer(2)/Integer(3), which then SymPy converts to Rational(2, 3). If False, no preparsing is done, and Python evaluates 2/3 to 0.6666667, which will then be converted by SymPy to a Float.

equation_config.integers_as_exact

False

e = Eqn(a + 1 / 2, b + 3 / 4)
e

$\displaystyle a + 0.5 = b + 0.75\qquad (e)$

equation_config.integers_as_exact = True

e = Eqn(a + 1 / 2, b + 3 / 4)
e

$\displaystyle a + \frac{1}{2} = b + \frac{3}{4}\qquad (e)$

It is reccommended to set this options to True only when executing purely symbolic computations, not when using other numerical libraries, such as Numpy, because it will create hard to debug situations. Consider executing this: np.cos(np.pi / 4). If integers_as_exact = True, this will raise an error because 4 is first replaced with sympy’s Integer(4), then np.pi / Integer(4) becomes a symbolic expression and np.cos is unable to evaluate it.

import numpy as np
np.cos(np.pi / 4)

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
AttributeError: 'Float' object has no attribute 'cos'

The above exception was the direct cause of the following exception:

TypeError                                 Traceback (most recent call last)
Cell In[62], line 3
      1 import numpy as np
----> 3 np .cos (np .pi /Integer (4 ))

TypeError: loop of ufunc does not support argument 0 of type Float which has no callable cos method

Example 1

Given the equations of the volume and mass of a cylindrical pressure vessel, rewrite the mass in terms of the volume and the ratio \(L/D\).

Parameters:

• V: volume;

• M: mass;

• D: diameter;

• p: pressure inside the vessel;

• sigma: stress level on the thin walls;

• rho: density of the material.

equation_config.show_label = False
var("D, L, V, M, rho, p, sigma")
Veq = Eqn(V, pi * D**3 / 6 + pi * D**2 * L / 4)
Meq = Eqn(M, pi * D**2 * rho * p / (2 * sigma) * (L + D / 2))
display(Veq, Meq)

$\displaystyle V = \frac{\pi D^{3}}{6} + \frac{\pi D^{2} L}{4}$

$\displaystyle M = \frac{\pi D^{2} p \rho \left(\frac{D}{2} + L\right)}{2 \sigma}$

Veq = (Veq / D**3).expand()
Veq

$\displaystyle \frac{V}{D^{3}} = \frac{\pi}{6} + \frac{\pi L}{4 D}$

Veq = Veq.collect(pi)
Veq

$\displaystyle \frac{V}{D^{3}} = \pi \left(\frac{1}{6} + \frac{L}{4 D}\right)$

Meq = (Meq / D**3).expand()
Meq

$\displaystyle \frac{M}{D^{3}} = \frac{\pi p \rho}{4 \sigma} + \frac{\pi L p \rho}{2 D \sigma}$

Meq = Meq.collect(pi * rho * p / sigma)
Meq

$\displaystyle \frac{M}{D^{3}} = \frac{\pi p \rho \left(\frac{1}{4} + \frac{L}{2 D}\right)}{\sigma}$

# isolate D**3
Veq = 1 / Veq * V
Veq

$\displaystyle D^{3} = \frac{V}{\pi \left(\frac{1}{6} + \frac{L}{4 D}\right)}$

eq = Meq.subs(Veq)
eq

$\displaystyle \frac{\pi M \left(\frac{1}{6} + \frac{L}{4 D}\right)}{V} = \frac{\pi p \rho \left(\frac{1}{4} + \frac{L}{2 D}\right)}{\sigma}$

eq = eq.cross_multiply()
eq

$\displaystyle \pi M \sigma \left(\frac{1}{6} + \frac{L}{4 D}\right) = \pi V p \rho \left(\frac{1}{4} + \frac{L}{2 D}\right)$

eq = eq / (pi * sigma)
eq

$\displaystyle M \left(\frac{1}{6} + \frac{L}{4 D}\right) = \frac{V p \rho \left(\frac{1}{4} + \frac{L}{2 D}\right)}{\sigma}$

res = eq / eq.lhs * M
res

$\displaystyle M = \frac{V p \rho \left(\frac{1}{4} + \frac{L}{2 D}\right)}{\sigma \left(\frac{1}{6} + \frac{L}{4 D}\right)}$