📄 best.cpp
字号:
#include <stdio.h>
#include <math.h>
#define n 9
#define m 9
#define m_j 3
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 A[],float B[]) //C[i]<--D[i]
{ int i; for(i=0;i<n;i++) A[i]=(-1)*B[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] ; }
int jieyue(float t) // jieyue function diserve the MAX way
{ if(t>0) return 1; else return 0; }
float F(float xc[]) //f(x)_the first part of F_ALL function
{ float f_sum=0;int i; for(i=0;i<n;i++) f_sum=f_sum+xc[i]; return f_sum; }
float F_II(float u,float v[],float xc[]) // the second part of F
{ int i;float t=0;
for(i=0;i<m;i++) t=t+pow(jieyue(v[i]-2*u*xc[i])*(v[i]-2*u*xc[i]),2)-pow(v[i],2);return t/(4*u); }
float F_III(float lamta[],float xc[])
{ return lamta[0]*(xc[0]+xc[1]+xc[2]-30)+lamta[1]*(xc[3]+xc[4]+xc[5]-300)+lamta[2]*(xc[6]+xc[7]+xc[8]-3000);}
float F_IV(float u,float xc[]) // the forth part of F
{ return u*(pow(xc[0]+xc[1]+xc[2]-30,2)+pow(xc[3]+xc[4]+xc[5]-300,2)+pow(xc[6]+xc[7]+xc[8]-3000,2));}
float F_ALL(float xc[],float u,float v[],float lamta[]) // the completed f_function
{ return F(xc)+F_II(u,v,xc)-F_III(lamta,xc)+F_IV(u,xc) ; }
void DF(float G[],float X[],float u,float v[],float lamta[]) //G[]=gradient(x[])
{ int i;
for(i=0;i<3;i++) G[i]=1+jieyue(v[i]-2*u*X[i])*(2*u*X[i]-v[i])-lamta[0]+2*u*(X[0]+X[1]+X[2]-30);
for(i=3;i<6;i++) G[i]=1+jieyue(v[i]-2*u*X[i])*(2*u*X[i]-v[i])-lamta[1]+2*u*(X[3]+X[4]+X[5]-300);
for(i=6;i<n;i++) G[i]=1+jieyue(v[i]-2*u*X[i])*(2*u*X[i]-v[i])-lamta[2]+2*u*(X[6]+X[7]+X[8]-3000); }
void ls(float MX[],float NP[],float RX[],float u,float v[],float lamta[]) //RX[]=ls.txt(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_ALL(MX,u,v,lamta) ;
DF(gg,MX,u,v,lamta); g1=product(gg,NP);
for(k=0;k<50;k++)
{ add(1,MX,t,NP,RX); //MX has been changed
f2=F_ALL(RX,u,v,lamta); DF(gg,RX,u,v,lamta); //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; }
}
float fai(float xc[],float u,float v[]) //////// Jugment Condition
{ float t=0; int i;
for(i=0;i<m;i++) t=t+pow(min(xc[i],v[i]/(2*u)),2);
t=t+pow(xc[0]+xc[1]+xc[2]-30,2)+pow(xc[3]+xc[4]+xc[5]-300,2)+pow(xc[6]+xc[7]+xc[8]-3000,2);
return sqrt(t);
}
void main()
{ int i=0,j=0,H=0,I=1,J=1,C=10,k;
float f,f0,t,sum=0,max,p[n],p0[n],g[n],g0[n],h1,h2,h3,rfa,V[n],LAMTA[m_j],U,Fai,Fai0;
float x0[n]={10,20,30,100,200,300,1000,2000,3000},x[n],e1,e2,e3,e=1e-6;
FILE *fp1,*fp2;
fp1=fopen("ls.txt","w+");
fprintf(fp1,"START=========================================================================WRITE\n",I++);
fclose(fp1); ++I;
for(i=0;i<n;i++) V[i]=30;
for(i=0;i<m_j;i++) LAMTA[i]=50;
U=10;
fp1=fopen("ls.txt","a");
fprintf(fp1,"\n=============================H order will be passed.===============================\n");
for(i=0;i<n;i++) { fprintf(fp1,"x%-d=%-17.5e",i,x0[i]);if((i+1)%3==0) fprintf(fp1,"\n");}
fprintf(fp1,"\nF=%-17.5eF_II=%-17.5eF_III=%-17.5eF_IV=%-17.5e\n",F(x0),F_II(U,V,x0),F_III(LAMTA,x0),F_IV(U,x0));
fprintf(fp1,"\nH=%-18dFx=%-17.5ek=%-18dj=%-17d\n",H,F_ALL(x0,U,V,LAMTA),k,j);
fclose(fp1);
for(i=0;i<15;i++)
{ // 最外层
Fai0=fai(x0,U,V);
eq(x,x0);
f=F_ALL(x0,U,V,LAMTA);
DF(g,x0,U,V,LAMTA);
//fp1=fopen("ls.txt","w+");
// for(i=0;i<n;i++) fprintf(fp1,"DX%-d=%-17.5e",i,g[i]);
// fclose(fp1); printf("Fai0=%-17.5e",Fai0);
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,U,V,LAMTA); //x0 has been put into x
fp1=fopen("ls.txt","a");
// fprintf(fp1,"ls.txt%d============================================================================ls.txt%d\n",J,J);
// for(i=0;i<n;i++)
// fprintf(fp1,"x%-d=%-17.5e",i,x[i]);
// fclose(fp1); ++J;
f=F_ALL(x,U,V,LAMTA);
DF(g,x,U,V,LAMTA);
e1=e2=1e-5;e3=1e-4; //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_ALL(x,U,V,LAMTA)-F_ALL(x0,U,V,LAMTA))/(fabs(F_ALL(x0,U,V,LAMTA))+1);
// fp1=fopen("ls.txt","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.txt","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]);if((i+1)%3==0) fprintf(fp1,"\n");}
fprintf(fp1,"\nH=%-18dFx=%-17.5ek=%-18dj=%-17d\n",H,F_ALL(x,U,V,LAMTA),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; }
// fp1=fopen("ls.txt","a");
// fprintf(fp1,"\n");
// fprintf(fp1,"H=%-18dFx=%-17.5ek=%-18dj=%-17d\n",H,F_ALL(x,U,V,LAMTA),k,j);
// fprintf(fp1,"\t\tJump over the H order and go down.\n");
// fclose(fp1);
if(k==n) { break; } //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; }
if(product(p,g)<-1*e)
eq(p0,p);
else uni_eq(p0,p);
} //end for
j++;}while(H==0&&j<100);
Fai=fai(x,U,V);
if(Fai<e) break;
if(Fai/Fai0>0.5) U=C*U;
for(i=0;i<n;i++) V[i]=jieyue(V[i]-2*U*x[i])*(V[i]-2*U*x[i]);
LAMTA[0]=LAMTA[0]-2*U*(x[0]+x[1]+x[2]-30);
LAMTA[1]=LAMTA[1]-2*U*(x[3]+x[4]+x[5]-300);
LAMTA[2]=LAMTA[2]-2*U*(x[6]+x[7]+x[8]-3000);
printf("Fai=%-17.5e",Fai);
eq(x0,x);H=0;
}////////////////最外层
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_ALL(x,U,V,LAMTA),k,j);
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -