Add ensemble-consistent snow masking for SWC-DA

docs/users_guide/running_tsmp_pdaf/input_enkfpf.md

@@ -590,17 +590,33 @@ are allowed.
   CLM5.0's `watmin` from `clm_varcon.F90` (current value `0.01`). If
   yes, set SM to `watmin`.
 
+(enkfpf:clm:swc_mask_snow)=
 ### CLM:swc_mask_snow ###
 
 `CLM:swc_mask_snow`: (integer) Switch for masking columns with snow
 cover from SWC updates.
 
-Snow covers larger than 1mm are switched off for the update.
+Columns with snow depth ≥ 1 mm are excluded from the update.
 
-Only takes effect if `CLM:update_swc``is switched on.
+Only takes effect if `CLM:update_swc` is switched on.
 
 Default setting is `0`: No masking of columns with snow cover.
 
+(enkfpf:clm:swc_mask_snow_ens)=
+### CLM:swc_mask_snow_ens ###
+
+`CLM:swc_mask_snow_ens`: (integer) Switch for ensemble-consistent snow
+masking. When set to `1`, a column is excluded from the SWC update if
+**any** ensemble member has snow depth ≥ 1 mm there, rather than each
+member masking independently based on its own snow state. Requires an
+MPI reduction across ensemble members via `COMM_couple`.
+
+Only takes effect if both `CLM:update_swc` and `CLM:swc_mask_snow` are
+switched on.
+
+Default setting is `0`: each member applies the snow mask based on its
+own snow depth (original behaviour).
+
 (enkfpf:cosmo)=
 ## [COSMO] ##
 
@@ -900,6 +916,8 @@ Default: 0, output turned off.
  |           | `statevec_max_layer`    | 25            |
  |           | `t_printensemble`       | -2            |
  |           | `watmin_switch`         | 0             |
+ |           | `swc_mask_snow`         | 0             |
+ |           | `swc_mask_snow_ens`     | 0             |
  | `[COSMO]` |                         |               |
  |           | `nprocs`                | 0             |
  |           | `dtmult`                | 0             |

interface/model/common/enkf.h

@@ -97,6 +97,7 @@ GLOBAL int clmstatevec_max_layer;
 GLOBAL int clmt_printensemble;
 GLOBAL int clmwatmin_switch;
 GLOBAL int clmswc_mask_snow;
+GLOBAL int clmswc_mask_snow_ens;
 GLOBAL int dtmult_cosmo;
 GLOBAL int pf_olfmasking;
 GLOBAL int pf_olfmasking_param;

interface/model/common/read_enkfpar.c

@@ -88,6 +88,7 @@ void read_enkfpar(char *parname)
   clmt_printensemble    = iniparser_getint(pardict,"CLM:t_printensemble",-2);
   clmwatmin_switch      = iniparser_getint(pardict,"CLM:watmin_switch",0);
   clmswc_mask_snow      = iniparser_getint(pardict,"CLM:swc_mask_snow",0);
+  clmswc_mask_snow_ens  = iniparser_getint(pardict,"CLM:swc_mask_snow_ens",0);
 
   /* get settings for COSMO */
   nproccosmo      = iniparser_getint(pardict,"COSMO:nprocs",0);

interface/model/eclm/enkf_clm_mod_5.F90

@@ -65,6 +65,7 @@ module enkf_clm_mod
   integer(c_int),bind(C,name="clmt_printensemble")       :: clmt_printensemble
   integer(c_int),bind(C,name="clmwatmin_switch")         :: clmwatmin_switch
   integer(c_int),bind(C,name="clmswc_mask_snow")            :: clmswc_mask_snow
+  integer(c_int),bind(C,name="clmswc_mask_snow_ens")        :: clmswc_mask_snow_ens
   real(c_double),bind(C,name="clmcrns_bd")      :: clmcrns_bd
 
   integer  :: nstep     ! time step index
@@ -631,6 +632,10 @@ subroutine update_clm_swc(tstartcycle, mype)
     use clm_varcon      , only : denh2o, denice, watmin
     use clm_varcon      , only : ispval
     use clm_varcon      , only : spval
+    use mpi, only: MPI_Allreduce
+    use mpi, only: MPI_IN_PLACE
+    use mpi, only: MPI_INTEGER
+    use mpi, only: MPI_MAX
 
     implicit none
 
@@ -649,6 +654,9 @@ subroutine update_clm_swc(tstartcycle, mype)
     real(r8)  :: watmin_set        ! minimum soil moisture for setting swc (mm)
     real(r8)  :: swc_update        ! updated SWC in loop
 
+    integer, allocatable :: snow_mask(:)  ! 1 = snow in any ensemble member, 0 = clear
+    integer :: MPIerr_swc
+
     integer :: i,j,cc
     character (len = 31) :: fn2    !TSMP-PDAF: function name for state vector outpu
     character (len = 32) :: fn3    !TSMP-PDAF: function name for state vector outpu
@@ -667,6 +675,24 @@ subroutine update_clm_swc(tstartcycle, mype)
 
     snow_depth => waterstate_inst%snow_depth_col ! snow height of snow covered area (m)
 
+    ! Build ensemble-consistent snow mask: if clmswc_mask_snow_ens==1,
+    ! reduce the per-column snow flag across all ensemble members via
+    ! COMM_couple_clm so that a column is masked whenever ANY member
+    ! has snow depth >= 1 mm.
+    allocate(snow_mask(clm_begc:clm_endc))
+    do j = clm_begc, clm_endc
+      if (snow_depth(j) >= 0.001_r8) then
+        snow_mask(j) = 1
+      else
+        snow_mask(j) = 0
+      end if
+    end do
+    if (clmswc_mask_snow == 1 .and. clmswc_mask_snow_ens == 1) then
+      call MPI_Allreduce(MPI_IN_PLACE, snow_mask(clm_begc), &
+                         clm_endc - clm_begc + 1, MPI_INTEGER, MPI_MAX, &
+                         COMM_couple_clm, MPIerr_swc)
+    end if
+
         ! Set minimum soil moisture for checking the state vector and
         ! for setting minimum swc for CLM
         if(clmwatmin_switch==3) then
@@ -690,8 +716,11 @@ subroutine update_clm_swc(tstartcycle, mype)
           ! do j=clm_begg,clm_endg
             do j=clm_begc,clm_endc
 
-              ! If snow is masked, update only, when snow depth is less than 1mm
-              if( (clmswc_mask_snow == 0) .or. snow_depth(j) < 0.001 ) then
+              ! If snow is masked, update only when no snow is present.
+              ! snow_mask reflects either the local member's snow state
+              ! (clmswc_mask_snow_ens=0, default) or the ensemble-wide
+              ! maximum (clmswc_mask_snow_ens=1).
+              if( (clmswc_mask_snow == 0) .or. snow_mask(j) == 0 ) then
               ! Update only those SWCs that are not excluded by ispval
               if(state_clm2pdaf_p(j,i) /= ispval) then
 
@@ -797,6 +826,8 @@ subroutine update_clm_swc(tstartcycle, mype)
         END IF
 #endif
 
+    deallocate(snow_mask)
+
   end subroutine update_clm_swc