数值分析讲稿

§4.

消去法的变形二、平方根法工程实际计算中，线性方程组的系数矩阵常常具有对称正定性，即其各阶顺序主子式及全部特征值均大于零。矩阵的这一特性使它的三角分解也有更简单的形式，从而导出一些特殊的解法，如平方根法与改进的平方根法。定理：设A是对称正定矩阵，则存在唯一的非奇异下三角阵L,使得A

LLT且L的对角元素皆为正.定理证明（1）111111

121u1nu12uuu22u22uu

u

u

u

11n

11

U

unn

unn

证：因A对称正定，其各阶顺序主子式均大于零，故有A

LU

其中L为单位下三角矩阵，U为上三角阵。令D

diag(u

,

,u

),

P

D1U

,则P为单位上三角阵。11

nn

DP

n1,n

un1,n1

A

PT

DP故

A

LU

LDP

AT

PT

DLT由LU分解的唯一性

PT

L定理证明（2）uiiu11

D1

0,

Di

/

Di1

0(i

2,

3, ,

n)由于A是对称正定的,其顺序主子式均大于零。故T令

D

diag

(

u11

,

,

unn

)

D

D

D则A

PT

DP

PT

DT

DP

(DP)T

DP

LLT其中L

(DP)T

为非奇异下三角阵,且对角元素皆为正数。唯一性：假定存在非奇异下三角阵G

L,其对角元素皆为正数，且使得

A

LLT

GGT

于是有LT

(GT

)1

L1LLT

(GT

)1

L1GGT

(GT

)1

L1G因LT

(GT

)1为上三角阵，L1G为下三角阵，故由上式得LT

(GT

)1

L1G

I即G

L,与假设

。平方根(Cholesky分解法）法2)2jj

jjjkl2122lllll0

llj

1k

10

0

l11ln1

l11l

0

21

22

0

n2

l

0

n1

n

2nn

nn

0

(i

1,

2,(

j

1,

2,

,

n),

(aij由A

LLT其中lii,n).由矩阵乘法及l

jk

0(当j

k时),1得0(i

j

1,

,

n);j

1k

1l

(a

lij

likljk

)

/

l

jj这里规定

0。计算顺序是按列进行，即k

1l11

li1

(i

2,3,

,n)

l22

li

2

(i

3,

,n)

。i1k

i1(i

1,

2, ,

n).(i

n,

n

1,

,1).yi

(bi

lik

yk

)

/

liik

1nxi

(

yi

lki

xk

)

/

lii当矩阵A完成平方根分解后，求解Ax

b，即求解两个三角形方程组(1)Ly

b,

求y;

（2）LT

x

y，求x.由于A的对称性，平方根法的乘除运算量为n3/6数量级，约是Gauss消去法的一半。上机计算时，所需单元也少，只要A的下三角部分和右端项b，计算中L存放在A的单元，y,x

在b的单元.但这种方法在求L时需作n次开方运算，这样又增加了计算量，为了避免开方，可使用改进的平方根方法.改进平方根法21dll

1lA

LDLT

21

d1

l1

d

n

ln1

n1

n

21

l21

1n

2

1

其中ljj

1,

l

jk

0

(

j

k

),由比较法得2i

iiik

kiiiiik

ikjil

d

;t

l

;,

n).i1k

1i1k

1i1k

1d

a

l

ji

(a

ji

)

/

d

(

j

i

1,

,

n).d

a

l

(a

)

/

d

(

j

i

1,ljk

dklik对于i

1,

2,

,

n,上式虽避免了开方运算，但增加了相乘因子，引进变量

对于i

1,

2,

,

n,

有i1tik

lik

dkji

tikljk

ik

1,

2,1.i

i

ik

kni

i

i

ki

ki1k

1k

i1i

2,

,

n;xn

yn

/

dn

;i

n

1,

y

b

l

y

,

x

y

/

d

l

x

,对称正定矩阵A按LDLT

分解和按LLT分解计算量差不多，但LDLT

分解不需要开方计算。求解Ly

b,

DLT

x

y计算公式

y1

b1;三、追赶法2222

2iiai

bi

ci

b1

c1

d1b

c

x1

a

x

d

x

d

an1

bn1

cn1

xn1

adn1

d

b

x

n

n

n

n

在数值计算中，如三次样条插值或用差分方法解常微分方程边值问题，常常会遇到求解以下形式的方程组简记Ax

d.此系数矩阵的非零元素集中分布在主对角线及其相邻两次对角线上，称为三对角矩阵。方程组称为三对角方程组。22i

i

in

nc1u

ccln

un

b1

c1

0

b

a

c

,

b

a

0l

21A

LU

l3

1n1

u1

1

aici

0(i

2,3,

,n

1)对角占优定理：设三对角方程组系数矩阵满足下列条件：则它可分解为

1其中ci

(i

1,

2, ,

n

1)为已给出的，且分解是唯一的定理证明（1）ii

i

i1i i1

i

ii

i

i1i i1

ia

l

u

(i

2,

3,

,

m)b

c

l

uu1

b1

l

a

/

u(i

2,

3,

,

m)u

b

c

l将上式右端按乘

则展开,并与A进行比较,得b1

u1如果ui

0

(i

1，2,,m)，则由上式可得定理证明（2）1

1(k

)nu

can1

bn1

cn1

a

bA

n

ck

1ak

/

uk

1

(k

2,

3,按Gause消去法步骤易得，经k

1次消元后，三对角方程的系数矩阵变为其中uk

bk,n)。kka

bcan

bnck

1b

cn1

n1n1

uk

ck

ak

1

bk

1

ck

1Ob1

c1O

O0

uk

1akO

OO定理证明（3）21

1

1

1u

b

bb

c1a2

b1b2

c1a2

b1b2

c1a2

b1c2u2

b2

c由于A满足定理所给条件，显然有

u1

b1

0.又因为

b1

c1

,

b2

a2

c2

,

于是b1b2

b1a2

b1c2

c1a2

b1c2故u2

0且矩阵A

仍满足定理条件。依此类推可得出(2)ui

0(i

1,2,L

,n)。因此由上面公式唯一确定了L和U。1

1n

n

b

c

0

ci

bi

ai

b

a

0aici

0(i

2,

3,,

n

1)从而有uk

bk

ck

1ak

/

uk

1

(k

2,

3, ,

n)追赶法的计算公式,

m);,

