高斯迭代法解线性方程组 下载本文

内容发布更新时间 : 2024/7/6 4:34:32星期一 下面是文章的全部内容请认真阅读。

一. 问题描述

用Gauss-Seidel迭代法求解线性方程组

由Jacobi迭代法中,每一次的迭代只用到前一次的迭代值。使用了两倍的存储空间,浪

x费了存储空间。若每一次迭代充分利用当前最新的迭代值,即在计算第i个分量i用最新分量x1(k?1)(k?1)时,

,x2(k?1)???xi-1(k?1)代替旧分量x1,x2(k)(k)???xi-1,可以起到节省存储

(k)空间的作用。这样就得到所谓解方程组的Gauss-Seidel迭代法。

二. 算法设计

将A分解成A?L?D?U,则Ax?b等价于(L?D?U)x?b 则Gauss-Seidel迭代过程

Dx(k?1)?b?Lx(k?1)?Ux(k)故

(D?L)x(k?1)?b?Ux(k)

若设(D?L)存在,则

?1x(k?1)?(D?L)?1Ux(k)?(D?L)?1b

G?(D?L)?1U,f?(D?L)?1b

则Gauss-Seidel迭代公式的矩阵形式为

x(k?1)?Gx(k)?f

其迭代格式为

(0)(0)Tx(0)?(x1(0),x2,???,xn) (初始向量),

x或者

(i?1)ii?1i?11(k?1)1,2,???;i?0,1,2,???n) ?(bi??aijxj??aijx(jk)) (k?0,aiij?1j?i?1?xi(k?1)?xi(k)??xi(k?k?0,1,2,???;i?0,1,2,???n)?i?1i?11 ? (i?1)(k?1)xi?(bi??aijxj??aijx(jk))?aiij?1j?i?1?

三. 程序框图

开始 读入数据n,初始向量,增广矩阵 1?k(bi??axj?1i?1(k?1)ijj??aijx(jk))/aii?yij?i?1i?1 maxxi?yi???k?1?kyi?xii?1,2,3???,n?? k=N? ??输出输出迭代失败标志 y1,y2???yn 结束

四. 结果显示

TestBench

利用Gauss-Seidel迭代法求解下列方程组

?5x1?2x2?x3??12??(0)??x1?4x2?2x3?20, 其中取x?1。 ?2x?3x?10x?323?1运行程序

依次输入:

1.方阵维数

2.增广矩阵系数 3.初始向量

得到:

迭代12次后算出 x[1] = -4.0 x[2] = 3.0 x[3] = 2.0

五. 程序

#include #include #include #include

#define MAX_n

100

#define PRECISION 0.0000001 #define MAX_Number 1000

void VectorInput(float x[],int n) //输入初始向量 { }

void MatrixInput(float A[][MAX_n],int m,int n) //输入增广矩阵 { }

void VectorOutput(float x[],int n) //输出向量 { }

int IsSatisfyPricision(float x1[],float x2[],int n) //判断是否在规定精度内 {

int i; int i;

for(i=1;i<=n;++i)

printf(\,i,x[i]); int i, j;

printf(\输入系数矩阵:\\n\); for(i=1;i<=m;++i) { }

printf(\增广矩阵行数%d : \,i); for(j=1;j<=n;++j)

scanf(\,&A[i][j]);

for(i=1;i<=n;++i) { }

printf(\,i); scanf(\,&x[i]); int i;

}

int Jacobi_(float A[][MAX_n],float x[],int n) //具体计算 { do{ }

int main() //主函数 {

if(k>=MAX_Number) { }

printf(\迭代次数为%d 次\,k); return 0; return 1; else

for(i=1;i<=n;++i) { }

printf(\); for(i=1;i<=n;++i) { } ++k;

x[i]=A[i][n+1]; for(j=1;j<=n;++j)

if(j!=i) x[i]-=A[i][j]*x[j]; x[i]/=A[i][i]; return 1;

if(fabs(A[i][i])>PRECISION)

printf(\,i,x[i]); x_former[i]=x[i];

k=0;

printf(\初始向量x0:\\n\); VectorInput(x,n); float x_former[MAX_n]; int i,j,k; for(i=1;i<=n;++i)

if(fabs(x1[i]-x2[i])>PRECISION) return 1; return 0;

else

}while(IsSatisfyPricision(x,x_former,n) && k