To make it clearer I want to plot the solutions of the 2nd order differential equation of the damping oscillation for a pendulum. Link to wiki about the equations used : https://en.wikipedia.org/wiki/Harmonic_oscillator
from sympy.interactive import printing
printing.init_printing(use_latex=True)
import numpy as np
import scipy as sp
from sympy import*
mport sympy as syp
`from scipy.integrate import odeint
import matplotlib.pyplot as plt
t=syp.Symbol('t')
x=syp.Function('x')(t)
m=2.0
k=5.0
a=0.5
z=a/(2.0*np.sqrt(m*k))
w=np.sqrt(k/m)
eq=x.diff(t,t)+2.0*z*w*x.diff(t)+w**2.0*x
dsolve(eq,t,0,ics={eq(1.0):0,eq(2.0):5})
You're not constructing the ics
argument as intended:
In [6]: dsolve(eq, ics={x.subs(t, 1.0): 0, x.subs(t, 2.0): 5})
Out[6]:
-0.125⋅t
x(t) = (-0.0346285740992263⋅sin(1.57619002661481⋅t) - 6.42012708343871⋅cos(1.57619002661481⋅t))⋅ℯ
The answer comes out nicer (subjectively) if you don't use floats. Also I find it more natural to keep the variable x
as the function x
rather than the applied function x(t)
e.g.:
In [15]: x = Function('x')
In [16]: x
Out[16]: x
In [17]: x(t)
Out[17]: x(t)
In [18]: eq = x(t).diff(t, 2) + x(t).diff(t)/4 + 5*x(t)/2
In [19]: eq
Out[19]:
d
──(x(t)) 2
5⋅x(t) dt d
────── + ──────── + ───(x(t))
2 4 2
dt
In [20]: dsolve(eq, x(t), ics={x(1): 0, x(2): 5})
Out[20]:
⎛ 1/4 ⎛√159⋅t⎞ ⎞ -t
⎜5⋅ℯ ⋅sin⎜──────⎟ ⎟ ───
⎜ ⎝ 8 ⎠ 1/4 ⎛√159⋅t⎞⎟ 8
x(t) = ⎜────────────────── - 5⋅ℯ ⋅cos⎜──────⎟⎟⋅ℯ
⎜ ⎛√159⎞ ⎝ 8 ⎠⎟
⎜ tan⎜────⎟ ⎟
⎝ ⎝ 8 ⎠ ⎠