地球物理计算方法-208秋-习题课

地球物理计算方法地球物理与学院课堂情况反馈复习问题（线性方程组）数值方法（三角分解法）误差分析（ 方程组与条件数）定理5：设三对角系数矩阵A对角占优矩阵，则它可以唯一地分解成矩阵L和U的乘积，A=LU（LU分解）其中，U为单位上二对角矩阵，L为下二对角矩阵。1

u1

12d122anadu1

ua

dn1n1L

0

dn

0

,

U

0n1

1

0

复习Ux

y复习LU分解方法求解方程组（25）令A=LU，方程组（25）Ax=f可表示成LUx=f，进而可以分解成两个方程：Ly

f由Ly=f

得，2d12andn

yn

fnada

dn1n1

0

f

yn1

n1

0

y1

f1

解方程组得，

y1

f1

/

d1i

i

i1

y

(

f

a

y

ii)

/

d

,

i

2,

3,...,

n.（33）追的过程：消去复习由Ux=y得，21

u

yn

x

y

n1

n1

n1

1

xn0

x1

y1

1

u1

1

u0其解为：

i

i

i i

+1x

y

u

x

，i

n

1,

n

2,...,1xn

yn赶的过程：回代（30）复习复习总结起来，方程组（25）的求解可归纳为如下三个步骤：

un1

dn步骤1，按式（32）计算d1

u1

d2

步骤2，方程组（33）顺序求y1

y2

yn步骤3，方程组（30）逆序求xn

xn1

x1.LU分解追的过程：消去赶的过程：回代定理6

设a32anna21an1

an

2A

a31ll32

l33为正定矩阵，则有如下三角阵l11l

21

22ln1

ln

2

ln30

lnn

L

l31使成立，A

L

LT若限定L的主对角线元素取正值，则这种分解是唯一的。（34）复习复习平方根法求解对称正定方程组Ax=b设对称正定矩阵A可以分解为A=LLT，则方程组Ax=b可以归结为两个方程组来求解，TLy

bL

x

y定理7：对称正定矩阵A可分解成A=LDLT

的形式，其中，d12dd

n

0

D

0为对角阵，而l32ln,n1l

21ln1

ln

2L

l31是单位下三角阵。复习复习令A=LDLT，方程组Ax=b可以重新写为：L(DLTx)=b，方程组的求解可归纳为求解如下两个方程组：T1Ly

bL

x

D

y其计算公式分别为：i1

lik

ykk

1,

i

1,

2,...,

n.yi

bi和xi

yi

/

di,

i

n,

n

1,...,1.n

lki

xkk

i1改进平方根法（Cholesky方法）优点：不用进行开方运算（运算量~n3/6）。小结：各方法的适用条件Gauss消元：所有顺序主子式不为0，即可得唯一解（充要条件）；A对角元占优，则Gauss消元可进行（充分条件）；LU分解：A对角元占优，则方程组可解；平方根法：对称正定矩阵，限定对角元为正，分解唯一；复习可以发现当系数矩阵和右端项有微小的扰动，其解也会发生变化。解对系数矩阵和右端项的扰动比较敏感，称这类方程组为方程组。几何学解释：方程组的解为两条接近平行的直线的交点，当其中一条直线的斜率稍微有变化时，新的交点和原交点相差比较远。复习方程组程度分析:右端项b的扰动系数矩阵A的扰动1

A

Ax

xb

b

AA11

A1

Ax

x复习x

b

x

b

cond

(

A)cond

(

A)Ax

xA

A

A1

cond

(

A)和上述两个式子表明：右端项b和系数矩阵A的扰动对解的影响与条件数cond(A)的大小有关，cond(A)越大，扰动对解的影响就越大。因而条件数的值刻画了方程组的程度。记cond(A)

A1

A误差估计式可分别表示为：称为矩阵A的条件数。复习复习近似解：x余量：r

b

Ax定理8：设

x

是方程组Ax=b的一个近似解，其精确解为

x*

r为，x的余量，r

b

Ax

则有x*x*

x

cond

(

A)rb对于 方程组来说，有可能余量r很小而条件数很大，此时不能单独由r的大小判断近似解的精度。习题课3第五章线性方程组迭代法+第六章线性方程组直接法线性方程组迭代法迭代公式迭代过程收敛性向量与矩阵范数雅可比迭代公式高斯-塞德尔迭代公式超松弛法雅可比迭代公式的矩阵表示？高斯-塞德尔迭代公式的矩阵表示？超松弛法迭代公式的矩阵表示？一般的线性方程组

Ax

b

，设将系数矩阵A

