anchor.diff:    IMPLICIT NONE
   
     ! This interface block is needed to pass dynamically allocated arrays
  +  ! Useful comments will be found near the actual subroutines.
     interface
        subroutine inhist(xfn, firstn, number, nextn)
          double precision :: xfn(:,:)
  
positive.code:   !     Start variable declaration
    !**********************************************************************
  
    ! lx_ing is the log_2(ngx_i)
    ! ngx_i is the number of fluid grid nodes in the i-direction
    ! nfsize is the number of solid nodes
    ! nfsize2 is the number of elements
    integer lxng,lyng,lzng,ngx,ngy,ngz,nfsize,nfsize2
    ! ngx_im1 is the number of grid nodes in the i-direction, minus one.
    ! This is only really useful for the z-direction.
    ! nbx_i is ngx_i+2, which is only really useful in the x and y directions.
    INTEGER NGXM1,NGYM1,NGZM1,NBX,NBY,NBZ
    ! These are just the ngx_i in double precision
    double precision :: FLNGX,FLNGY,FLNGZ
    PARAMETER(LXNG=$lngx$,LYNG=LXNG,LZNG=LXNG)
    PARAMETER(NGX=2**LXNG,NGY=2**LYNG,NGZ=2**LZNG)
    PARAMETER(NBX=NGX+2,NBY=NGY+2,NBZ=NGZ+2)
    PARAMETER(NGXM1=NGX-1,NGYM1=NGY-1,NGZM1=NGZ-1)
    PARAMETER(FLNGX=NGX,FLNGY=NGY,FLNGZ=NGZ)
    integer nrects
    integer nrect
    double precision fnrect
    INTEGER KLOK,KLOK1,KLOKEND,NSTEP
    double precision :: T,H,h64,TD,VSC,TIME,RHO,PI,RADX,FOSTAR
    double precision, allocatable :: rad(:)
  
negative.code:   ALLOCATE(NEXTN(1:NFSIZE),XFN(1:3,1:NFSIZE),FRC(1:3,1:NFSIZE),Foptical(1:3,1:NFSIZE))
    ALLOCATE(elmnew(1:3,1:NFSIZE2))
    ALLOCATE(shpint(1:3,1:NFSIZE2),shpfs(1:7,1:NFSIZE2))
    allocate(xpi(3,nrectnodes))
  
    pi = 3.14159265358979323846d0 ! Taken from Wikipedia; 20 digits
    !     Physical parameters -- using cgs system
  