PROGRAM d10r5(input,output);
(* driver for routine AMOEBA *)
CONST
   np=3;
   mp=4;
   ftol=1.0e-6;
TYPE
   glmpnp = ARRAY [1..mp,1..np] OF real;
   glmp = ARRAY [1..mp] OF real;
   glnp = ARRAY [1..np] OF real;
VAR
   i,iter,j,ndim : integer;
   x : glnp;
   y : glmp;
   p : glmpnp;

(*$I MODFILE.PAS *)
(*$I BESSJ0.PAS *)

FUNCTION func(x: glnp): real;
(* calling function must define type
TYPE
   glnp = ARRAY [1..np] OF real;
where np is the physical dimension of the argument x. *)
BEGIN
   func := 0.6-bessj0(sqr(x[1]-0.5)+sqr(x[2]-0.6)+sqr(x[3]-0.7))
END;

(*$I AMOEBA.PAS *)

BEGIN
   p[1,1] := 0.0; p[1,2] := 0.0; p[1,3] := 0.0;
   p[2,1] := 1.0; p[2,2] := 0.0; p[2,3] := 0.0;
   p[3,1] := 0.0; p[3,2] := 1.0; p[3,3] := 0.0;
   p[4,1] := 0.0; p[4,2] := 0.0; p[4,3] := 1.0;
   ndim := np;
   FOR i := 1 to mp DO BEGIN
      FOR j := 1 to np DO BEGIN
         x[j] := p[i,j]
      END;
      y[i] := func(x)
   END;
   amoeba(p,y,ndim,ftol,iter);
   writeln;
   writeln('Iterations: ',iter:3);
   writeln('Vertices of final 3-d simplex and');
   writeln('function values at the vertices:');
   writeln;
   writeln('i':3,
      'x[i]':10,'y[i]':12,'z[i]':12,'function':14);
   writeln;
   FOR i := 1 to mp DO BEGIN
      write(i:3);
      FOR j := 1 to np DO write(p[i,j]:12:6);
      writeln(y[i]:12:6)
   END;
   writeln;
   writeln('True minimum is at (0.5,0.6,0.7)')
END.