严水仙
【期刊名称】《赣南师范学院学报》
【年(卷),期】2018(039)003
【总页数】5页(P10-14)
【关键词】常系数线性微分方程;Matlab;矩阵指数
【作 者】严水仙
【作者单位】赣南师范大学数学与计算机科学学院,江西赣州341000
【正文语种】中 文
【中图分类】O175
微分方程课程是高校不少理工科专业(如数学、力学、控制等) 的重要基础理论课程.常微分方程是描述自然科学、工程技术和社会科学中的运动、演化和变化规律的重要连续型模型. 物理、化学、材料、医学、经济学等领域中的许多原理和规律都可以描述成相应的微分方程, 如生物种中的生态平衡、流行病存在的阈值定理、化学反应中的稳定性、遗传基因变异、股票的涨幅趋势、利率的浮动、市场均衡价格的变化等.描述、认识和分析其中的规律可以通过研究相应的微分方程数学模型来实现.[1]
在微分方程的理论中，线性微分方程组是非常值得重视的一部分内容，它是了解并掌握非线性微分方程、非线性动力系统、非线性控制等课程的基础. 常系数线性微分方程组的求解是线性微分方程组理论中最简单、最直观的部分，熟悉并掌握常系数线性微分方程的求解将有利于更好的理解线性系统的基本理论.
matlab求导Matlab是由美国的Cleve Moler博士等[2-3]于1980年提出的以矩阵运算为基础，把计算、程
序设计等融合到了一个简单易用的交互式工作环境中.可实现工程计算、算法研究、符号运算、建模和仿真、原型开发、数据分析及可视化、科学和工程绘图、应用程序设计等功能.Matlab强大的运算功能和图形使其成为目前世界上应用最为广泛的科学计算软件之一, 在教学中能快速的计算方程的解并描绘直观的几何图形.[4-6]鉴于此,本文主要介绍借助于Matlab来求解常系数线性微分方程组，通过利用Matlab命令，计算系数矩阵的特征值、特征向量、矩阵指数求解线性微分方程组.
1 常系数线性微分方程的基本理论[1]
定理1[1] 如果A(t)是n×n阶矩阵函数，f(t)是n维列向量函数.它们都在区间a≤t≤b上连续，则对区间a≤t≤b上的任意t0∈[a,b]及任一常数n维列向量η,方程组
x′=A(t)x+f(t)(1)
存在唯一解φ(t)，定义于整个区间a≤t≤b上，且满足初值条件φ(t0)=η.
定理2[1] 齐次线性微分方程组x′=A(t)x一定存在n个线性无关的解x1(t),x2(t),…,xn(t).
定理 3[1] 齐次线性微分方程组x′=A(t)x一定存在一个基解矩阵Φ(t).如果ψ(t)是方程组的任意解，那么
ψ(t)=Φ(t)c,(2)
这里c是确定的n维常数列向量.
定理4[1] 如果Φ(t)，ψ(t)在区间a≤t≤b上是x′=A(t)x的两个基解矩阵，那么，存在一个非奇异n×n常数矩阵C，使得在a≤t≤b区间上ψ(t)=Φ(t)C.
定理5[1] 设Φ(t)是齐次线性微分方程组x′=A(t)x的基解矩阵，是非齐次线性微分方程组x′=A(t)x+f(t)的某一个解，则方程组x′=A(t)x+f(t)的任一解φ(t)都可表示为
(3)
这里c是确定的n维常数列向量.
矩阵指数exp A的定义和性质： 如果A是一个n×n常数矩阵，则定义为矩阵指数，其中E为n阶单位矩阵，且规定A0=E,0!=1,对于所有元均为0的矩阵，易知expO=E.根据矩阵指数的定
义及级数的收敛性，易知对一切矩阵A都是绝对收敛，在t的任何有限区间上是一致收敛.矩阵指数expA有如下性质：
(Ⅰ)如果矩阵A,B是可交换的矩阵，即AB=BA,则exp(A+B)=exp A exp B.
(Ⅱ)对于任何矩阵A,(exp A)-1存在，且(exp A)-1=exp(-A).
(Ⅲ)如果T是非奇异矩阵(可逆矩阵)，则exp(T-1AT)=T-1(exp A)T.
设常系数线性微分方程组为
x′=Ax(4)
其中A是n×n阶常数矩阵,x=(x1,x2,…,xn)T.
定理6 矩阵
Φ(t)=exp At(5)
为方程组(4)的基解矩阵，且Φ(0)=E.
证明 由矩阵指数的定义易得Φ(0)=E.(1.5)式对t求导， 我们得到Φ′=(exp A t) exp At=AΦ(t),所以Φ(t)=exp At是方程(4)的解矩阵.又因为det Φ(0)=det E=1,因此Φ(t)是(4)的基解矩阵.
定理6已经给出了常系数线性微分方程组(4)的基解矩阵，但是exp At是由At的矩阵级数定义的，矩阵中的每个元是什么没有具体给出，下面我们讨论exp At的计算方法.
2 常系数线性微分方程组基解矩阵的求解
定理7 如果矩阵A具有n个线性无关的特征向量v1,v2,…,vn，它们对应的特征值分别为λ1,λ2,…,λn(不必各不相同)，那么矩阵ψ(t)=[eλ1tv1,eλ2tv2,…,eλntvn],(-∞≤t≤+∞)是常系数线性微分方程组(4)的一个基解矩阵.
证明 因为λ1,λ2,…,λn为矩阵A的特征值，对应的特征向量为v1,v2,…,vn，有det(λiE-A)=0,且(λiE-A)vi=0,i=1,2,…,n.
又因为eλit≠0，则λieλitvi=Aeλitvi，即φi(t)=eλitvi是方程组(4)的一个解.因此，矩阵ψ(t)=[eλ1tv1,eλ2tv2,…,eλntvn]是(4)的一个解矩阵.
v1,v2,…,vn是矩阵A线性无关的特征向量组，所以det ψ(0)=det[v1,v2,…,vn]≠0，因此可知，ψ(t)=[eλ1tv1,eλ2tv2,…,eλntvn]是(4)的一个基解矩阵.
由此可知，Φ(t)=exp At和ψ(t)=[eλ1tv1,eλ2tv2,…,eλntvn]均为方程组(4)的基解矩阵，根据定理4的结论可知，存在一个非奇异的常数矩阵C,使得Φ(t)=exp At=ψ(t)C.令t=0，我们得到C=ψ-1(0)，因此exp At=ψ(t)ψ-1(0).
所以，常系数线性微分方程的求解转化为求系数矩阵的特征值和特征向量.
例1 求下列方程组的基解矩阵
解 系数矩阵为 的特征方程为系数矩阵的特征值为λ1=1,λ2=2,λ3=3.
设λ1=1时对应的特征向量v1=(a,b,c)T，则(λ1E-A)v1=0，解得特征向量为其中c1为参数.
设λ2=2时对应的特征向量v2=(a,b,c)T，则(λ2E-A)v2=0，解得其中c2为参数.
设λ3=3时对应的特征向量v3=(a,b,c)T，则(λ3E-A)v3=0，解得其中c3为参数.
由定理7知，方程组的基解矩阵为：
根据定理3，方程组的通解为φ(t)=ψ(t)C,其中C=(c1,c2,c3)T.
常系数线性微分方程的求解实际就是求系数矩阵的特征值和特征向量，或者通过拉普拉斯变换法求解以及直接求解系数矩阵的矩阵指数exp(At). 然而，不管矩阵的的特征值、特征向量、拉普拉斯变换、还是矩阵指数，直接求解均比较复杂，特别是系数矩阵的特征值出现重根、复根等情况时，求解通解就显得更为困难.下面介绍使用Matlab数学软件求解常系数线性方程组.
3 Matlab在求解常系数微分方程中的应用
根据线性代数、高等代数的理论，矩阵对角化过程是一个复杂的计算过程，特别是系数矩阵的特征值出现重根、复根等情况时，特征向量计算比较麻烦. 学习使用Matlab软件在求解数学问题的原理将能够让学生更好的理解数学思想，减少重复性计算等. 下面将介绍几个利用Matlab计算基解矩阵的例子.
例2 求下列线性微分方程组满足初始条件φ(0)=η的解
其中
解 在Matlab软件中直接输入如下命令，
>> A=[0 1 0;0 0 1;-6 -11 -6];
>> syms t;
>> expm(A*t) %运行得线性方程组的基解矩阵exp(At)
ans =
[3*exp(-t) - 3*exp(-2*t) + exp(-3*t), (5*exp(-t))/2 - 4*exp(-2*t) + (3*exp(-3*t))/2, exp(-t)/2 - exp(-2*t) + exp(-3*t)/2]
发布评论