# Symbolic Computation
Symbolic computation deals with symbols, representing them exactly, instead of numerical approximations (floating point). 

We will start with the following [borrowed](https://docs.sympy.org/latest/tutorial/intro.html) tutorial to introduce the concepts of SymPy. Devito uses SymPy heavily and builds upon it in its DSL. 


```python
import math

math.sqrt(3)
```




    1.7320508075688772




```python
math.sqrt(8)
```




    2.8284271247461903



$\sqrt(8) = 2\sqrt(2)$, but it's hard to see that here


```python
import sympy
sympy.sqrt(3)
```




$\displaystyle \sqrt{3}$



SymPy can even simplify symbolic computations


```python
sympy.sqrt(8)
```




$\displaystyle 2 \sqrt{2}$




```python
sympy.sqrt(20)
```




$\displaystyle 2 \sqrt{5}$




```python
from sympy import symbols
x, y, z = symbols('x y z')
expr = x + 2*y
expr
```




$\displaystyle x + 2 y$



Note that simply adding two symbols creates an expression. Now let's play around with it. 


```python
expr + 1
```




$\displaystyle x + 2 y + 1$




```python
expr - x
```




$\displaystyle 2 y$



Note that `expr - x` was not `x + 2y -x`


```python
x*expr
```




$\displaystyle x \left(x + 2 y\right)$




```python
from sympy import expand, factor
expanded_expr = expand(x*expr)
expanded_expr
```




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




```python
factor(expanded_expr)
```




$\displaystyle x \left(x + 2 y\right)$




```python
from sympy import diff, sin, exp

diff(sin(x)*exp(x), x)
```




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




```python
from sympy import limit

limit(sin(x)/x, x, 0)
```




$\displaystyle 1$




```python
from sympy import symbols
equation = x**2 + z*x*y + z**3
equation
```




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




```python
x * equation
```




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




```python
expanded_equation = expand(x * equation)
expanded_equation
```




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




```python
factor(expanded_equation)
```




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




```python
diff(expanded_equation,x)
```




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



### Exercise

Solve $x^2 - 2 = 0$ using sympy.solve


```python
# Type solution here
from sympy import solve
solve(x**2 - 2, x)
```




$\displaystyle \left[ - \sqrt{2}, \  \sqrt{2}\right]$



## Pretty printing


```python
from sympy import init_printing, Integral, sqrt

init_printing(use_latex='mathjax')
```


```python
Integral(sqrt(1/x), x)
```




$\displaystyle \int \sqrt{\frac{1}{x}}\, dx$




```python
from sympy import latex

latex(Integral(sqrt(1/x), x))
```




    '\\int \\sqrt{\\frac{1}{x}}\\, dx'



More symbols.
Exercise: fix the following piece of code


```python
# NBVAL_SKIP
# The following piece of code is supposed to fail as it is
# The exercise is to fix the code
expr2 = x + 2*y +3*z
```

### Exercise 

Solve $x + 2*y + 3*z$ for $x$


```python
# Solution here
from sympy import solve
expression = x + 2 * y + 3 * z
solve(expression,x)
```




$\displaystyle \left[ - 2 y - 3 z\right]$



Difference between symbol name and python variable name


```python
x, y = symbols("y z")
```


```python
x
```




$\displaystyle y$




```python
y
```




$\displaystyle z$




```python
# NBVAL_SKIP
# The following code will error until the code in cell 16 above is
# fixed
z
```




$\displaystyle z$



Symbol names can be more than one character long


```python
crazy = symbols('unrelated')

crazy + 1
```




$\displaystyle unrelated + 1$




```python
x = symbols("x")
expr = x + 1
x = 2
```

What happens when I print expr now? Does it print 3?


```python
print(expr)
```

    x + 1


How do we get 3?


```python
x = symbols("x")
expr = x + 1
expr.subs(x, 2)
```




$\displaystyle 3$



## Equalities


```python
x + 1 == 4
```




    False




```python
from sympy import Eq

Eq(x + 1, 4)
```




$\displaystyle x + 1 = 4$



Suppose we want to ask whether $(x + 1)^2 = x^2 + 2x + 1$


```python
(x + 1)**2 == x**2 + 2*x + 1
```




    False




```python
from sympy import simplify

a = (x + 1)**2
b = x**2 + 2*x + 1

simplify(a-b)
```




$\displaystyle 0$



### Exercise 
Write a function that takes two expressions as input, and returns a tuple of two booleans. The first if they are equal symbolically, and the second if they are equal mathematically.

## More operations


```python
z = symbols("z")
expr = x**3 + 4*x*y - z
expr.subs([(x, 2), (y, 4), (z, 0)])
```




$\displaystyle 36$




```python
from sympy import sympify

str_expr = "x**2 + 3*x - 1/2"
expr = sympify(str_expr)
expr
```




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




```python
expr.subs(x, 2)
```




$\displaystyle \frac{19}{2}$




```python
expr = sqrt(8)
```


```python
expr
```




$\displaystyle 2 \sqrt{2}$




```python
expr.evalf()
```




$\displaystyle 2.82842712474619$




```python
from sympy import pi

pi.evalf(100)
```




$\displaystyle 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068$




```python
pi.evalf(1000)
```




$\displaystyle 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587006606315588174881520920962829254091715364367892590360011330530548820466521384146951941511609433057270365759591953092186117381932611793105118548074462379962749567351885752724891227938183011949129833673362440656643086021394946395224737190702179860943702770539217176293176752384674818467669405132000568127145263560827785771342757789609173637178721468440901224953430146549585371050792279689258923542019956112129021960864034418159813629774771309960518707211349999998372978049951059731732816096318595024459455346908302642522308253344685035261931188171010003137838752886587533208381420617177669147303598253490428755468731159562863882353787593751957781857780532171226806613001927876611195909216420199$




```python
from sympy import cos

expr = cos(2*x)
expr.evalf(subs={x: 2.4})
```




$\displaystyle 0.0874989834394464$



### Exercise



```python
from IPython.core.display import Image 
Image(filename='figures/comic.png')
```

Write a function that takes a symbolic expression (like pi), and determines the first place where 789 appears.
Tip: Use the string representation of the number. Python starts counting at 0, but the decimal point offsets this

## Solving an ODE


```python
from sympy import Function

f, g = symbols('f g', cls=Function)
f(x)
```




$\displaystyle f{\left(y \right)}$




```python
f(x).diff()
```




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




```python
diffeq = Eq(f(x).diff(x, x) - 2*f(x).diff(x) + f(x), sin(x))
diffeq
```




$\displaystyle f{\left(y \right)} - 2 \frac{d}{d y} f{\left(y \right)} + \frac{d^{2}}{d y^{2}} f{\left(y \right)} = \sin{\left(y \right)}$




```python
from sympy import dsolve

dsolve(diffeq, f(x))
```




$\displaystyle f{\left(y \right)} = \left(C_{1} + C_{2} y\right) e^{y} + \frac{\cos{\left(y \right)}}{2}$



## Finite Differences


```python

f = Function('f')
dfdx = f(x).diff(x)
dfdx.as_finite_difference()
```




$\displaystyle - f{\left(y - \frac{1}{2} \right)} + f{\left(y + \frac{1}{2} \right)}$




```python
from sympy import Symbol

d2fdx2 = f(x).diff(x, 2)
h = Symbol('h')
d2fdx2.as_finite_difference(h)
```




$\displaystyle - \frac{2 f{\left(y \right)}}{h^{2}} + \frac{f{\left(- h + y \right)}}{h^{2}} + \frac{f{\left(h + y \right)}}{h^{2}}$



Now that we have seen some relevant features of vanilla SymPy, let's move on to Devito, which could be seen as SymPy finite differences on steroids!

40c65c0eed5c9ae91349ddd95d4cd39352f86cde
ipynb
Jupyter Notebook
tinkered with examples/userapi/00_sympy.ipynb
ofmla/Devito-playbox
f3547c7c1bfd82cc32b51179c178f685ecf12e84
MIT
2022-02-08T20:09:00.000Z
2022-02-08T20:09:00.000Z
tinkered with examples/userapi/00_sympy.ipynb
ofmla/Devito-playbox
f3547c7c1bfd82cc32b51179c178f685ecf12e84
MIT
tinkered with examples/userapi/00_sympy.ipynb
ofmla/Devito-playbox
f3547c7c1bfd82cc32b51179c178f685ecf12e84
MIT
Qwen/Qwen-72B
1. YES
2. YES