风花雪月 发表于 2006-10-11 16:06

Powell法优化子程序

本程序包含5个C文件
mpowell.c
powell.c
funct.c
jtf.c
hjfgf.c


本程序由heizi友情提供

风花雪月 发表于 2006-10-11 16:08

mpowell.c代码如下：

#include "powell.c"
main()
{double p[]={1,2};
double ff,x;
ff=powell(p,0.3,0.001,0.0001,2,x);
printf("x=%f,x=%f,ff=%f\n",x,x,ff);
getch();
}


powell.c代码如下：
#include "hjfgf.c"
double oneoptim(double x0[],double s[],double h0,double epsg,int n,double x[])
{double *a,*b,ff;
a=(double *)malloc(n*sizeof(double));
b=(double *)malloc(n*sizeof(double));
jtf(x0,h0,s,n,a,b);
ff=gold(a,b,epsg,n,x);
free(a);
free(b);
return (ff);
}
double powell(double p[],double h0,double eps,double epsg,int n,double x[])
{int i,j,m;
double *xx,*ss,*s;
double f,f0,f1,f2,f3,fx,dlt,df,sdx,q,d;
ss=(double *)malloc(n*(n+1)*sizeof(double));
s=(double *)malloc(n*sizeof(double));
for(i=0;i<n;i++)
{for(j=0;j<=n;j++)
   *(ss+i*(n+1)+j)=0;
   *(ss+i*(n+1)+i)=1;
}
for(i=0;i<4;i++)
xx=(double *)malloc(n*sizeof(double));
for(i=0;i<n;i++)
*(xx+i)=p;
for(;;)
{for(i=0;i<n;i++)
    {*(xx+i)=*(xx+i);
   x=*(xx+i);
    }
   f0=f1=objf(x);
   dlt=-1;
   for(j=0;j<n;j++)
    {for(i=0;i<n;i++)
      {*(xx+i)=x;
       *(s+i)=*(ss+i*(n+1)+j);
      }
   f=oneoptim(xx,s,h0,epsg,n,x);
   df=f0-f;
   if(df>dlt)
      {dlt=df;
       m=j;
      }
    }
    sdx=0;
    for(i=0;i<n;i++)
   sdx=sdx+fabs(x-(*(xx+i)));
    if(sdx<eps)
   {free(ss);
      free(s);
      for(i=0;i<4;i++)
      free(xx);
      return(f);
   }
    for(i=0;i<n;i++)
   *(xx+i)=x;
    f2=f;
    for(i=0;i<n;i++)
    {*(xx+i)=2*(*(xx+i)-(*(xx+i)));
   x=*(xx+i);
    }
    fx=objf(x);
    f3=fx;
    q=(f1-2*f2+f3)*(f1-f2-dlt)*(f1-f2-dlt);
    d=0.5*dlt*(f1-f3)*(f1-f3);
    if((f3<f1)||(q<d))
   {if(f2<=f3)
       for(i=0;i<n;i++)
      *(xx+i)=*(xx+i);
      else
       for(i=0;i<n;i++)
      *(xx+i)=*(xx+i);
   }
    else
    {for(i=0;i<n;i++)
   {*(ss+(i+1)*(n+1))=x-(*(xx+i));
      *(s+i)=*(ss+(i+1)*(n+1));
   }
   f=oneoptim(xx,s,h0,epsg,n,x);
   for(i=0;i<n;i++)
    *(xx+i)=x;
   for(j=m+1;j<=n;j++)
   for(i=0;i<n;i++)
    *(ss+i*(n+1)+j-1)=*(ss+i*(n+1)+j);
    }
}
}

风花雪月 发表于 2006-10-11 16:10

funct.c代码如下：
#include "stdio.h"
#include "stdlib.h"
#include "math.h"
double objf(double x[])
{double ff;
ff=x*x+x*x-x*x-10*x-4*x+60;
return(ff);
}


jtf.c代码如下：
#include "funct.c"
void jtf(double x0[],double h0,double s[],int n,double a[],double b[])
{int i;
double *x,h,f1,f2,f3;
for(i=0;i<3;i++)
x=(double *)malloc(n*sizeof(double));
h=h0;
for(i=0;i<n;i++)
*(x+i)=x0;
f1=objf(x);
for(i=0;i<n;i++)
*(x+i)=*(x+i)+h*s;
f2=objf(x);
if(f2>=f1)
{h=-h0;
    for(i=0;i<n;i++)
    *(x+i)=*(x+i);
   f3=f1;
    for(i=0;i<n;i++)
    {*(x+i)=*(x+i);
   *(x+i)=*(x+i);
    }
   f1=f2;
   f2=f3;
   }
   for(;;)
   {h=2*h;
   for(i=0;i<n;i++)
   *(x+i)=*(x+i)+h*s;
   f3=objf(x);
   if(f2<f3) break;
   else
    { for(i=0;i<n;i++)
       {*(x+i)=*(x+i);
      *(x+i)=*(x+i);
       }
      f1=f2;
      f2=f3;
    }
   }
   if(h<0)
    for(i=0;i<n;i++)
    {a=*(x+i);
   b=*(x+i);
    }
   else
    for(i=0;i<n;i++)
    {a=*(x+i);
   b=*(x+i);
   }
   for(i=0;i<3;i++)
   free(x);
}


hjfgf.c代码如下：
#include "jtf.c"
double gold(double a[],double b[],double eps,int n,double xx[])
{int i;
double f1,f2,*x,ff,q,w;
for(i=0;i<2;i++)
x=(double *)malloc(n*sizeof(double));
for(i=0;i<n;i++)
{*(x+i)=a+0.618*(b-a);
   *(x+i)=a+0.382*(b-a);
}
f1=objf(x);
f2=objf(x);
do
   {if(f1>f2)
   {for(i=0;i<n;i++)
      {b=*(x+i);
       *(x+i)=*(x+i);
       }
   f1=f2;
   for(i=0;i<n;i++)
      *(x+i)=a+0.382*(b-a);
   f2=objf(x);
   }
    else
   { for(i=0;i<n;i++)
       {a=*(x+i);
       *(x+i)=*(x+i);}
   f2=f1;
    for(i=0;i<n;i++)
   *(x+i)=a+0.618*(b-a);
    f1=objf(x);
   }
q=0;
for(i=0;i<n;i++)
   q=q+(b-a)*(b-a);
w=sqrt(q);
}while(w>eps);
for(i=0;i<n;i++)
   xx=0.5*(a+b);
ff=objf(xx);
for(i=0;i<2;i++)
free(x);
return(ff);
}

lzjld502 发表于 2006-11-14 07:08

xiexiele hehe:@)

zhoyi 发表于 2007-6-8 09:56

注释到没有，很难看明白。要是看明白了，也会编写了 。建议上传时把注释也写上。不管怎样，还是谢谢louzhu。

hnsqzrj 发表于 2007-6-8 20:36

怎么不能运行

怎么不能运行

yangwenli 发表于 2008-12-14 23:00

感谢楼主提供的程序:@D

to楼上:编译在一起就可运行:)

crywolf911 发表于 2009-11-10 20:27

问题是怎么编译在一起啊，我是新手啊~！

modest_xiu 发表于 2012-5-7 17:37

谢谢~~~楼主辛苦
页: [1]
查看完整版本: Powell法优化子程序