python - 将 matlab ODE 求解器转换为 python

Question

我正在尝试将 MATLAB 代码转换为 python，但得到的答案完全不同。我尝试使用 scipy.ode、solve_ivp 和 odeint。运行代码时，我得到的值范围从 1 到 0.2，但在 MATLAB 中，它们的范围从 30 到 70。
MATLAB 代码:

function dydt = onegroup(t,y,tsi,rho)
%Point Reactor Kinetics equation parameters
Lambda      = 10^-4;
beta        = 0.0065;
lambda      = 0.08;
%Reactivity
rho         = interp1(tsi,rho,t);

dydt = zeros(2,1);
%One Group-Delayed Precursor
dydt(1) = -lambda*y(1)+beta*y(2);
%Power
dydt(2) = ((rho-beta)/Lambda)*y(2)+(lambda*y(1))/Lambda;
end

输入文件如下:

%%
clear;
clc;
tsi=linspace(0,20,21);
rho=ones(1,21)*0.0025;
y0= [1; 0];
ts=[0 20];

ode23(@(t,y) onegroup(t,y,tsi,rho),ts,y0)

我的python代码如下:

from scipy.integrate import ode                            
from numpy import arange,vstack,array,sqrt,ones                   
from pylab import plot,close,show,xlabel,ylabel,title,grid,pi,clf     
from scipy.interpolate import interp1d                                                    

# Function defining derivates of the positions and velocities.                                                            
def dydt(t,y,tsi,rho):

    Lambda      = 10^-4
    beta        = 0.0065
    lambda2      = 0.08

    rho = interp1d(tsi,rho, fill_value = 'extrapolate')


    #one group delayed precursor
    dydt1 = (-lambda2*y[0] + beta*y[1])
    #power
    dydt2 = (((rho(t)-beta)/Lambda)*y[1]+(lambda2*y[0])/Lambda)
    return array([dydt1, dydt2])

'''
#Final time
x = 21
#Time steps
dt = 1
tsi=np.arange(0,x,dt)
j0 = 0
times = np.arange(0,x,dt)

dt = 1

rho=np.ones(x)*0.0025
y0= [1,0]
t0 = [0,x-1]
'''

# Initial Conditions
y0, t0 = [1.,0.], 0.

# Model parameters 
k = arange(0,21) 
m = ones(21)*0.0025      

# CREATE ODE OBJECT
i = ode(dydt) # Create an ode object and bind the rhs function.
i.set_integrator('dopri5')  # Which integrator to use. 
i.set_initial_value(y0,t0)  # The initial values
# Define parameters for the derivatives function. 
# These will be passed  to the function at each time.
i.set_f_params(k,m)         

tf = 21        # Final time
dt = 1                  # Output interval
yf=[y0]                    # List for storing the output
times = arange(t0,tf,dt)    # Times to evaluate a solution. 

# Main loop for the integration
for t in times[1:]:
    i.integrate(i.t+dt)
    yf.append(i.y)

yf = array(yf)

Answer 1

^在python中是按位logical and .

用

Lambda = 1e-4

对于那个参数。