PDAF-related: GRACE-DA changes

src/clm5/biogeophys/SoilHydrologyType.F90

@@ -53,6 +53,11 @@ Module SoilHydrologyType
      real(r8), pointer :: i_0_col           (:)     ! col VIC average saturation in top soil layers 
      real(r8), pointer :: ice_col           (:,:)   ! col VIC soil ice (kg/m2) for VIC soil layers
 
+#ifdef USE_PDAF
+     ! Yorck
+     real(r8), pointer :: wa_col_inc        (:)     ! increment col water in the unconfined aquifer (mm)
+     real(r8), pointer :: wa_col_mean       (:)     ! mean col water in the unconfined aquifer (mm)
+#endif
    contains
 
      ! Public routines
@@ -117,6 +122,11 @@ subroutine InitAllocate(this, bounds)
     allocate(this%zwts_col          (begc:endc))                 ; this%zwts_col          (:)     = nan
 
     allocate(this%wa_col            (begc:endc))                 ; this%wa_col            (:)     = nan
+#ifdef USE_PDAF
+    ! Yorck
+    allocate(this%wa_col_inc        (begc:endc))                 ; this%wa_col_inc        (:)     = nan
+    allocate(this%wa_col_mean       (begc:endc))                 ; this%wa_col_mean       (:)     = nan
+#endif
     allocate(this%qcharge_col       (begc:endc))                 ; this%qcharge_col       (:)     = nan
     allocate(this%fracice_col       (begc:endc,nlevgrnd))        ; this%fracice_col       (:,:)   = nan
     allocate(this%icefrac_col       (begc:endc,nlevgrnd))        ; this%icefrac_col       (:,:)   = nan

src/clm5/biogeophys/WaterStateType.F90

@@ -54,6 +54,57 @@ module WaterstateType
      real(r8), pointer :: ice1_grc               (:)   ! grc initial gridcell total h2o ice content
      real(r8), pointer :: ice2_grc               (:)   ! grc post land cover change total ice content
      real(r8), pointer :: tws_grc                (:)   ! grc total water storage (mm H2O)
+#ifdef USE_PDAF
+     ! Yorck additions, variables for getting monthly means which can be compared to a GRACE measurements, other variables will just
+     ! provide instantaneous values
+
+     real(r8), pointer :: h2osno_col_mean             (:)   ! col snow water (mm H2O)
+     real(r8), pointer :: h2osoi_liq_col_mean         (:,:) ! col liquid water (kg/m2) (new) (-nlevsno+1:nlevgrnd)
+     real(r8), pointer :: h2osoi_ice_col_mean         (:,:) ! col ice lens (kg/m2) (new) (-nlevsno+1:nlevgrnd)
+     real(r8), pointer :: h2ocan_patch_mean           (:)   ! patch canopy water (mm H2O)
+     real(r8), pointer :: snocan_patch_mean           (:)   ! patch canopy water (mm H2O)
+     real(r8), pointer :: h2osfc_col_mean             (:)   ! col surface water (mm H2O)
+     real(r8), pointer :: total_plant_stored_h2o_col_mean(:) ! col water that is bound in plants, including roots, sapwood, leaves, etc
+                                                             ! in most cases, the vegetation scheme does not have a dynamic
+                                                             ! water storage in plants, and thus 0.0 is a suitable for the trivial case.
+                                                             ! When FATES is coupled in with plant hydraulics turned on, this storage
+                                                             ! term is set to non-zero. (kg/m2 H2O)
+
+
+     real(r8), pointer :: tws_hactive                 (:)    ! TWS for hydrological active columns
+     real(r8), pointer :: tws_hactive_mean            (:)    ! TWS for hydrological active columns
+
+
+     ! also increments for all variables are added
+     real(r8), pointer :: h2osno_col_inc              (:)   ! col snow water (mm H2O)
+     real(r8), pointer :: h2osoi_liq_col_inc          (:,:) ! col liquid water (kg/m2) (new) (-nlevsno+1:nlevgrnd)
+     real(r8), pointer :: h2osoi_ice_col_inc          (:,:) ! col ice lens (kg/m2) (new) (-nlevsno+1:nlevgrnd)
+     real(r8), pointer :: h2ocan_patch_inc            (:)   ! patch canopy water (mm H2O)
+     real(r8), pointer :: h2osfc_col_inc              (:)   ! col surface water (mm H2O)
+     real(r8), pointer :: total_plant_stored_h2o_col_inc (:) ! col water that is bound in plants, including roots, sapwood, leaves, etc
+                                                             ! in most cases, the vegetation scheme does not have a dynamic
+                                                             ! water storage in plants, and thus 0.0 is a suitable for the trivial case.
+                                                             ! When FATES is coupled in with plant hydraulics turned on, this storage
+                                                             ! term is set to non-zero. (kg/m2 H2O)
+     real(r8), pointer :: tws_grc_inc                  (:)   ! grc total water storage (mm H2O)
+
+     ! for analysis of the statevector, also variables for putting the state inside are initialized
+     real(r8), pointer :: tws_state_before             (:)    ! TWS state
+     real(r8), pointer :: h2osoi_liq_state_before      (:,:)  ! soil liq state
+     real(r8), pointer :: h2osoi_ice_state_before      (:,:)  ! soil ice state
+     real(r8), pointer :: h2osfc_state_before          (:)    ! surface water state
+     real(r8), pointer :: h2osno_state_before          (:)    ! snow state
+     real(r8), pointer :: h2ocan_state_before          (:)    ! canopy state
+
+     real(r8), pointer :: tws_state_after              (:)    ! TWS state
+     real(r8), pointer :: h2osoi_liq_state_after       (:,:)  ! soil liq state
+     real(r8), pointer :: h2osoi_ice_state_after       (:,:)  ! soil ice state
+     real(r8), pointer :: h2osfc_state_after           (:)    ! surface water state
+     real(r8), pointer :: h2osno_state_after           (:)    ! snow state
+     real(r8), pointer :: h2ocan_state_after           (:)    ! canopy state
+
+     ! END Yorck
+#endif
 #ifdef COUP_OAS_PFL
      real(r8), pointer :: pfl_psi_col            (:,:) ! ParFlow pressure head   COUP_OAS_PFL
      real(r8), pointer :: pfl_h2osoi_liq_col     (:,:) ! ParFlow soil liquid     COUP_OAS_PFL
