1004.txt

数值分析的c语言几种算法实现

int SOR(int n, double a[][MAX_SIZE], double b[], double x[], double w, double TOL, int MAXN)
{
	
	
   double S, S1;
   int i,j,k;
   double dX;
   double min,max;
   for(i=0;i<n;i++) x[i]=0.0;
   for(i=0;i<n;i++)
   {
      min=fabs(a[i][i]);
      k=i;     
	  for(j=i;j<n;j++)
      {
			if(fabs(a[j][i])>min) 
			{
				k=j;
				min=fabs(a[j][i]);
			 }
      }
      if(min>1e-10)
	  {
		    for(j=0;j<n;j++)
			{ max=a[i][j];a[i][j]=a[k][j];a[k][j]=max;}
			max=b[i];b[i]=b[k];b[k]=max;//b[i]+=b[k];
      }
     else if(min<1e-10)
	 {
		 min=0;
         k=i;
         for(j=i-1;j>=0;j--)
			if(fabs(a[j][i])>min) 
			{
				k=j;
				min=fabs(a[j][i]);
			 }

			if(min>1e-10)
			{
			  for(j=0;j<n;j++)
					a[i][j]+=a[k][j];
			  b[i]+=b[k];
			}
			else 
				return -1;
    }
  }
  /*for(i=0;i<n;i++)
	{  
		for(j=0;j<n;j++)	  
		  printf("%lf ",a[i][j]);
		  printf("%lf\n",b[i]);
    }*/
   for (k=0; k<MAXN; k++)
   {
     S = 0.0;
     for (i=0; i<n; i++)
     {
       double sum = 0.0;       
	   for ( j=0; j<n; j++)
			 sum += a[i][j]*x[j];
       
	   if (fabs(a[i][i]) < 1e-10)
			return -1;

       dX = w * (b[i] - sum ) /a[i][i];
       S1 = dX;
	   // printf("%lf\n",dX[i]);
       if (fabs(S1) > fabs(S)) S = S1;
	   //printf("%lf\n",S);
       x[i] = x[i] + S1;
	  // printf("%lf\n",pow((double)2,(double)127));
       if(fabs(x[i])>pow((double)2,(double)127)) 
	   {
		 //  printf("%lf\n",S);
		 //  printf("%d\n",k);
		  // cout << k << endl;
		   return -2;
	   }
     }
     if (fabs(S)<=TOL)
		 return 1;
   }
   return 0;
}