用matlab实现线性常系数差分方程的求解

模型建立 记第k年沙丘鹤的数量为xk,年均增长率为 r，则第k+1年鹤的数量为
xk+1=(1+r)xk k=0,1,2······
已知x0=100, 在较好，中等和较差的自然 环境下 r=0.0194, -0.0324,和-0.0382 我们利用 Matlab编程，递推20年后观察沙丘鹤的 数量变化情况
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计算了
K=(0:20)’; Y1=zwfz(100,21,0.18); Y2=zwfz(100,21,0.19); Y3=zwfz(100,21,0,20); Round([k,y1’,y2’,y3’]) Plot(k,y1,k,y2,’:’,k,y3,’o’), Gtext(‘b=0.18’),gtext(‘b=0.19’),gtext(‘b=0.20’)
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）作，而不是必须编一个函数

matlab差分方程MATLAB是一种广泛使用的计算机辅助工具，其中包含了许多实用算法和解决方案。

差分方程是MATLAB中非常重要的一种工具，可以用于模拟和解决各种差分方程问题。

下面将介绍如何使用MATLAB来解决差分方程问题。

首先在MATLAB窗口中打开一个新的脚本文件（Ctrl+N），左侧显示脚本编辑器的窗口。

在窗口中输入以下内容：function dy = diffeq(t,y)dy = zeros(2,1);dy(1) = y(2);dy(2) = -0.1*y(2) - y(1) - 10*(y(1)^3);在这个脚本中，我们定义了一个名为“diffeq”的函数，它有两个参数（t和y）。

该函数返回一个长度为2的dy向量，dy是y的导数（dy/dt）。

在本例中，我们使用了系统描述的常见方法：x'=f（x，t），y是系统状态向量。

换句话说，我们将梯度设置为我们想要模拟的方程。

一旦函数被定义，我们现在可以开始运行模拟。

接下来，我们将使用MATLAB的ODE求解器来解决我们的差分方程问题。

我们可以这样编写代码：[t,y] = ode45(@diffeq,[0 30],[1 0]);在这里，ode45是MATLAB中用于解决常微分方程的函数，它需要三个参数。

第一个参数是定义我们的方程的函数（即我们之前声明的diffeq函数），第二个参数是我们期望的时间范围（从0到30，单位为秒），第三个参数是初值（在这个例子中，我们使用y(0)=1和y'(0)=0作为初值）。

运行后，MATLAB会将结果存储在两个向量t和y中，我们可以使用下面的代码来显示不同时间点t的y值：plot(t,y(:,1),'-')在此代码中，我们使用plot函数来绘制y的前一个元素（我们的状态向量正在被建模）以及时间t之间的关系。

结果应该是一个类似于与时间的函数y（t）的曲线。

这个值可以根据不同的初值和系统变量被改变。

Matlab中的差分方程求解差分方程是微分方程的离散形式，常常在计算机科学、物理学、经济学等领域中广泛应用。

而Matlab作为一种数值计算软件，提供了强大的工具和函数用于差分方程的求解和分析。

本文将介绍Matlab中差分方程求解的基本方法和常见应用。

一、差分方程的基本概念差分方程是一种通过递归关系描述变量之间关系的数学表达式。

与微分方程不同，差分方程是以离散时间点为基础的，适用于描述离散系统的动态行为。

一般来说，差分方程可以分为线性差分方程和非线性差分方程两类。

线性差分方程的一般形式为：y[n] = a*y[n-1] + b*y[n-2] + ... + c*x[n]，其中y[n]为方程的解，x[n]为给定的输入，a、b、c为系数。

而非线性差分方程则没有这种简单的表达形式，通常需要通过迭代或数值方法求解。

二、在Matlab中，可以利用函数和工具箱来求解差分方程。

下面将介绍几种常见的求解方法。

1. 符号计算方法Matlab的符号计算工具箱提供了一系列用于求解差分方程的函数，例如dsolve()函数。

这些函数可以根据给定的差分方程自动进行符号运算，得到方程的解析解。

符号计算方法适用于简单的线性差分方程，对于复杂的非线性差分方程则很难求解。

