Я пытаюсь решить систему Лоренца , используя метод Рунге Кутты 4-го порядка, где
dx/dt=a*(y-x)
dy/dt=x(b-z)-y
dx/dt=x*y-c*z
Поскольку эта система не зависит от времени, возможно, эту часть в итерации игнорировать, поэтому я просто
йХ = Р (х, у, г)
def func(x0):
a=10
b=38.63
c=8/3
fx=a*(x0[1]-x0[0])
fy=x0[0]*(b-x0[2])-x0[1]
fz=x0[0]*x0[1]-c*x0[2]
return np.array([fx,fy,fz])
def kcontants(f,h,x0):
k0=h*f(x0)
k1=h*f(f(x0)+k0/2)
k2=h*f(f(x0)+k1/2)
k3=h*f(f(x0)+k2)
#note returned K is a matrix
return np.array([k0,k1,k2,k3])
x0=np.array([-8,8,27])
h=0.001
t=np.arange(0,50,h)
result=np.zeros([len(t),3])
for time in range(len(t)):
if time==0:
k=kcontants(func,h,x0)
result[time]=func(x0)+(1/6)*(k[0]+2*k[1]+2*k[2]+k[3])
else:
k=kcontants(func,h,result[time-1])
result[time]=result[time-1]+(1/6)*(k[0]+2*k[1]+2*k[2]+k[3])
В результате должны получиться аттракторы Лоренца, однако мой код расходится вокруг пятой итерации, и это потому, что контуры, которые я создаю в kconstants
, делают, однако я проверил, и я почти уверен, что имплементация Rge Kutta не виновата ... (по крайней мере, я так думаю)
редактирование:
Нашел похожий пост , но пока не могу понять, что я делаю не так