@@ -211,6 +262,42 @@ subroutine InitAllocate(this, bounds)
     allocate(this%ice2_grc               (begg:endg))                     ; this%ice2_grc               (:)   = nan
     allocate(this%tws_grc                (begg:endg))                     ; this%tws_grc                (:)   = nan
 
+#ifdef USE_PDAF
+    ! Yorck additions (see above)
+    allocate(this%h2osoi_ice_col_mean    (begc:endc,-nlevsno+1:nlevgrnd)) ; this%h2osoi_ice_col_mean    (:,:) = nan
+    allocate(this%h2osoi_liq_col_mean    (begc:endc,-nlevsno+1:nlevgrnd)) ; this%h2osoi_liq_col_mean    (:,:) = nan
+    allocate(this%h2osno_col_mean        (begc:endc))                     ; this%h2osno_col_mean        (:)   = nan
+    allocate(this%h2ocan_patch_mean      (begp:endp))                     ; this%h2ocan_patch_mean      (:)   = nan
+    allocate(this%snocan_patch_mean      (begp:endp))                     ; this%snocan_patch_mean      (:)   = nan
+    allocate(this%h2osfc_col_mean        (begc:endc))                     ; this%h2osfc_col_mean        (:)   = nan
+    allocate(this%total_plant_stored_h2o_col_mean(begc:endc))             ; this%total_plant_stored_h2o_col_mean(:) = nan
+
+    allocate(this%tws_hactive            (begg:endg))                     ; this%tws_hactive            (:)   = nan
+    allocate(this%tws_hactive_mean       (begg:endg))                     ; this%tws_hactive_mean       (:)   = nan
+
+    allocate(this%h2osoi_ice_col_inc     (begc:endc,1:nlevsoi))           ; this%h2osoi_ice_col_inc     (:,:) = nan
+    allocate(this%h2osoi_liq_col_inc     (begc:endc,1:nlevsoi))           ; this%h2osoi_liq_col_inc     (:,:) = nan
+    allocate(this%h2osno_col_inc         (begc:endc))                     ; this%h2osno_col_inc         (:)   = nan
+    allocate(this%h2ocan_patch_inc       (begp:endp))                     ; this%h2ocan_patch_inc       (:)   = nan
+    allocate(this%h2osfc_col_inc         (begc:endc))                     ; this%h2osfc_col_inc         (:)   = nan
+    allocate(this%total_plant_stored_h2o_col_inc (begc:endc))             ; this%total_plant_stored_h2o_col_inc (:) = nan
+    allocate(this%tws_grc_inc            (begg:endg))                     ; this%tws_grc_inc            (:)   = nan
+
+    allocate(this%tws_state_before       (begg:endg))                     ; this%tws_state_before       (:)   = nan
+    allocate(this%h2osoi_liq_state_before(begg:endg,1:nlevsoi))           ; this%h2osoi_liq_state_before(:,:) = nan
+    allocate(this%h2osoi_ice_state_before(begg:endg,1:nlevsoi))           ; this%h2osoi_ice_state_before(:,:) = nan
+    allocate(this%h2osno_state_before    (begg:endg))                     ; this%h2osno_state_before    (:)   = nan
+    allocate(this%h2osfc_state_before    (begg:endg))                     ; this%h2osfc_state_before    (:)   = nan
+    allocate(this%h2ocan_state_before    (begg:endg))                     ; this%h2ocan_state_before    (:)   = nan
+
+    allocate(this%tws_state_after        (begg:endg))                     ; this%tws_state_after        (:)   = nan
+    allocate(this%h2osoi_liq_state_after (begg:endg,1:nlevsoi))           ; this%h2osoi_liq_state_after (:,:) = nan
+    allocate(this%h2osoi_ice_state_after (begg:endg,1:nlevsoi))           ; this%h2osoi_ice_state_after (:,:) = nan
+    allocate(this%h2osno_state_after     (begg:endg))                     ; this%h2osno_state_after     (:)   = nan
+    allocate(this%h2osfc_state_after     (begg:endg))                     ; this%h2osfc_state_after     (:)   = nan
+    allocate(this%h2ocan_state_after     (begg:endg))                     ; this%h2ocan_state_after     (:)   = nan
+    ! END Yorck
+#endif
     allocate(this%total_plant_stored_h2o_col(begc:endc))                  ; this%total_plant_stored_h2o_col(:) = nan
 
     allocate(this%snw_rds_col            (begc:endc,-nlevsno+1:0))        ; this%snw_rds_col            (:,:) = nan
