powell.c

优化设计常用程序

#include"stdio.h"
#include"conio.h"
#include"math.h"
#include"f:\sq.h"
#include"f:\golden.h"
#define  N 5
powell(float,int);
tidu();
float f(float);
float f1(float *);
static float d1[5],d[5][5],x0[5],x1[5][5],g0[5],g1[5],h0[5][5],yk[5],sk[5];
int n;
main()
{
  float E;
  int i;
  clrscr();
  printf("\nPlease input weishu n and the E:\n");
  scanf("%d%f",&n,&E);
  printf("\nPlease input the first point:\n");
  for(i=0;i<n;i++)   scanf("%f",&x0[i]);
  powell(E,n);
  getch();
}
powell(float E,int n)
{
  float A,B,delta[5],F0,F2,F3,deltam;
  double a,b,EE=0.0;
  int i,j,l,maxm=0;
  int k;
  for(i=0;i<n;i++) for(j=0;j<n;j++) if(i==j) d[i][j]=1;else d[i][j]=0;
FIRST:
  for(i=0;i<n;i++) x1[0][i]=x0[i];
  for(l=1;l<=n;l++)
  {  for(i=0;i<n;i++) d1[i]=d[l-1][i];
     for(i=0;i<n;i++) x0[i]=x1[l-1][i];
     sq(&a,&b,j);
     A=golden(a,b,E);
     for(i=0;i<n;i++)
     { x1[l][i]=x1[l-1][i]+A*d1[i];
       printf("\nx1[%d]=%f",i,x1[l][i]);
     }
     delta[l]=f1(x1[l-1])-f1(x1[l]);
  }
  for(i=0;i<n;i++)
  {  d[n][i]=x1[n][i]-x1[0][i];
     x1[n+1][i]=2*x1[n][i]-x1[0][i];
  }
  F0=f1(x1[0]);F2=f1(x1[n]);F3=f1(x1[n+1]);
  deltam=-1000;
  for(i=1;i<=n;i++)
     if(delta[i]>deltam) {deltam=delta[i];maxm=i;}
  if((F3<F2)&&((F0-2*F2+F3)*pow(F0-F2-deltam,2)<0.5*deltam*pow(F0-F3,2)))
    {
       for(i=0;i<n;i++) d1[i]=d[n][i];
       for(i=0;i<n;i++) x0[i]=x1[n][i];
       sq(&a,&b,j);
       A=golden(a,b,E);
       for(i=0;i<n;i++)
	 { x0[i]=x1[n][i]+A*d1[i];
	 }
       for(i=maxm-1;i<n;i++)
	  for(j=0;j<n;j++)  d[i][j]=d[i+1][j];
    }
  else if(F2<F3)  for(i=0;i<n;i++) x0[i]=x1[n][i];
  else for(i=0;i<n;i++) x0[i]=x1[n+1][i];
  EE=0.0;
  for(i=0;i<n;i++)
    EE=EE+pow(x1[n][i]-x1[0][i],2);
  if(fabs(EE)<E) goto  EXIT;
  else goto  FIRST;
EXIT:  for(i=0;i<n;i++)
    printf("\nX%d=%f",i,x0[i]);
}
 float f(float a)
 {
   float f=0;
     f=10*pow((x0[0]+a*d1[0]+x0[1]+a*d1[1]-5),2)+pow((x0[0]+a*d1[0]-(x0[1]+a*d1[1])),2);
   return(f);
 }
float f1(float *xx)
  {
     return(10*pow(*xx+*(xx+1)-5,2)+pow(*xx-*(xx+1),2));

  }