创新实验之脉冲微分方程的数值解实验报告

1、编号： 数学与计算科学学院创新性实验结题报告书实验题目： 脉冲微分方程的数值解 实验成绩（教师填写）： 学 院： 数学与计算科学学院 专 业： 数学与应用数学 学生姓名： 小组成员： 指导教师： 2014年 10月 17日 摘 要：脉冲微分方程在生态学、环境科学、电力工程及自动控制系统的实时仿真等科学与工程应用领域中有广泛的应用。但常见的脉冲微分方程中,只有极少数能够获得理论解的解析表达式,而许多实际问题不需要求解脉冲微分方程的解析解,只需要数值解,因此利用数值解法求解实际问题具有十分重要的理论意义及实用价值。关键字：脉冲微分方程、数值解一、实验目的及主要创新性点：本次实验我们要研究的是脉冲时

2、滞微分方程。脉冲时滞微分方程在生态学、环境科学、电力工程及自动控制系统的实时仿真等科学与工程应用领域中有广泛的应用。但常见的脉冲时滞微分方程中，只有极少数能够获得理论解的解析表达式，而许多实际问题不需要求解脉冲时滞微分方程的解析解，只需要数值解，因此利用数值解法求解实际问题具有十分重要的理论意义及实用价值。2、 实验要求：1、了解脉冲微分方程的原理及其应用；2、使用数学软件来求解脉冲时滞微分方程的数值解。三、实验原理：基本概念：(1)考虑微分系统， ，其中，为开集，为n维欧式空间,;2）对任意存在两个集合；3）算子任意。设是系统,满足的解。它具有如下特点：动点 起始于并沿曲线运动，在处，与集合

3、相遇，这时算子将点作用到，其中。而点从出发继续沿解曲线运动，直到时刻再次遇到集合,点被作用到,其中。同前面一样，沿着系统的解曲线继续 运动，如果系统的解存在，则一直运动下去。将具有上述运动过程的综合起来称为脉冲微分方程系统，有所描述的曲线称为积分曲线，积分曲线所代表的函数称为脉冲微分系统的解，与集合相遇的时刻称为脉冲时刻。假设脉冲微分系统的解在脉冲时刻是左连续的，即上述的脉冲微分方程的解有以下几个特点：（a) 微分方程的解是连续函数，若积分弧度与没有交集或者交其不动点；（b) 微分方程的解是分段连续函数并且有有限个第一类间断点，积分弧和相交同时不是的不动点；（c) 微分方程的解是分段连续函数并

4、且第一类间断点的个数是可数的，积分弧与相交且这些点不是的不动点。1.脉冲时滞微分方程初值问题的数值解法 考虑如下脉冲时滞微分方程初值问题: ， （1） 求初值问题(1)的数值解,就是寻求准确解在一系列离散节点上的近似值其中:脉冲时刻,满足,时,；是脉冲算子;,是正整数的集合；是时滞；是初始值函数.显然,只有当初值问题(1)的解存在且唯一时,使用数值解法才有意义.当是连续且满足李普希茨条件时,李普希茨常数为,即对任意恒有时初值问题的解存在并且唯一.2.基于Runge-Kutta法的脉冲时滞微分方程初值问题的数值解龙格库塔法(Runge-Kutta)是一种微分方程近似求解的数值方法.龙格库塔法通过

5、计算不同点上的函数值,并对这些函数值作线性组合,构造近似公式,再将近似公式与精确解的泰勒展开式进行比较,使前面的若干项相同,从而使近似公式的局部截断误差阶数尽可能地高.用个函数的值表示的Runge-Kutta方法称为级Runge-Kutta方法.截断误差为的称为阶Runge-Kutta方法。四级四阶龙格库塔法即其截断误差为,每步都要计算4次的值,其标准公式即古典Runge-Kutta法为: (3)其中表示计算过程中选取的步长,通常采用等距节点.为了得到系统的精确一般步长选的尽可能小.若,则将区间由右向左以步长h进行分割.设 则代表网域,其中元素称为节点。(3)式中的表示点处的斜率；表示利用求得

