Решение дифференциального уравнение и построение спектра мощности (Фурье) python

есть дифференциальное уравнение с затуханием X''+2*sigmaХ'+ sinx= Acos(omega*t) для нелинейного гармонического осциллятора. я его решил, но как построить спектр мощности(Фурье) мой код

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import math
def f(y,t, A):
    omega = 1.01
    gamma = 0.01
    wj = 1
    y1, y2 = y #вводим имена новых функций
    return [y2,A*math.cos(omega*t)-gamma*y2-wj**2*math.sin(y1)] 
t = np.linspace(0,10,200)
y0 = [0.1,0.2] #начальные условия
A = 0.0
w1 = odeint(f, y0, t, args=(A,))
A = 0.1
w2 = odeint(f, y0, t, args=(A,))
A = 0.2
w3 = odeint(f, y0, t, args=(A,))
A = 3
w4 = odeint(f, y0, t, args=(A,))
A = 4
w5 = odeint(f, y0, t, args=(A,))
A = 5
w6 = odeint(f, y0, t, args=(A,))
y1 = w1[:,0] #вектор значений u
y2 = w1[:,1] #вектор значение u'
y3 = w2[:,0] #вектор значений u
y4 = w2[:,1] #вектор значение u'
y5 = w3[:,0] #вектор значений u
y6 = w3[:,1] #вектор значение u'
plt.plot(t,y1,'r-',linewidth=2, label = "A = 0")
plt.plot(t,y3,'b-',linewidth=2, label = "A = 0.1")
plt.plot(t,y5,'k-',linewidth=2, label = "A = 0.2")
plt.ylabel("u(t)")
plt.xlabel("t")
plt.legend()
plt.show()
plt.plot(t,y2,'r:',linewidth=2, label = "A = 0")
plt.plot(t,y4,'b:',linewidth=2, label = "A = 0.1")
plt.plot(t,y6,'k:',linewidth=2, label = "A = 0.2")
plt.ylabel("u'(t)")
plt.xlabel("t")
plt.legend()
plt.show()
#график фазовой траектории
t = np.linspace (0,10,150)
A = 0.0
[y1,y2] = odeint(f, y0, t, args=(A,), full_output = False).T
A = 0.1
[y3,y4] = odeint(f, y0, t, args=(A,), full_output = False).T
A = 0.2
[y5,y6] = odeint(f, y0, t, args=(A,), full_output = False).T
#A = 0.3
#[y1,y2] = odeint(f, y0, t, args=(A,), full_output = False).T
plt.plot(y1,y2,'r-',linewidth=2, label = "A = 0")
plt.plot(y3,y4,'b-',linewidth=2, label = "A = 0.1")
plt.plot(y5,y6,'k-',linewidth=2, label = "A = 0.2")
plt.ylabel("u'(t)")
plt.xlabel("u(t)")
plt.legend()
plt.show()

Ответы (0 шт):