ffnewuoa_edp.edp

/*  README NEWUOA:      This is the Fortran version of NEWUOA. Its purpose is to seek the least value of a function F of several variables, when derivatives are not available, where F is specified by the user through a subroutine called CALFUN. The algorithm is intended to change the variables to values that are close to a local minimum of F. The user, however, should assume responsibility for finding out if the calculations are satisfactory, by considering carefully the values of F that occur. The method is described in the report "The NEWUOA software for unconstrained optimization without derivatives", which is available on the web at www.damtp.cam.ac.uk, where you have to click on Numerical Analysis and then on Reports, the number of the report being NA2004/08. Let N be the number of variables. The main new feature of the method is that quadratic models are updated using only about NPT=2N+1 interpolation conditions, the remaining freedom being taken up by minimizing the Frobenius norm of the change to the second derivative matrix of the model.      The new software was developed from UOBYQA, which also forms quadratic models from interpolation conditions. That method requires NPT=(N+1)(N+2)/2 conditions, however, because they have to define all the parameters of the model. The least Frobenius norm updating procedure with NPT=2N+1 is usually much more efficient when N is large, because the work of each iteration is much less than before, and in some experiments the number of calculations of the objective function seems to be only of magnitude N.      The attachments in sequence are a suitable Makefile, followed by a main program and a CALFUN routine for the Chebyquad problems, in order to provide an example for testing. Then NEWUOA and its five auxiliary routines, namely NEWUOB, BIGDEN, BIGLAG, TRSAPP and UPDATE, are given. Finally, the computed output that the author obtained for the Chebyquad problems is listed.      The way of calling NEWUOA should be clear from the Chebyquad example and from the comments of that subroutine. It is hoped that the software will be helpful to much future research and to many applications. There are no restrictions on or charges for its use. If you wish to refer to it, please cite the DAMTP report that is mentioned above, which has been submitted for publication in the proceedings of the 40th Workshop on Large Scale Nonlinear Optimization (Erice, Italy, 2004). December 16th, 2004                    M.J.D. Powell (mjdp@cam.ac.uk) */ 

load "ffnewuoa" 
/*       SUBROUTINE CALFUN (N,X,F,IWF)       IMPLICIT REAL*8 (A-H,O-Z)       DIMENSION X(*),Y(10,10),IWF(*)       DO 10 J=1,N       Y(1,J)=1.0D0    10 Y(2,J)=2.0D0*X(J)-1.0D0       DO 20 I=2,N       DO 20 J=1,N    20 Y(I+1,J)=2.0D0*Y(2,J)*Y(I,J)-Y(I-1,J)       F=0.0D0       NP=N+1       IW=1       DO 40 I=1,NP       SUM=0.0D0       DO 30 J=1,N    30 SUM=SUM+Y(I,J)       SUM=SUM/DFLOAT(N)       IF (IW .GT. 0) SUM=SUM+1.0D0/DFLOAT(I*I-2*I)       IW=-IW    40 F=F+SUM*SUM       RETURN       END */ 
int N=2; 
real[int,int] Y(N+1,N); 
real[int] X(N); 
func real J(real[int] & X) 
{ 
  Y(0,:)=1.; 
  Y(1,:)=2.*X; 
  for(int i=2;i<=N;++i) 
    for(int j=0;j<N;++j) 
      Y(i,j) = Y(1,j) * Y(i-1,j)- Y(i-2,j); 

  real f=0; 
  int np=N+1,iw=1; 
  for(int i=0;i<=N;i++) 
    { 
      real s=Y(i,:).sum/N; 
      int i1=i+1; 
      if(iw>0) s+= 1./(i1*i1-2.*i1); 
      iw=-iw; 
      f += s*s; 
    } 
  return 2.14+f; 
} 
real delta=1; 
mesh Th=square(30,30,[(x-0.5)*delta,(y-0.5)*delta]); 
fespace Vh(Th,P1); 
Vh u; 
for(int i=0;i<Th.nv;++i) 
  {X[0]=Th(i).x;X[1]=Th(i).y; 
    u[][i]=J(X); 
  } 
plot(u,wait=1); 

for(int i=0;i<N;++i) 
  X[i]=(i+1)/(N+1.); 
real mincost=newuoa(J,X,rhobeg=2*X[0],rhoend=1e-6,npt=2*N+1); 
cout << " min " << mincost << "  at :" << X << endl;