example3_14.java

清华大学2002年出版的《科学与工程数值计算算法Java》配套源码

import java.applet.*;
import java.awt.*;

public class Example3_14 extends Applet
{
 public static void  QR(double[][] a,double[][] Q,int m,int n )
  {int i,j,k,l=Math.max(m,n);
   double[] t=new double[l];
   double alpha,beta;
   for(i=0;i<m;i++)
    for(j=0;j<m;j++)
     {if (i==j) Q[i][j]=1;
      else Q[i][j]=0; }

  for(k=0;k<Math.min(n,m)-1;k++)
  { if(a[k][k]==0)alpha=0;
    else
   {alpha=0;
    for(l=k;l<m;l++)alpha+=a[l][k]*a[l][k];
    alpha=Math.sqrt(alpha);
    alpha*=a[k][k]/Math.abs(a[k][k]);}
    beta=(alpha+a[k][k])*alpha;
    for(i=0;i<m;i++)
     {t[i]=Q[i][k]*alpha;
      for(l=k;l<m;l++)t[i]+=Q[i][l]*a[l][k];
     }

    for(i=0;i<m;i++)
    { Q[i][k]-=t[i]/alpha;
	  for(j=k+1;j<m;j++)Q[i][j]-=(t[i]*a[j][k]/beta);}
    for(j=0;j<n;j++)
     {t[j]=alpha*a[k][j];
      for(l=k;l<m;l++)t[j]+=a[l][k]*a[l][j];
     }

   for(j=k+1;j<n;j++)a[k][j]-=t[j]/alpha;
   for(i=k+1;i<m;i++)
     for(j=k+1;j<n;j++)a[i][j]-=(a[i][k]*t[j]/beta);
   for(l=k+1;l<m;l++)a[l][k]=0;
   a[k][k]=-alpha;
   }
 }

 public static void QRSolve(double[][] a,double[] b,int m,int n)
 {double[][] Q=new double[m][m];
  double[] X=new  double[m];
  int i,j;
  QR(a,Q,m,n);
  for(i=0;i<n;i++)
    { X[i]=0;
	  for(j=0;j<m;j++)
	  X[i]+=Q[j][i]*b[j];
    }
  X[n-1]/=a[n-1][n-1];
  for(i=n-2;i>=0;i--)
    {for(j=i+1;j<n;j++) X[i]-=a[i][j]*X[j];
     X[i]/=a[i][i];}
  for(j=0;j<n;j++)b[j]=X[j];
}


  public void paint(Graphics g)
  { double[][] a={{1,1,-1},{2,1,0},{1,-1,0},{-1,2,1}};
   double[] b={2,-3,1,4};
   int m=4,n=3;
   QRSolve(a,b,m,n);
   for(m=0;m<n;m++)
   g.drawString("X"+m+"="+b[m],10,20*m+20);
  }

}