Gauss-Seidel迭代法求方程组的解

注意c语言中,abs()为求整型的绝对值的函数。对于float型需用fabs(),代码被这个不起眼的错误坑了好久,最后得出结论不熟悉的系统函数一定要慎用,只有在确切知道该函数的顶以后才取使用它,否则自己写函数都比使用不明函数来的好。

c语言版

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <ctime>
const int maxn=1001,maxl=100;
double  x[maxn],y[maxn],b[maxn],A[maxn][maxn];
int flag;
int gdsl(int n){
	int k=0;
	for(int i=0;i<n;i++){
		x[i]=0;
	}
	while(1){
		for(int i=0;i<n;i++)
			y[i]=x[i];
		for(int i=0;i<n;i++){
			double z=b[i];
			for(int j=0;j<n;j++){
				if(j!=i)z-=A[i][j]*x[j];
			}
			if(abs(A[i][i])<1e-10||k==maxl){
				flag=0;return 0;
			}
			z/=A[i][i];x[i]=z;
			printf("%lf ",x[i]);
		}
		puts("");
		for(int i=0;i<n;i++)
            printf("%lf ",y[i]);
		puts("@");
		double norm=0;
		for(int i=0;i<n;i++)
			if(fabs(y[i]-x[i])>norm)norm=fabs(y[i]-x[i]);
			printf("%lf#\n",norm);
		if(norm<1e-5)break;
		++k;
	}
	return k;
}
int main(){
	int n;
	scanf("%d",&n);
	for(int i=0;i<n;i++){
		for(int j=0;j<n;j++)
			scanf("%lf",&A[i][j]);
	}
	for(int i=0;i<n;i++)
		scanf("%lf",&b[i]);
	int tmp=gdsl(n);
	printf("%d\n",tmp);
	for(int i=0;i<n;i++)
		printf("%lf\n",x[i]);
	return 0;
}

matlab版

function[x,k,flag,nor]=gdsl(A,b)
% A:方程组的系数矩阵
% b:方程组的右边项
% delata:精度1e-5
% maxl:最大迭代次数
n=length(A);k=0;
x=zeros(n,1);y=zeros(n,1);flag='OK!';
delta=1e-5;maxl=100;
while 1
    y=x;
    for i=1:n
        z=b(i);
        for j=1:n
            if j~=i
                z=z-A(i,j)*x(j);
            end
        end
        if abs(A(i,i))<1e-10|k==maxl
            flag='Fail!';return;
        end
        
        z=z/A(i,i);x(i)=z;
    end
    if norm(y-x,inf)<delta
        break;
    end
    k=k+1;
end

直接调用gdsl函数即可。


评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值