anchor.diff:          read(10) array
           close(10) 
           ! wrap array around for periodicity
  -        array(0,:,) = array(nx,:,)
  -        array(nx+1,:,) = array(1,:,)
  -        array(:,0,) = array(:,ny,)
  -        array(:,ny+1,) = array(:,1,)
  +        array(0,:) = array(nx,:)
  +        array(nx+1,:) = array(1,:)
  +        array(:,0) = array(:,ny)
  +        array(:,ny+1) = array(:,1)
       else
           array = default
       endif
  
positive.code:     end do
  
  ! wrap fields around for periodic simulations
      wrap_fields_3D(unew,nx,ny,layers)
      wrap_fields_3D(vnew,nx,ny,layers)
      wrap_fields_3D(hnew,nx,ny,layers)
      wrap_fields_2D(etanew,nx,ny)
  
  
  !    Accumulate average fields
  
negative.code:     integer nx,ny, i,j, maxits, n
      double precision dt
      double precision rjac, eps
      double precision rhs(nx,ny) 
      double precision res(nx,ny)
      double precision norm, norm0, norm_old
      double precision relax_param
  