<jupyter_start><jupyter_text>Python - Symbolic Mathematics (`sympy`)<jupyter_code>%matplotlib inline

import sympy as sp
import numpy as np
import matplotlib.pyplot as plt
sp.init_printing()<jupyter_output><empty_output><jupyter_text>`sympy` treats stuff fundementally different than `numpy`<jupyter_code>np.sqrt(8)
sp.sqrt(8)
np.pi
sp.pi<jupyter_output><empty_output><jupyter_text>sympy has its own way to handle rational numbers<jupyter_code>sp.Rational(3,5)
sp.Rational(3,5) + sp.Rational(1,2)<jupyter_output><empty_output><jupyter_text>Least Common Multiple<jupyter_code>sp.lcm(2,5)<jupyter_output><empty_output><jupyter_text>Adding `.n()` to the end of a sympy expression will `evaluate` expression<jupyter_code>sp.pi.n()
sp.pi.n(100)<jupyter_output><empty_output><jupyter_text>`nsimplify()` will sort-of do the reverse<jupyter_code>sp.nsimplify(0.125)
sp.nsimplify(4.242640687119286)
sp.nsimplify(sp.pi, tolerance=1e-2)
sp.nsimplify(sp.pi, tolerance=1e-5)
sp.nsimplify(sp.pi, tolerance=1e-6)<jupyter_output><empty_output><jupyter_text>... to $\infty$ and beyond<jupyter_code>sp.oo
sp.oo + 3
1e9 < sp.oo<jupyter_output><empty_output><jupyter_text>Primes<jupyter_code>list(sp.primerange(0,100))
sp.nextprime(2019)
sp.factorint(2272019)<jupyter_output><empty_output><jupyter_text>Symbolic You have to explicitly tell `SymPy` what symbols you want to use.<jupyter_code>x, y, z = sp.symbols('x y z')
a, b, c = sp.symbols('a b c')
mu, rho = sp.symbols('mu rho')<jupyter_output><empty_output><jupyter_text>Expressions are then able use these symbols<jupyter_code>my_equation = 2*x + y

my_equation
my_equation + 3
my_equation - x
my_equation / x
my_greek_equation = mu**2 / rho * (a + b)

my_greek_equation<jupyter_output><empty_output><jupyter_text>`SymPy` has all sorts of ways to manipulates symbolic equations<jupyter_code>sp.simplify(my_equation / x)
another_equation = (x + 2) * (x - 3)

another_equation
sp.expand(another_equation)
long_equation = 2*y*x**3 + 12*x**2 - x + 3 - 8*x**2 + 4*x + x**3 + 5 + 2*y*x**2 + x*y

long_equation
sp.collect(long_equation,x)
sp.collect(long_equation,y)<jupyter_output><empty_output><jupyter_text>You can evaluate equations for specific values<jupyter_code>trig_equation = a*sp.sin(2*x + y) + b*sp.cos(x + 2*y)

trig_equation
trig_equation.subs({a:2, b:3, x:4, y:5})
trig_equation.subs({a:2, b:3, x:4, y:5}).n()
sp.expand(trig_equation, trig=True)
sp.collect(sp.expand(trig_equation, trig=True),sp.cos(x))<jupyter_output><empty_output><jupyter_text>You can evaluate/simplify equations sybolically<jupyter_code>my_equation_xyz = sp.sqrt((x * (y - 4*x)) / (z * (y - 3*x)))

my_equation_xyz
my_equation_x = (3 * a * y) / (9 * a - y)

my_equation_x
my_new_xyz = my_equation_xyz.subs(x, my_equation_x)

my_new_xyz
sp.simplify(my_new_xyz)<jupyter_output><empty_output><jupyter_text>System of equations $$\begin{array}{c}9x - 2y = 5 \\-2x + 6y = 10 \\\end{array}\hspace{3cm}\left[\begin{array}{cc}9 & -2 \\-2 & 6 \\\end{array}\right]\left[\begin{array}{c}x\\y\end{array}\right]=\left[\begin{array}{c}5\\10\end{array}\right]$$<jupyter_code>a_matrix = sp.Matrix([[9, -2],
                      [-2, 6]])

b_matrix = sp.Matrix([[5],
                      [10]])
