Gauss-Seidel迭代法求解线性方程组
- 格式:docx
- 大小:139.44 KB
- 文档页数:6
一. 问题描述
用Gauss-Seidel 迭代法求解线性方程组
由Jacobi 迭代法中,每一次的迭代只用到前一次的迭代值。使用了两倍的存储空间,浪
费了存储空间。若每一次迭代充分利用当前最新的迭代值,即在计算第i 个分量)
1(+k i
x 时,
用最新分量)
1(1
+k x ,⋅⋅⋅+)
1(2
k x )
1(1
-+k i x 代替旧分量)
(1
k x ,⋅⋅⋅)
(2
k x )
(1
-k i x ,可以起到节省存储
空间的作用。这样就得到所谓解方程组的Gauss-Seidel 迭代法。
二. 算法设计
将A 分解成U D L A --=,则b x =A 等价于b x =--U)D (L 则Gauss-Seidel 迭代过程
)
()1()1(k k k Ux Lx b Dx ++=++
故
)()1()(k k Ux b x L D +=-+
若设1
)(--L D 存在,则
b L D Ux L D x k k 1)(1)1()()(--+-+-=
令
b L D f U L D G 11)()(---=-=,
则Gauss-Seidel 迭代公式的矩阵形式为
f Gx x k k +=+)()1(
其迭代格式为
T
n x x x x )
()0()0(2)0(1)0(,,,⋅⋅⋅= (初始向量),
)(1
1
1
1
1
)()1()1(∑∑-=-+=++--=i j i i j k j ij
k j ij i ii i i
x a
x a b a x
)210i 210(n k ⋅⋅⋅=⋅⋅⋅=,,,;,,,
或者
⎪⎩
⎪⎨⎧--=⋅⋅⋅=⋅⋅⋅==∆+=∑∑-=-+=+++)
(1)210i 210(111
1)()
1()1()()1(i j i i j k j ij k j ij i ii i i i k i k i x a x a b a x n k k x x x ,,,;,,,
三. 程序框图
四. 结果显示
TestBench
利用Gauss-Seidel 迭代法求解下列方程组
⎪⎩⎪⎨⎧=+-=++--=++3
1032202412
25321
321321x x x x x x x x x , 其中取→
=1)
0(x 。 运行程序
依次输入:
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) //输入初始向量
{
int i;
for(i=1;i<=n;++i)
{
printf("x[%d]=",i);
scanf("%f",&x[i]);
}
}
void MatrixInput(float A[][MAX_n],int m,int n) //输入增广矩阵
{
int i, j;
printf("\n输入系数矩阵:\n");
for(i=1;i<=m;++i)
{
printf("增广矩阵行数%d : ",i);
for(j=1;j<=n;++j)
scanf("%f",&A[i][j]);
}
}
void VectorOutput(float x[],int n) //输出向量
{
int i;
for(i=1;i<=n;++i)
printf("\nx[%d]=%f",i,x[i]);
}
int IsSatisfyPricision(float x1[],float x2[],int n) //判断是否在规定精度内{
int i;
for(i=1;i<=n;++i)
if(fabs(x1[i]-x2[i])>PRECISION) return 1;
return 0;
}
int Jacobi_(float A[][MAX_n],float x[],int n) //具体计算
{
float x_former[MAX_n];
int i,j,k;
printf("\n初始向量x0:\n");
VectorInput(x,n);
k=0;
do{
for(i=1;i<=n;++i)
{
printf("\nx[%d]=%f",i,x[i]);
x_former[i]=x[i];
}
printf("\n");
for(i=1;i<=n;++i)
{
x[i]=A[i][n+1];
for(j=1;j<=n;++j)
if(j!=i) x[i]-=A[i][j]*x[j]; if(fabs(A[i][i])>PRECISION)
x[i]/=A[i][i];
else
return 1;
}
++k;
}while(IsSatisfyPricision(x,x_former,n) && k if(k>=MAX_Number) return 1; else { printf("\nGauss-Seidel迭代次数为%d 次",k); return 0; } } int main() //主函数 {