eulerode.m

通信程序源代码和模型通信程序源代码和模型通信程序源代码和模型

function [t_out, x_out]=eulerode(odefuncname,tspan, x0, options, par)
% 欧拉法求状态方程odefuncname的数值解
% odefuncname为描述状态方程函数名称（也即文件名）
% tspan为求解的时间序列（要求是列向量）
% x0为系统初始状态列向量
% options 此处不使用
% par 是传入odefuncname的参数，与odefuncname 函数的参数par设计相匹配
% 返回值：t_out 计算的时间点序列, 等于输入的tspan
% x_out是对应于t_out的状态变量解矩阵
t_out=tspan;
x_out=zeros(length(t_out),length(x0));% 初始化输出解矩阵, 以加快程序计算速度
                                % x_out的第k行状态值对应着t_out中的第k个时刻
for k=1:length(t_out)
    x_out(k,:)=x0';             % 将第k步计算结果存入输出解矩阵(注意要转置)
    if k<length(t_out) 
        h=t_out(k+1)-t_out(k);  % 到下一时刻的计算步长
    else
        h=0;                    % 最后一个步长
    end
    x0=x0+h*eval([odefuncname, '(t_out(k),x0,[], par)']);% 欧拉算法
    % 用eval调用状态方程函数, 其参数flag不使用, 置为空矩阵或任意数值
end