anchor.diff:      ! Allocate domains and prepare comms
       call md%setup()
   
  +    ! Apply repeated smooths. 
  +    grid_smooth_iterations = nint(10 * mesh_refine**2)
  +
       ! Set initial conditions
       do j = 1, size(md%domains)
  -        call set_initial_conditions(md%domains(j))
  +        call set_initial_conditions(md%domains(j), grid_smooth_iterations)
       end do
   
       ! These steps are important in complex nested models, but not really needed here.
  
positive.code:     use logging_mod, only: log_output_unit
      implicit none
  
      contains 
  
      subroutine set_initial_conditions(domain)            
          class(domain_type), intent(inout):: domain
          integer(ip):: i, j
          real(dp), allocatable:: x(:), y(:)
          type(multi_raster_type) :: elevation_rast
  
  
negative.code:     ! grid size (number of x/y cells)
      integer(ip), parameter :: global_nx(2) = [945_ip, 615_ip] 
  
      integer(ip):: j, grid_smooth_iterations
  
      ! Set the model name
      md%output_basedir = './OUTPUTS/'
  