数值分析_lec4_解线性方程组_LU分解法_1013
- 格式:pdf
- 大小:924.46 KB
- 文档页数:25
lu分解迭代求精法
LU分解是一种常用的矩阵分解方法,用于求解线性方程组或矩阵的逆。
LU分解迭代求精法是指在已经进行了LU分解的基础上,通过迭代逐步求取精确解的方法。
下面以迭代求解线性方程组为例来说明LU分解迭代求精法的步骤:
1. 对系数矩阵A进行LU分解,得到下三角矩阵L和上三角矩阵U。
2. 将线性方程组Ax=b转化为LUx=b。
3. 设定初始解向量x^(0)的初始值。
4. 迭代计算x^(k+1) = U^(-1)(b - Lx^(k)),其中k为迭代次数。
5. 判断精度是否满足要求,如果满足则停止迭代,否则返回第4步。
迭代求解的精度取决于迭代次数和初始解向量的选择,一般可以使用收敛准则来判断是否满足要求。
常用的收敛准则包括残差准则和误差准则。
需要注意的是,LU分解迭代求精法虽然可以提高解的精度,但在某些情况下可能会导致迭代过程发散或慢速收敛。
因此,在实际应用中需要根据具体问题选择合适的迭代方法和初始解向量,并进行收敛性分析以保证求解结果的可靠性和精度。
数值分析实验报告(1)学院:信息学院班级:计算机0903班姓名:***学号:********课题一A.问题提出给定下列几个不同类型的线性方程组,请用适当的方法求解线性方程组1、设线性方程组⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--------------------------1368243810041202913726422123417911101610352431205362177586832337616244911315120130123122400105635680000121324⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡10987654321x x x x x x x x x x =⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡-2119381346323125 x *= ( 1, -1, 0, 1, 2, 0, 3, 1, -1, 2 )T2、设对称正定阵系数阵线方程组⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡----------------------19243360021411035204111443343104221812334161206538114140231212200420424⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡87654321x x x x x x x x = ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡---4515229232060 x * = ( 1, -1, 0, 2, 1, -1, 0, 2 )T3、三对角形线性方程组⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡------------------4100000000141000000001410000000014100000000141000000001410000000014100000000141000000001410000000014⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡10987654321x x x x x x x x x x = ⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡----5541412621357 x *= ( 2, 1, -3, 0, 1, -2, 3, 0, 1, -1 )TB.(1)对上述三个方程组分别用Gauss 顺序消去法与Gauss 列主元消去法;平方根 与改进平方根法;追赶法求解(选择其一) (2)编写算法通用程序(3)在应用Gauss 消去时,尽可能利用相应程序输出系数矩阵的三角分解式C.(1)通过该课题的程序编制,掌握模块化结构程序设计方法 (2)掌握求解各类线性方程组的直接方法,了解各种方法的特点 (3)体会高斯消去法选主元的必要性 实验步骤:(高斯消去法,列主元,LU )1顺序高斯消去法2.LU 分解法3.列主元高斯消去法(如下图)(1)高斯消去法运行结果如下(2)对方程的系数矩阵进行LU分解并求出方程组的解(3)列主元高斯消去法实验体会总结:利用gauss消去法解线性方程组的时候,如果没有经过选主元,可能会出现数值不稳定的现象,使得方程组的解偏离精确解。
LU 分解法解线性方程组一. 实验目的1.熟练运用已学计算方法求解方程组2.加深对计算方法技巧,选择正确的计算方法来求解各种方程组3.培养使用电子计算机进行科学计算和解决问题的能力二.实验环境VC++6.0 语言:C++三.实验内容编写一个解线性方程组的L U -算法。
算法包括:1. 矩阵的输入输出2. 判断矩阵是否奇异3. 进行分解4. 输出L 、U5. 回带解出y 、x四.实验原理1.若一个线性方程组系数矩阵为n 阶方阵A 且各阶顺序主子式均不为零,则A 的L U -分解存在且唯一。
2.在满足1的条件下可推导得出以下公式 (1):11()j ij ij ik ij ijk l a l u u -==-∑11j j i j i j kk ik u a lu -==-∑(2): 11i i i ik kk y b l y -==-∑1()ni i ik k iik i x y u x u =+=-∑3.公式(1)用于求解矩阵L 、U ,公式(2)用于回带求解y 、x 。
从公式中可以看出:L 对角线上元素为1,U 第一行与A 第一行相同。
4.L U-分解的具体过程和计算顺序如下:(1)第一步分解:1111u a =(2)第二步分解:212111l a u =1212u a=22222112u a l u =-(3)第三步分解:313111l a =3232311222()l a l u =-1313u a = 23232113u a l u =-333331133223u a l u l u =--(n).第n 步分解:依次计算:1.1..n n n in nnl l u u -⎧⎪⎪⎪⎪⎨⎪⎪⎪⎪⎩1112121222121.1.......1n n n n nn u u u l uu l l u ⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⨯⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦五.实验流程图1.程序总的流程图2.考虑主要目的这里只给出L U 分解的流程图:六.VC++源程序#include<iostream>using namespace std;double Fun(int n1,double a1[11][11]); //声明求矩阵行列式的递归函数函数double Fun(int n1,double a1[11][11]){int i_1,j_1,c; //c为数组b的行double b[11][11]; //用于存放余子式int p=0,q=0;int sum=0;if(n1==1)return a1[1][1];for(i_1=1;i_1<=n1;i_1++){for(c=1;c<=n1-1;c++){if(c<i_1)p=0;elsep=1;for(j_1=1;j_1<=n1-1;j_1++){b[c][j_1]=a1[c+p][j_1+1];}}if((i_1-1)%2==0)q=1;elseq=(-1);sum=sum+a1[i_1][1]*q*Fun(n1-1,b);}return sum;}//主程序,用于矩阵输入输出,判断,分解,回带void main(){int n,i,j,k;double temp=1.0;double a[11][11]={0};double L[11][11]={0};double b[11]={1};double x[11]={1};double m=0.0;double y[11]={1};double temp1=0.0;double temp2=0.0;//完成系数矩阵与右端项的输入与输出cout<<"!!!请输入待解方程组未知数的个数: n "<<endl; cin>>n;cout<<"!!!请输入待解方程组的系数"<<endl;for(i=1;i<=n;i++){for(int j=1;j<=n;j++){cin>>a[i][j];}}cout<<"!!!请输入待解方程组的系数矩阵如下:"<<endl;for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<a[i][j]<<" ";}cout<<endl;}cout<<"!!!请输入待解方程组的右端项"<<endl;//初始化LUfor(int p=1;p<=n;p++){L[p][p]=1;}for(p=1;p<=n;p++){}cout<<endl;for(i=1;i<=n;i++){cin>>b[i];}cout<<"b="<<endl;for(i=1;i<=n;i++){cout<<b[i]<<endl;}for(i=2;i<=n;i++){temp=Fun(i,a);if(temp==0){cout<<endl<<endl;cout<<endl<<"output:temp ="<<temp<<endl;cout<<"矩阵奇异";break;}}cout<<endl;//LU分解if(temp!=0){for(i=1;i<=n;i++){for(j=1;j<=i-1;j++){if(a[j][j]==0)cout<<"请重新输入系数矩阵与右端项:"<<endl;elsem=0.0;for(k=1;k<=j-1;k++){m=m+L[i][k]*U[k][j];}L[i][j]=(a[i][j]-m)/U[j][j];}for(j=1;j<=i;j++){m=0.0;for(k=1;k<=j-1;k++){m=m+L[j][k]*U[k][i];}U[j][i] = a[j][i] - m;}}}if(temp!=0){cout<<"L="<<endl;for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<L[i][j]<<" ";}cout<<endl<<endl;}cout<<endl<<endl<<"U="<<endl; for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<U[i][j]<<" ";}cout<<endl<<endl;}//回带求解y[]cout<<"y="<<endl;for(int s=1;s<=n;s++){for(int t=1;t<=s-1;t++){temp1=temp1+L[s][t]*y[t];}y[s]=b[s]-temp1;cout<<" "<<y[s]<<" "<<endl;temp1=0.0;}cout<<endl;//回带求解x[]x[1]=y[1];for(s=n;s>=1;s--){for(int t=s+1;t<=n;t++){temp2=temp2+U[s][t]*x[t];}x[s]=(y[s]-temp2)/U[s][s];temp2=0.0;}cout<<"x="<<endl;for(i=1;i<=n;i++)cout<<" "<<x[i]<<endl;}}七.程序执行结果1.系数矩阵主子式有为0时,随便输入一个矩阵结果如下:2. 系数矩阵各阶主子式不为0时,输入课本70页习题12所示矩阵系数5791068109710895765A⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦与右端项2618229b⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦,顺序执行程序,结果经matlab验证无误,结果如下:11 / 12验证:七.实验优缺点此算法可以一定程度上判断是否可解并解出其解,但由于未添加主元选取的步骤,导致对于有的方程组会判断失误并认为方程无解。