实验二 线性方程组的数值解法 （一）利用矩阵LU分解求解线性方程组
2.1 实验目的
①掌握利用矩阵LU分解求解线性方程组的解的基本思路和步骤； ②培养编程与上机调试能力。 2.2 算法描述
2.2.1 矩阵LU分解求解线性方程组的基本思路
?LU， 则Ax?b若A
(LU)x?bLy?b (下三角方程组）Ux?y (上三角方程组）
这样就将一般的线性方程组化简为三角方程组，从而简化了线性方程组的求
解过程。
2.2.2 矩阵LU分解求解线性方程组的计算步骤 步骤1： 利用矩阵LU分解得到矩阵L和U 步骤2： 求解三角方程组 2.3 实验内容
用矩阵LU分解求解如下线性方程组：
?123??x1??1??223??x???2? ???2?????333????x3????3??2.4 程序源代码
void shuru (n,a)
int n;
float a[5][6]; {int i,j; float b;
printf(\ for (i=0;i<=n-1;i++) {for (j=0;j<=n;j++) {scanf(\ printf(\ a[i][j]=b; }
printf(\ } }
void LUfenjie(n,a)
int n;
float a[5][6]; {int i,j,k;
for (i=1;i<=n-1;i++)
a[i][0]=a[i][0]/a[0][0]; for(j=1;j<=n;j++)
a[1][j]=a[1][j]-a[1][0]*a[0][j]; for (i=2;i<=n-1;i++) {for (j=1;j<=i-1;j++) {for (k=0;k<=j-1;k++)
a[i][j]=a[i][j]-a[1][k]*a[k][j]; a[i][j]=a[i][j]/a[j][j]; }
for(j=i;j<=n;j++) for(k=0;k<=i-1;k++)
a[i][j]=a[i][j]-a[1][k]*a[k][j]; } }
void huidai(n,a) int n;
float a[5][6]; {int i,j;
a[n-1][n]=a[n-1][n]/a[n-1][n-1]; for(i=n-2;i>=0;i--)
{ for(j=n-1;j>=i+1;j--)
a[i][n]=a[i][n]-a[i][j]*a[j][n]; a[i][n]=a[i][n]/a[i][i]; } } main() {int i,n;
float a[5][6];
printf(\ scanf(\ printf(\ shuru(n,a); LUfenjie(n,a); huidai(n,a);
for(i=0;i<=n-1;i++)
printf(\ printf(\ }
（二）高斯列主元消去法
2.1 实验目的
① 掌握高斯列主元消去法的基本思路和步骤； ② 培养编程与上机调试能力。 2.2 算法描述
2.2.1 高斯消去法基本思路。
设有方程组Ax?b，设A是可逆矩阵。高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵B??A?b?，将其中的A变换成一个上三角矩阵，然后求解这个三角形方程组。
2.2.2列主元高斯消去法计算步骤
将方程组用增广矩阵B??A?b???aij?n?(n?1)表示。
步骤1：消元过程，对k?1,2,?,n?1 （1） 选主元，找ik??k,k?1,?,n?使得
aik,k?maxaik
k?i?n（2） 如果aik,k?0，则矩阵A奇异，程序结束；否则执行（3）。
（3） 如果ik?k，则交换第k行与第ik行对应元素位置，akj?akij，
j?k,?,n?1。
（4） 消元，对i?k,?,n，计算lik?aik/akk,对j?k?1,?,n?1，计算
aij?aij?likakj.
步骤 2：回代过程：
（1） 若ann?0,则矩阵奇异，程序结束；否则执行（2）。
n??（2） xn?an,n?1/ann;对i?n?1,?,2,1，计算xi??ai,n?1??aijxj?/aii
j?i?1??2.3实验内容
2x1?x2?2x3?5解线性方程组 5x1?x2?x3?8
x1?3x2?4x3??42.4 程序源代码
#include
void zhuyuan (k,n,a) int k,n;
float a[5][6];
{int t,i,j; float p,q;
p=fabs(a[k][k]); t=k;
for(i=k+1;i<=n;i++) if (fabs(a[i][k])>p) {p=fabs(a[i][k]>p); t=i; }
for(j=k;j<=n+1;j++) {q=a[k][j];
a[k][j]=a[t][j]; a[t][j]=q; } }
void shuru(n,a) int n;
float a[5][6]; {int i,j;
printf(\ for(i=0;i<=n+1;i++) {for(j=0;j<=n+1;j++)
{scanf(\ printf(\ }
printf(\ } }
void xiaoyuan (n,a) int n;
float a[5][6]; {int k,i,j;
for(k=0;k a[i][j]=a[i][j]-a[i][k]*a[k][j]/a[k][k]; } } void huidai(n,a,x) int n; float a[5][6],x[5]; {int k,j; x[n]=a[n][n+1]/a[n][n]; for (k=n+1;k>=0;k--) {x[k]=a[k][n+1]; for(j=k+1;j<=n;j++) x[k]=x[k]-a[k][j]*x[j]; x[k]=x[k]/a[k][k]; } } main() {int n,i; float a[5][6],x[5]; printf(\ scanf(\ printf(\ n=n-1; shuru(n,a); xiaoyuan(n,a); huidai(n,a,x); for(i=0;i<=n;i++) printf(\ printf(\ }