a_matrix, b_matrix
a_matrix.inv()
a_matrix.inv() * a_matrix
a_matrix.inv() * b_matrix<jupyter_output><empty_output><jupyter_text>Solving equations - `solve`$$\begin{array}{c}9x - 2y = 5 \\-2x + 6y = 10 \\\end{array}$$<jupyter_code>equation_a = 9*x - 2*y - 5
equation_b = -2*x + 6*y - 10
sp.solve([equation_a, equation_b], [x,y])
yet_another_equation = x**3 + x + 10

yet_another_equation
sp.solve(yet_another_equation,x)<jupyter_output><empty_output><jupyter_text>... complex numbers<jupyter_code>sp.I
a_complex_number = 2 + 3 * sp.I

a_complex_number
sp.re(a_complex_number), sp.im(a_complex_number)<jupyter_output><empty_output><jupyter_text>... solving more symbolically<jupyter_code>symbolic_equation = a*x**2 + b*x +c

symbolic_equation
sp.solve(symbolic_equation, x)<jupyter_output><empty_output><jupyter_text>Calculus<jupyter_code>symbolic_equation
sp.diff(symbolic_equation,x)
sp.diff(symbolic_equation,x,2)
sp.integrate(symbolic_equation,x)
sp.integrate(symbolic_equation,(x,0,5))   # limits x = 0 to 5
sp.integrate(symbolic_equation,(x,0,5)).subs({a:2, b:7, c:3}).n()
trig_equation
sp.diff(trig_equation,x)
sp.integrate(trig_equation,x)<jupyter_output><empty_output><jupyter_text>Ordinary differential equation - `dsolve`<jupyter_code>f = sp.Function('f')
f(x)
sp.Derivative(f(x),x,x)
equation_ode = sp.Derivative(f(x), x, x) + 9*f(x)

equation_ode
sp.dsolve(equation_ode, f(x))<jupyter_output><empty_output><jupyter_text>Limits<jupyter_code>limit_equation = (1 + (1 / x)) ** x

limit_equation<jupyter_output><empty_output><jupyter_text>$$\lim _{x\to 5 }\left(1+{\frac {1}{x}}\right)^{x}$$<jupyter_code>sp.limit(limit_equation, x, 5)
sp.limit(limit_equation, x, 5).n()<jupyter_output><empty_output><jupyter_text>$$\lim _{x\to \infty }\left(1+{\frac {1}{x}}\right)^{x}$$<jupyter_code>sp.limit(limit_equation, x, sp.oo)
sp.limit(limit_equation, x, sp.oo).n()<jupyter_output><empty_output><jupyter_text>Summation $$ \sum{\frac {x^{a}}{a!}} $$<jupyter_code>sum_equation = x**a / sp.factorial(a)

sum_equation<jupyter_output><empty_output><jupyter_text>$$ \sum _{a=0}^{3}{\frac {x^{a}}{a!}} $$<jupyter_code>sp.summation(sum_equation, [a, 0, 3])
sp.summation(sum_equation.subs({x:1}), [a, 0, 3])
sp.summation(sum_equation.subs({x:1}), [a, 0, 3]).n()<jupyter_output><empty_output><jupyter_text>$$ \sum _{a=0}^{10}{\frac {x^{a}}{a!}} $$<jupyter_code>sp.summation(sum_equation.subs({x:1}), [a, 0, 10]).n()<jupyter_output><empty_output><jupyter_text>$$ \sum _{a=0}^{\infty}{\frac {x^{a}}{a!}} $$<jupyter_code>sp.summation(sum_equation, [a, 0, sp.oo])<jupyter_output><empty_output><jupyter_text>Let's do some graphing stuff ... $$\large y_1 = \frac{x^3}{4} - 3x$$<jupyter_code>my_np_x = np.linspace(-2*np.pi,2*np.pi,200)
my_np_y1 = my_np_x ** 3 / 4 - 3 * my_np_x
fig,ax = plt.subplots(1,1)
fig.set_size_inches(10,4)

fig.tight_layout()

ax.set_ylim(-7,7)
ax.set_xlim(-np.pi,np.pi)

ax.set_xlabel("This is X")
ax.set_ylabel("This is Y")

ax.plot(my_np_x, my_np_y1, color='r', marker='None', linestyle='-', linewidth=4);<jupyter_output><empty_output><jupyter_text>Fourier Series<jupyter_code>my_sp_y1 = x ** 3 / 4 - 3 * x

my_sp_y1
my_fourier = sp.fourier_series(my_sp_y1, (x, -sp.pi, sp.pi))