@@ -402,6 +489,13 @@ subroutine InitHistory(this, bounds)
          avgflag='A', long_name='total water storage', &
          ptr_lnd=this%tws_grc)
 
+#ifdef USE_PDAF
+    ! Yorck: add also TWS_hactive in history outputs
+    this%tws_hactive(begg:endg) = spval
+    call hist_addfld1d (fname='TWS_hactive',  units='mm',  &
+         avgflag='A', long_name='total water storage of hydrological active cells', &
+         ptr_lnd=this%tws_hactive)
+#endif
     ! (rgk 02-02-2017) There is intentionally no entry  here for stored plant water
     !                  I think that since the value is zero in all cases except
     !                  for FATES plant hydraulics, it will be confusing for users

src/clm5/main/atm2lndType.F90

@@ -147,6 +147,11 @@ module atm2lndType
      real(r8) , pointer :: t_mo_patch                   (:)   => null() ! patch 30-day average temperature (Kelvin)
      real(r8) , pointer :: t_mo_min_patch               (:)   => null() ! patch annual min of t_mo (Kelvin)
 
+#ifdef USE_PDAF
+     ! Yorck
+     real(r8), pointer :: volr_grc_inc                  (:)   => null() ! rof volr total volume increment (m3)
+     real(r8), pointer :: volr_grc_mean                 (:)   => null() ! rof volr total volume mean (m3)
+#endif
    contains
 
      procedure, public  :: Init
@@ -552,6 +557,13 @@ subroutine InitAllocate(this, bounds)
     allocate(this%pfl_psi_grc                   (begg:endg,1:nlevgrnd)); this%pfl_psi_grc            (:,:) = ival
     allocate(this%pfl_h2osoi_liq_grc            (begg:endg,1:nlevgrnd)); this%pfl_h2osoi_liq_grc     (:,:) = ival
 #endif
+
+#ifdef USE_PDAF
+     ! Yorck
+     allocate(this%volr_grc_inc                  (begg:endg))        ; this%volr_grc_inc                  (:)   = ival
+     allocate(this%volr_grc_mean                 (begg:endg))        ; this%volr_grc_mean                 (:)   = ival
+#endif
+
     ! anomaly forcing
     allocate(this%bc_precip_grc                 (begg:endg))        ; this%bc_precip_grc                 (:)   = ival
     allocate(this%af_precip_grc                 (begg:endg))        ; this%af_precip_grc                 (:)   = ival

src/clm5/main/clm_driver.F90

@@ -1068,6 +1068,9 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro
          solarabs_inst, drydepvel_inst,       &
          vocemis_inst, fireemis_inst, dust_inst, ch4_inst, glc_behavior, &
          lnd2atm_inst, &
+#ifdef USE_PDAF
+         soilhydrology_inst, soilstate_inst, &
+#endif
          net_carbon_exchange_grc = net_carbon_exchange_grc(bounds_proc%begg:bounds_proc%endg))
     deallocate(net_carbon_exchange_grc)
     call t_stopf('lnd2atm')

src/clm5/main/clm_varcon.F90

@@ -106,6 +106,14 @@ module clm_varcon
   ! Keep this negative to avoid conflicts with possible valid values
   integer , public, parameter :: ispval = -9999             ! special value for int data
 
+#ifdef USE_PDAF
+  ! ------------------------------------------------------------------------
+  ! Yorck variables
+  ! ------------------------------------------------------------------------
+
+  integer :: averaging_var = 0 ! averaging integer for averaging TWS variables according to GRACE interval
+  integer :: set_averaging_to_zero = ispval ! integer indicating at which point in time the averageTemp has to be set to zero due to missing observation values
+#endif
   ! ------------------------------------------------------------------------
   ! These are tunable constants from clm2_3
   ! ------------------------------------------------------------------------