freuroth.cpp

共轭梯度的几种方法对标准函数的测试结果比较 附比较结果

//cF1                        Extended Freudenstein & Roth
//*
//*                       Initial Point: [0.5, -2, ...,0.5, -2].
#include<iostream.h>
#include<math.h>
#define N 2
void CALF(int n,double x[N],double &f){
      double t1,t2;     
      f=0.0;
      int n1=n/2; 
      for(int i=0;i<n1;i++){
        t1=-13.0+x[2*i-1]+5.0*x[2*i]*x[2*i]-pow(x[2*i],3)-2.0*x[2*i]; 
        t2=-29.0+x[2*i-1]+pow(x[2*i],3)+pow(x[2*i],2)-14.0*x[2*i];   
        
        f=f+t1*t1+t2*t2;}
}
  
//C
void CALG(int n,double x[N],double g[N]){  
      double t1,t2;
      int n1=n/2; 
      int j=0;   
      for(int i=0;i<n1;i++) 
        t1=-13.0+x[2*i-1]+5.0*x[2*i]*x[2*i]-pow(x[2*i],3)-2.0*x[2*i]; 
        t2=-29.0+x[2*i-1]+pow(x[2*i],3)+pow(x[2*i],2)-14.0*x[2*i];  
        
        g[j]=2.0*(t1+t2);   
        g[j+1]=2.0*t1*(10.0*x[2*i]-3.0*x[2*i]*x[2*i]-2.0)+2.0*t2*(3.0*x[2*i]*x[2*i]+2.0*x[2*i]-14.0);   
        j=j+2; 
}

double DOIT(double m1[N],double m2[N]){
	double sum=0;
	for(int i=0;i<N;i++)
		sum+=m1[i]*m2[i];
	return sum;
}//DOIT
int wolf1(double x[],double work[],double g[],double alpha,double c1,double c2){
	double f1,f2;
	double temp[N];
	CALF(N,x,f1);
	for(int i=0;i<N;i++)
		temp[i]=x[i]+alpha*work[i];
	CALF(N,temp,f2);
	double temp1=DOIT(g,work);
	double temp2=-c1*alpha*temp1;
	if((f1-f2)>=temp2)return 1;
	else return 0;
}//wolf1
int wolf2(double x[],double work[],double g[],double alpha,double c1,double c2){
	double ng[N];
	double temp[N];
	for(int i=0;i<N;i++)
		temp[i]=x[i]+alpha*work[i];
	CALG(N,temp,ng);
	double temp1=DOIT(ng,work);
	double temp2=c2*DOIT(g,work);
	if(temp1>=temp2)return 1;
	else return 0;
}//wolf2
void main()
{	
	double x[N],G[N],work[N];
	double c1=0.1,c2=0.5;
	double alpha=1;
	double a=0,b=1000;
	int k=0;
	for(int i=0;i<N;i++)
		if(i%2==0)x[i]=0.5;
			else x[i]=-2;
   CALG(N,x,G);
	for(i=0;i<N;i++)
		work[i]=-G[i];
	wolf1(x,work,G,alpha,c1,c2);
step2:
	if(wolf1(x,work,G,alpha,c1,c2)&&wolf2(x,work,G,alpha,c1,c2)){
	for(i=0;i<N;i++)
		   x[i]=x[i]+alpha*work[i];
	goto step5;}
	if(wolf1(x,work,G,alpha,c1,c2)==0){
		k++;
		goto step3;}
	if(wolf2(x,work,G,alpha,c1,c2)==0){
		k++;
		goto step4;}
step3:
	b=alpha;
	alpha=(a+alpha)/2;
	goto step2;
step4:
	a=alpha; 
	alpha=(2*alpha)>(alpha+b)/2?(alpha+b)/2:(2*alpha);
	goto step2;
step5:
	cout<<"迭带步数"<<k<<'\n';
	cout<<"alpha:"<<alpha<<'\n';

}