为A=L+D+U;其中D为对角阵，L和U分别为严格下三角和严格上三角阵。如果将所给方程L

D

U

x

b

改写为Dx

L据此建立的迭代公式即为雅可比迭代公式，其中，L+U=A-D，xk

1

DL

+

D1b或：x(k

1)(I-D1A)x(k

)22(k

0,1,2,)————————Jacobi迭代法式4于是有迭代公式xxxnnx(

k

)n2

2x(

k

)n1

1x(

k

)

](

k

1)n2n

n(

k

1)22

21

1

2x(

k

)

]1n

nx(

k

)12

2111a111a221ann1(

k

1)

..........

0

x(

k

)

]

a

[b

a

..........

ax(

k

)

0

x(

k

)

[b

a

..........

a

a

0

x(

k

)

[bxk

1

DL

+式7迭代方法称为与Jacobi迭代法(式4)对应的Seidel方法。23xxx(

k

)n

n1

1

n2

2

n(

k

)2n

n(

k

1)n(

k

1)2221

12(

k

)1n

nx(

k

)12

2111a111a221ann1(

k

1)

0

x(

k

)

b

a

a

0

x(

k

)

bx(

k

1)

x(

k

1)x(

k

1)

b

a

a

0

x

a

x

a

xxk

1

D1

Lxk

1

Uxk

D1b如果改写为据此建立的迭代公式xk

1

D1

Lx由此解出xk1

，则有x(k

1)

(D

L)1

Ux(k

)

(D

L)1b即为高斯—塞德尔迭代公式。D

L

x(

k

)i

1,

2,...,

nn(

k

)a

x

),

ij

jj

i

1iij

jiii

1(

k

1)j

1(bi

a

x(

k

1)xi

(1

)x

aSOR方法的计算表达式为,x(k

1)

(1)x(k

)

D1(b

Lx(k

1)

Ux(k

)

)整理之后的超松弛法的矩阵形式为，x(k

1)

(D

L)1[(1)D

U]x(k

)

(D

L)1bP171，习题8xk

1

D1

Lx(k

1)

1 (k

)

1令A=L+D+U，则对应于上式的矩阵形式为，

x

(D

L)

Ux

(D

L)

b线性方程组迭代法迭代公式迭代过程收敛性向量与矩阵范数雅可比迭代公式高斯-塞德尔迭代公式超松弛法为了研究迭代过程的收敛性，需要对向量“大小”引入某种，其范数记为

，它是一个实数，1

2Tn,

x度量。任给向量x

x

,x

,x,

当且仅当x

0

时,(2)

对任意实数

及任意向量

x

,(3)

对于任意向量x

与y

,有其中性质(3)被称为向量范数的三角不等式。且满足：(1)对任意向量x

,

x

0x

0

（非负性）（齐次性）x

xx

y

x

y根据向量范数的定义，常用的范数有：2-范数（向量长度）1-范数范数的通用表达式：np1/

pp

ix

(

x

)i1无穷-范数1122nnix2x

x

2

2

x2

x2

i1

n1

1

2

n

ix

x

x

x

xi1x

max

x1

,

x2

, ,

x对给定的n阶方阵A，比值为A的范数，记为

。Ax

x

0

的上确界称Ax由定义知，对任意向量x，有

Ax

Ax

。矩阵范数具有以下基本性质：(1)

对任意方阵A

，A

0

，当且仅A

0

当时对任意实数

和任意方阵A

，有

A

对任意两个同阶方阵A

和B

，有AA

B

A

BAB

A

B(相容性条件)A

0(正定性)(齐次性)(三角不等式)xAxA

maxx常用的矩阵范数：1A

max1

jnni1ijaA

(

AT

A)2(3)无穷范数（行范数）(2)2-范数（谱范数）(1)1-范数（列范数）nj1ijaA

max1in(4)F-范数(Frobenious

范数)1n

n2

2

i1

j

1A

aij

F

x

A

x

A

(

A)

Ax

)，Ax

A

x(

Ax

x

定理：(

A)

A

，A

为A的任意矩阵范数(A)

max

i

为A的谱半径。1ininn定义：设A

R

的特征值为

(i

1,2,,n)称定理3

对给定的方阵

G，若为非奇异。G

1，则矩阵I

1、迭代收敛的充分条件设将方程组Ax

b

改写成x

Gxxk

1

Gxk

d的形式，据此建立迭代公式这类迭代公式的收敛性与矩阵G的性态有关；定理4

若迭代矩阵

G满足均收敛。（充分条件）

1

，则迭代公式对任意初值x0G若线性方程组的系数矩阵为对角占优阵，则称这个线性方程组为对角占优方程组。i

是对角占优的，如果其主对角线元其它元素绝对值之和：定理5

若A为对角占优矩阵，则它是非奇异的。定理6

对角占优方程组的雅可比迭代公式和高斯—塞德尔迭代公式均收敛。（充分条件）称n

阶方阵A

素的绝对值大于nj

1,

j

i

aij

aiii

1,

2, ,

n2、对角占优方程组3、充要条件作业点评：P170-171，习题五：2，5，6，8，11定理4的推论：

（事后误差估计法）若迭代矩阵G满足||G||＜1，则有满足精度要求的最大迭代次数（事先误差估计法）||G||或L

越小收敛越快x

*

x(k

1)x(k

1)

x(k

)G1

GkGx

*

x(k

)

x(1)

x(0)1

G1k

1x

*

xk

x1

Lk

xkLkx

*

x

x1

x01

L对比：非线性方程求根压缩映像原理：P171，习题11迭代结束的条件对比：非线性方程求根控制误差ε的方法(1)

先计算满足误差要求的迭代次数n，再迭代。由可得ln

|

x

x

|n

1

0

(1

L)ln

L(2)

事后误差估计法。由于因而可用|xn+1-xn|≤ε来控制迭代过程。1n1

nx

*

xn

x

x

1

Lx

*

xnLn

x

x

1

L

1

0(1)若使‖x(k)–x*‖<

,只需即可以事先估计达到某一精度需要迭代多少步；

x(0)1

Bx(1)B

kx(1)

x(0)k

ln(

ε(1

B

)

)

/

ln

B对比：线性方程组迭代法控制误差ε的方法即通过 来控制迭代过程。(2)若使‖x(k+1)

–x*‖<

,只需x

*

x(k

1)x(k

1)

x(k

)B1

B

x(k

1)

x(k

)例：已知方程组用J迭代法求方程组的解，取x(0)=(0,0,0)T，若使误差x(k)-x*<10-5，问需要迭代多少次？于是有，x(1)-x(0)=1.4，B=0.5

。

3

1000

0解

由题知，x(1)=(1.4,0.5,1.4)T，

B

1015

1

310110

3)

/

ln1.40.5105k应满足

k

ln(故取k=19，即需要迭代19次。例题选讲5.1

迭代公式的设计

0.2x(k

)

0.721

2

3

1

2

32

1

33

1

2x(k

1)(k

1)

0.1x

0.2x(k

)

0.83x(k

1)

0.2x(k

1)

0.2x(k

1)21

33

1

2

0.84

x(k

1)

0.2x(k

)

0.2x(k

)

0.84

0.72

0.83x(k

1)x(k

1)

0.1x(k

)

x(k

1)

0.1x(k

1)

0.2x(k

1)

0.1x(k

)

0.2x(k

1)（P168，题1-2，高斯-塞德尔预报-校正系统）D

L

x

D

U

x例题选讲5.2

迭代过程的收敛性（P169-170，题1-4）题4线性方程组直接法消去法误差分析： 问题与条件数三角分解法约当消去法高斯消去法高斯选主元法追赶

U分解平方根法改进平方根法例题选讲:6.2.三角分解的两种模式定义：设A为n阶矩阵(n≥2).称A=LU为矩阵A的三角分解，其中L是下三角矩阵，U是上三角矩阵。为了保证分解的唯一性，引入如下定义：定义：如果L是单位下三角阵，U是上三角阵，称三角分解A=LU为Doolittle分解；如果L是下三角阵，U是单位上三角阵，则称A=LU为Crout分解。定理：如果n阶(n≥2)矩阵A的前n-1个顺序主子式不为零，则A有惟一的Doolittle分解和惟一的Crout分解。22231u

uu1l32u33u

1u1n

lu11

u12

u13

21

L

l31

,U

n1n

ln1

ln

2

ln31unn

0

0231u12

u131

u1ll32

l33l11l1

21

22uuu1n

2n

n1n

2n

L

l31

,U

ln1

ln

2

ln3lnn

0

0例题选讲:6.3.对称的乔里斯基分解题3：设a

aA

A1

Ta1

1

nn

为n阶对称正定矩阵，又1

1，其中L

为单位下三角阵，D

为具有正对角元素的对角阵，试确定c1与dn使成立，1 1

1

1TA

L

D

L11

10T11

0

dcT

L

0

D

0

LA

cT

1

n

1

作业点评：P197-199，习题六：1，2，3,6，7，9，10迭代法：用某种极限过程逐步近线性方程组精确解的方法。具有存贮量少、程序计算简单和原始系数矩阵在计算过程中始终不变等优点，但要求方程组的系数矩阵具有某种特殊性质，以保证迭代过程的收敛性。解大型稀疏矩阵方程组常用。（代数本质：矩阵的加法分解）直接法：经有限步算术运算求得线性方程组精确解的方法。实际计算中由于舍入误差的存在，这种方法求得的仍是线性方程组的近似解。系数矩阵进行递推和更新，未知量保持不变。（代数本质：矩阵的乘法分解）线性方程组迭代法与直接法的比较矩阵分解是一种重要的数值方法设计思想，将复杂的矩阵分解为简单矩阵的求和/积的形式，从而简化矩阵运算（如求逆矩阵、求解线性方程组、求特征值等）。加法分解：线性方程组迭代法（雅可比、高斯-塞德尔）乘法分解：线性方程组直接法、矩阵特征值（LU分解、Cholesky分解、QR分解、SVD分解）迭代-校正-

近（二分法、方程求根+线性方程组迭代法）非线性问题线性化（利用泰勒展开取到一次项）以直代曲（分段线性插值、欧拉折

）差商/差分代微分（数值微分、常微分方程、有限差分法）连续问题离散化（数值积分/微分、常微分方程差分法）（6）平均（龙贝格加速、龙格-库塔+亚当姆斯法、牛顿下山法、SOR法）（7）矩阵分解（加法分解、乘法分解）一些重要的数值算法设计思想总复习Ch0：绪论Ch1：函数插值与拟合方法Ch2：数值积分方法Ch3：常微分方程的数值解法Ch4：方程求根的迭代法Ch5：线性代数方程组的迭代解法

Ch6：线性代数方程组的直接解法误差理论数值

近方程数值解数值代数绪论算法设计计算方法及地球物理数值计算误差分析秦九韶算法二分法误差来源与分类避免误差危害的原则有效数字、误差与误差限（绝对/相对）误差

与估计重点：重要概念：计算方法的两个方面、地球物理正反演重要的计算方法：秦九韶算法、二分法误差的来源与分类（四种、其中两种与数值计算相关）误差与误差限（绝对/相对）有效数字、绝对误差限、相对误差限三者之间的关系误差的

与估计避免误差危害的原则（4条）重要习题：2，8，9，1052插值与拟合函数插值函数

近最小二乘拟合泰勒插值拉格朗日插值牛顿插值三次样条插值埃尔米特插值分段插值重点：构造拉格朗日/牛顿插值多项式、余项定理（估计截断误差）拉格朗日插值基函数的性质差商的定义、性质（三条）与计算（差商表）埃尔米特插值（计算题必考！）龙格现象以及如何克服（分段低次插值）分段线性插值、分段三次埃尔米特插值及其截断误差估计三次样条函数的概念、边界条件插值问题与拟合问题的异同点最小二乘拟合（线性拟合、非线性模型线性化）重要习题：5,7,9,12,16,17,22,23,31,33,355455Hermite插值问题解法：待定系数法基于承袭性基于基函数重节点差商数值积分数值微分数值积分数值微分机械求积公式牛顿-柯特斯公式复化求积公式龙贝格加速算法高斯积分重点：机械求积公式、代数精度、插值型求积公式Newton-Cotes公式（梯形公式、Simpson公式及其余项）复化梯形公式、复化Simpson公式及其余项梯形递推算法（变步长梯形算法）Romberg加速算法（梯形-Simpson-Cotes-Romberg）高斯积分公式（两点、三点）及其代数精度、区间变换求积公式的设计（充分利用对称性）数值微分：中点格式及其加速算法（类比Romberg加速算法）插值型求导公式及其余项（两点公式、三点公式、高阶导数公式）求导公式的设计57重要习题：2,3,4,10,17,20,P86题8一但给定一系列求积节点，则求积系数可通过以下两种方法来确定：求解线性方程组（代数精度条件）利用插值型求积公式根据充要条件或方程组解的存在唯一性，二者是等价的。设计求导公式同样可以利用上述两种方法。构造求积公式的2种情况：求解线性方程组即可

牛顿-科特斯型在此基础上复化求积法龙贝格加速算法（1）若给定[a,b]上的求积节点xi,这时只有

Ai

未知,（2）若不给定

xi

,这时要求解关于

xi

Ai

的非线性方程组

高斯型换元法,定积分区间

[a,b]→[-1,1]常微分方程差分法初值问题边值问题欧拉方法龙格-库塔方法亚当姆斯方法收敛性稳定性方程组与高阶方程方法理论分析重点：Euler格式、Euler隐式公式、Euler两步格式、梯形格式、改进Euler格式Euler格式的分类（显示/隐式、单步/多步）局部截断误差、精度分析（T