anchor.diff:              !! 'nesting_boundary_flux_*')
   
           !
  -        ! Output file content
  +        ! Output file content -- bespoke binary format (that is no longer properly supported, consider removing)
           !
           character(len=charlen):: metadata_ascii_filename
               !! Name of ascii file where we output metadata
  
positive.code:         !
          ! Setup for Chezy or Manning friction
          !
          if(domain%friction_type == 'manning') then
              ! For efficiency we ensure the friction power is known at compile time
              ! Previously it was set here
  
          else if(domain%friction_type == 'chezy') then
              ! For efficiency we ensure the friction power is known at compile time
              ! Previously it was set here
  
              if(domain%timestepping_method == 'cliffs') then
                  write(log_output_unit, *) "ERROR: chezy_friction has not yet been implemented for our cliffs solver"
                  call generic_stop
              end if
  
          else
              write(log_output_unit, *) ' Error: Unrecognized value of domain%friction_type '
              call generic_stop
          end if
  
  #ifdef DEBUG_ARRAY
          ! Optional for debugging -- this will be written to the netcdf file
  
negative.code:             !! U holds main flow variables -- STG, UH, VH, ELV
              !! First 2 dimensions = space, 3rd = number of quantities
          real(dp), allocatable :: max_U(:,:,:)
              !! Store various flow maxima. FIXME: Consider storing in reduced precision.
          !
          ! Multi-dimensional arrays below are only required for nonlinear solver. Would be possible
          ! to further reduce memory usage, but not completely trivial. Not all of these are required
  