求助好心人：老版本matlab的Runge-Kutte rk45怎么转为ode45

youbindi

clear;
clc;
A_1 = 0.11*pi/180;
A_2 = 0.08*pi/180;
t0 = 0;
tfinal = 20;
tspan = [t0, tfinal];
x0 = [0,0]';
[t,x] = rk45('sfile45', tspan , x0); %老版本rk45怎么转换为ode45格式呢？？  
x_err = (A_1*sin(2*x(:,1))+A_2*sin(4*x(:,1)))*180/pi;
x_out = (x(:,1)/N)+(x_err*pi/180);
plot(t,x_err,'r');
hold on;
plot(t,x_out,'--b');
grid;


% function sfile45          % sfile45.m
function[sys,x0] = sfile45(t,x,u,flag)   %关键是怎么函数怎么改写呢，谢谢好心人!
N = 50;
u_d = pi*67.5/180;
A_1 = 0.11*pi/180;
A_2 = 0.08*pi/180;
B_in = 1.83e-4;  
B_L = 1.0e-3;   
J_in = 4.5e-3;  
J_L = 5.0e-2;   
K_P = 1.17;      
K_D = 0.3225;   
u = u_d;
if(flag == 0)   
    sys = [2,0,2,1,0,0];
    x0 = [0,0];
elseif(abs(flag) == 1)
    sys = zeros(1,2);
    sys(1) = x(2);
    sys(2) = (-(B_in+B_L*(2*A_1*cos(2*x(1))+4*A_2*cos(4*x(1))+...
        1/N)^2)*x(2)+(J_L*4*(A_1*sin(2*x(1)+4*A_2*sin(4*x(1))))*...
        (2*A_1*cos(2*x(1))+4*A_2*cos(4*x(1))+1/N)*(x(2)^2))-...
        K_D*(2*A_1*cos(2*x(1))+4*A_2*cos(4*x(1))+(x(2)/N))-K_P*...
        (A_1*sin(2*x(1))+A_2*sin(4*x(1))+(x(1)/N))+...
        K_P*u)/(J_in+J_L*(2*A_1*cos(2*x(1))+4*A_2*cos(4*x(1))+1/N));
elseif(abs(flag) == 3)
    sys = zeros(1,2);
    sys(1) = x(1);
    sys(2) = x(2);
else
    sys = [];
end;

youbindi

给自己顶一下。。。支持一下自己:lol :lol

sogooda

1. 貌似新版本里不用flag了。
2. u是附加参数吗？看rk45的调用好像没有传递u啊，如果不需要就去掉u
3. [t,x] = rk45('sfile45', tspan , x0); 可以改成[t,x] = rk45(@sfile45, tspan , x0);

youbindi

u是控制参数，flag是判断微分还是积分还是雅可比计算。。。但是就不知道怎么改成符合ode45调用的函数！！！

sogooda

1. 
   
2. [t,x] = rk45(@sfile45, tspan , x0,,u);
   
3. ...
   
4. function x0= sfile45(t,x,u)
   

复制代码