选主元的三角分解法

2012-2013（1）专业课程实践论文选主元（列）的三角分解法
范俊，0818180124，R数学08-1班
毛勇，0818180117，R数学08-1班
一、算法理论
从直接三角分解公式可看出，当时计算将中断或者当绝对值很小时，按分解公式计算可能要求舍入误差的累积，但如果非奇异，就可通过交换的行实现矩阵的分解，因此可采用与列主元素消去法类似的方法
（可以证明下述方法与列主元素消去法等价），将直接三角分解法修改为（部分）选主元的三角分解法。

设第步分解已完成，这时有
第步分解需要用到式
及式
为了避免用小的数作除数，引进量
于是有
,,
用作为，交换的行与行元素(将位置的新元素仍记作及
)，于是有。

由此在进行第步分解计算。

该程序即利用选主元的三角分解法解方程求方程的根。

先选择列主元，再构造LU矩阵，通过求解LY = B和UX = Y得出需要的解。

注：方程维数在程序中需按题意自定义。

否
源代码源代码：
LU_Decomposition.cpp
#include<iostream.h>
#include<math.h>
#define N 4 //矩阵维数，可自定义
static double A[N][N]; //系数矩阵
static double B[N]; //右端项
static double Y[N]; //中间项
static double X[N]; //输出
static double S[N]; //选取列主元的比较器
int i,j,k; //计数器
void main()
{
cout << "请输入线性方程组（ai1，ai2，ai3......ain, yi）："<<endl;
for ( i = 0; i < N ;i++)
{
for (int j = 0; j< N ;j++ )
cin >> A[i][j];
cin >>B[i];
}
for (k = 0 ; k < N; k++)
{
//选列主元
int index =k;
for (i = k ; i < N; i++)
{
double temp = 0;
for (int m = 0 ; m < k ;m++)
{
temp = temp + A[i][m]*A[m][k];
}
S[i] = A[i][k] - temp;
if(S[index] < S[i])
{
index = i;
}
} //交换行
double temp;
for( i = k ; i < N ;i++ )
{
temp = A[index][i];
A[index][i] = A[k][i];
A[k][i] = temp;
}
temp = B[index];
B[index] = B[k];
B[k] = temp; // 构造L、U矩阵
for (j = k ; j < N; j++)
{
double temp = 0;
for (int m = 0 ;m < k; m++ )
{
temp = temp + A[k][m] * A[m][j];
}
A[k][j] = A[k][j] - temp; //先构造U一行的向量}
for( i = k+1; i < N; i++)
{
double temp = 0;
for (int m =0 ; m < k; m++ )
{
temp = temp + A[i][m] * A[m][k];
}
A[i][k] = (A[i][k] - temp)/A[k][k]; //再构造L一列的向量}
}
//求解L Y = B
Y[0] = B[0];
for (i = 1; i < N ; i++)
{
double temp = 0;
for (int j =0 ; j < i; j++ )
{
temp = temp + A[i][j] * Y[j];
}
Y[i] = B[i] - temp;
}
//求解UX = Y
X[N-1] = Y[N-1]/A[N-1][N-1];
for (i = N-2 ; i >= 0 ; i-- )
{
double temp = 0;
for (int j =i+1 ; j < N; j++ )
{
temp = temp + A[i][j] * X[j];
}
X[i] = (Y[i] - temp)/A[i][i];
}
//打印X
cout << "线性方程组的解（X1，X2，X3......Xn）为："<<endl; for( i = 0; i < N ;i++)
{
cout << X[i] <<" ";
}
}
四、算法实现
例1. 用选主元（列）的三角分解法解
解：利用软件的解答方法为：
（1）先将程序中的维数N改为3。

（2）输入1 2 3 14回车，输入2 5 2 18 回车，输入3 1 5 20回车。

（3）显示结果 1 2 3
例2. 用选主元（列）的三角分解法解
解：利用软件的解答方法为：
（1）先将程序中的维数N改为3。

（2）输入0 3 4 1回车，输入1 -1 1 2回车，输入2 1 2 3回车。

（3）显示结果 1.27273 -1.18182 0.818182。