my_fourier
my_fourier.truncate(3).n(2)
my_np_1term = -4.1 * np.sin(my_np_x)
my_np_2term = -4.1 * np.sin(my_np_x) + 0.91 * np.sin(2*my_np_x)
my_np_3term = -4.1 * np.sin(my_np_x) + 0.91 * np.sin(2*my_np_x) - 0.47 * np.sin(3*my_np_x)
fig,ax = plt.subplots(1,1)
fig.set_size_inches(10,4)

fig.tight_layout()

ax.set_ylim(-7,7)
ax.set_xlim(-np.pi,np.pi)

ax.set_xlabel("This is X")
ax.set_ylabel("This is Y")

ax.plot(my_np_x, my_np_y1, color='r', marker='None', linestyle='-', linewidth=8)

ax.plot(my_np_x, my_np_1term, color='b', marker='None', linestyle='--', label="1-term")
ax.plot(my_np_x, my_np_2term, color='g', marker='None', linestyle='--', label="2-term")
ax.plot(my_np_x, my_np_3term, color='k', marker='None', linestyle='--', label="3-term")

ax.legend(loc = 0);<jupyter_output><empty_output><jupyter_text>Another Function$$\large y_2 = 2\,\sin(5x) \ e^{-x}$$<jupyter_code>my_np_y2 = 2 * np.sin(5 * my_np_x) * np.exp(-my_np_x)
fig,ax = plt.subplots(1,1)
fig.set_size_inches(10,4)

fig.tight_layout()

ax.set_ylim(-10,10)
ax.set_xlim(-np.pi,np.pi)

ax.set_xlabel("This is X")
ax.set_ylabel("This is Y")

ax.plot(my_np_x, my_np_y2, color='r', marker='None', linestyle='-', linewidth=4);<jupyter_output><empty_output><jupyter_text>Taylor Expansions<jupyter_code>my_sp_y2 = 2 * sp.sin(5 * x) * sp.exp(-x)
my_sp_y2
my_taylor = sp.series(my_sp_y2, x, x0 = 0)

my_taylor<jupyter_output><empty_output><jupyter_text>if you want more terms* n = magnitude of the highest term* n = 8 means all terms up to x$^{8}$ or $\mathcal{O}(8)$<jupyter_code>my_taylor = sp.series(my_sp_y2, x, x0 = 0, n=8)

my_taylor
my_taylor.removeO()
my_taylor.removeO().n(2)<jupyter_output><empty_output><jupyter_text>General Equation Solving - `nsolve` $$\large y_1 = \frac{x^3}{4} - 3x\\\large y_2 = 2\,\sin(5x) \ e^{-x}$$  Where do they cross? - The graph<jupyter_code>fig,ax = plt.subplots(1,1)
fig.set_size_inches(10,4)

fig.tight_layout()

ax.set_ylim(-5,5)
ax.set_xlim(-np.pi,4)

ax.set_xlabel("This is X")
ax.set_ylabel("This is Y")

ax.plot(my_np_x, my_np_y1, color='b', marker='None', linestyle='--', linewidth = 4)
ax.plot(my_np_x, my_np_y2, color='r', marker='None', linestyle='-', linewidth = 4);<jupyter_output><empty_output><jupyter_text>Where do they cross? - The `sympy` solution<jupyter_code>my_sp_y1, my_sp_y2
my_guess = 3.3

sp.nsolve(my_sp_y1 - my_sp_y2, x, my_guess)
all_guesses = (3.3, 0, -0.75)

for val in all_guesses:
    result = sp.nsolve(my_sp_y1 - my_sp_y2, x, val)
    print(result)<jupyter_output><empty_output><jupyter_text>Your guess has to be (somewhat) close or the solution will not converge:<jupyter_code>my_guess = -40

sp.nsolve(my_sp_y1 - my_sp_y2, x, my_guess)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18263835405435142, 'fraction_numerical': 0.03366707008471779, 'mean_word_length': 6.619446772841576, 'pattern_counts': {'":': 0, '<': 138, '<?xml version=': 0, '>': 137, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '34', 'dir': 'jupyter-structured-clean-dedup', 'id': '230050', 'license': 'MIT', 'n_tokens_mistral': 4252, 'n_tokens_neox': 4053, 'n_words': 606, 'path': 'Python_SymPy.ipynb', 'repo_name': 'UWashington-Astro300/Astro300-Wtr19', 'size': '9089'}