Modsimpmethod.rar

#include<stdio.h> #include<math.h> #include<iostream.h> int i,j; double d[3][100]={{0,1,0,0},{0,0,1,0},{0,0,0,1}},f[100];//d[][]为单纯形的顶点，本算例中未知数个数为3，则顶点个数为4 double g,h,l,q,s=1,t=2,u=0.5,v=0.0001,y=0;//s为反射系数，t为扩展系数，u为压缩系数，v为允许误差 int o,F,r,D,e,lj=0,N=4;//N为顶点的个数，o为最大值点的位置，F为最小值点的位置，r为次大值点的位置 void function1(int e)//求解函数f[e] { f[e]=(d[0][e]-3)*(d[0][e]-3)+2*(d[1][e]+2)*(d[1][e]+2)+(d[2][e]-4)*(d[2][e]-4);//函数为f=(x1-3)^2+2(x2+2)^2+(x3-4)^2,求其最小值 } void function2() { while((++lj)<100)//最大迭代次数 { for(i=0,g=f[i];i<N-1;i++)//求最大值点 { if(g<f[i+1]) { g=f[i+1];o=i+1; } else if(i==0) o=i; } for(i=0,h=f[i];i<N-1;i++)//求最小值点 { if(h>f[i+1]) { h=f[i+1];F=i+1; } else if(i==0) F=i; } for(i=0,l=f[i];i<N-1;i++)//求次大值点 { if(i==o&&i==0) { l=f[i+1]; r=i+1; continue; } if(l<f[i+1]&&(i+1)!=o) { l=f[i+1];r=i+1; } else if(i==0) r=i; } for(i=0;i<N-1;i++)//求除最大值点，其余点的形心 { d[i][N]=0; for(j=0;j<N;j++) { if(j!=o) { d[i][N]=d[i][N]+d[i][j]; } } d[i][N]=d[i][N]/(N-1);//平均值 } function1(N); for(i=0;i<N-1;i++) { d[i][N+1]=d[i][N]+(d[i][N]-d[i][o])*s;//反射 } function1(N+1); if(f[N+1]<f[F]) { for(i=0;i<N-1;i++) { d[i][N+2]=d[i][N]+(d[i][N+1]-d[i][N])*t;//扩展 } function1(N+2); if(f[N+2]<f[N+1]) { for(i=0;i<N-1;i++) { d[i][o]=d[i][N+2]; } f[o]=f[N+2]; y=0; for(i=0;i<N;i++) { y=y+(f[i]-f[N])*(f[i]-f[N]); } y=sqrt(y/(N*1.0)); if(y<v) { for(i=0;i<N-1;i++) { printf("d[%d][%d]=%f",i,F,d[i][F]); printf("\n"); } printf("f[%d]=%f",F,f[F]); break; } else { function2(); break; } } else { for(i=0;i<N-1;i++) { d[i][o]=d[i][N+1]; } f[o]=f[N+1]; y=0; for(i=0;i<N;i++) { y=y+(f[i]-f[N])*(f[i]-f[N]); } y=sqrt(y/(N*1.0)); if(y<v) { for(i=0;i<N-1;i++) { printf("d[%d][%d]=%f",i,F,d[i][F]); printf("\n"); } printf("f[%d]=%f",F,f[F]); break; } else { function2(); break; } } } if(f[N+1]>=f[F]&&f[N+1]<=f[r]) { for(i=0;i<N-1;i++) { d[i][o]=d[i][N+1]; } f[o]=f[N+1]; y=0; for(i=0;i<N;i++) { y=y+(f[i]-f[N])*(f[i]-f[N]); } y=sqrt(y/(N*1.0)); if(y<v) { for(i=0;i<N-1;i++) { printf("d[%d][%d]=%f",i,F,d[i][F]); printf("\n"); } printf("f[%d]=%f",F,f[F]); break; } else { function2(); break; } } if(f[N+1]>f[r]) { if(f[N+1]>f[o])//压缩 { D=o; } else { D=N+1; } for(i=0;i<N-1;i++) { d[i][N+3]=d[i][N]+(d[i][D]-d[i][N])*u; } function1(N+3); if(f[N+3]<=f[D]) { for(i=0;i<N-1;i++) { d[i][o]=d[i][N+3]; } f[o]=f[N+3]; y=0; for(i=0;i<N;i++) { y=y+(f[i]-f[N])*(f[i]-f[N]); } y=sqrt(y/(N*1.0)); if(y<v) { for(i=0;i<N-1;i++) { printf("d[%d][%d]=%f",i,F,d[i][F]); printf("\n"); } printf("f[%d]=%f",F,f[F]); break; } else { function2(); break; } } else { for(i=0;i<N-1;i++) { for(j=0;j<N;j++) { d[i][j]=d[i][j]+0.5*(d[i][F]-d[i][j]);//收缩 } } for(i=0;i<N;i++) { function1(i); } y=0; for(i=0;i<N;i++) { y=y+(f[i]-f[N])*(f[i]-f[N]); } y=sqrt(y/(N*1.0)); if(y<v) { for(i=0;i<N-1;i++) { printf("d[%d][%d]=%f",i,F,d[i][F]); printf("\n"); } printf("f[%d]=%f",F,f[F]); break; } else { function2(); break; } } } } } void main() { for(i=0;i<N;i++) { function1(i); } function2(); }