2. 数值迭代方法对于非线性差分方程，常常采用数值迭代的方法来求解。

Matlab提供了多种迭代函数，例如fsolve()函数和fminsearch()函数。

这些函数可以根据给定的差分方程和初始值，通过迭代计算得到方程的数值解。

数值迭代方法适用于各种类型的差分方程，但需要注意选择合适的初始值和迭代算法以确保收敛。

3. 差分方程求解函数除了符号计算和数值迭代方法，Matlab还提供了一些专门用于求解差分方程的函数，例如ode23()函数和ode45()函数。

这些函数可以根据给定的差分方程和初始条件，通过数值方法求解方程的数值解。

相比于数值迭代方法，差分方程求解函数更加高效和稳定，适用于大规模的复杂差分方程。

• title('(a)'); xlabel('(n)');ylabel('(n)');%ep141.m: 调用filter解差分方程y(n)-ay(n-1)=x(n)
• a=0.8; ys=1;
%设差分方程系数a=0.8，初始状态：y(-1)=1
• xn=[1,zeros(1,30)] %x(n)=单位脉冲序列，长度N=31
• B=1; A=[1,-a]; %差分方程系数
• xi=filtic(B,A,ys); %由初始条件计算等效初始条件的输入序列xi
• yn=filter(B,A,xn,xi); %调用filter解差分方程，求系统输出信号y(n)
图形表现
批注本地保存成功开通会员云端永久保存去开通
用matlab解差分方程
本例子截取书上第19页，
matlab求解程序
• %ep141.m: 调用filter解差分方程y(n)-ay(n-1)=x(n)
• a=0.8; ys=1;
%设差分方程系数a=0.8，初始状态：y(-1)=1
• xn=[1,zeros(1,30)] %x(n)=单位脉冲序列，长度N=31
• B=1; A=[1,-a]; %差分方程系数
• xi=filtic(B,A,ys); %由初始条件计算等效初始条件的输入; %调用filter解差分方程，求系统输出信号y(n)
• n=0:length(yn)-1;
• subplot(3,2,1);stem(n,yn,'.')

用Matlab求解差分方程
采用filter函数实现线性常系数差分方程的递推 求解 yn=filter(B,A,xn) 计算输入信号xn的零状态响 应yn yn=filter(B,A,xn,xi) 计算输入信号xn的全响应 yn，xi为等效初始条件的输入序列 xi=filtic(B,A,ys,xs) 由初始条件计算xi的函数
程序结果
(a) 2 1.5
y(n)
1 0.5 0
0
5
10
15 n (b)
20
25
30
1
ห้องสมุดไป่ตู้
y(n)
0.5
0
0
5
10
15 n
20
25
30
xn=[1 zeros(1,20)] B=[2,3] A=[1,0.5,0.06] ys=[1,2] xi=filtic(B,A,ys) yn1=filter(B,A,xn) yn2=filter(B,A,xn,xi) subplot(2,1,1) n1=0:length(yn1)-1 stem(n1,yn1,'.') axis([0,21,-3,3]) subplot(2,1,2) n2=0:length(yn2)-1 stem(n2,yn2,'.')
%初始状态: y(-1)=1
xn=［1, zeros(1, 30)］; %x(n)=单位脉冲序列, 长度N=31 B=1; A=［1, -a］; xi=filtic(B, A, ys); % %由初始条件计算等效初始条件 的输入序列xi
yn=filter(B, A, xn, xi); %调用filter解差分方程， 求
系统输出信号y(n)。
n=0:length(yn)-1; subplot(3, 2, 1); stem(n, yn, '.') title('(a)'); xlabel('n'); ylabel('y(n)') 程序中取差分方程系数a=0.8时，得到系统输出y(n)

数字信号处理课程设计题目：试实现线性常系数差分方程的求解学院：专业：班级：学号：组员：指导教师：题目：用Matlab 实现线性常系数差分方程求解一． 设计要求1． 掌握线性常系数差分方程的求解2． 熟练掌握Matlab 基本操作和各类函数调用 3．结合Matlab 实现线性常系数差分方程的求解二．设计原理1．差分与差分方程与连续时间信号的微分及积分运算相对应，离散时间信号有差分及序列求和运算。

