source/maxlik/src_MD_T_maxlik-NEWCORR-PMF/testamebsa.f

   1 c      implicit real*8 (a-h,o-z)
   2       parameter (ndim=2)
   3       real p(ndim+1,ndim),y(ndim+1),pb(ndim),pom(ndim),yb,ftol,temptr
   4       real f
   5       external f
   6       temptr=0.1d0
   7       iter=100
   8       ftol=1.0d-4
   9       yb=1.0d6
  10       do i=1,ndim
  11         do j=1,ndim+1
  12           p(j,i)=1.0d1
  13         enddo
  14       enddo
  15
  16       p(1,1)=-100.0d0
  17       p(1,2)=-50.0d0
  18       p(2,1)=-80.0d0
  19       p(2,2)=6.0d1
  20       p(3,1)=4.0d0
  21       p(3,2)=20.0d0
  22       do i=1,ndim+1
  23         do j=1,ndim
  24           pom(j)=p(i,j)
  25         enddo
  26         y(i)=f(pom)
  27       enddo
  28
  29       call amebsa(p,y,ndim+1,ndim,ndim,pb,yb,ftol,f,iter,temptr)
  30
  31       print *,"pb",pb
  32       print *,"yb",yb
  33       print *,"iter",iter
  34
  35       stop
  36       end
  37
  38       real function f(x)
  39       real x(2)
  40       f=0.5*((x(1)-1)**2+(x(2)+1)**2)+10*(1-cos(x(1)-1))*(1-cos(x(2)+1))
  41       return
  42       end
  43
  44 CcccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccC
  45       REAL FUNCTION ran1(idum)
  46       INTEGER idum,IA,IM,IQ,IR,NTAB,NDIV
  47       REAL AM,EPS,RNMX
  48       PARAMETER (IA=16807,IM=2147483647,AM=1./IM,IQ=127773,IR=2836,
  49      *     NTAB=32,NDIV=1+(IM-1)/NTAB,EPS=1.2e-7,RNMX=1.-EPS)
  50 C     "Minimal" random number generator of Park and Miller with Bays-Durham shuffle and
  51 C     added safeguards. Returns a uniform random deviate between 0.0 and 1.0 (exclusive of
  52 C     the endpoint values). Call with idum a negative integer to initialize; thereafter, do not
  53 C     alter idum between successive deviates in a sequence. RNMX should approximate the largest
  54 C     floating value that is less than 1. Recommended for seqences smaller than 100e6 (i.e. 5%
  55 C     of ran1's period.
  56 C     (from Press etal, Numerical Recipes in Fortran 77)
  57 CcccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccC
  58       INTEGER j,k,iv(NTAB),iy
  59       SAVE iv,iy
  60       DATA iv /NTAB*0/, iy /0/
  61       if (idum.le.0.or.iy.eq.0) then  ! initialize
  62          idum=max(-idum,1)            ! load the shuffle table (after 8 warmups)
  63          do j=NTAB+8,1,-1
  64             k=idum/IQ
  65             idum=IA*(idum-k*IQ)-IR*k
  66             if (idum.lt.0) idum=idum+IM
  67             if (j.le.NTAB) iv(j)=idum
  68          enddo
  69          iy=iv(1)
  70       endif
  71       k=idum/IQ                       ! start here when not initializing
  72       idum=IA*(idum-k*IQ)-IR*k
  73       if (idum.lt.0) idum=idum+IM
  74       j=1+iy/NDIV
  75       iy=iv(j)
  76       iv(j)=idum
  77       ran1=min(AM*iy,RNMX)
  78       return
  79       END