线性方程组的古典迭代算法

#include <math.h>

 //雅克比迭代

int Jacobi(double **a,double *b,double *x,int dim,

double eps=0.0001,int CNT=100000)

{

    double *y=new double [dim];

    double err=1.0;

    int i,j,count=0;

    for(i=0;i<dim;i++)

        y[i]=0;

    while(err>eps&&count<CNT)

    {

        err=0;

        count++;

        for(i=0;i<dim;i++)

        {

            x[i]=0;

            for(j=0;j<i;j++)

                x[i]+=a[i][j]*y[j];

            for(j=i+1;j<dim;j++)

                x[i]+=a[i][j]*y[j];

            x[i]=(b[i]-x[i])/a[i][i];

            if(err<fabs(x[i]-y[i]))

                err=fabs(x[i]-y[i]);

        }

        for(i=0;i<dim;i++)

            y[i]=x[i];

    }

    return count;

}

//高斯——赛德尔迭代

 

int G_S(double **a,double *b,double *x,int dim,double eps=0.0001,int CNT=100000)

{

    double y;

    double err=1.0;

    int i,j,count=0;

    for(i=0;i<dim;i++)

        x[i]=0;

    while(err>eps&&count<CNT)

    {

        err=0;

        count++;

        for(i=0;i<dim;i++)

        {

            y=x[i];x[i]=0;

            for(j=0;j<i;j++)

                x[i]+=a[i][j]*x[j];

            for(j=i+1;j<dim;j++)

                x[i]+=a[i][j]*x[j];

            x[i]=(b[i]-x[i])/a[i][i];

            if(err<fabs(y-x[i]))

                err=fabs(x[i]-y);

        }

    }

    return count;

}

//超松弛迭代

 

int SOR(double **a,double *b,double *x,int dim,double w=1,double eps=0.0001,int CNT=100000)

{

    double y;

    double err=1.0;

    int i,j,count=0;

    for(i=0;i<dim;i++)

        x[i]=0;

    while(err>eps&&count<CNT)

    {

        err=0;

        count++;

        for(i=0;i<dim;i++)

        {

            y=x[i];

            x[i]=0;

            for(j=0;j<i;j++)

                x[i]+=a[i][j]*x[j];

            for(j=i+1;j<dim;j++)

                x[i]+=a[i][j]*x[j];

            x[i]=(b[i]-x[i])/a[i][i];

            x[i]=(1-w)*y+w*x[i];

            if(err<fabs(y-x[i]))

                err=fabs(x[i]-y);

        }

    }

    return count;

}