shekel.cc

在linux下面实现的单纯性算法的源代码

// SHEKEL FAMILY
//   coded by Erin Parker
//   modified by Chris Siefert
//   tested, comments added by P.L. Shepherd
// Shekel Family on [0,10]^4
//-----------------------------

//if m==5
//f(4,4,4,4) = -10.1531

/* actually, I get f(1,1,1,1) = -5.0551956 and f(8,8,8,8)
 * = 5.100757.   --pls 1/2001
 */
//f(1,1,1,1) = -5.1007
//f(8,8,8,8) = -5.0551
//f(6,6,6,6) = -2.6828
//f(3,7,3,7) = -2.6303

//if m==7
//f(4,4,4,4) = -10.4029
//f(8,8,8,8) = -5.1288
//f(1,1,1,1) = -5.0876

/* For this one, I get -3.7227518
 *  --pls 1/2001
 */
//f(5,5,3,3) = -3.7243
//f(3,7,3,7) = -2.7659
//f(6,6,6,6) = -2.7519

/* For this one, I get -1.8370824
 *  --pls 1/2001
 */   
//f(2,9,2,9) = -1.83759

//if m==10
//f(4,4,4,4) = -10.5364           // -10.536284   pls 1/2001
//f(8,8,8,8) = -5.1756
//f(1,1,1,1) = -5.1285
//f(5,5,3,3) = -3.8354            //  -3.833635   pls 1/2001
//f(6,6,6,6) = -2.8711            //  -2.8709955  pls 1/2001
//f(3,7,3,7) = -2.8066
//f(7,3.6,7,3.6) = -2.42734       //  -2.4265188  pls 1/2001
//f(6,2,6,2) = -2.4217            //  -2.4208314  pls 1/2001
//f(2,9,2,9) = -1.8584            //  -1.8589173  pls 1/2001      
//f(8,1,8,1) = -1.6765            //  -1.6752526  pls 1/2001 


#include <math.h>
#include <iostream.h>
#include <iomanip.h>
#include <malloc.h>
#include <stdlib.h>

#define n 4 //do not change this... this is the number of dimensions
#define m 5 //change for shekel5, shekel7 or shekel10

#define LBOUND 0.0
#define UBOUND 10.0

int main()

{
  double * x=(double*)malloc(n*sizeof(double));
  double f;
  int flag;
  cin>>x[0]>>x[1]>>x[2]>>x[3];
  //deal with constrained minimization
    if((x[0]<0||x[0]>10)||(x[1]<0||x[1]>10)||(x[2]<0||x[2]>10)||(x[3]<0||x[3]>10)){
      cerr<<"SHEKEL FAILED!!!!\n";
      flag = 1; 
      f = 10000000;
      free(x);
      return 1;
    }
    else {
  //initialize matrix a
  double a[10][n];
  for(int k=0; k<n; k++)
     a[0][k]=4.0;
  for(int k=0; k<n; k++)
     a[1][k]=1.0;
  for(int k=0; k<n; k++)
     a[2][k]=8.0;
  for(int k=0; k<n; k++)
     a[3][k]=6.0;
  a[4][0]=3.0;
  a[4][1]=7.0;
  a[4][2]=3.0;
  a[4][3]=7.0;
  a[5][0]=2.0;
  a[5][1]=9.0;
  a[5][2]=2.0;
  a[5][3]=9.0;
  a[6][0]=5.0;
  a[6][1]=5.0;
  a[6][2]=3.0;
  a[6][3]=3.0;
  a[7][0]=8.0;
  a[7][1]=1.0;
  a[7][2]=8.0;
  a[7][3]=1.0;
  a[8][0]=6.0;
  a[8][1]=2.0;
  a[8][2]=6.0;
  a[8][3]=2.0;
  a[9][0]=7.0;
  a[9][1]=3.6;
  a[9][2]=7.0;
  a[9][3]=3.6;

  //initiate vector c
  double c[10];
  c[0]=0.1;
  c[1]=0.2;
  c[2]=0.2;
  c[3]=0.4;
  c[4]=0.4;
  c[5]=0.6;
  c[6]=0.3;
  c[7]=0.7;
  c[8]=0.5;
  c[9]=0.5;

  f=0.0;
  double sum=0.0;
  double temp;
  for(int i=0; i<m; i++)
    {
      for(int j=0; j<n; j++)
	{
          temp=x[j]-a[i][j];
          sum+=temp*temp;
	}
      f-=(1.0/(sum+c[i]));
      sum=0.0;
    }
  }
  flag=0;//cs
  free(x);
  cout<<setprecision(8)<<f<<endl;
  return flag;
}