simulation - Dynamical System in Python: Getting around Stiffness -
i have been learning how use python (still beginner), , trying plot four-dimensional dynamical system (and it's stochastic equivalent). parameters, observe dynamical behavior expect. however, parameters (such below in code), encountering described "stiffness". plots seem stuck @ fixed value , sit there time (analytically, should cycle between 0 , 1). wondering if there way around this. knowledge python limited. right i'm integrating using fourth order runge-kutta though encountering same problem using odeint, , stochastic version integrated using euler-maruyama (this 1 i'm familiar , seems common 1 used). don't know if there's easy way fix this? i've been trying read online, honestly, lot of things i'm reading way beyond level of understanding of python.
import scipy sp import pylab plt import matplotlib.pyplot import numpy np scipy.integrate import odeint import scipy.integrate spi #constants c13 = 6.2 c14 = 4.2 c21 = 7.3 c32 = 2.4 c34 = 12.7 c42 = 5.7 c43 = 5.0 e12 = 0.5 e23 = 2.5 e24 = 2.0 e31 = 2.0 e41 = 4.8 b = np.array([1, 1, 1,1]) n=len(b) = = np.array([[1, 1+c21, 1-e31, 1-e41], [1-e12, 1, 1+c32, 1+c42], [c13+1, 1-e23, 1, 1+c43], [c14+1, 1-e24, 1+c34, 1]]) #time t = 1000 dt = 0.01 t = sp.arange(0.0, t, dt) n = t.size sqrtdt = np.sqrt(dt) #initial condition ic = [0.5,0.5,0.6,0.7] def model(m,t): dmdt = np.dot(np.diag(m), b.t-np.dot(a,np.dot(np.diag(m),m))) return dmdt y = odeint(model, ic, t) y1 = y[:,0] y2 = y[:,1] y3 = y[:,2] y4 = y[:,3] plt.subplot(4, 1, 1) plt.plot(t, y1, 'b') plt.ylim([0,1.2]) plt.subplot(4, 1, 2) plt.plot(t, y2, 'r') plt.ylim([0,1.2]) plt.subplot(4, 1, 3) plt.plot(t, y3, 'g') plt.ylim([0,1.2]) plt.subplot(4, 1, 4) plt.plot(t, y4, 'purple') plt.xlim([0,t]) plt.ylim([0,1.2]) plt.xlabel('time') plt.show()
Comments
Post a Comment