6、的点处的斜率；表示利用求得的点处的斜率；表示利用求得的点处的斜率。（3）式是常微分方程的数值解的龙格库塔法.而所谓时滞,就是一种时间差,是时间延迟,是控制系统普遍存在的一种现象,与始点和终点都有关,即某一时刻系统(1)的数值解与系统(1)在时刻及时刻的数值解都有关。将(3)式应用于方程(1),得 (4) 即当计算的值时,需要用到的值,由于的值也在网格点上,因此将其代替时滞项的值进行迭代运算。所谓脉冲,就是在短时间内突变,随后又迅速返回其初始值的物理量,它可以是周期性重复的,也可以是非周期性的或单次的.即在某一个脉冲时刻,将系统的解修改为所给定的脉冲值,而在该脉冲时刻和下一个脉冲时刻之间系统被重

7、新积分求脉冲时滞微分方程初值问题的算法具体描述如下:()设定,的初始值,并设定存储空间为，()设,应用(3)式分步计算区间的历史函数值；()设脉冲时刻为，,，···，···，则在每一个脉冲时刻,将该时刻的函数值修改为脉冲值,即将该时刻对应的存储单元的值修改为，再应用(4)式分步计算区间的系统(1)的数值解。求脉冲微分方程数值解的相应算法描述如下:()设定,的初始值,并设定存储空间为，()设,脉冲时刻为，,，···，···，则在每一个脉冲时刻,将该时刻的函数值修改为脉冲值,即将

8、该时刻对应的存储单元的值修改为，再应用(3)式分步计算区间的脉冲微分方程的数值解。求时滞微分方程数值解的相应算法描述如下:()设定,的初始值,并设定存储空间为，()设,应用(3)式分步计算区间的历史函数值；()在区域应用(4)式分步计算时滞微分方程的数值解。4、 实验内容与步骤：4x”+35.78x+2000x=F(t)F(t)=100 0<t<=0.1F(t)=-1000(t-0.2) 0.1<t<=0.2要求用ord45解出x'和x(数值解）t=0时x=0,x'=0function hhx0=10 5;T1,X1 = ode45(odefun1,0:

9、0.01:0.1,x0)T2,X2 = ode45(odefun2,0.1:0.01:0.2,X1(end,1),X1(end,2)plot(T1,X1,T2,X2)function dx=odefun1(t,x)x1=x(1);x2=x(2);dx1=x1;dx2=(1000-35.78*x1-2000*x2)/4; dx=dx1;dx2;function dx=odefun2(t,x)x1=x(1);x2=x(2);dx1=x1;dx2=(1000*(t-0.2)-35.78*x1-2000*x2)/4; dx=dx1;dx2;五、实验过程原始记录（数据，图表，计算等）：新建一个函数脚本文

10、件t0=0; tN=15; y0=5;0; h=0.01;t=t0:h:tN; N=length(t); j=1;for i=1:Nt1=t0+h;K1=R_K(t0,y0);K2=R_K(t0+h/2,y0+h*K1/2); K3=R_K(t0+h/2,y0+h*K2/2);K4=R_K(t0+h,y0+h*K3);y1=y0+(h/6)*(K1+2*K2+2*K3+K4);Y1(j)=y1(1); Y2(j)=y1(2);t0=t1; y0=y1; j=j+1;end plot (t,Y1,t,Y2); grid;新建一个函数M-K.m文件Function ydot=f(t,y)ydot=zeros(size(y);ydot(1)=y(2);_ ydot(2)=-y(2)-y(1); end用脚本文件调用函数文件求解，初值，下图为求解图形，针对不同的微分方程适当修改函数M文件即可例1 求如下脉冲微分方程初值的数值解： , , , , , 应用上述方法计算所得的数值解仿真结果如图1 图1例2 求如下脉冲时滞微分方程初值的数值解： , , , , , 应用上述方法计算所得的数值解仿真结果如图2 图2六，实验结果分析或总结：寻找所求脉冲微分方程问题与已学的微分方程问题的联系，利用已知的龙哥库塔法，加