anchor.diff:              !! but sometimes we want several grids (e.g. east/north velocity)
           type(c_ptr) :: grid_update_context = c_null_ptr
               !! Used to store any data required by the get_grid_at_time function
  +        logical :: even_x_spacing, even_y_spacing
  +            !! If x/y are evenly spaced then we can quickly find the nearest index
   
           procedure(get_grid_at_time), pointer, nopass :: get_grid_at_time => NULL()
               !! User-provided subroutine to update the g1 and g2 values
  
positive.code:             ng_local = 1_ip
          end if
  
          do i = 1, (size(x, kind=ip)-1)
              if(x(i+1) <= x(i)) then
                  write(log_output_unit, *) "x must be increasing", x
                  call generic_stop
              end if
          end do
          do i = 1, (size(y, kind=ip)-1)
              if(y(i+1) <= y(i)) then
                  write(log_output_unit, *) "y must be increasing", y
                  call generic_stop
              end if
          end do
  
          ! Grid size
  
negative.code:             !print*, xs(i), gsi%x(i0), gsi%x(i1), i0, i1, w0_x, w1_x
  
              ! Get y-indices for bilinear interpolation: j0 < j1
              if(gsi%even_y_spacing) then
                  j0 = nint( (ys(i) - gsi%y(1)) / ((gsi%y(gsi%nx(2)) - gsi%y(1))/(gsi%nx(2)-1.0_dp)) ) + 1
                  j0 = min(max(j0, 1), gsi%nx(2))
              else
                  call nearest_index_sorted(gsi%nx(2), gsi%y, ys(i), j0)
              end if
              if(ys(i) > gsi%y(j0)) then
                  j1 = min(j0+1, gsi%nx(2))
              else
  