如何用Matlab数值方法绘制一阶非线性微分方程的解x(t)

1. 定义微分方程匿名函数 #

Matlab 常微分方程数值求解器要求输入的微分方程函数接收t、x两个输入参数（即使表达式中不含t也需要保留该参数位），将你给出的 dx/dt 表达式改写为支持元素级运算的匿名函数即可：

% 注意乘方、除法前加. 做元素级运算，适配求解器的数组输入逻辑
dxdt = @(t,x) -(1.0*(6.84e+45*x.^2 + 5.24e+32*x - 2.49e+42))./(2.47e+39*x + 7.12e+37);

2. 配置求解参数 #

• 初始条件：对应你给出的 t=0 时 x=0，设置为 x0 = 0;
• 求解时间区间：根据你需要观察的解的变化范围调整，示例设置为从0到1e-5（该方程系数量级很大，解的变化速度会非常快，可以根据实际输出调整区间大小）：tspan = [0 1e-5];
• 可选精度配置：如果需要更高的求解精度，可以通过odeset自定义求解器参数：

options = odeset('RelTol',1e-8,'AbsTol',1e-10);

3. 调用求解器计算数值解 #

普通场景直接用自适应步长求解器ode45即可，如果求解时报错提示方程刚性强，可以换成ode15s求解：

% 不自定义精度用这个
[t,x] = ode45(dxdt, tspan, x0);
% 要自定义精度就用下面这句替换上面的
% [t,x] = ode45(dxdt, tspan, x0, options);

4. 绘制x(t)曲线 #

plot(t,x,'LineWidth',1.5)
xlabel('t')
ylabel('x(t)')
grid on
title('一阶非线性微分方程数值解曲线')

注意事项 #

求解时如果出现发散、奇点报错，大概率是你的求解区间内分母2.47e+39*x + 7.12e+37趋近于0，可尝试缩小求解时间区间，或者先验证x的合理取值范围，避免进入分母为0的无意义区域。

内容的提问来源于stack exchange，提问作者Andrew