设有序列f(k),则称…,f(k+2),f(k+1),…,f(k －1),f(k －2),…为f(k)的移位序列。

序列的差分可以分为前向差分和后向差分。

一阶前向差分定义为()(1)()f k f k f k ∆=+- （3.1—1）一阶后向差分定义为()()(1)f k f k f k ∆=-- （3.1—2）式中Δ和Δ称为差分算子。

由式（3.1—1）和式（3.1—2）可见，前向差分与后向差分的关系为()(1)f k f k ∆=∆- （3.1—3）二者仅移位不同，没有原则上的差别，因而它们的性质也相同。

此处主要采用后向差分，并简称其为差分。

由查分的定义，若有序列1()f k 、2()f k 和常数1a ，2a 则1122112211221112221122[()()][()()][(1)(1)][()(1)][()(1)]()()a f k a f k a f k a f k a f k a f k a f k f k a f k f k a f k a f k ∆+=+--+-=--+--=∆+∆ （3.1—4）这表明差分运算具有线性性质。

二阶差分可定义为2()[()][()(1)]()(1)()2(1)(2)f k f k f k f k f k f k f k f k f k ∆=∆∆=∆--=∆-∆-=--+- （3.1—5）类似的，可定义三阶、四阶、…、n 阶差分。

一般地，n 阶差分10()[()](1)()nn n jj n f k f k f k j j -=⎛⎫∆=∆∆=-- ⎪⎝⎭∑ （3.1—6）式中!,0,1,2,,()!!n n j n j n j j ⎛⎫== ⎪-⎝⎭ （3.1—7）为二项式系数序列f(k)的求和运算为()ki f i =-∝∑（3.1—8）差分方程是包含关于变量k 的未知序列y(k)及其各阶差分的方程式，它的一般形式可写为,(),(),,()0n F k y k y k y k ⎡⎤∆∆=⎣⎦ （3.1—9a ）式中差分的最高阶为n 阶，称为n 阶差分方程。

高阶线性常系数差分方程
如果第k+1时段变量Xk+1不仅取决 于第k时段变量Xk，而且与以前时段变 量有关，就要用高阶差分方程来描述
一年生植物的繁殖
一年生植物春季发芽，夏天开花，秋季 产种，没有腐烂，风干，被人为掠取的 那些种子可以活过冬天，其中一部分能 在第2年春季发芽，然后开花，产种，其 中的另一部分虽未能发芽，但如又能活 过一个冬天，则其中一部分可在第三年 春季发芽，然后开花，产种，如此继续， 一年生植物只能活1年，而近似的认为， 种子最多可以活过两个冬天，试建立数 学模型研究这种植物数量变化的规律， 及它能一直繁殖下去的条件。
1, 2

5 10 2
b
植物能一直繁殖下去的条件是b>0.191
线性常系数差分方程组
汽车租赁公司的运营
一家汽车租赁公司在3个相邻的城市运营，为方便顾客起见公司 承诺，在一个城市租赁的汽车可以在任意一个城市归还。根据经
验估计和市场调查，一个租赁期内在A市租赁的汽车 在A,B,C市归还的比例分别为0.6,0.3,0.1;在B市 租赁的汽车归还比例0.2,0.7,0.1;C市租赁的归还 比例分别为0.1,0.3,0.6。若公司开业时将600辆 汽车平均分配到3个城市，建立运营过程中汽 车数量在3个城市间转移的模型，并讨论时间 充分长以后的变化趋势。

x2
(k
1)



0.3
0.7
0.3


x2
(
k
)

x3 (k 1) 0.1 0.1 0.6 x3 (k )
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

MATLAB中的差分方程建模与求解方法引言差分方程是数学中常见的一种方程类型，是一种离散形式的微分方程。

在实际问题中，差分方程能够提供对系统的离散描述，对于动态模型的建立和求解具有重要作用。

MATLAB作为一种功能强大的数值计算软件，其内置了丰富的工具箱和函数，为差分方程的建模和求解提供了便利。

