anchor.diff:              ! Relax snow height to maximum (eg, 5 m)
               snow_to_ice  = dmax1(0.d0,hsnow-hsmax) 
               hsnow        = hsnow - snow_to_ice
  -            massbal_ice  = snow_to_ice/tstic - melt_ice + refrozen_rain   ! Use to force ice sheet model
  +            massbal_ice  = snow_to_ice/tstic - melt_ice + refr   ! Use to force ice sheet model
               new_ice      = massbal_ice*tstic ! update new ice budget: remove or add ice
   
               if (.not. bnd%massbal) then
                   where (land_ice_ocean == 2)
  -                    massbal = massbal_snow + massbal_ice
  +                    massbal = massbal_snow + massbal_ice - snow_to_ice/tstic
                   elsewhere
  -                    massbal = massbal_snow + dmax1(0.0_dp,massbal_ice)
  +                    massbal = massbal_snow + dmax1(0.0_dp,massbal_ice - snow_to_ice/tstic)
                   end where
               end if
   
  
positive.code:                 ! potential refeezing snow
                  refrozen_snow = dmax1(refr-refrozen_rain,0.0_dp)
                  ! actual refreezing snow
                  refrozen_snow = f_rz*dmin1(refrozen_snow,melt_snow)
                  ! actual refreezing
                  refr = refrozen_rain + refrozen_snow
              end if 
  
              ! 3) runoff
              runoff = melt + rf - refrozen_rain
  
              ! 4) accumulation: sum of all incoming solid water
              if (.not. bnd%acc) then
                  acc = sf - subl/rhow + refrozen_rain
                  !acc = sf - subl/rhow + refr
              end if
  
              
              ! 5) surface mass balance
              massbal_snow = sf - subl/rhow - melt_snow
  
negative.code:                 melt_ice  = melt-melt_snow
  
                  ! actual melt is sum of melted snow and ice (melted snow over land)
                  where (land_ice_ocean == 2)
                      melt = melt_snow + melt_ice
                  elsewhere
                      melt = melt_snow
  