第6章 微分方程系统求解的伪谱方法 (1)
谱方法解微分方程
录
1. 谱方法解微分方程 ............................................................................................................................ 1 1.1 解微分方程的加权余量法(METHODS OF WEIGHT RESIDUALS) ................................................... 1 1.1.1 基本思想.............................................................................................................................. 1 1.1.2 MWR 的基本方法 .................................................................................................................. 2 1.2 正交函数系与谱近似 .................................................................................................................. 3 1.2.1 正交函数系与正交多项式 .................................................................................................. 3 1.2.2 谱近似...............................
谱方法介绍
摘要:近些年来，无限维动力系统得到了很大的发展.随着对它研究的深入和计算能力的迅速提高，使得与之相关的数值研究越来越被人们关注.谱方法作为一种数值求解偏微分方程的方法，它具有无穷阶收敛性.因此，谱方法也就引起人们更多的关注.关键词：谱方法；偏微分；收敛；逼近；1偏微分方程及其谱方法的介绍偏微分方程主要借助于未知函数及其导数来刻画客观世界的物理量的一般变化规律。
理论上，对偏微分方程解法的研究已经有很长的历史了。
最初的研究工作主要集中在物理，力学，几何学等方面的具体问题，其经典代表是波动方程，热传导方程和位势方程（调和方程）。
通过对这些问题的研究，形成了至今仍然使用的有效方法，例如，分离变量法，fourier变换法等。
早期的偏微分方程研究主要集中在理论上，而在实际操作中其研究方法和研究结果都难以得到广泛的应用。
求解的主要方法为：有限差分法，有限元法，谱方法。
谱方法起源于Ritz-Galerkin方法，它是以正交多项式（三角多项式，切比雪夫多项式，勒让得多项式等）作为基函数的Galerkin方法、Tau方法或配置法，它们分别称为谱方法、Tau方法或拟谱方法（配点法），通称为谱方法。
谱方法是以正交函数或固有函数为近似函数的计算方法。
从函数近似角度看．谱方法可分为Fourier方法．Chebyshev或Legendre方法。
前者适用于周期性问题，后两者适用于非周期性问题。
而这些方法的基础就是建立空间基函数。
下面介绍几种正交多项式各种节点的取值方法及权重。
1) Chebyshev-Gauss:2) Chebyshev-Gauss-Radau: x0 =1,3) Chebyshev-Gauss-Lobatto: x0 =1, xN =1,4)Legendre-Gauss: xj 是的零点且5）Legendre-Gauss-Radau: xj 是的N+1个零点且6）Legendre-Gauss-Lobatto: x0=-1,xN=1其它N-1个点是的零点且下面介绍谱方法中最重要的Jacobi正交多项式其迭代公式为：其中：Jacobi正交多项式满足正交性：而Chebyshev多项式是令时Jacobi多项式的特殊形式，另外Legendre多项式是令时Jacobi多项式的特殊形式。
LGL伪谱算法
Lagrange-Gauss-Lobatto （LGL ）伪谱算法一般的非线性优化控制问题可以描述为00min [()](,)d s.t.(,),(0),()ft f f f LUL UJ x t L tf t ϕ=+===≤≤≤≤⎰x u x x u x x x x x x xu u u (1.1)其中：x 为运动状态，u 为控制力，0x 为初始状态，fx 为目标状态，L Ux x 、为运动状态的约束界，L U u u 、为控制力的约束界，f t为末端时刻，J 为性能指标。
引入新的时间1t τα=-以及时间因子2/f t α=，将时间区间[0,]f t 变换为[1,1]-。
令()N L τ为定义在[1,1]-上的N 阶Legendre 多项式。
令j τ为()N L τ'的零点，11j N ≤≤-，其中撇号为对τ的导数。
令01τ=-，1N τ=，则x n 维运动状态()τx ，u n 维控制力()τu 的多项式近似分别为00()()()()Nj j j N j j j τϕττϕτ====∑∑xx uu (1.2)其中()j j τ=x x ，()j j τ=u u 为运动状态和控制力在插值点处的值，基函数()j ϕτ为2(1)()1()(1)()N j N j jL N N L ττϕττττ'-=+- (1.3)不仿构造合适的0{}N j j =x 和0{}Nj j =u ，使得()τ'x 在k τ处的估计误差为零，即有()()()(,)0Nk k k j j k k k j f ττϕτα='''∆=-=-=∑q x x x u(1.4)令0011[,,,,,,]T T T T T T TN N =Z x u x u x u ，则可将优化问题(1.1)转化为min ()s.t.()0L UJ =Φ∆=≤≤Z Z Z Z Z (1.5)从而实现了对连续控制问题转化为非线性规划问题，可以采用遗传算法、神经网络、单纯形法等求解，得到最优控制率*()τu 。
解析微分方程的常见近似解法与稳定性分析
解析微分方程的常见近似解法与稳定性分析微分方程是数学中的重要概念，广泛应用于自然科学、工程技术等领域。
解析微分方程通常是一项艰巨的任务，但常见的近似解法可以在某些情况下提供有效的近似解，并且对解的稳定性进行分析。
本文将介绍几种常见的近似解法，并探讨它们的稳定性。
一、欧拉法欧拉法是最简单的近似解法之一，适用于一阶常微分方程。
它基于差分近似的思想，将微分方程转化为差分方程。
具体步骤如下：1. 将微分方程表示为dy/dx=f(x,y)，其中f(x,y)为已知函数。
2. 将x的区间[a,b]分成n个小区间，每个区间的长度为h=(b-a)/n。
3. 定义x的序列x0,x1,...,xn，其中xi=a+i*h。
4. 利用差分近似，得到y的递推公式：yi+1 = yi + h*f(xi,yi)。
欧拉法的稳定性分析较为简单，通常通过步长h来评估。
当步长h较小时，欧拉法的近似解较为准确，并且稳定性较好。
然而，当步长h过大时，欧拉法的误差会较大，并且可能导致解的不稳定性。
二、改进的欧拉法改进的欧拉法是对欧拉法的一种改进，主要通过引入中点来提高近似解的准确性。
具体步骤如下：1. 将微分方程表示为dy/dx=f(x,y)，其中f(x,y)为已知函数。
2. 将x的区间[a,b]分成n个小区间，每个区间的长度为h=(b-a)/n。
3. 定义x的序列x0,x1,...,xn，其中xi=a+i*h。
4. 利用差分近似，得到y的递推公式：yi+1 = yi + h*f(xi+0.5h, yi+0.5h*f(xi,yi))。
改进的欧拉法相比于欧拉法，具有更高的精度和稳定性。
通过引入中点，它能够更好地逼近真实解，并减小近似误差。
三、龙格-库塔法龙格-库塔法是一类常见的高阶近似解法，包括二阶和四阶龙格-库塔法。
它们通过计算多个函数值来提高近似解的准确性。
以四阶龙格-库塔法为例，具体步骤如下：1. 将微分方程表示为dy/dx=f(x,y)，其中f(x,y)为已知函数。
伪谱的定义、性质和计算
近代数学选讲——伪谱引言：二十世纪九十年代以前，研究矩阵的传统工具是特征值（谱），它们可以揭示线性和非线性系统的特征，包括稳定性、共振、矩阵迭代的可行性等，因此它们是数学学科的一个重要的标准工具。
在计算数学方面，该问题的理论和数值计算也取得了很多成果。
然而，在科学和工程应用中，人们经常遇到这样的现象：根据特征值或谱的性质所作的判断与许多观察的现象或数值结果不相匹配。
究其原因，主要是这些问题所包含的矩阵往往是非正规的，甚至是高度非正规的。
所以，特征值（谱）对分析非正规矩阵是一个不完美的工具。
作为谱的自然延伸，伪谱是一个针对非正规系统的新工具。
摘要：本文首先介绍了伪谱的定义及性质，然后介绍了经典的FOV方法来粗略地给出了伪谱范围的矩形界定，之后介绍了伪谱计算的两种方法，即随机扰动法和SVD方法，最后给出了伪谱的一个应用。
关键字：伪谱定义及性质、矩形界定、伪谱计算记号及说明：文中所有矩阵均为定义在复数域上的n 阶方阵，H A 表示矩阵A 的共轭转置，I 表示相应阶的单位阵，)Re(z ，)Im(z 分别表示复数z 的实部，虚部数值，),(εz B 表示以z 为中心，ε为半径的闭圆域。
正文：伪谱的定义及性质伪谱的定义：假定有矩阵n n C A ⨯∈，A 的谱是指矩阵A 的特征值的全体，可表示如下：}0)det(:{)(=-∈=ΛA zI C z A我们知道，当)(A z Λ∈时，1)(--A zI 是没有意义的。
如果我们定义+∞=--||)(||1A zI 。
那么当||)(||1--A zI 有限而且非常大时，又会如何？这就导致了伪谱最初的一个定义。
给定ε，矩阵A 的-ε伪谱（)(A εΛ）定义：}||)(:||{)(11--≥-∈=ΛεεA zI C z A ，（1） 等价地，伪谱也可以用扰动矩阵的特征值来定义：}||||),(:{)(εε≤+Λ∈∈=ΛE E A z C z A ，（2） 也就是说，A 的-ε伪谱是A 的任何一个-ε扰动矩阵的特征值全体。
求解偏微分方程三种数值方法
数值模拟偏微分方程的三种方法介绍（有限差分方法、有限元方法、有限体积方法）I．三者简介有限差分方法（Finite Difference Methods）是数值模拟偏微分方程最早采用的方法，至今仍被广泛使用。
该方法包括区域剖分和差商代替导数两个步骤。
首先将求解区域划分为差分网格，用有限个网格节点代替连续的求解区域。
其次，利用Taylor级数展开等方法将偏微分方程中的导数项在网格节点上用函数值的差商代替进行离散，从而建立以网格节点上的值为未知量的代数方程组。
该方法是一种直接将微分问题变为代数问题的近似数值解法，数学概念直观，表达简单，是发展较早且十分成熟的数值方法。
差商代替导数后的格式称为有限差分格式，从格式的精度来考虑，有一阶格式、二阶格式和高阶格式。
从差分的空间离散形式来考虑，有中心格式和迎风格式。
对于瞬态方程，考虑时间方向的离散，有显格式、隐格式、交替显隐格式等。
目前常见的差分格式，主要是以上几种格式的组合，不同的组合构成不同的差分格式。
差分方法主要适用于结构网格，网格的大小一般根据问题模型和Courant 稳定条件来决定。
有限元方法（Finite Element Methods）的基础是虚位移原理和分片多项式插值。
该方法的构造过程包括以下三个步骤。
首先，利用虚位移原理得到偏微分方程的弱形式，将计算区域划分为有限个互不重叠的单元（三角形、四边形、四面体、六面体等），在每个单元上选择合适的节点作为求解函数的插值点，将偏微分方程中的变量改写成由各变量或其导数的节点值与所选用的分片插值基函数组成的线性表达式，得到微分方程的离散形式。
利用插值函数的局部支集性质及数值积分可以得到未知量的代数方程组。
有限元方法有较完善的理论基础，具有求解区域灵活（复杂区域）、单元类型灵活（适于结构网格和非结构网格）、程序代码通用（数值模拟软件多数基于有限元方法）等特点。
有限元方法最早应用于结构力学，随着计算机的发展已经渗透到计算物理、流体力学与电磁学等各个数值模拟领域。
常微分方程的常见解法
Euler折线法
近似导数
y(x0)
y(x1) h
y( x0 )
记为
y( x1 ) y( x0 ) hy( x0 ) y0 h f ( x0 , y0 )
解:设t时刻雪球的体积为
，表面积为 ，
由题得
球体与表面积的关系为
引入新常数
再利用题中的条件得
分离变量积分得方程得通解为
再利用条件 确定出常数C和r代入关系式得 t的取值在 之间。
方程为全微分方程的充要条件
定理2.1 设函数
和
在一个矩形区域
中连续且有连续的一阶偏导数，则
是全微分方程的充要条件为:
（2.3.3)
nan (x x0 )n1
f
x,
an
(
x
x0
)n
n1
n0
展开后比较两端同次幂的系数确定
an ,
y
y0
N n1
cn1 (x n
x0 )n
例：用待定系数法求
dy x2 y2 ,
的近似解。
dx
y(0) 1
解： 令 y a n (x x0 )n, 由 y (0) 1 得 a0 1 n0
([diff(y(x),x)=-y(x)],y(x),
# 定义微分方程
x=-2..2,
# 指定x范围
[[y(-2)=2],[y(-2)=1],[y(-2)=-2]], # 给出3个初始值
dirgrid=[17,17],
第六章 量子力学微扰理论与近似方法
102第六章 近似计算方法§6.1 微扰理论 一、非简并定态微扰论 1、定态微扰论的主要思想在量子力学中，当体系的哈密顿算符不显含时间时，属于定态问题，通过解其基本方程：ˆn n nH E Ψ=Ψ 可以求出Hˆ的本征值和本征函数。
如果H ˆ比较复杂，但是如果H ˆ可以写成两部分： H H H ˆˆˆ0′+= （0ˆH 和H ′ˆ都不显含时间），而且满足下列条件：（1）0ˆH 的本征方程：(0)(0)(0)0ˆnn n H E ψψ= 可以精确求解，即n ε和n Φ是已知的。
（2）0ˆH 和H ′ˆ的差别很大，或者说H ′ˆ很小，可以看作0ˆH 的基础上加一个小的微扰H ′ˆ，故H′ˆ称为微扰项。
这样，我们就可以通过微扰理论来近似求解。
(0)(1)(2)n n n n E E E E =+++ (0)(1)(2)n n n n ψψψψ=+++2、定态微扰计算假设微扰时体系的能量是哈密顿算符0ˆH 的第n 个本征值(0)nE ，这个本征值无简并，即体系于定态(0)n ψ。
当体系受到一个与时间无关的微扰H ˆ′作用时，它将处于一个新的能级nE 和状态n Ψ。
n E 和n Ψ是H H H ˆˆˆ0′+=的本征值和本征函数.即满足： ˆn n nH E Ψ=Ψ 微扰论的主要思想：H ˆ′代表一个微小的扰动，那么我们就有理由认为n E 和(0)n E 相差不多，nΨ和(0)n ψ也十分接近。
（1）、非简并能量的一级修正在非简并微扰情况下，由一级微扰确定一级近似波函数和一级能量修正103010010ˆˆn n n nE E H H Ψ′+Ψ=Ψ′+Ψ 两边左乘()*0n Ψ，并对整个空间积分得：()()()()()()()()τττd H d E d E H n n n n n n n n ∫∫∫Ψ′Ψ−ΨΨ=Ψ−Ψ0*00*01100*0ˆˆ 注意到0ˆH 是厄密算符，所以有： ()()()()()()[]0*ˆˆ0001100*0=Ψ−Ψ=Ψ−Ψ∫∫ττd E H d E H n n n n n n 从而得到()()()τd H E nn n 0*01ˆΨ′Ψ=∫ 即()n H n E n′=1 （2）、非简并能量的二级修正令()()()001l ll n a Ψ=Ψ∑得：000ˆˆn n n n nE E E H H Ψ′′+Ψ′′+Ψ′′=Ψ′′+Ψ′′ ()()()()()()()001010010ˆnn n ll l n l l llH E a E a EΨ′−Ψ=Ψ−Ψ∑∑ 将()()n m m ≠Ψ*0左乘上式两边后，对整个空间积分，所以有()()()()mn n m ml ll n ml l lH d H a E a H ′−=Ψ′Ψ−=−∫∑∑τδδ0*01010ˆˆ 其中()()ml l m d δτ=ΨΨ∫0*()()mnm l n H a E E ′=−100 ()01mn mnm E E H a −′=()()0001m mn mn n E E H Ψ−′=Ψ∑左乘()*0n Ψ，并对整个空间积分得104()()()()()()()2111200*0ˆn nl ll n nl ll n n n E a E H a d E H ++′−=Ψ−Ψ∑∑∫δτ 当n l ≠时，利用0ˆH 的厄密性可得 ()∑∑−′=′=ll n nlnlll n E E H H a E 022即()∑−′′=ll n n E E l H n l Hn E 02ˆ（3）、非简并波函数的一级修正(1)'(0)(0)(0)mn n m mn mH E E ψψ′=−∑ 二、简并定态微扰论 1、简并的处理 （1）问题假设(0)n E 是k 度简并的，0ˆH 属于本征值(0)n E 的本征函数有k 个: k φφφ,,,21 ，且它们已经是相互正交的。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
2016/9/25
航空航天中的计算方法
Page 21
6.4 最佳配点分布 Legendre-Gauss-Lobatto：
2016/9/25
航空航天中的计算方法
Page 22
6.5 微分矩阵与两点边值问题求解 6.5 微分矩阵与两点边值问题求解 6.5.1 微分矩阵的概念
伪谱法将微分方程近似解用Lagrange插值表示：
2016/9/25 航空航天中的计算方法 Page 4
6.1 谱方法及伪谱法的概念 6.1 谱方法及伪谱法的概念 以N+1个全局基函数的加权和近似某一连续函数：
( x ) kk ( x ) y( x ) y
k 0 N
k ( x ) e ikx
Fourier谱方法
其中：k ( x) 为多项式或三角函数。 残差函数： R( x; 0 , 1 , , N ) 例，二阶微分方程求解 y( x ) f ( x, y( x ), y( x )) ( x) f ( x, y ( x), y ( x)) R( x) y 残差为 某种准则下使残差最小，确定系数。
航空航天中的计算方法
授课教师：陈琪锋 中南大学航空航天学院
第二部分 边值问题求解方法
第6章 微分方程求解的伪谱法
பைடு நூலகம்
内容提要 6.1 6.2 6.3 6.4 6.5 谱方法及伪谱法的概念 谱方法与Lagrange插值 正交多项式 最优配点分布 微分矩阵与两点边值问题求解
[1] John P. Boyd, Chebyshev and Fourier Spectral Methods (Second Edition), DOVER Publications, Inc., 2000.Chap.1,3-6 [2] Shen, J., and Tang, T., Spectral and High-Order Methods with Applications（谱方法和高精度算法及其应用）, Science Press, Beijing, 2006, Chap.（1.1-1.3；2.1,2.4）.
x0 1, x N 1,
(a , b) ( 1,1)
( x ) (1 i N 1) xi zeros of LN
i
2 1 (0 i N ) 2 N ( N 1) [ LN ( xi )]
Legendre-Gauss-Lobatto点没有显式表达式，需数值求解
p1 ( x ) x 1
0 p1 , p0 ( x ) xdx 1 ( x )dx
a a b b
pn1 ( x ) ( x n1 ) pn ( x ) n 1 pn 1 ( x ), n 1
(a , b) ( 1,1) 时，得到Legendre多项式 Ln ( x) 当 ( x ) 1 ， (a , b) ( 1,1) 时，得到Chebyshev多项 当 ( x ) (1 x 2 ) ， 式 Tn ( x )
Chebyshev多项式 Tn ( x ) ：
T0 ( x ) 1, T1 ( x ) x , T2 ( x ) 2 x 2 1 Tn 1 ( x ) 2 xTn ( x ) Tn1 ( x ), n 1
2016/9/25
航空航天中的计算方法
Page 14
6.3 正交多项式 正交多项式曲线图：
i 0 N
N
航空航天中的计算方法
Page 9
6.2 谱方法与Lagrange插值 6.2.2 Runge现象 对任意光滑函数f(x)，根据均匀分布的N+1个插值点的函
数值，构造N次Lagrange插值近似，误差随N增大趋于0？ 例： 1 插值点随均匀分布时， f ( x) , x 5, 5 1 x2 误差随点数增多不收敛
与任何不高于n次的多项式正交。
若多项式序列 pn ( x )n 0 是正交的，则多项式 pn1 ( x ) 的零点是互不相同的实数，且位于开区间 (a , b) 内。
2016/9/25
航空航天中的计算方法
Page 12
6.3 正交多项式 6.3.2 正交多项式的生成 根据正交多项式的定义（首一情况为例） p0 ( x ) 1
谱方法应用具有高阶次的全局基函数在整个计算域上
2016/9/25
航空航天中的计算方法
Page 7
6.1 谱方法及伪谱法的概念 伪谱方法精度高、收敛快、存贮省，适用于问题的几何 特征平滑和规则时 伪谱法的问题： 如何选择最优的基函数？
如何选择最优的配点？
2016/9/25
航空航天中的计算方法
Ln ( x )
2016/9/25 航空航天中的计算方法
Tn ( x )
Page 15
6.4 最佳配点分布 6.4 最佳配点分布 6.4.1 Gauss求积与Lagrange插值
将积分表示为被积函数在若干点处的函数值加权和：
1 1
f ( x )dx i f ( xi )
i 0
N
若适当选取 i 和 xi ，可使公式对次数≤ 2N+1的多项式被积 函数均精确成立，节点 xi ( i 0,1, , N ) 称为高斯点。 等价于将函数 f 用Lagrange插值近似为插值多项式，然后求 积分。若选用Gauss点插值，能实现最高精度。 最佳配点（插值点）为Gauss点
2016/9/25 航空航天中的计算方法
谱方法
Page 5
6.1 谱方法及伪谱法的概念
在与未知量个数相对的特定点处令残差为零：配点法
R( i ) 0, i 1, 2, , N a 1 2 N b
采用最佳配点的谱方法，即伪谱法。 加权残差为零：加权残差法
2016/9/25 航空航天中的计算方法 Page 11
6.3 正交多项式 任意n次多项式q(x)均可表示为正交多项式 p0 , p1 , , pn 的线性加权和： q( x ) bn pn bn 1 pn 1 b0 p0 若多项式序列 pn ( x )n 0 是正交的，则多项式 pn1 ( x )
两端点附近的误差大
2016/9/25
端点附近插值点增多，中间可减少
航空航天中的计算方法 Page 10
6.3 正交多项式 6.3 正交多项式 6.3.1 函数正交性与正交多项式
函数f(x)与g(x) 在加权Sobolev空间 L (a , b) 上正交，是指
f , g : f , g ( x ) f ( x ) g( x )dx 0
2016/9/25 航空航天中的计算方法 Page 17
6.4 最佳配点分布 Gauss-Radau求积点 q( x ) p N 1 ( x ) p N ( x ) 定义： q( a ) 0 若采用 x0 a，以及多项式 q( x ) ( x a )的零点 x1 , x2 , , x N 作 为求积点，称为Gauss-Radau求积点。 由方程组： a pk ( x ) ( x )dx i pk ( xi ), 0 k N
Gauss-Lobatto求积点包括端点a和b，适用于两点边值问题
航空航天中的计算方法 Page 19
6.4 最佳配点分布 6.4.3 常用正交多项式的Gauss点 Chebyshev多项式的Gauss点
( x ) (1 x 2 )
1 2
(a , b) ( 1,1)
Chebyshev-Gauss-Lobatto：
i 0
其中，插值基函数：
Ci ( x )
k 0, k i
N
x xk xi xk
C i ( x j ) ij
任意N次多项式
PN ( x ) ii ( x )
i 0
2016/9/25
等价
Lagrange插值形式
PN ( x ) PN ( xi )C i ( x )
2016/9/25 航空航天中的计算方法 Page 16
6.4 最佳配点分布 6.4.2 几类Gauss点 Gauss求积点
对于带权函数的Gauss求积：
b
b
a
f ( x ) ( x )dx i f ( xi )
i 0
N
其中Gauss点为 xi 正交多项式 pn1 的零点。
由方程组：
由方程组：
2016/9/25
b
a
pk ( x ) ( x )dx i pk ( xi ), 0 k N
i 0 N
N
可唯一解出 i ( i 0,1, , N )，并且
b a
q( x ) ( x )dx i q( xi ),
i 0
for all
q P2 N 1
Page 8
6.2 谱方法与Lagrange插值 6.2 谱方法与Lagrange插值 6.2.1 Lagrange插值
对函数f(x)，根据N+1个插值点的函数值，构造N次插值 N 多项式近似： PN ( x ) f ( xi )C i ( x ) PN ( xi ) f ( xi )
b
a
w i ( x ) R( x )dx 0,
i 1, 2, , N ,
wi ( x ) 为权函数
Galerkin法：wi ( x ) i ( x ) 。
2016/9/25
航空航天中的计算方法
Page 6
6.1 谱方法及伪谱法的概念 谱方法、有限单元法、有限差分法的区别： 有限单元法将区间分成一些子区间，在子区间选择局部 多项式基函数 有限差分是局部计算