n);i

i1k

k

k

k

1k

k

k

1

ku1

b1

iu

b

c

l

i

i i1

i

y1

d1

y

d

l

y(k

2,

3,x

(

y

c

x

)

/

u(k

n

1,

n

2,

,1).

kA

LU分解公式：

l

a

/

u

(i

2,

3,解Ly

d得：再解Ux

y得:

xn

yn

/

un追赶法的基本思想与Gauss消去法及三角分解法相同,只是由于系数中出现了大量的零,可使计算公式简化,减少了计算量。可证,当系数矩阵为严格对角占优时,此方法具有良好的数值稳定性。§5.向量和矩阵的范数一、向量范数向量范数定义:设对任意向量x

Rn

,

按一定的规则有一实数与之对应,记为

x

,

若

x

满足1,x

0

,且x0当且仅当

x

0;

(正定)2,x

x,为任意实数

(齐次)3,x

yxy,

对任意x,y

Rn(三角不等式)则称

x

为

向量x的范数.向量范数例2122nix

)i1x

x2

x2

(1

n1nx

x1

xn

xii

11inx

max{

x1

,,

xn

}

max{

xi

}npip1/

px

i

1

x1in

1in1in

1in

max

xi

max

yi

x

y

可验证上面范数均满足范数定义的条件。以-范数为例:满足条件1,2显然。由于x,

y

Rn为向量，而其分量x

,

y

(i

1,

,

n)i

i为实数，故有x

y

max

xi

yi

max

xi

yi

例：计算向量x

(1,2,3)T

的各种范数。解：

x

6,

x

3,

x

14.1

2如果Rn中两个范数

和

'

，存在实数m,M

0，使得对任意n维向量x,都有m

x

x

'

M

x

，则称这两个范数是等价的。对两个等价范数而言，同一向量序列有相同的极限。21

2i

1x2xx2nnxn1in1in

2

x2

x2

x

.2

n

max

xi

.

x2

x22

n

2

x

x不难证明，1－范数，2－范数和－范数是等价的。例：

x

max

x

设

x

x

j则

x

x

j

2－范数和－范数等价。如不作说明，今后

是指任意一种向量范数。二、矩阵的范数定义：对任意n阶方阵A，按一定的规则由一实数与之对应，记为

A

。若

A

满足1,

A

0

,

且

A

0当且仅当

A

0;

(正定)2,

A

A

,

为任意实数

(齐次)3,4，A

B

A

B

,

对任意A,B两个n阶方阵(三角)AB

A

B

（相容性条件）则称A

为矩阵A的范数。xxxxxx

1Ax

A

xx

xA

max

Ax

max

Ax定理：设A为n阶方阵，

是Rn中的向量范数，则A

max

Ax

是一种矩阵范数，称其为由向量范数

诱导出的矩阵范数。证：设A

(aij

)为任意n阶方阵，x为任意n维非零向量。因为x

为范数是1的单位向量，故x

xxx

y

1

x

1max

Ax

max

A

x

max

Ay

max

Axxxxx

1x

1Ax

0

Ax

A

0.

max

Ax

A

.x

11，显然

A

0.若A

0,则

A

max

Ax

0.反之，若

A

0

2，

A

max

Ax

max

Axx

1

x

13,对任意两个n阶方阵A和B，A

B

max (

A

B)x

max

Ax

Bxx

1

x

1

max(

Ax

Bx

)

max

Ax

max

Bxx

1

x

1

x

1

A

B

.xxx

1A

max

Ax

max

AxAxx4，对任意n维非零向量x,有

A

即

Ax

A

x

.故有

AB

max

(

AB)x

max

A(Bx)x

1

x

1

max

A

Bx

max

A

B

xx

1

x

1

A

B5，对任意n维向量x，都有

Ax

A

x

。这一性质称为矩阵范数与向量范数的相容性。可由三种常用的向量范数诱导出矩阵范数。xxA

max

Ax21

12

2Tx

1A

max

Ax

其中

是

A

A,

的最大特征值。又称为谱范数。设A

(aij

)为n阶方阵。矩阵范数例与前述三种向量范数相容的三种矩阵范数：1

1nx

1

1i1A

max

Ax

max1

jnaij

,

为矩阵的列向量的1－范数的最大值称为矩阵的列和范数。1nj

1A

max

Axx

1

max1inaij

,

为矩阵的行向量的1－范数的最大值称为矩阵的行和范数2

1

2如果将矩阵范数看作Rn

空间上的向量范数，则由向量范数的等价性可得矩阵范数的等价性。例：A

34

,计算A的各种范数。

解：

A

6,

A

7,

A

5.46.1

2矩阵的误差可用矩阵范数表示：设A*是A的近似矩阵，A

A*

、A

A*

/

A

分别称为A*的关于范数

的绝对误差与相对误差。矩阵A的谱半径i1in(i

1,

2,

(

A)

max

i定义：设A

Rnn的特征值为,n)

称(

Ax

x

x

,

Ax

A

x

x

A

x

A

(

A)

A

)定理：

(

A)

A为A的谱半径。,

A

为

A

的任意矩阵范数.§5.误差分析一、矩阵的条件数x

1.00001x

2

x

0

1

2

2

x1

1.x

1.00001x

2.00001x

1

1

2

2一个实际问题化为数学问题，初始数据往往会有误差，即有扰动，从而使计算结果产生误差。例：方程组

x1

x2

2

x1

2

.而方程组

x1

x2

2比较这两个方程组可以看出，他们只是右端项有微小的差别，最大相对误差为1

105

,但它们的解却大不相同，解分量21的相对误差至少为。2定义：如果矩阵A或常数项b的微小变化，引起方程组Ax

b解的巨大变化，则称此方程组为“”方程组，矩阵A称为“”矩阵（相对于方程组而言）。否则称方程组为“良态”方程组，A称为“良态”矩阵。矩阵的“

”性质是矩阵本身的特性。为了定量刻划方程组的“

”程度，下面对方程组Ax

b就系数矩阵或右端项分别有扰动的两种情形进行

。右端项b的扰动对解的影响Ax

bA

AAbxb

b

b

xA1A1A1A1A

x

b

x

A1

b

x

A1

b

b设b有扰动

b，相应解x的扰动记为

x,即A(x

x)

b

b由Ax

b，两边取范数又因为

x

A此式表明,当右端项有扰动时,解的相对误差不超过右端项的相对误差的

A

倍。系数矩阵A的扰动对解的影响AAxA

A

x

A

A

x

AA1A1A1A11

x

A1

A(x

x)

A

(

x

x

)1

A1

A1

A

A如果右端项无扰动，系数矩阵A有扰动

A，相应的解x的扰动仍记为

x,则(

A

A)(x

x)

b

A

x

A(x

x)

0如果

A充分小，使得

A1

x

(1

A1

A

)

A

1,则由上式得A1A1上式表明，当系数矩阵有扰动时，解的扰动仍与A有关。一般地，A

越大，解的扰动也越大。条件数的定义max222min.vA-1

(

AT

A)

(

AAT

)A1A1(1)

cond

(

A)

cond

(

A)

A

综上分析可知,量

A-1 A实际上刻划了解对原始数据变化的灵敏程度,即刻划了方程组的“

”程度。定义：设A为非奇异阵，称数cond

(A)v

A

(v

1,2或v）为矩阵A的条件数。常用的条件数，有(2)

A的谱条件数2,n

1当A为对称矩阵时，cond

(A)其中1，n为A的绝对值最大和绝对值最小的特征值。Au

u

u

A1u

A1u

1

uA

条件数的性质A

A1

A

I

1.v

v1、对任何非奇异矩阵A，都有cond

(A)v

1.1由定义

cond

(

A)v

Av2、设A为非奇异矩阵且c

（0

常数），则cond

(cA)v

cond

(

A)v3、如果A为正交矩阵，则cond

(A)2＝1；

如果A为非奇异矩阵，R为正交矩阵,则cond

(RA)2

cond

(

AR)2

cond

(

A)2

.max22min.

(

AT

A)

(

AT

A)1cond

(

A)2

A

A

1121111n

1

1

3Hn

2

n

1

1

n

1

n2n

1例：Hilbert矩阵计算H3的条件数。n(1)

H3(

)同

可计算一般阵当0.500

0.3333

0.

5简记为(H3

H其解为x

xH3

b

x

(0.0895,

0.5120,

0.4910)T

H3

0.18103

0.06%

b

0.182%,

x

0.5120

51.2%.x

1由于x

x

(1.0895,0.4880,1.491)T

,x

(1,1,1)T这表明H3与b相对误差不超过0.2%，而引起解的相对误差超过50％.”方程的经验判断。根”的。计算条件数需要求矩阵的逆，因而比较据数值经验，在下列情况下，方程组常是“（1）在用主元素法时出现小主元；如果A的最大特征值和最小特征值之比（