matlab解微分方程组（带参数）

Question

function dy=eq2(t,y)
dy(1)=-(a*y(2))/(4*exp(a*t/4));
dy(2)=-(a/4)*(exp(a*t/4))*(y(1)+0.5)+(a/4)*y(2)-y(3)*((exp(a*t/4))^2);
dy(3)=4*y(2);
方程是上面的,a是可以给定的参数（结束的时间），y的初值是[1.4 0.1 0.1]，画出a取值在0~10所对应的y(1)结束的值，t的范围是0~a.

Answer 1

可以借助于嵌套函数或匿名函数实现附加参数的传递，例如

function main
y0 = [1.4; 0.1; 0.1];
A = linspace(eps, 10, 20);
Y = A * NaN;
for ii = length(A)
a = A(ii);
y = ode45(@eq2, [0 a], y0);
Y(ii) = y(end, 1);
end
plot(A, Y)

function dy=eq2(t,y)
dy = y*0;
dy(1)=-(a*y(2))/(4*exp(a*t/4));
dy(2)=-(a/4)*(exp(a*t/4))*(y(1)+0.5)+(a/4)*y(2)-y(3)*((exp(a*t/4))^2);
dy(3)=4*y(2);
end
end

但微分方程组似乎是刚性的，不过换用ode15s、ode23s等适合刚性系统的算法效果也不理想（可以调用ode*函数时不返回参数，观察求解的过程）。追问

那画图呢，怎么画出来

追答

上面的程序带画图语句啊。
只不过方程组是刚性的，求解时间长，我没等结果出来，你可以试试看。