一、差分方程的建模差分方程的建模是将实际问题转化为数学方程的过程。

在MATLAB中，差分方程的建模可以通过定义离散系统的状态和状态转移方程来实现。

下面以一个简单的例子说明差分方程的建模过程。

假设有一个人口增长模型，人口数在每年增加10%，则差分方程可以表示为：P(n+1) = P(n) + 0.1 * P(n)，其中P(n)表示第n年的人口数，P(n+1)表示第n+1年的人口数。

在MATLAB中，可以通过定义一个函数来描述差分方程的状态转移方程，代码如下：```matlabfunction Pn = population_growth(Pn_minus_1)growth_rate = 0.1;Pn = Pn_minus_1 + growth_rate * Pn_minus_1;end```上述代码定义了一个名为"population_growth"的函数，该函数的输入参数为上一年的人口数"Pn_minus_1"，输出为当前年的人口数"Pn"。

其中，growth_rate表示人口增长率，根据差分方程的定义，将上一年的人口数乘以增长率再加上本身，即可得到当前年的人口数。

二、差分方程的求解方法在MATLAB中，差分方程的求解可以通过多种方法实现。

下面介绍两种常用的差分方程求解方法：欧拉法和四阶龙格-库塔法。

1. 欧拉法（Euler's method）欧拉法是差分方程求解中最简单直观的一种方法。

其基本思想是通过离散化的方式逐步逼近连续函数的解。

具体步骤如下：1) 将时间段分割成若干个小区间；2) 根据差分方程的状态转移方程，在每个小区间内进行计算；3) 迭代计算直到达到指定的时间点。

x1 ( k + 1) = 0.6 x1 ( k ) + 0.2 x2 ( k ) + 0.1x3 ( k ) x2 ( k + 1) = 0.3 x1 ( k ) + 0.7 x2 ( k ) + 0.3 x3 ( k ) x ( k + 1) = 0.1x ( k ) + 0.1x ( k ) + 0.6 x ( k ) 1 2 3 3
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计算了
500
400
300
200
100
0
0
1
2
3
4
5
6
7
8
9
10
按年龄分组的种群增长
野生或饲养的动物因繁殖而增加，因自然死亡 和人为屠杀而减少，不同年龄动物的繁殖率， 死亡率有较大差别，因此在研究某一种群数量 的变化时，需要考虑年龄分组的种群增长。 将种群按年龄等间隔的分成若干个年龄组，时 间也离散化为时段，给定各年龄组种群的繁殖 率和死亡率，建立按年龄分组的种群增长模型， 预测未来各年龄组的种群数量，并讨论时间充 分长以后的变化趋势。
Matlab求解差分方程问题 用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计算了
Xk= a1bcXk-1 + a2b(1-a1)bcXk-2
• • • • • • • • Function x=zwfz(x0,n,b) C=10;a1=0.5;a2=0.25; p=a1*b*c;q=a2*b*(1-a1)*b*c; X1=x0; X2=p*(x1); for k=3:n X(k)=p*(xk-1)+q*(xk-2); end
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市，可以看到 变化时间充分长以后，各城市汽车数量 趋于稳定，与初始值无关
直接输入x（:,1）的值即可
x(:,1)=[600,0,0]; round(x'); plot(k,x),grid

差分方程的解法分析及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),% 计算并显示零状态响应该程序的运行过程与手算过程对应,显示在命令窗的运行结果与手算结果相同.。

