第九讲 差分方程 Matlab语言程序设计 教学课件
差分方程及matlab求解
称为n阶差分方程的通解.在通解中给任意常数 C1,C2,…,Cn 以确定的值所得的解,称为n阶差分方程的 特解.
例如,函数yt=at+C(a为已知常数,C为任意常数)是差
分方程yt+1-yt=a的通解.而函数yt=at,yt=at-1,…均是这个 差分方程的特解. 由差分方程的通解来确定它的特解,需要给出确 定特解的定解条件. n阶差分方程F(t,yt,yt+1,…，yt+n)=0常见的定解条件为 初始条件.
方程的通解
1 t 1 1 t -1 2 t 1 t 1 t +1 yt = A( ) + ( ) ( 2 - 1) = A ( ) + 2 2 3 2 2 3 2 A = A - 为任意常数 . 3
2.待定系数法求特解
情形Ⅰ f(t)为常数． 方程变为yt+1+ayt=b, a,b均为非零常数． 试以 yt = (为待定常数)形式的特解代入方程得 +a =(1+a) =b． 当a≠-1时,可求得特解
三、 差分方程的解
定义4 如果将已知函数yt=j(t)代入方程F(t,yt,yt+1,…， yt+n)=0, 使 其 对 t=…,-2,-1,0,1,2,… 成 为 恒 等 式 , 则 称 yt=j(t)为方程的解.含有n个任意(独立)常数C1,C2,…,Cn 的解
yt=j(t,C1,C2,…,Cn)
解,那么,非齐次线性差分方程的通解为： y(t)=yA(t)+ y (t)
即
y(t)=A1y1(t)+A2y2(t)+…+Anyn(t)+ y(t)， 这里A1，A2,…,An为n个任意(独立)常数．
数学建模差分方程PPT课件
的方程都是差分方程。 方程中所含未知函数角标的最大值与最小值的差数称为差分
方程的阶。 若一个函数代入差分方程后，方程两端恒等，则称此函数为
差分方程的解。如果解中所含相互独立的任意常数的个数等于方 程的阶数，则称该解为差分方程的通解。满足初始条件的解称为 特解。
• 第一阶段: w(k)每周减1千克, c(k)减至下限10000千卡
w (k)w (k1)1 w ( k 1 ) w ( k ) c ( k 1 ) w ( k )
c(k1) 1[w(k)1] w (k)w (0)k
c(k1) w (0) 1(1k)
1 8000
0.025
120 200 k 00Cm 10000 k 10
2 x k 2 x k 1 x k 2 ( 1 ) x 0 , k 1 , 2 ,
二阶线性常系数差分方程
x0为平衡点 研究平衡点稳定，即k, xkx0的条件
模型的推广 2 x k 2 x k 1 x k 2 ( 1 ) x 0
方程通解
xk
c1
k 1
c2
k 2
(c1, c2由初始条件确定)
相当于70千克的人每天消耗2000千卡 ~ 3200千卡；
3）运动引起的体重减少正比于体重，且与运动 形式有关；
4）为了安全与健康，每周体重减少不宜超过1.5 千克，每周吸收热量不要小于10000千卡。
减肥计划
某甲体重100千克，目前每周吸收20000千卡热量， 体重维持不变。现欲减肥至75千克。
1 2 k 是（3）的 k 重根，则只要将 Y1 (i),Y2 (i),,Yk (i) 换为
差分方程的解法分析及MATLAB实现(程序)
差分方程的解法分析及MATLAB 实现（程序）摘自：张登奇,彭仕玉.差分方程的解法分析及其MATLAB 实现[J]. 湖南理工学院学报.2014(03) 引言线性常系数差分方程是描述线性时不变离散时间系统的数学模型,求解差分方程是分析离散时间系统的重要内容.在《信号与系统》课程中介绍的求解方法主要有迭代法、时域经典法、双零法和变换域法[1].1 迭代法例1 已知离散系统的差分方程为)1(31)()2(81)1(43)(-+=-+--n x n x n y n y n y ,激励信号为)()43()(n u n x n =,初始状态为21)2(4)1(=-=-y y ，.求系统响应. 根据激励信号和初始状态,手工依次迭代可算出2459)1(,25)0(==y y . 利用MATLAB 中的filter 函数实现迭代过程的m 程序如下:clc;clear;format compact;a=[1,-3/4,1/8],b=[1,1/3,0], %输入差分方程系数向量,不足补0对齐n=0:10;xn=(3/4).^n, %输入激励信号zx=[0,0],zy=[4,12], %输入初始状态zi=filtic(b,a,zy,zx),%计算等效初始条件[yn,zf]=filter(b,a,xn,zi),%迭代计算输出和后段等效初始条件2 时域经典法用时域经典法求解差分方程:先求齐次解；再将激励信号代入方程右端化简得自由项,根据自由项形式求特解；然后根据边界条件求完全解[3].用时域经典法求解例1的基本步骤如下.（1）求齐次解.特征方程为081432=+-αα,可算出41 , 2121==αα.高阶特征根可用MATLAB 的roots 函数计算.齐次解为. 0 , )41()21()(21≥+=n C C n y n n h （2）求方程的特解.将)()43()(n u n x n =代入差分方程右端得自由项为 ⎪⎩⎪⎨⎧≥⋅==-⋅+-1,)43(9130 ,1)1()43(31)()43(1n n n u n u n n n 当1≥n 时,特解可设为n p D n y )43()(=,代入差分方程求得213=D . （3）利用边界条件求完全解.当n =0时迭代求出25)0(=y ,当n ≥1时,完全解的形式为 ,)43(213 )41()21()(21n n n C C n y ⋅++=选择求完全解系数的边界条件可参考文[4]选)1(),0(-y y .根据边界条件求得35,31721=-=C C .注意完全解的表达式只适于特解成立的n 取值范围,其他点要用)(n δ及其延迟表示,如果其值符合表达式则可合并处理.差分方程的完全解为)(])43(213 )41(35)21(317[)1(])43(213 )41(35)21(317[)(25)(n u n u n n y n n n n n n ⋅+⋅+⋅-=-⋅+⋅+⋅-+=δ MATLAB 没有专用的差分方程求解函数,但可调用maple 符号运算工具箱中的rsolve 函数实现[5],格式为y=maple('rsolve({equs, inis},y(n))'),其中:equs 为差分方程表达式, inis 为边界条件,y(n)为差分方程中的输出函数式.rsolve 的其他格式可通过mhelp rsolve 命令了解.在MATLAB 中用时域经典法求解例1中的全响应和单位样值响应的程序如下.clc;clear;format compact;yn=maple('rsolve({y(n)-3/4*y(n-1)+1/8*y(n-2)=(3/4)^n+1/3*(3/4)^(n-1),y(0)=5/2,y(-1)=4},y(n))'),hn=maple('rsolve({y(n)-3/4*y(n-1)+1/8*y(n-2)=0,y(0)=1,y(1)=13/12},y(n))'),3 双零法根据双零响应的定义,按时域经典法的求解步骤可分别求出零输入响应和零状态响应.理解了双零法的求解原理和步骤,实际计算可调用rsolve 函数实现.yzi=maple('rsolve({y(n)-3/4*y(n-1)+1/8*y(n-2)=0,y(-1)=4, y(-2)=12},y(n))'),yzs=maple('rsolve({y(n)-3/4*y(n-1)+1/8*y(n-2)=(3/4)^n+1/3*(3/4)^(n-1),y(0)=1,y(-1)=0},y(n))'),4 变换域法设差分方程的一般形式为)()(00r n x b k n y a r Mr k N k -=-∑∑==.对差分方程两边取单边z 变换,并利用z 变换的位移公式得])()([])()([1010m r m r r M r l k l k k N k z m x z X z b z l y z Y z a ---=-=---=-=∑∑∑∑+=+整理成)()()()()()(00z X z X z B z Y z Y z A +=+形式有. )(, )(110110M M N N z b z b b z B z a z a a z A ----+++=+++=. )()(, )()(110110∑∑∑∑=--=--=--=--==M r r m m r r N k k l l k k z m x b s X zl y a s Y可以看出,由差分方程可直接写出 )(z A 和 )(z B ,系统函数)(/)()(z A z B z H =,将系统函数进行逆z 变换可得单位样值响应.由差分方程的初始状态可算出 )(0z Y ,由激励信号的初始状态可算出 )(0z X ,将激励信号进行z 变换可得 )(z X ,求解z 域代数方程可得输出信号的象函数 , )()()()()()(00z A z Y z X z X z B z Y -+= 对输出象函数进行逆z 变换可得输出信号的原函数)(n y .利用z 变换求解差分方程各响应的步骤可归纳如下:（1）根据差分方程直接写出 )(z A 、 )(z B 和)(z H ,)(z H 的逆变换即为单位样值响应；（2）根据激励信号算出 )(z X ,如激励不是因果序列则还要算出前M 个初始状态值；（3）根据差分方程的初始状态 )(, ),2( ),1(N y y y -⋅⋅⋅--和激励信号的初始状态 )(, ),2( ),1(M x x x -⋅⋅⋅--算出 )(0z Y 和 )(0z X ；（4）在z 域求解代数方程)()()()()()(00z X z X z B z Y z Y z A +=+得输出象函数 )(z Y , )(z Y 的逆变换即为全响应；（5）分析响应象函数的极点来源及在z 平面中的位置,确定自由响应与强迫响应,或瞬态响应与稳态响应；（6）根据零输入响应和零状态响应的定义,在z 域求解双零响应的象函数,对双零响应的象函数进行逆z 变换,得零输入响应和零状态响应.用变换域法求解例1的基本过程如下. 根据差分方程直接写出2181431 )(--+-=z z z A ,1311 )(-+=z z B .系统函数的极点为41,21. 对激励信号进行z 变换得)43/( )(-=z z z X .激励象函数的极点为3/4. 根据差分方程的初始状态算出102123 )(-+-=z z Y .根据激励信号的初始状态算出 0)(0=z X . 对z 域代数方程求解,得全响应的象函数)323161123/()83243125( )(2323-+-+-=z z z z z z z Y . 进行逆z 变换得全响应为)(])43(213 )41(35)21(317[)(n u n y n n n ⋅+⋅+⋅-= 其中,与系统函数的极点对应的是自由响应;与激励象函数的极点对应的是强迫响应. )(z Y 的极点都在z 平面的单位圆内故都是瞬态响应.零输入响应和零状态响应可按定义参照求解.上述求解过程可借助MATLAB 的符号运算编程实现.实现变换域法求解差分方程的m 程序如下: clc;clear;format compact;syms z n %定义符号对象% 输入差分方程、初始状态和激励信号%a=[1,-3/4,1/8],b=[1,1/3], %输入差分方程系数向量y0=[4,12],x0=[0], %输入初始状态,长度分别比a 、b 短1,长度为0时用[]xn=(3/4)^n, %输入激励信号,自动单边处理,u(n)可用1^n 表示% 下面是变换域法求解差分方程的通用程序，极点为有理数时有解析式输出 %N=length(a)-1;M=length(b)-1;%计算长度Az=poly2sym(a,'z')/z^N;Bz=poly2sym(b,'z')/z^M;%计算A(z)和B(z)Hz=Bz/Az;disp('系统函数H(z):'),sys=filt(b,a),%计算并显示系统函数hn=iztrans(Hz);disp('单位样值响应h(n)='),pretty(hn),%计算并显示单位样值响应Hzp=roots(a);disp('系统极点:');Hzp,%计算并显示系统极点Xz=ztrans(xn);disp('激励象函数X(z)='),pretty(Xz),%激励信号的单边z 变换Y0z=0;%初始化Y0(z),求Y0(z)注意系数标号与变量下标的关系for k=1:N;for l=-k:-1;Y0z = Y0z+a(k+1)*y0(-l)*z^(-k-l);endenddisp('初始Y0(z)'),Y0z,%系统初始状态的z 变换X0z=0;%初始化X0(z),求X0(z)注意系数标号与变量下标的关系for r=1:M;for m=-r:-1;X0z = X0z+b(r+1)*x0(-m)*z^(-r-m);endenddisp('初始X0(z)'),X0z,%激励信号起始状态的z 变换Yz=(Bz*Xz+X0z-Y0z)/Az;disp('全响应的z 变换Y(z)'),pretty(simple(Yz)),yn=iztrans(Yz);disp('全响应y(n)='),pretty(yn),% 计算并显示全响应Yziz=-Y0z/Az;disp('零输入象函数Yzi(z)='),pretty(Yziz),%零激励响应的z 变换yzin=iztrans(Yziz);disp('零输入响应yzi(n)='),pretty(yzin),% 计算并显示零输入响应 Yzsz=(Bz*Xz+X0z)/Az;disp('零状态象函数Yzs(z)='),pretty(Yzsz),%零状态响应的z 变换yzsn=iztrans(Yzsz);disp('零状态响应yzs(n)='),pretty(yzsn),% 计算并显示零状态响应该程序的运行过程与手算过程对应,显示在命令窗的运行结果与手算结果相同.。
Matlab解差分方程
一阶线性常系数差分方程的解、平衡点及其稳定性
x a x b k 1 k
• 自然环境下,b=0 • 人工孵化条件下
xk a k x0
k k 1 xa xba ( 1 a ) k 0
1 ak a x0 b 1 a
k
• 令xk=xk时，xk→x,称平衡点是稳定的
• • • • • • •
A=[0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6]; >> n=10; >> for k=1:n x(:,1)=[200,200,200]'; x(:,k+1)=A*x(:,k); end >> round(x)
作图观察数量变化趋势
300
• • • • •
用Matlab求解差分方程问题
一阶线性常系数差分方程
高阶线性常系数差分方程
线性常系数差分方程组
差分方程是在离散时段上描述现 实世界中变化过程的数学模型
• 例1、 某种货币1年期存款的年利率是r ， 现存入M元，问年后的本金与利息之和 是多少？
• Xk+1=(1+r)xk , k = 0 , 1 , 2 · · · · ·
利用plot 绘图观察数量变化趋势
• 可以用不同线型和颜色绘图 • r g b c m y k w 分别表示 红绿兰兰绿洋红黄黑白色 ： + o * . X s d 表示不同的线型
• plot(k,y1,k,y2,k,y3) 在同一坐标系下画图 plot(k,y2,':')
>> plot(k,y2,'--')
0.1 • A 0.6 0.7 0.3 B C
差分方程的解法分析及MATLAB实现
差分方程的解法分析及MATLAB实现差分方程是描述离散时序系统行为的数学工具。
在离散时间点上，系统的行为由差分方程给出，这是一个递归方程，其中当前时间点的状态取决于之前的状态和其他外部因素。
解差分方程的方法可以分为两类：直接解法和转化为代数方程的解法。
直接解法通过求解差分方程的递归形式来得到解析或数值解。
转化为代数方程的解法则将差分方程转化为代数方程进行求解。
一、直接解法的步骤如下：1.将差分方程表示为递归形式，即将当前时间点的状态表示为之前时间点的状态和其他外部因素的函数。
2.根据初始条件，确定初始时间点的状态。
3.根据递归形式，计算出后续时间点的状态。
以下是一个简单的差分方程的例子：y(n)=2y(n-1)+1，其中n为时间点。
按照上述步骤求解该差分方程：1.将差分方程表示为递归形式：y(n)=2y(n-1)+12.根据初始条件，假设y(0)=1，确定初始时间点的状态。
3.根据递归形式，计算出后续时间点的状态：y(1)=2y(0)+1=2*1+1=3y(2)=2y(1)+1=2*3+1=7y(3)=2y(2)+1=2*7+1=15...依此类推计算出所有时间点的状态。
二、转化为代数方程的解法的步骤如下：1.假设差分方程的解具有指数形式，即y=r^n，其中r为待定参数。
2.将差分方程代入上述假设中，得到r的方程。
3.解得r的值后，再根据初始条件求解出常数值。
4.得到差分方程的解析解。
以下是一个复杂一些的差分方程的例子：y(n)=2y(n-1)+3y(n-2)，其中y(0)=1，y(1)=2按照上述步骤求解该差分方程：1.假设差分方程的解具有指数形式：y=r^n。
2.代入差分方程得到：r^n=2r^(n-1)+3r^(n-2)。
3.整理得到：r^2-2r-3=0。
4.解得r的值为：r1=-1，r2=35.根据初始条件求解出常数值：y(0)=c1+c2=1，y(1)=c1-c2=2、解得c1=1.5，c2=-0.56.得到差分方程的解析解：y(n)=1.5*(-1)^n+-0.5*3^n。
Matlab解差分方程
一阶线性常系数差分方程
高阶线性常系数差分方程
线性常系数差分方程组
差分方程是在离散时段上描述现 实世界中变化过程的数学模型
• 例1、 某种货币1年期存款的年利率是r ， 现存入M元，问年后的本金与利息之和 是多少？
• Xk+1=(1+r)xk , k = 0 , 1 , 2 · · · · ·
x ( k 1 ) 0 . 6 x ( k ) 0 . 2 x ( k ) 0 . 1 x ( k ) 1 1 2 3 x ( k 1 ) 0 . 3 x ( k ) 0 . 7 x ( k ) 0 . 3 x ( k ) 2 1 2 3 x ( k 1 ) 0 . 1 x ( k ) 0 . 1 x ( k ) 0 . 6 x ( k ) 3 1 2 3
k=0:10; plot(k,x) ，grid gtext('x1(k)'), gtext('x2(k)'), gtext('x3(k)')
280 260 240 220 200 180 160 140 120
x2(k)
x1(k)
x3(k) 0 1 2 3 4 5 6 7 8 9 10
• 可以看到时间充分长以后3个城市汽车数 量趋于180，300，120 • 可以考察这个结果与初始条件是否有关 • 若最开始600辆汽车都在A市，可以看到 变化时间充分长以后，各城市汽车数量 趋于稳定，与初始值无关
Matlab实现
• • • • • • • 首先建立一个关于变量n ，r的函数 function x=sqh(n,r) a=1+r; x=100; for k=1:n x(k+1)=a*x(k); end
Matlab解差分方程
• • • • • •
function x=czqc(n) A=[0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6]; x(:,1)=[200,200,200]'; for k=1:n x(:,k+1)=A*x(:,k); end
如果直接看10年或者20年发展趋势，可以直接在命令窗 口（commond window）作，而不是必须编一个函数
c 1 q )c n ( 0
n
以cn=c0/2代入即求解。
一阶线性常系数差分方程
• 濒危物种的自然演变和人工孵化 • 问题 Florida沙丘鹤属于濒危物种，它在较好
自然环境下，年均增长率仅为1.94%，而在中 等和较差环境下年均增长率分别为 -3.24% 和 -3.82%，如果在某自然保护区内开始有100只 鹤，建立描述其数量变化规律的模型，并作 数值计算。
• • • • • • •
A=[0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6]; >> n=10; >> for k=1:n x(:,1)=[200,200,200]'; x(:,k+1)=A*x(:,k); end >> round(x)
作图观察数量变化趋势
300
• • • • •
直接输入x（:,1）的值即可
x(:,1)=[600,0,0]; round(x'); plot(k,x),grid
600 500
400
300
200
100
0
0
1
2
3
4
5
6
7
8
9
10
按年龄分组的种群增长
matlab差分方程积分
matlab差分方程积分Matlab中进行差分方程积分可以采用以下步骤：1. 首先，定义差分方程。
假设我们要解决的差分方程是\(y(n+1) - y(n) = r*(y(n) - y(n-1))\)，其中\(r\) 是常数。
2. 在Matlab中，我们可以使用内置的函数来求解差分方程。
例如，我们可以使用`diff` 函数来计算差分，`poly` 函数来找到多项式的根，以及`roots` 函数来找到多项式的零点。
3. 在定义了差分方程之后，我们可以使用`for` 循环来迭代计算每个\(y(n)\) 的值。
在每次迭代中，我们可以使用前一个\(y(n)\) 和\(y(n-1)\) 的值来计算新的\(y(n)\) 值。
4. 最后，我们可以将所有的\(y(n)\) 值绘制成图表，以直观地观察解的动态行为。
可以使用Matlab 中的`plot` 函数来绘制图表。
下面是一个简单的示例代码：```matlab% 定义差分方程r = 0.5;y = zeros(1, 10); % 初始化y 值y(1) = 1; % 初始条件y(2) = r*y(1);for n = 3:10y(n) = r*(y(n-1) - y(n-2)) + y(n-1);end% 绘制解的动态行为figure;plot(y);xlabel('n');ylabel('y(n)');title('Solution of the Difference Equation');```这个代码将计算并绘制差分方程的解的动态行为。
你可以根据需要修改代码以适应你的特定问题。
1。
差分方程模型matlab
差分方程模型matlab差分方程模型在数学和工程领域中具有重要的应用。
它是描述动态系统行为的一种数学模型，通常由一系列离散时刻的状态变量和状态转移方程组成。
MATLAB作为一种功能强大的数值计算软件，为差分方程模型的建模和求解提供了便捷的工具和环境。
本文将介绍差分方程模型在MATLAB中的使用方法和应用场景。
首先，我们将探讨差分方程模型的基本原理和概念，然后详细介绍在MATLAB中的建模步骤和求解技巧。
最后，我们会给出一些在实际问题中使用差分方程模型的案例，并展示其在系统分析、控制和优化等方面的优势。
差分方程模型是描述离散系统行为的数学模型，常用于描述在给定时间步长下变量之间的关系。
它与连续时间的微分方程模型相对应，但在很多情况下，离散系统更符合实际情况。
差分方程模型可以描述许多系统，例如电路、金融市场、人口增长等。
在MATLAB中建立差分方程模型需要以下步骤：1. 定义变量：首先需要确定模型涉及的状态变量，然后在MATLAB 中声明这些变量。
可以使用向量或矩阵表示多个变量。
2. 构建状态转移方程：差分方程模型通过状态转移方程描述系统变量在不同时间步长之间的变化规律。
在MATLAB中，可以使用循环或矩阵运算构建状态转移方程。
3. 设定初值条件：差分方程模型通常需要给定初始条件，即在 t=0 时刻各个变量的值。
在MATLAB中，可以使用向量或矩阵存储初始条件。
4. 求解差分方程：在MATLAB中可以使用函数或求解器来求解差分方程模型。
常用的函数包括 `solve`、`ode45`、`ode15s`等，它们可以根据模型的具体特点选择合适的求解方法。
在实际应用中，差分方程模型在系统分析、控制和优化等方面具有广泛的应用。
例如，在系统分析中，可以通过建立差分方程模型来预测系统的行为和变化趋势。
在控制问题中，差分方程模型可以描述系统动态行为，从而设计和优化控制策略。
在优化问题中，差分方程模型可以作为约束条件或目标函数进行求解。
《数学建模与数学实验》电子课件-赵静、但琦 第9讲 差分方程
例1
a
n a1
an1 1, a2
na n2 2
其中， a1 1, a2 2 叫初始值.
例 2 设第一月初有雌雄各一的一对小兔.假定两月后长成成兔，
同时（即第三月）开始每月初产雌雄各一的一对小兔，
新增小兔也按此规律繁殖.设第 n 月末共有 Fn 对兔子， 试建立关于 Fn 的差分方程.
(1 2
5 )n
由初始条件 F1 1, F2 1 ，得
c1
(1 2
5
)
c2
(1
2
5)
1
联立解出
c1
1 (
2
5 )2
c2
1 (
2
5 )2
1
c1 1 , c2 1
5
5
故
Fn
1 [(1 5 )n (1 5 ) n ]
52
2
定理 2 (重根)差分方程 an b1an1 b2an2 bk ank 0 ,
况下每个人的行进轨迹。
D
B
O
C
图1
求解过程:
1. 建立平面直角坐标系:A(x1,y1), B(x2,y2), C(x3,y3), D(x4,y4)。
2. 取时间间隔为 Δt,计算每一点在各个时刻的坐标。 设某点在 t 时刻的坐标为: (xi,yi) 则在 t+Δt 时刻的坐标为:(xi+vΔtcosα,yi+vΔtsinα) 其中
有如下问题需要解决: 1) 分析如何实现可持续捕获(即每年开始捕捞时渔场中各年龄组鱼 群不变)，并在此前提下得到最高收获量； 2) 合同要求某渔业公司在5年合同期满后鱼群的生产能力不能受到 太大的破坏,承包时各年龄组鱼群数量为122，，，3.29(×109条),在 固定努力量的捕捞方式下,问该公司应采取作怎样的捕捞策略,才能 使总收获量最高.
Matlab解差分方程
Xk-1决定的部分是 a1bcXk-1,由Xk-2决定的部分是 a2b(1-a1)bcXk-2
Xk= a1bcXk-1 + a2b(1-a1)bcXk-2
Xk= a1bcXk-1 + a2b(1-a1)bcXk-2
• 实际上，就是Xk= pXk-1 + qXk-2 我们需 要知道x0,a1,a2,c, 考察b不同时，种子繁 殖的情况。在这里假设 • X0=100,a1=0.5,a2=0.25,c=10,b=0.18~0.20 • 这样可以用matlab计算了
模型及其求解
• 记一棵植物春季产种的平均数为c,种子能 活过一个冬天的(1岁种子)比例为b,活过 一个冬天没有发芽又活过一个冬天的（2 岁种子）比例仍为b,1岁种子发芽率a1，2 岁种子发芽率a2。
• 设c,a1,a2固定，b是变量，考察能一直繁殖的条件 • 记第k年植物数量为Xk，显然Xk与Xk-1,Xk-2有关，由
Xk+1=aXk +b ,a=1+r
• • • • • •
function x=fhsqh(n,r,b) a=1+r; X=100; For k=1:n X(k+1)=a*x(k)+b; end
• k=(0:20) ; %一个行向量 • y1=(20,-0.0324,5); 也是一个行向量 • round( [ k ’, y 1 ’] ) 对k,y1四舍五入，但 是 不改变变量的值 plot( k , y1) k y1 是行向量列向量都可以 也可以观察200年的发展趋势，以及在较差 条件下的发展趋势，也可以考察每年孵 化数量变化的影响。
方程(1)的解可以表为 C1,c2 由初始条件x0,x1确定。
1 ,x ( k ) 1 ,x 0 ( k ) 1 , 2 k 1 , 2 k
MATLAB与差分方程
0.7 0.1
x(k)=[x1(k), x2(k), x3(k)]T
0.6 0.2 0.1 A 0 . 3 0 . 7 0 . 3 0 .1 0 .1 0 .6
x(k 1) Ax(k ), k 0,1,2,
例 汽车租赁公司的运营
模型及其求解
开始时600辆汽车平均分配到3个城市
开始时600辆汽车全分配给A市
例 按年龄分组的种群增长
• 动物因自然或人工繁殖而增加，因自然死亡 和人为屠杀而减少; • 不同年龄动物的繁殖率、死亡率有较大差别; • 将种群按年龄等间隔地分成若干个年龄组， 时间离散化为时段; • 在稳定环境下假定各年龄组种群的繁殖率和 死亡率与时段无关; • 建立按年龄分组的种群增长模型; • 预测未来各年龄组的种群数量;
如果每年人工孵化5只鹤放入该保护区， 在中等自然环境下鹤的数量将如何变化?
例1 濒危物种(Florida 沙丘鹤)的自然演变 和人工孵化
模型及其求解 记第k年沙丘鹤的数量为xk， 自然环境下年平均增长率为r
xk 1 axk , a 1 r , k 0,1,2,
160 140 120 100 r=-0.0324 80 60 40 r=-0.0382 r=0.0194
高阶线性常系数差分方程的平衡点及其稳定性
a0 xk n a1 xk n1 an1 xk 1 an xk b,
特征方程 特征根 平衡点
k 1,2,
a0 n a1n1 a n1 a n 0
1, 2, …n
x b /(a0 a1 an1 an )
x1 (k 1) bi xi (k ) (设至少1个bi>0)
怎样用Matlab求解差分方程题解读
Xk+1=aXk +b ,a=1+r
function x=fhsqh(n,r,b) a=1+r; X=100; For k=1:n X(k+1)=a*x(k)+b; end
k=(0:20) ; %一个行向量 y1=(20,-0.0324,5); 也是一个行向量 round( [ k ’, y 1 ’] ) 对k,y1四舍五入，但 是 不改变变量的值 plot( k , y1) k y1 是行向量列向量都可以 也可以观察200年的发展趋势，以及在较差 条件下的发展趋势，也可以考察每年孵 化数量变化的影响。
高阶线性常系数差分方程
如果第k+1时段变量Xk+1不仅取决 于第k时段变量Xk，而且与以前时段变 量有关，就要用高阶差分方程来描述
一年生植物的繁殖
一年生植物春季发芽，夏天开花，秋季 产种，没有腐烂，风干，被人为掠取的 那些种子可以活过冬天，其中一部分能 在第2年春季发芽，然后开花，产种，其 中的另一部分虽未能发芽，但如又能活 过一个冬天，则其中一部分可在第三年 春季发芽，然后开花，产种，如此继续， 一年生植物只能活1年，而近似的认为， 种子最多可以活过两个冬天，试建立数 学模型研究这种植物数量变化的规律， 及它能一直繁殖下去的条件。
用矩阵表示
x1 (k 1) 0.6 0.2 0.1 x1 ( k ) x ( k 1) 0.3 0.7 0.3 x ( k ) 2 2 x (k 1) 0.1 0.1 0.6 x ( k ) 3 3
function x=czqc(n) A=[0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6]; x(:,1)=[200,200,200]'; for k=1:n x(:,k+1)=A*x(:,k); end
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化 ❖ 人工孵化是挽救濒危物种的措施之一，如果每
年孵化5只鹤放入保护区，观察在中等自然条 件下沙丘鹤的数量如何变化 ▪ Xk+1=aXk+5 ,a=1+r ▪ 如果我们想考察每年孵化多少只比较合适，
Matlab语言程序设计
理学院 徐 屹
东北电力大学
Northeast Dianli College
第九讲 差分方程
❖ 污水处理厂每天可将处理池的污水浓度降低一 个固定比例q，问多长时间才能将污水浓度降 低一半？
❖ 记第k天的污水浓度为ck,则第k+1天的污水浓 度为 ck+1=(1-q)ck,k=0,1,2,···· 从k=0开始递推n次得
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化 模型建立：
❖ 记第k年沙丘鹤的数量为xk,年均增长率为r，则 第k+1年鹤的数量为
❖
xk+1=(1+r)xk k=0,1,2······
❖已知x0=100, 在较好，中等和较差的自然环境
下 r=0.0194, -0.0324,和-0.0382 我们利
❖ end
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化 ❖ 在command窗口里调用sqh函数
k=(0:20)';
>> y1=sqh(20,0.0194); >> y2=sqh(20,-0.0324); >> y3=sqh(20,-0.0382); >> round([k,y1',y2',y3'])
可以令 ▪ Xk+1=aXk+b ,a=1+r
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
❖function x=fhsqh(n,r,b)
❖ a=1+r;
❖ x=100;
❖for k=1:n
❖ x(k+1)=a*x(k)+b;
❖ end
fhsqh
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
二、高阶线性常系数差分方程
例2：一年生植物的繁殖
exno93,exno94
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
一阶线性常系数差分方程的解、平衡点及其稳定性
❖ 自然环境下,b=0
xk a k x0
❖ 人工孵化条件下
x k a kx 0 b (1 a L a k 1 )
❖ 令xk=xk+1=x得 ❖ 差分方程的平衡点
cn (1q)nc0
以cn=c0/2代入即求解。
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
❖ 濒危物种的自然演变和人工孵化 ❖问题 Florida沙丘鹤属于濒危物种，它在较好
自然环境下，年均增长率仅为1.94%，而在中 等和较差环境下年均增长率分别为 -3.24% 和 -3.82%，如果在某自然保护区内开始有100只 鹤，建立描述其数量变化规律的模型，并作 数值计算。
用Matlab编程，递推20年后观察沙丘鹤的数量
变化情况
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
Matlab实现：
❖首先建立一个关于变量n ，r的函数
❖function x=sqh(n,r)
❖ a=1+r;
❖ x=100;
❖for k=1:n
❖ x(k+1)=a*x(k);
sqh
❖ plot(k,y1,k,y2,k,y3) 在同一坐标系下画图 >> subplot(2,2,1),plot(k,y1,'r:') >> subplot(2,2,2),plot(k,y2,'k-') >> subplot(2,2,3),plot(k,y3,'g*') >> subplot(2,2,4),plot(k,y1,'r:',k,y2,'k',k,y3,'g*')
Xk= a1bcXk-1 + a2b(1-a1)bcXk-2
二、高阶线性常系数差分方程
例2：一年生植物的繁殖
Xk= a1bcXk-1 + a2b(1-a1)bcXk-2 实际上，就是Xk= pXk-1 + qXk-2 我们需要知道x0,a1,a2,c, 考 察b不同时，种子繁殖的情况。在这里假设 X0=100,a1=0.5,a2=0.25,c=10,b=0.18、0.19、0.20 这样可以用matlab计算了
二、高阶线性常系数差分方程
二、高阶线性常系数差分方程
例2：一年生植物的繁殖
记一棵植物春季产种的平均数为c,种子能活过一个冬天的(1岁 种子)比例为b,活过一个冬天没有发芽又活过一个冬天的（2岁 种子）比例仍为b,1岁种子发芽率a1，2岁种子发芽率a2。 设c,a1,a2固定，b是变量，考察能一直繁殖的条件 记第k年植物数量为Xk，显然Xk与Xk-1,Xk-2有关，由 Xk-1决定 的部分是 a1bcXk-1,由Xk-2决定的部分是 a2b(1-a1)bcXk-2
ak x0
b1 ak 1 a
x
b
1 a
❖ k→∞时，xk→x,称平衡点是稳定的
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
二、高阶线性常系数差分方程 如果第k+1时段变量Xk+1不仅取决于第k时段 变量Xk，而且与以前时段变量有关，就要用 高阶差分方程来描述.
exno91
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化 利用plot 绘图观察数量变化趋势
❖ 可以用不同线型和颜色绘图 ❖rgb c m ykw ❖ 红 绿 兰 兰绿 洋红 黄 黑 白色 ❖ ： + o * . X s d 表示不同的线型
一、一阶线性常系数差分方程
例1：濒危物种的自然演化和人工孵化
❖k=(0:20); %一个行向量 ❖y1=fhsqh(20,-0.0324,5); % 也是一个行向量 ❖round([ k’,y1’])% 对k,y1四舍五入，但是不改变变
量的值 ❖plot(k,y1) %k,y1是行向量列向量都可以 ❖ 也可以观察200年的发展趋势，以及在较差条件下的
发展趋势，也可以考察每年孵化数量变化的影响。