Решение системы ODE с медленным замедлением - PullRequest
Новая
       21

Решение системы ODE с медленным замедлением

0 голосов
/ 05 июня 2019

Я пытаюсь решить систему ODE с scipy.integrate.ode, используя метод Адама или bdf. Scipy обеспечивает хорошие и быстрые результаты, когда некоторые параметры в ODE устанавливаются на ноль, но в противном случае кажется, что поиск решения очень сильно замедляется. Мне просто интересно, что вызывает это замедление? 

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import ode
def func(t,y,p):
    wk,w1,w2,B,t0,u,vi = p
    x0,x1,x2,x3,x4,x5=y
    dE = np.exp(1j*(w1-w2)/2*t)
    dEc=np.exp(-1j*(w1-w2)/2*t)
    mw1 = B/2*np.exp(1j*vi)*np.exp(1j*(wk-w1)*t/2)
    mw2 = B/2*np.exp(1j*vi)*np.exp(1j*(wk-w2)*t/2)
    mw1c = B/2*np.exp(-1j*vi)*np.exp(-1j*(wk-w1)*t/2)
    mw2c = B/2*np.exp(-1j*vi)*np.exp(-1j*(wk-w2)*t/2)
    f = [mw1*x2+x1*mw2,
         t0*x4*dE+t0*x5*dE+mw1*x3+mw2c*x0,
         -t0*x4*dEc-t0*x5*dEc+mw1c*x0+mw2*x3,
         mw1c*x1+mw2c*x2,
         t0*x1*dEc-t0*x2*dE+x4*u,
         t0*x1*dEc-t0*x2*dE+x5*u]
    g = [-1j*2*np.pi*i for i in f]
    return g
w1=2e9
w2=1e9
t0=0
B=5e6
vi=0
wk=2e9
u=1e12
p = [wk,w1,w2,B,t0,u,vi]
y0 = [1,0,0,0,0,0]
t = np.linspace(0,100e-9,1001)
dt = t[1]-t[0]

r = ode(func).set_integrator('zvode',atol=1e-8,rtol=1e-6
       ,nsteps=200000,method='bdf')
r.set_initial_value(y0,0).set_f_params(p)
psi = []
while r.successful() and r.t < t[-1]:
    psi.append(r.integrate(r.t+dt))
for i in range(6):
    plt.plot(abs(np.array(psi)[:,i])**2)
plt.show()

Данный код приведет к синусоидальному решению, которое является правильным и решается очень быстро. Однако, когда я устанавливаю параметр t0 в конечное значение, например, t0 = 5 (предпочтительно 4,3e9). Тогда кажется невероятное замедление поиска решения.

...