Imitate ode45 function from MATLAB in Python
Solution 1
As @LutzL mentioned, you can use the newer API, solve_ivp
.
results = solve_ivp(obj_func, t_span, y0, t_eval = time_series)
If t_eval
is not specified, then you won't have one record per one timestamp, which is mostly the cases I assume.
Another side note is that for odeint
and often other integrators, the output array is a ndarray
of a shape of [len(time), len(states)]
, however for solve_ivp
, the output is a list(length of state vector)
of 1-dimension ndarray(which length is equal to t_eval
).
So you have to merge it if you want the same order. You can do so by:
Y =results
merged = np.hstack([i.reshape(-1,1) for i in Y.y])
First you need to reshape to make it a [n,1]
array, and merge it horizontally.
Hope this helps!
Solution 2
The function scipy.integrate.solve_ivp uses the method RK45 by default, similar the method used by Matlab's function ODE45 as both use the Dormand-Pierce formulas with fourth-order method accuracy.
vdp1 = @(T,Y) [Y(2); (1 - Y(1)^2) * Y(2) - Y(1)];
[T,Y] = ode45 (vdp1, [0, 20], [2, 0]);
from scipy.integrate import solve_ivp
vdp1 = lambda T,Y: [Y[1], (1 - Y[0]**2) * Y[1] - Y[0]]
sol = solve_ivp (vdp1, [0, 20], [2, 0])
T = sol.t
Y = sol.y
Ordinary Differential Equations (solve_ivp)
Migui Mag
Updated on March 24, 2021Comments
-
Migui Mag about 3 years
I am wondering how to export MATLAB function ode45 to python. According to the documentation is should be as follows:
MATLAB: [t,y]=ode45(@vdp1,[0 20],[2 0]); Python: import numpy as np def vdp1(t,y): dydt= np.array([y[1], (1-y[0]**2)*y[1]-y[0]]) return dydt import scipy integrate l=scipy.integrate.ode(vdp1([0,20],[2,0])).set_integrator("dopri5")
The results are completely different, Matlab returns different dimensions than Python.
-
extraymond about 6 yearsThx for cleanup.
-
crobar over 5 yearsDoes this do variable step sizes like ode45?
-
Lin over 5 yearsdon't you hate when someone come, ask a question and do not pick one as correct? (also when someone resurrect an old post)
-
Lutz Lehmann over 3 yearsNote that Matlab's ode45 has a "Refine" option that defaults to 4. That is, per computed point it adds 3 additional interpolated points from the segment to the return vectors. There is no such option in python solve_ivp, so that with the same tolerances leading to approximately the same work ode45 will return more points that give smoother looking plots.