anchor.diff:          do j = 1, domain%nx(2)
               y = domain%y(j)
   
  +            domain%U(:,j,STG) = exp(sin(2*pi*x/10.0_dp))*cos(2*pi*y/10.0_dp)
  +
               call elevation_rast%get_xy(x, y, domain%U(:,j,ELV), domain%nx(1), bilinear=1_ip)
  -
  -            domain%U(:,j,STG) = exp(sin(2*pi*x/10.0_dp))*cos(2*pi*y/10.0_dp)
           end do
   
           deallocate(x,y)
          
  -        ! Wet everywhere 
  +        ! Wet everywhere and not too shallow
           domain%U(:,:,ELV) = min(domain%U(:,:,ELV), -100.0_dp)
   
           ! Ensure topography is smooth
  
positive.code: module local_routines 
      !! Test the convergence of the code in a periodic domain
      use global_mod, only: dp, ip, wall_elevation, pi, charlen
      use domain_mod, only: domain_type, STG, UH, VH, ELV
      use read_raster_mod, only: multi_raster_type
  
negative.code:     do j = 1, size(md%domains)
          call set_initial_conditions(md%domains(j))
      end do
  
      ! These steps are important in complex nested models, but not really needed here.
      call md%make_initial_conditions_consistent()
      call md%set_null_regions_to_dry()
  
      ! Fixed timestep
  