powell优化算法的过程

Powell优化算法是利用仪器测井理建立误差函数（非相关函数），借助Powell方向加速法求出非相关函数达到最小时的解，对于气，水两相流动，从预设的气，水流量初始值出发，碰手沿不同的广向进行搜索，可求出气，水两相流动中可能碧含最大产量。与目前常用的生产测井解释方法相比，文中提出的方法具有精度高，实用性强等优悔吵笑点，在测井曲线有缺陷时，仍有可能得到较好的结果。

powell.c代码如下：

CODE:

＃i nclude "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[4],*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=0i<ni++)

{for(j=0j<=nj++)

*(ss+i*(n+1)+j)=0

*(ss+i*(n+1)+i)=1

}

for(i=0i<4i++)

xx[i]=(double *)malloc(n*sizeof(double))

for(i=0i<ni++)

*(xx[0]+i)=p[i]

for()

{for(i=0i<ni++)

{*(xx[1]+i)=*(xx[0]+i)

x[i]=*(xx[1]+i)

}

f0=f1=objf(x)

dlt=-1

for(j=0j<nj++)

{for(i=0i<ni++)

{*(xx[0]+i)=x[i]

*(s+i)=*(ss+i*(n+1)+j)

}

f=oneoptim(xx[0],s,h0,epsg,n,x)

df=f0-f

if(df>dlt)

{dlt=df

m=j

}

}

sdx=0

for(i=0i<ni++)

sdx=sdx+fabs(x[i]-(*(xx[1]+i)))

if(sdx<eps)

{free(ss)

free(s)

for(i=0i<4i++)

free(xx[i])

return(f)

}

for(i=0i<ni++)

*(xx[2]+i)=x[i]

f2=f

for(i=0i<ni++)

{*(xx[3]+i)=2*(*(xx[2]+i)-(*(xx[1]+i)))

x[i]=*(xx[3]+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=0i<ni++)

*(xx[0]+i)=*(xx[2]+i)

else

for(i=0i<ni++)

*(xx[0]+i)=*(xx[3]+i)

}

else

{for(i=0i<ni++)

{*(ss+(i+1)*(n+1))=x[i]-(*(xx[1]+i))

*(s+i)=*(ss+(i+1)*(n+1))

}

f=oneoptim(xx[0],s,h0,epsg,n,x)

for(i=0i<ni++)

*(xx[0]+i)=x[i]

for(j=m+1j<=nj++)

for(i=0i<ni++)

*(ss+i*(n+1)+j-1)=*(ss+i*(n+1)+j)

}

}

}

或者%powell算法，用于寻找无约束最优值点

function powel=powell(x0,n,q,r,h,a)

d=eye(n)%n个线性无关的初始搜索方向

k=1

kk=1

xx(1,1:n)=x0

while (kk)

y(1,1:n)=xx(k,1:n)

for j=1:n

s(j)=HJ(y(j,1:n),d(j,1:n),q) %调用0.618算法

y(j+1,1:n)=y(j,1:n)+s(j).*d(j,1:n)

end

d(n+1,1:n)=y(n+1,1:n)-y(1,1:n)

if (norm(d(n+1,1:n),2)<r)

kk=0

break

else

ww=0

m=1

for i=1:n

gg=ff(y(i,1:n),q)-ff(y(i+1,1:n),q)

if (gg>=ww)

m=i

end

end

cha=ff(y(1,1:n),q)-2*ff(y(n+1,1:n),q)+ff(2*y(n+1,1:n)-y(1,1:n),q)

cha1=2*(ff(y(m,1:n),q)-ff(y(m+1,1:n),q))

if (cha<cha1)

s(n+1)=HJ(y(n+1,1:n),h,a,d(n+1,1:n),q)

xx(k+1,1:2)=y(n+1,1:n)+s(n+1).*d(n+1,1:n)

for j=m+1:n

d(j,1:n)=d(j+1,1:n)

end

k=k+1

else

xx(k+1,1:n)=y(n+1,1:n)

k=k+1

end

end

end

powel=y(n+1,1:n)

function w=HJ(x0,h,d,dd,q) %0.618算法

[a,b]=JTF(x0,h,d,dd,q) %调用进退法算法，确定范围

r=0.618

r1=a+(1-r)*(b-a)

r2=a+r*(b-a)

y1=ff(x0+r1.*dd,q)

y2=ff(x0+r2.*dd,q)

k=1

while (abs(r1-r2)>=0.1)

if y1<y2

b=r2

r2=r1

y2=y1

r1=a+(1-r)*(b-a)

y1=ff(x0+r1.*dd,q)

else

a=r1

r1=r2

y1=y2

r2=a+r*(b-a)

y2=ff(x0+r2.*dd,q)

end

end

w=(r1+r2)/2

%进退法

function [a,b]=JTF(x0,h,d,dd,q)

r0=0

y0=ff(x0+r0.*dd,q)

k=0

l=1

while (l)

r1=r0+h

y1=ff(x0+r1.*dd,q)

if y1<y0

h=d*h

r=r0

r0=r1

y0=y1

else

if k==0

h=-h

r=r0

else

l=0

break

end

end

k=k+1

end

a=min(r,r1)

b=max(r,r1)

Powell是好胡人名，鲍威尔。

以下文章首次提出了枯肢Powell算法：

Powell, M. J. D. (1964). "An efficient method for finding the minimum of a function of several variables without calculating derivatives."友败拦 The Computer Journal 7(2): 155.

---