#166 Solution of a differential equation using Bubnov-Galerkin method with Sympy package
The problem and solution in this pdf file:
ex166
import sympy
from matplotlib import pyplot as plt
import seaborn as sns
import numpy as np
from sympy.utilities.lambdify import lambdify
from sympy import simplify
from sympy import *
import os
os.chdir(r'D:\projects\wordpress\ex166')
a, b, x = sympy.symbols("a, b, x")
C0, C1, C2, C3 = sympy.symbols("C0, C1, C2, C3")
u0 = sympy.Function("u0")
u1 = sympy.Function("u1")
u2 = sympy.Function("u2")
u3 = sympy.Function("u3")
y = sympy.Function("y")
u0=1 - x
u1=x*(1 - x)
u2=x**2*(1 - x)
u3=x**3*(1 - x)
y = u0 + C1*u1 + C2*u2 + C3*u3
dy = sympy.diff(y, x)
d2y = sympy.diff(dy, x)
R=d2y+x*dy+y-2*x
a=0
b=1
i1 = sympy.integrate(u1*R, (x, a, b))
i2 = sympy.integrate(u2*R, (x, a, b))
i3 = sympy.integrate(u3*R, (x, a, b))
sols=sympy.solve([i1, i2, i3], (C1,C2,C3), dict=True)
C1_sol= [sol[C1] for sol in sols][0]
C2_sol= [sol[C2] for sol in sols][0]
C3_sol= [sol[C3] for sol in sols][0]
C1_sol_float = "{:.4f}".format(float(C1_sol))
C2_sol_float = "{:.4f}".format(float(C2_sol))
C3_sol_float = "{:.4f}".format(float(C3_sol))
print("C1 =",C1_sol_float)
print("C2 =",C2_sol_float)
print("C3 =",C3_sol_float)
y = u0 + C1_sol*u1 + C2_sol*u2 + C3_sol*u3
y_eq=nfloat(simplify(y),5)
print(y_eq)
func = lambdify(x, y,'numpy')
xvals = np.arange(a,b,.01)
yvals = func(xvals)
sns.set_theme()
plt.rcParams["figure.figsize"] = (8,8)
fig, ax = plt.subplots(1,1,subplot_kw=dict(aspect='equal'))
ax.plot(xvals, yvals)
ax.set_xlabel('x')
ax.set_ylabel('y(x)')
plt.text(0,0,'y = {}'.format(y_eq))
plt.savefig("figure_ex166.pdf")
plt.show()
Recent Comments