best.bak

另带完整程序实验说明

#include <stdio.h>
#include <math.h>
#define n 5

FILE *fp1,*fp2;

float min(float a,float b)
      { if(a>b) return b; else return a; }
void eq(float A[],float B[])                                  //A[i]<--B[i]
      { int i; for(i=0;i<n;i++) A[i]=B[i];}
void uni_eq(float C[],float D[])                              //C[i]<--D[i]
      { int i; for(i=0;i<n;i++) C[i]=(-1)*D[i];}
float product(float E[],float L[])                            //return E[]*L[]
      { float j=0;int i;
	for(i=0;i<n;i++) j=j+E[i]*L[i]; return j; }
void  add(float c1,float M[],float c2,float N[],float R[])    //R[]=c1*M[]+c2*N[]
      { int i; for(i=0;i<n;i++)  R[i]=c1*M[i]+c2*N[i] ;    }
float F(float xc[])                                           //f(x)
      { return pow(xc[0],2)*pow(xc[1],2)+pow(xc[0],2)+pow(xc[1],2)+pow(xc[2],2)+pow(xc[3],2)+pow(xc[4]-10,2); }
void DF(float G[],float X[])                                  //G[]=gradient(x[])
      { G[0]=2*X[0]*pow(X[1],2)+2*X[0];G[1]=2*X[1];G[2]=2*X[2];G[3]=2*X[3];G[4]=2*pow((X[4]-10),1); }
void ls(float MX[],float NP[],float RX[])                     //RX[]=ls(MX[],NP[])
      { float p=0.1,B=2,v=0.4,a=0,b=1e+8,t=1,f1,f2,g1=0,g2=0,gg[n];
	int i,k=0;
	f1=F(MX);
	DF(gg,MX);  g1=product(gg,NP);
	for(k=0;k<100;k++)
	   { add(1,MX,t,NP,RX);                               //MX has been changed
	     f2=F(RX); DF(gg,RX);                             //g=DF(x)
	     g2=product(gg,NP);
	     if(g2<v*g1)  { a=t;t=min(0.5*(a+b),2*a);continue;}
	     else  if(f1-f2<(-1)*p*t*g1)
		   {b=t;t=0.5*(a+b);continue; }
	     break;      }

	}
  void main()
{  int i=0,j=0,H=0,I=1,J=1,k,row,colum;
   float f,f0,t,sum=0,max,p[n],p0[n],g[n],g0[n],h1,h2,h3,rfa;
   char C;
   float x0[n]={1000,2000,3000,4000,5000},x[n],e1,e2,e3,e=1e-8;
   FILE *fp1,*fp2;
   fp1=fopen("ls","w+");
       fprintf(fp1,"START=========================================================================WRITE\n",I++);
   fclose(fp1); ++I;

   eq(x,x0);
   f=F(x0);
   DF(g,x0);
   j=0;
   do
   {    k=0;  uni_eq(p0,g);
       for(;H!=1&&k<300;k++)
       { eq(x0,x);
	   f0=f;
	   eq(g0,g);
	   ls(x0,p0,x);		                     //x0 has been put into x
	   fp1=fopen("ls","a");
       fprintf(fp1,"ls%d============================================================================ls%d\n",J,J);
       for(i=0;i<n;i++)
       fprintf(fp1,"x%-d=%-17.5e",i,x[i]);
       fclose(fp1); ++J;
	   f=F(x);
	   DF(g,x);

	   e1=e2=1e-15;e3=1e-14;                    //H order
	   h1=h2=h3=0;
	      for(i=0;i<n;i++)
	      { h1=h1+fabs(g[i]);
		h2=h2+fabs(x[i]-x0[i]);
		h3=h3+fabs(x0[i]);  }
		h2=h2/(h3+1);
		h3=fabs(F(x)-F(x0))/(fabs(F(x0))+1);
	     fp1=fopen("ls","a");
	     fprintf(fp1,"\n\nh1=%-17.eh2=%-17.eh3=%-17.e",h1,h2,h3);
	     fprintf(fp1,"\ne3=%-17.ee1=%-17.ee2=%-17.e\n",e3,e1,e2);
	     fclose(fp1);
		if(h1<e3&&h2<=e1&&h3<=e2)
		{ H=1;
	     fp1=fopen("ls","a");
	     fprintf(fp1,"\n=============================H order has been passed.===============================\n");
	     for(i=0;i<n;i++)    fprintf(fp1,"x%-d=%-17.5e",i,x[i]);
	     fprintf(fp1,"\nH=%-18dFx=%-17.5ek=%-18dj=%-17d\n",H,F(x),k,j);
	     fprintf(fp1,"\nh1=%-17.eh2=%-17.eh3=%-17.e\n",h1,h2,h3);
	     fprintf(fp1,"\ne3=%-17.ee1=%-17.ee2=%-17.e\n",e3,e1,e2);
	     fclose(fp1);
		  break;  break; }
		fp1=fopen("ls","a");
		fprintf(fp1,"\n");
	     fprintf(fp1,"H=%-18dFx=%-17.5ek=%-18dj=%-17d\n",H,F(x),k,j);
	     fprintf(fp1,"\t\tJump over the H order and go down.\n");
	     fclose(fp1);
	   if(k==n) { break; continue;}         //break from 'for' and continue 'do while'
	      rfa=product(g,g)/product(g0,g0);
	      add(-1,g,rfa,p0,p);               //p=-g+a*p0
	   if(fabs(product(p,g))<=e)  { break; continue; }
	   if(product(p,g)<-1*e)
	      eq(p0,p);
	      else  uni_eq(p0,p);
	      }                                 //end for
   j++;}while(H==0&&j<100);

   for(i=0;i<n;i++)                             //output on the square-eye
   printf("x0%d=%f\t",i+1,x[i]);
   printf("\nFx=%.5f\tk=%d\tj=%d\n",F(x),k,j);

   }