Моя цель состоит в том, чтобы построить следующий набор связанных ОДУ:Python черчения с помощью dopri5
Метод odeint дает мне проблемы (например, некоторые взятия X_i на отрицательные значения, которые аналитически не делает), поэтому я решил использовать решатель Runge-Kutta четвертого порядка. Я следовал примеры здесь, чтобы использовать метод dopri5 (Using adaptive step sizes with scipy.integrate.ode):
import scipy as sp
import pylab as plt
import numpy as np
import scipy.integrate as spi
#Constants
c13 = 6.2
c14 = 1.0
c21 = 7.3
c32 = 2.4
c34 = 12.7
c42 = 5.7
c43 = 5.0
e12 = 1.5
e23 = 2.5
e24 = 2.0
e31 = 3.0
e41 = 4.8
#Time
t_end = 700
t_start = 0
t_step = 1
t_interval = sp.arange(t_start, t_end, t_step)
#Initial Condition
ic = [0.2,0.3,0.3,0.5]
def model(t,ic):
Eqs= np.zeros((4))
Eqs[0] = (ic[0]*(1-ic[0]*ic[0]-ic[1]*ic[1]-ic[2]*ic[2]-ic[3]*ic[3])-c21*((ic[1]*ic[1])*ic[0])+e31*((ic[2]*ic[2])*ic[0])+e41*((ic[3]*ic[3])*ic[0]))
Eqs[1] = (ic[1]*(1-ic[0]*ic[0]-ic[1]*ic[1]-ic[2]*ic[2]-ic[3]*ic[3])+e12*((ic[0]*ic[0])*ic[1])-c32*((ic[2]*ic[2])*ic[1])-c42*((ic[3]*ic[3])*ic[1]))
Eqs[2] = (ic[2]*(1-ic[0]*ic[0]-ic[1]*ic[1]-ic[2]*ic[2]-ic[3]*ic[3])-c13*((ic[0]*ic[0])*ic[2])+e23*((ic[1]*ic[1])*ic[2])-c43*((ic[3]*ic[3])*ic[2]))
Eqs[3] = (ic[3]*(1-ic[0]*ic[0]-ic[1]*ic[1]-ic[2]*ic[2]-ic[3]*ic[3])-c14*((ic[0]*ic[0])*ic[3])+e24*((ic[1]*ic[1])*ic[3])-c34*((ic[2]*ic[2])*ic[3]))
return Eqs
ode = spi.ode(model)
ode.set_integrator('dopri5')
ode.set_initial_value(ic,t_start)
ts = []
ys = []
while ode.successful() and ode.t < t_end:
ode.integrate(ode.t + t_step)
ts.append(ode.t)
ys.append(ode.y)
t = np.vstack(ts)
x1,x2,x3,x4 = np.vstack(ys).T
plt.subplot(1, 1, 1)
plt.plot(t, x1, 'r', label = 'x1')
plt.plot(t, x2, 'b', label = 'x2')
plt.plot(t, x3, 'g', label = 'x3')
plt.plot(t, x4, 'purple', label = 'x4')
plt.xlim([0,t_end])
plt.legend()
plt.ylim([-0.2,1.3])
plt.show()
, который производит следующий сюжет:
Однако мои участки имеют что-то странное них - x1 Значение (1, 0, 0) является фиксированной точкой динамической системы, для меня не имеет большого значения, что она превышает один (возможно, была причина, если шипы показали какой-то повторяющийся узор, но они выглядят случайными в сюжете, поэтому мне интересно, может ли это иметь больше общего с численная интеграция, а не фактическая динамика). Изменение параметров также меняет этот всплеск (некоторые значения параметров, которые я пробовал, все еще имеют этот пик, но это существенно меньше). Поэтому у меня есть 2 вопроса:
1) Что вызывает появление этих шипов? Моя оригинальная мысль заключалась в том, что это был «t_step», большой, поэтому я попытался поиграть с ним (что приводит меня к моему второму вопросу)
2) Я попытался изменить размер шага от 1 до 0,1. Однако он произвел этот сюжет, который на самом деле выглядит совсем иначе, чем если бы размер шага составлял 1 (что я считал, что более низкий размер шага даст более точные графики?). На этом графике кажется, что время, которое каждый x_i проводит около 1, больше, чем если бы размер шага был 1 вместо 0,1, а высоты спайков равны - как я узнаю, какой сюжет более точным из истинной системы, учитывая, что графики имеют существенные отличия?
Реальная суть проблемы, я работаю смотрит взаимодействие между x3 и x4, но я хочу, чтобы убедиться, что моделирование правильно настроены так, что я могу получить точные результаты о x3 и x4!