数值分析——三角分解法（LU分解法）C++

从例题着手,求解方程组：

10x1+20x2+20x3+50x4=200

50x1+40x2+10x3+10x4=250

30x1+10x2+40x3+20x4=210

20x1+40x2+40x3+40x4=340

根据方程组可以得到系数矩阵Aij[4][4], 通过矩阵A的元素进行LU分解；由Ax=b，从题干得到b[4]矩阵。现在已知了A、b矩阵，再通过LU分解法的递推公式计算出L U矩阵各个元素。

注意：（选用了最笨的办法）通过LU分解法的递推过程我们知道，要把系数矩阵A分解为L和U，需要先计算出U的第一行和L的第一列，接着计算出U的第二行L的第二列，继续U的第三行L的第三列......以此类推。整个递推的过程不是先计算出U的整个矩阵元素，再计算L的；两者是穿插进行的，直到所有行和列都陆续求出，得到最后完整的L[4][4], U[4][4]。

公式如下：

(1）U第一行，L第一列：u1i=a1i (i=0,1,2...n), li1=ai1/u11

(2)计算U的第r行，L的第r列元素：

**Uri = Ari - $\sum_{k=1}^{r-1}$ Lrk*Uki, (i=r,r+1,...n);**

**Lir = ( Air - $\sum_{k=1}^{r-1}$ Lik*Ukr)/Urr, i=r+1,....n 且r != n；**

(3)求解Ly=b， Ux=y 的计算公式：

y1=b1;

yi = bi - $\sum_{k=1}^{i-1}$ Lik*yk ,i=2,3,4...n;

xn=yn/Unn;

xi = ( yi - $\sum_{k=i+1}^{n}$ Uik*xk)/Uii , i=n-1,n-2,...1;

(根据上面公式编写代码)例题代码如下：


double Ai[4][4] = { 10,20,20,50,50,40,10,10,30,10,40,20,20,40,40,40 };  //Ai系数矩阵
double L[4][4] = { 0 }; //L初始化
double U[4][4] = { 0 }; //U初始化
double bi[4] = { 200,250,210,340 };  //b矩阵
double y[4] = { 0 };  //存放Ly=b 的结果y
double x[4] = { 0 }; //存放Ux=y 的最终结果xvoid U1i_and_Li1( ) //计算U的第一行U1i，L的第一列Li1
{for (int i = 1; i <= 4; i++){U[0][i-1] = Ai[0][i-1];L[i-1][0] = Ai[i-1][0] / U[0][0];}
}
//根据LU分解法递推公式分别计算U的某一行，L的某一列
void Uri(int m )  //形参m表示计算U的某一行
{double temp = 0.0;for (int r = 2; r <=m; r++){for (int j = r; j <= 4; j++){for (int k = 1; k <= r - 1; k++){temp = temp + L[r - 1][k - 1] * U[k - 1][j - 1];}U[r-1][j-1] = Ai[r-1][j-1] - temp;temp = 0.0;}}
}
void Lir( int m) //m：L的某一列
{double temp = 0.0;for (int r = 2; r <= m; r++){for (int i = r; i <= 4;i++){for (int k = 1; k <= r - 1; k++){temp = temp + L[i-1][k-1] * U[k-1][r-1];}L[i-1][r-1] = (Ai[i-1][r-1] - temp) / U[r-1][r-1];temp = 0.0;}}
}
void print(double a[4][4])  //打印结果
{for (int i = 0; i < 4; i++){for (int j = 0; j < 4; j++){cout << a[i][j] << " ";if ((i == 0 || i == 1 || i == 2 || i == 3) && j == 3)cout << endl;}}
}
void yi()  //计算yi
{double temp = 0;for (int i = 1; i <=4; i++){if (i == 1)y[0] = bi[0];else {for (int j = 1; j <=i - 1; j++){temp = temp + L[i-1][j-1] * y[j-1];}y[i - 1] = bi[i - 1] - temp;temp = 0;}           }
}
void xi()  //计算最终结果xi
{double temp = 0;x[3] = y[3] / U[3][3];for (int i = 3; i >=1; i--){for (int k = i + 1; k <=4; k++){temp = temp + U[i - 1][k - 1] * x[k - 1];}x[i - 1] = (y[i - 1] - temp) / U[i - 1][i - 1];temp = 0;}
}

main函数和结果

int main()
{U1i_and_Li1();cout << "U 的第一行为：" << U[0][0] << " " << U[0][1] << " " << U[0][2] << " " << U[0][3] << endl;cout << "L 的第一列为：" << L[0][0] << " " << L[1][0] << " " << L[2][0] << " " << L[3][0] << endl;Uri(2);Lir(2);cout << "U 的第二行为：" << U[1][0] << " " << U[1][1] << " " << U[1][2] << " " << U[1][3] << endl;cout << "L 的第二列为：" << L[0][1] << " " << L[1][1] << " " << L[2][1] << " " << L[3][1] << endl;//l12,l22,l32,l42Uri(3);Lir(3);cout << "U 的第三行为：" << U[2][0] << " " << U[2][1] << " " << U[2][2] << " " << U[2][3] << endl;cout << "L 的第三列为：" << L[0][2] << " " << L[1][2] << " " << L[2][2] << " " << L[3][2] << endl;Uri(4);Lir(4);cout << "U 的第四行为：" << U[3][0] << " " << U[3][1] << " " << U[3][2] << " " << U[3][3] << endl;cout << "L 的第四列为：" << L[0][3] << " " << L[1][3] << " " << L[2][3] << " " << L[3][3] << endl;cout << "U矩阵表示为：" << endl;print(U);cout << "L矩阵表示为：" << endl;print(L);yi();cout << "由Ly=bi计算出y[4]为：" << y[0] << " " << y[1] << " " << y[2] << " " << y[3] << endl;xi();cout<< "由Ux=y计算出x[4]为：" << x[0] << " " << x[1] << " " << x[2] << " " << x[3] << endl;return 0;
}

数值分析——三角分解法（LU分解法）C++相关推荐

Python02 雅克比迭代法 Gauss-Seidel迭代法列选主元法 LU分解法(附代码)
1. 实验结果 (1)在定义的矩阵类中设置需要求解的方程为: (2)在 test.py 中选择雅克比迭代法求解: 输入:最大容许迭代次数和精度要求: 输出:根据谱半径判断方法是否收敛,收敛时得到满足精 ...
数值分析3-解线性方程组的高斯消去法、LU分解法及列主元消去法的matlab程序和调试方法
对于形如Ax=b的线性方程组,在线性代数中是通过求逆的方式求解的,即x=A-1b,而在数值分析中,解线性方程组的方法是通过直接法或者迭代法来实现的,今天写的三个程序都属于直接法,分别为高斯消去法.LU ...
紧凑存储的杜利特尔分解法Doolittle(LU分解法)_解线性方程组的直接解法
紧凑存储的杜利特尔分解法Doolittle(LU分解法)_解线性方程组的直接解法标签:计算方法实验 /* 紧凑存储的杜利特尔分解法Doolittle:如果初始矩阵不要求保留的话,可以紧凑存储.因为每 ...
Doolittle分解法（LU分解法）的Python实现
在解一般的非奇异矩阵线性方程组的时候,或者在迭代改善算法中,需要使用LU分解法. 对于一个一般的非奇异矩阵A=(a11, a12,-,a1n,a21,-ann),可分解为一个下三角矩阵L和一个上三角矩 ...
计算方法：列主元消去法，LU分解法，雅可比迭代法，高斯塞德尔迭代法解线性方程（C++）
Matrix.h包括矩阵类Matrix的定义,Matrix.cpp包括该类成员函数的实现,LinearEqu.h包括线性方程类LinearEqu的定义,继承自Matrix类,其中solve()方法为列 ...
LU分解法 | matlab
% LU分解法 % M为输入的增广矩阵 % precision为输入的精度要求,如不输入或输入有误,则默认为10位if nargin == 2trydigits(precision);catchdis ...
计算机数值方法——LU分解法（C++\Python代码实现）
算法流程首先,LU分解法就是在高斯消元法的基础上,把矩阵AAA分解为一个上三角矩阵UUU与一个单位下三角矩阵LLL的乘积. 懒得敲LaTeX公式了,书上由具体的推导过程,这里我们重点介绍代码吧主要 ...
解线性方程组的直接方法:LU分解法及其C语言算法实现
在上一篇博客里面,笔者介绍了解线性方程组的列主元Guass消元法,这篇将介绍LU分解法及其算法实现. 什么是LU分解? 对于一个线性方程组Ax=b,其中A是非奇异系数矩阵,b是线性方程组右端项,在列主 ...
Matlab | Lab4——用LU 分解法、 Jacobi 迭代、 Gauss-Seidel 迭代解线性病态方程组（系数矩阵为Hilbert矩阵）
1.要求考虑线性方程组Hx=b,其中H为n阶Hilbert矩阵,即通过先给定解(例如取x的各个分量为1),再计算出右端向量b的办法给出一个精确解已知的问题． (1)分别编写Doolittle LU ...
matlab lu解线性方程,MATLAB报告用LU分解法求解线性方程组.doc
MATLAB报告用LU分解法求解线性方程组 <MATLAB>期末大报告报告内容:用LU分解法求解线性方程组学院(系): 专业: 班级: 学号: 学生姓名: 2014 年 6 月 ...

数值分析——三角分解法（LU分解法）C++

(1）U第一行，L第一列：u1i=a1i (i=0,1,2...n), li1=ai1/u11

**Uri = Ari - $\sum_{k=1}^{r-1}$ Lrk*Uki, (i=r,r+1,...n);**

**Lir = ( Air - $\sum_{k=1}^{r-1}$ Lik*Ukr)/Urr, i=r+1,....n 且r != n；**

y1=b1;

yi = bi - $\sum_{k=1}^{i-1}$ Lik*yk ,i=2,3,4...n;

xn=yn/Unn;

xi = ( yi - $\sum_{k=i+1}^{n}$ Uik*xk)/Uii , i=n-1,n-2,...1;

数值分析——三角分解法（LU分解法）C++相关推荐

最新文章

热门文章

数值分析——三角分解法（LU分解法）C++

(1）U第一行，L第一列：u1i=a1i (i=0,1,2...n), li1=ai1/u11

Uri = Ari - Lrk*Uki, (i=r,r+1,...n);

Lir = ( Air - Lik*Ukr)/Urr, i=r+1,....n 且r != n；

y1=b1;

yi = bi - Lik*yk ,i=2,3,4...n;

xn=yn/Unn;

xi = ( yi - Uik*xk)/Uii , i=n-1,n-2,...1;

数值分析——三角分解法（LU分解法）C++相关推荐

最新文章

热门文章

**Uri = Ari - $\sum_{k=1}^{r-1}$ Lrk*Uki, (i=r,r+1,...n);**

**Lir = ( Air - $\sum_{k=1}^{r-1}$ Lik*Ukr)/Urr, i=r+1,....n 且r != n；**

yi = bi - $\sum_{k=1}^{i-1}$ Lik*yk ,i=2,3,4...n;

xi = ( yi - $\sum_{k=i+1}^{n}$ Uik*xk)/Uii , i=n-1,n-2,...1;