FinX
FinX

Reputation: 13

Call Matlab's ODE45 from python

I would like to call Matlab's ode45 function from python just like MATLAB's isprime() function is called in the following python code

import matlab.engine
eng = matlab.engine.start_matlab()
tf = eng.isprime(37)
print(tf)

but I don't really know how to implement the code in python. I've tried

import matlab.engine

def  dydt(t,y):
              dydt= 2*y
              return dydt
          
eng = matlab.engine.start_matlab()
T,Y=eng.ode45(dydt,[0, 20],[2, 0])

and it returned :

line 73, in call out=_stdout, err=_stderr)

TypeError: unsupported Python data type: function.

I've already installed the MATLAB Engine API for Python according to the following link https://www.mathworks.com/help/matlab/matlab_external/install-the-matlab-engine-for-python.html and tested it successfully using examples in the following link https://www.mathworks.com/help/matlab/matlab_external/call-matlab-functions-from-python.html .

I found the following thread that might be relevant Imitate ode45 function from MATLAB in Python but from what I see they just used python's libraries to solve the ODEs and don't call the MATLAB's ODE45 from python. Maybe I just don't understand that correctly.

Could anyone help me make progress with that issue? Thanks!

Upvotes: 1

Views: 1467

Answers (2)

Paolo
Paolo

Reputation: 26034

Reproducing the example from the doc:

import matlab.engine
eng = matlab.engine.start_matlab()
eng.eval('tspan = [0 5];',nargout = 0)
eng.eval('y0 = 0;',nargout = 0)
t,y=eng.eval('ode45(@(t,y) 2*t, tspan, y0)',nargout = 2)

Now you can plot the results too:

eng.plot(t,y,'-o')

Upvotes: 1

norok2
norok2

Reputation: 26886

You could be using scipy.integrate.RK45(), or, possibly preferably, scipy.integrate.solve_ivp():

from scipy.integrate import solve_ivp
def exponential_decay(t, y): return -0.5 * y
sol = solve_ivp(exponential_decay, [0, 10], [2, 4, 8])
print(sol.t)
# [  0.           0.11487653   1.26364188   3.06061781   4.85759374
#    6.65456967   8.4515456   10.        ]
print(sol.y)
# [[2.         1.88836035 1.06327177 0.43319312 0.17648948 0.0719045
#   0.02929499 0.01350938]
#  [4.         3.7767207  2.12654355 0.86638624 0.35297895 0.143809
#   0.05858998 0.02701876]
#  [8.         7.5534414  4.25308709 1.73277247 0.7059579  0.287618
#   0.11717996 0.05403753]]

(example from the official doc.)

Upvotes: 1

Related Questions