时滞微分方程的matlab解法

有位小伙伴在matlab编程爱好者群中问道有关时滞微分方程的matlab解法，问题是选自由清华大学出版社出版、薛定宇著的《高等应用数学问题的MATLAB求解 (第四版)》的课后习题，问题的如下：

显然这是时滞固定的时滞微分方程，采用dde23函数即可对其进行求解，在给出解法之前先来看看dde23到底该怎么使用，其常用的调用格式如下：

sol = dde23(ddefun,lags,history,tspan);

以上语句中：

ddefun —— 待求解的时滞微分方程，要求：函数必须为以下形式：

dydt = ddefun(t,y,Z);

lags —— 固定正时滞向量 τ1, ..., τk;

history —— 按以下三种方式之一指定 history：

• ① 一个 t 函数，要求 y = history(t) 能够将 t ≤ t0 的解 y(t) 以列向量的形式返回
• ② 一个固定列向量（如果 y(t) 为常量）
• ③ 来自之前积分的解 sol（如果此调用继续该积分）
• tspan —— 从 t0=tspan(1) 到 tf=tspan(end) 的时间积分区间，其中 t0 < tf。
• sol —— 返回的计算结果，为结构体类型数据，通过sol.x来读取时间，sol.y来读取计算结果。
• 了解了dde23函数的使用方法，接下来就是如何实现的问题了。

在ddefun中，t对应时间t，y对应导数x'(t)、y'(t)，Z对应各个历史x(t-0.1)、x(t-0.2)、x(t-0.5)、y(t-0.1)、y(t-0.2)、y(t-0.5)，为了求解所有的量都需要定义成列向量的形式，即有：

y = [x'(t);y'(t)],

Z = [x(t-0.1) x(t-0.2) x(t-0.5); y(t-0.1) y(t-0.2) y(t-0.5)],

=[Z(1,1) Z(1,2) Z(1,3); Z(2,1) Z(2,2) Z(2,3)],

所以有：dydt = ddefun(t,y,Z).

lags = [0.1 0.2 0.5].

history = [t; exp(t)].

tspan = [0 0.8].

这里，tspan是自行定义的。

到这儿，其实代码已经完成了，照葫芦画瓢填进去即可，下面就一起来看看本题的源代码吧：

% 作者：巴山
% 欢迎关注matlab爱好者公众号
% Q群：531421022
% 知识星球：资源分享园地
% 抖音号：matlab. (有个点哦)
% B站号：matlabaihaozhe

clc;clear;close all
% 定义微分方程 —— ddefun
ddefun = @(t,x,Z)[Z(1,2)^2+Z(2,2)^2-6*Z(1,3)-8*Z(2,1);...
x(1)*(2*Z(2,2) - x(1)+5 - 2*Z(1,1)^2)];
% 定义时滞向量 —— lags
lags = [0.1 0.2 0.5];
% 定义历史 —— history
history =@(t)[t;exp(t)];
% 定义时间积分区间 —— tspan
tspan = [0 0.8];
sol=dde23(ddefun,lags,history,tspan);
figure('Color','w');
% plot(sol.x,sol.y)
plot(sol.x,sol.y(1,:),'r','LineWidth',2.0);
hold on
plot(sol.x,sol.y(2,:),'b--','LineWidth',2.0)
hold off
title('时滞微分方程组');
xlabel('时间 t');
ylabel('结果 y');
legend('x(t)','y(t)');

本文只起抛砖引玉之作用，今后若有小伙伴在matlab编程爱好者Q群中提有意思的matlab编程问题，同样会与大家分享，指不定哪天就能派上用场。

matlab中求解时滞微分方程相关的函数了dde23之外，还有ddensd 、ddesd、ddeset、ddeget、deval等着大家学习，这里就不赘述了。

更多关于微分方程、含积分方程的求解与拟合问题，敬待后续推送！

参考资料：

[1] 薛定宇. 《高等应用数学问题的MATLAB求解 (第四版)》. 清华大学出版社.

[2] https://ww2.mathworks.cn/help/matlab/ref/dde23.html