anchor.diff:                write(52,1005) j,soln(j,1),soln(j,2)
     515     continue
             write(52,*)
  -C         Save solution in best models if in top NSOLN obj values
  +
  +C         Find the worst saved solution
  +          test = 0.0d0
             do 517 j = 1,nsoln
  -              if (obj.lt.obj0(j)) then
  -                  do 516 jj = 1,nflt
  -                      best0(jj,1,j) = soln(jj,1)
  -                      best0(jj,2,j) = soln(jj,2)
  -  516             continue
  -                  obj0(j) = obj
  -                  goto 518
  +              if (obj0(j).gt.test) then
  +                  test = obj0(j) 
  +                  jj = j
                 endif
     517     continue
  -  518     continue
  -C         If solution is best, save it
  +          j = jj
  +C         Save solution in best models array if in top NSOLN objective values
  +          if (obj.lt.obj0(j)) then
  +              do 516 jj = 1,nflt
  +                  best0(jj,1,j) = soln(jj,1)
  +                  best0(jj,2,j) = soln(jj,2)
  +  516         continue
  +              obj0(j) = obj
  +          endif
  +C         If solution is better than previous, save it
             if (obj.lt.objmin) then
                 do 519 j = 1,nflt
                     best(j,1) = soln(j,1)
                     best(j,2) = soln(j,2)
                     objmin = obj
     519         continue
  +C         When using Metropolis-Hastings, allow solution chance to get worse
  +          elseif (algrthm.eq.'metropolis') then
  +              test = dexp((objmin-obj)/temp)
  +              random = ran2(idum)
  +              if (random.lt.test) then
  +                  do 527 j = 1,nflt
  +                      best(j,1) = soln(j,1)
  +                      best(j,2) = soln(j,2)
  +                      objmin = obj
  +  527             continue
  +              endif
             endif
   C         Cool off
  -          if (temp.gt.1e-4) then
  +          if (mod(i,reset).lt.tconst) then
  +              temp = temp
  +          elseif (temp.gt.1e-4) then
                 temp = temp*cool
             else
                 temp = temp*0.995d0
  +          endif
  +          if (mod(i,reset).eq.0) then
  +              temp = 1
             endif
     520 continue
         write(51,*)
  
positive.code: C----------------------------------------------------------------------C
  
        SUBROUTINE readannf(inityp,slip0,rake0,slip,rake,temp0,cool,
       1                    nitmax,annf,nflt,FLTMAX,nsoln)
  C----
  C Read annealing control file if it exists
  C----
        IMPLICIT none
        INTEGER nflt,FLTMAX
        CHARACTER*80 annf
        INTEGER inityp,nitmax
        INTEGER nsoln
        REAL*8 slip0,rake0,slip(FLTMAX,3),rake(FLTMAX,3),temp0,cool
        LOGICAL ex
        CHARACTER*120 line
        CHARACTER*1 dum
        INTEGER i,j,ns,nr
  C If annealing file not specified, use some default values
    551 if (annf.eq.'none') then
            inityp = 0
  
negative.code: C----
  C Write nsoln best solutions to file best_nsoln.dat
  C----
  C      write(nsolnc,'(I0.2)') nsoln
  C      open(unit=53,file='best_'//nsolnc//'.dat',status='unknown')
        open(unit=53,file='best0.dat',status='unknown')
        write(53,'(50F12.6)') (obj0(j),j=1,nsoln)
        do 521 i = 1,nflt
            write(53,'(50F12.6)') (best0(i,1,j),j=1,nsoln)
    521 continue
  