时滞微分方程的Matlab解法之一dde23

延迟微分方程

matlab提供了dde23求解非中性微分方程。dde23的调用格式如下：
sol = dde23(ddefun,lags,history,tspan)
lags是延迟量，比如方程中包含y1(t-0.2)和y2(t-0.3)则可以使用lags=[0.2,0.3]。
这里的ddefun必须采用如下的定义方式：
dydt = ddefun(t,y,Z)
其中的Z(:,1)就是y(t-lags(1)),Z(:,2)就是y(t-lags(2)）...
下面是使用dde23求解延迟微分方程的两个例子。

第一个例子：

代码如下：

ddex1dez = @(t,y,Z) [y(1)*(1 + 0.1*sin(t)-0.1*Z(1,1) - y(2)/(1+y(1)) );
    y(2)*( (2+sin(t))*10^(-5) + 9*Z(1,2)/(1+Z(1,2)) - Z(2,1) )];

 %y(1)表示x_1(t)，因为dde求解的结果中sol会有个x，为了区别用y(1)表示x_1(t)；Z(1,1)表示时滞项x_1(t-0.1)；Z(1,2)表示时滞项x_1(t-0.3)

sol = dde23(ddex1dez,[0.1, 0.3],[2 2],[0, 50]);%dde23(@....,tau,history,tspan);

 %[0.1, 0.3]是时滞，[2 2]是初值，[0, 50]是时间范围

figure;
% plot(sol.x,sol.y)
plot(sol.x,sol.y(1,:) )
hold on
plot(sol.x,sol.y(2,:),'-.' )
hold off

title('时滞微分方程组');
xlabel('time t');
ylabel('solution y');
legend('x1','x2');

第二个例子：

This example shows how to use dde23 to solve a system of DDEs with constant delays.

The differential equations are:

are solved on [0,5] with history:

for t ≤ 0.

1. Create a new program file in the editor. This file will contain a main function and two local functions.

2. Define the first-order DDE as a local function.

   function dydt = ddex1de(t,y,Z) ylag1 = Z(:,1); ylag2 = Z(:,2); dydt = [ylag1(1); ylag1(1)+ylag2(2); y(2)]; end

3. Define the solution history as a local function.

   function S = ddex1hist(t) S = ones(3,1); end

4. Define the delays, τ₁,…,τ_k in the main function.

   lags = [1,0.2];

5. Solve the DDE by calling dde23 in the main function. Pass the DDE function, the delays, the solution history, and interval of integration, [0,5], as inputs.

   sol = dde23(@ddex1de,lags,@ddex1hist,[0,5]);

   The dde23 function produces a continuous solution over the whole interval of integration [t₀,t_f].

6. Plot the solution returned by dde23. Add this code to your main function.

   plot(sol.x,sol.y); title('An example of Wille and Baker'); xlabel('time t'); ylabel('solution y'); legend('y_1','y_2','y_3','Location','NorthWest');

7. Evaluate the solution at 10 equally spaced points over the interval of integration. Then plot the the results on the same axes as sol.y. Add this code to the main function.

   tint = linspace(0,5,10); Sint = deval_r(sol,tint) hold on plot(tint,Sint,'o');

8. Run your program to generate and plot the results.

代码如下：

ddex1dez = @(t,y,Z) [Z(1,1);Z(1,1)+Z(2,2);y(2)];

lags = [1,0.2];

sol = dde23(ddex1dez,lags,[1 1 1],[0,5]);

 

plot(sol.x,sol.y);
title('An example of Wille and Baker');
xlabel('time t'); ylabel('solution y');
legend('y_1','y_2','y_3','Location','NorthWest');

 

tint = linspace(0,5,10);
Sint = deval(sol,tint)

hold on
plot(tint,Sint,'o');

或者 按如下代码执行：

clear;clc
lags=[1,0.2];
history=[1;1;1];
tspan=[0,5];
sol = dde23(@myddefun,lags,history,tspan)

plot(sol.x,sol.y) 

function dy = myddefun(t,y,Z)
dy=[
    Z(1,1);
    Z(1,1)+Z(2,2);
    y(2) ];

引用：

http://blog.sina.com.cn/s/blog_3ecbcc0701013bcd.html

http://wenku.baidu.com/view/31c21f54cc1755270722088b.html