anchor.diff:          
       end subroutine geo_LUT
       
  -    ! interpolate boundaries of fieldout to fieldin using geolut
  +    !>------------------------------------------------------------
  +    !!  Interpolate boundaries of fieldout to fieldin using geolut
  +    !!
  +    !!------------------------------------------------------------
       subroutine boundary_interpolate(fieldout, fieldin, geolut)
           implicit none
           real,intent(inout)::fieldout(:,:,:)
  
positive.code: 
      end subroutine boundary_interpolate
      
      ! interpolate fieldout to fieldin using geolut.  
      ! if boundary_only is true, call boundary_interpolate instead
      ! loops over y,z,x but geolut is only defined over x,y (for now)
      subroutine geo_interp(fieldout,fieldin,geolut,boundary_only)
          implicit none
          real,intent(inout)::fieldout(:,:,:)
  
negative.code:         endif
      end subroutine geo_interp
      
      !>------------------------------------------------------------
      !!  Interpolate fieldout to fieldin using geolut. 
      !!
      !!------------------------------------------------------------
      subroutine geo_interp2d(fieldout, fieldin, geolut)
          implicit none
          real, dimension(:,:),intent(inout) :: fieldout
  