程中出现．
定义3′ 含有两个或两个以上函数值yt,yt+1,…的函数方 程,称为(常)差分方程,出现在差分方程中未知函数下 标的最大差,称为差分方程的阶． n阶差分方程的一般形式为
F(t,yt,yt+1,…，yt+n)=0,
其中F为t,yt,yt+1,…，yt+n的已知函数 ,且yt和yt+n一定要 在差分方程中出现.
其中A1，A2，…，An为n个任意(独立)常数．
定理4(非齐次线性差分方程通解结构定理) 如果 y (t)是非齐次线性方程yt+n+a1(t)yt+n-1+a2(t)yt+n-2
+…+an-1(t)yt＋1+an(t)yt=f(t)的一个特解,yA(t)是其对应的齐
次线性方程 yt+n+a1yt+n-1 +a2yt+n-2 +…+an-1yt+1+anyt=0 的通
解,那么,非齐次线性差分方程的通解为： y(t)=yA(t)+ y (t)
即
y(t)=A1y1(t)+A2y2(t)+…+Anyn(t)+ y(t)， 这里A1，A2,…,An为n个任意(独立)常数．
第二节 一阶常系数线性差分方程
一阶常系数线性差分方程的一般形式为 yt+1+ayt=f(t) 和 yt+1+ayt=0, 其中 f(t) 为 t 的已知函数 ,a≠0 为常数 .分别称为一阶常 系数非齐次线性差分方程和其对应的齐次差分方程.
试以特解 yt =a+bt,(a,b为待定系数)代入方程得

解析法
通过数学变换和求解技巧，得到差分方程的解析解。
数值法
通过数值计算方法，如欧拉法、龙格-库塔法等，得到差分方程的近似解。
02
Matlab在解差分方程中的 应用
Matlab的符号计算功能
符号计算
Matlab提供了符号计算的功能， 可以用于解决差分方程的符号解。 通过符号计算，可以找到差分方 程的通解或特解，并分析解的性 质。
析解。
THANKS
要点二
递推法
通过递推公式求解高阶差分方程，适用于一阶高阶差分方 程。
04
差分方程的数值解法
欧拉方法
总结词
简单易行，但精度较低
详细描述
欧拉方法是一种简单的数值方法，用于求解差分方程。它基于差分方程的递推性质，通 过迭代的方式逐步逼近解。由于其简单性，欧拉方法在许多情况下是首选的数值方法。 然而，由于其精度较低，对于需要高精度解的问题，可能需要采用更高级的数值方法。
特征方程法
通过求解特征方程来求解线性差分方程，适用于多阶线性差分方程。
非线性差分方程的解析解法
迭代法
通过迭代公式求解非线性差分方程，适用于一阶非线性差分方程。
解析法
通过解析表达式求解非线性差分方程，适用于多阶非线性差分方程。
高阶差分方程的解析解法
要点一
降阶法
将高阶差分方程转化为低阶差分方程，再求解低阶差分方 程。
05
Matlab实现差分方程的解 法示例
解析解法的示例
差分方程
$y(n+1) - y(n) = 0$
解析解
$y(n) = y(0)$
解释
该差分方程表示每一项都与前一项相等，因此解为常 数。
数值解法的示例

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

差分方程求解matlab代码差分方程在数值分析中是很重要的应用，matlab可以帮助我们有效的求解差分方程。

首先，我们要明确待求解的差分方程。

比如研究一个某种物质随时间的变化，可以用分类微分方程去描述：dC/dt = k * C其中C表示物质的浓度，t表示时间，k表示一个常量，dC/dt表示物质的瞬时浓度。

要使用matlab求解上面这个差分方程，需要如下步骤：1. 将差分方程写成matlab可以求解的形式。

因为 matlab 提供了 ode45 函数来求解常微分方程，所以我们需要将上面的差分方程转换成对应的常微分方程：dC/dt = k * C可以得到：dC/C = k * dt又即有：dC/C = k * dt2. 使用matlab编写代码进行求解。

假设已知t0=0, C0=3，t = 10（s），k = 2（1/s），则matlab代码如下：t0=0;t1=10; %定义起始时间、终止时间C0=3; %定义起始浓度k=2; %定义常量f=inline('k*x','x','k'); %定义右端函数[t,c]=ode45(f,[t0,t1],C0); %调用ode45函数进行求解plot(t,c), grid on; xlabel('Time(s)'); ylabel('Concentration(M)'); %绘制图像最后，使用 plot 函数绘制出解曲线，看到在 t0=0 时，物质的浓度为3 ，在 t=10 时，物质的浓度为 24。

在这里，我们用matlab求解了简单的一维差分方程，可以看到matlab还提供很多用于求解差分方程的函数，比如 ode23， ode113 等，这些函数可以帮助我们在不同的问题中有效地求解差分方程。