anchor.diff:          ! While the result is good, there are many ways in which to question it.
           nd = 7
       else
  -        ! This case has a domain res ~ 1.5x1.5m around the monai inundation peak. 
  +        ! This case has a domain res ~ 1.5x1.5m around the monai inundation peak (with mesh_refine = 1.0_dp)
           !
           ! It leads to a peak of >28.0m (vs 31.7 obs, although the latter varies depending on which
           ! dataset is used). So the very_high_res_monai case does better.
  
positive.code:         real(dp), allocatable :: random_uniform(:)
          ! Use this to check the effect of very small perturbations of the initial condition
          ! Eventually  a tiny perturbation will lead to differences in the results. This seems comparable
          ! to the differences that eventually emerge in openmp/coarray versions, which seem
          ! to begin with differences of this order in elevation/stage, I suppose due to round-off
          ! scale differences in x/y etc that emerge from grid partitioning.
          real(dp), parameter :: random_perturbation_scale = 0.0e-10_dp
  
          ! Stage
  
negative.code:     ! A more detailed Monai domain 
      call md%domains(5)%match_geometry_to_parent(&
          parent_domain=md%domains(4), &
          lower_left  = [139.4150_dp, 42.0968_dp], &
          upper_right = [139.4264_dp, 42.1020_dp], &
          dx_refinement_factor = nest_ratio, &
  