矩阵doolittle分解算法

解线性方程组的Doolittle 分解

目的意义：

1.学习和掌握线性代数方程组的Doolittle 分解法。

2.运用Doolittle 分解法进行计算。

方法原理：

n 阶线性方程组的系数矩阵A 非奇异且有分解式A=LR ，其中L 为单位下三角矩阵，R 为上三角矩阵，即L=(l ij ),当ij 时，r ij =0，矩阵A 的这种分解方法为Doolittle 的分解。

⎥⎥⎥⎥⎦

⎤⎢⎢⎢⎢⎣

⎡⎥⎥⎥⎥⎦⎤⎢

⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡nn n n n n nn n n n n r r r r r r l l l a a a a a a a a a M O

Λ

Λ

ΛO M M Λ

ΛΛΛΛΛΛ

222112

11212121

22221112

11111

比较等号两边的第i 行和第j 列的元素，可知∑==

n

k kj

ik ij r

l

a 1

，因为

0,11,======++ij j j in i i r r l l ΛΛ，

所以∑==

n

k kj

ik ij r l

a 1

=ij kj i k ik r r l +∑-=1

1

，i<=j ，从而

.

,1,,1

1n i i j r l a r i k kj ik ij ij Λ+=-=∑-=当

n

i i j Λ,2,1++=时，

ii ji i

k i k kj jk kj jk ij r l r l r l a +==∑∑=-=1

1

1

，从而n i j r r l a l ii i k ki jk ji ji ,,1,/)(1

1

Λ+=-=∑-=，于是就得

到了计算LR 分解的一般计算公式。

算法描述：

Setp1：利用for 循环求出1

1

k ki

ki kp pi p r a l r -==-∑，1

1

()/k ik ik ip pk

kk p l a l

r r -==-

∑。

Step2：

1

1

i i i ik k k y b l y -==-∑,得出1

()/n

i i ik

k ii k i x y r

x r =+=-

∑

程序代码：

头文件：

#include typedef double Datatype;

class Matrix

{

private:

Datatype **ar;

int M;

int N;

public:

Matrix(int a=0,int b=0);

~Matrix();

void print();

void Init();

void Doolittle(Matrix &L,Matrix &R);

friend void Init_LR(Matrix &L,Matrix &R);

void get_Y(Matrix &b,Matrix &L);

void get_X(Matrix &Y,Matrix &R);

};

CPP文件：

#include

#include

#include

#include"Doolittle.h"

void Init_LR(Matrix &L,Matrix &R)

{

for(int i=0;i

{

for(int j=0;j

{

if(i<=j)

{

if(i!=j)

{

L.ar[i][j]=0;

}

else

{

L.ar[i][j]=1;

}

}

else

{

R.ar[i][j]=0;

}

}

}

}

Matrix::Matrix(int a,int b)///////////////////构造函数

{

M=a;

N=b;

ar=(Datatype **)malloc(sizeof(Datatype *)*M);

for(int i=0;i

{

ar[i]=(Datatype *)malloc(sizeof(Datatype)*N);

}

}

Matrix::~Matrix()////////////////////////////析构函数

{

for(int i=0;i

{

free(ar[i]);

}

free(ar);

}

void Matrix::print()////////////////////////打印函数

{

for(int i=0;i

{

for(int j=0;j

{

cout<

}

cout<

}

}

void Matrix::Init()/////////////////////////初始化函数

{

for(int i=0;i

{

cout<<"input the data of "<

for(int j=0;j

{

cin>>ar[i][j];

}

}