Implementation of CCPP timestep_init and timestep_final phases (ufs-c…

…ommunity#217)

- replace calls to CCPP step `time_vary` with `timestep_init` in `atmos_model.F90`, add call to CCPP step `timestep_final` 
- update `ccpp/driver/CCPP_driver.F90` with calls to CCPP `timestep_init` and `timestep_final`
- add `h2o_def.f` and `ozne_def.` to `ccpp_prebuild_config.py`
- update of `gfsphysics/GFS_layer/GFS_typedefs.F90`: cleanup work for o3 and h2o physics (required by the updates to the CCPP time vary physics)
- update metadata in `gfsphysics/GFS_layer/GFS_typedefs.meta` for the above changes, and clean up the index used for the surface wind enhancement due to convection in the `phy_f2d` array (use proper index, not just the last entry in the array)

atmos_model.F90

@@ -254,6 +254,7 @@ subroutine update_atmos_radiation_physics (Atmos)
   type (atmos_data_type), intent(in) :: Atmos
 !--- local variables---
     integer :: nb, jdat(8), rc
+
     procedure(IPD_func0d_proc), pointer :: Func0d => NULL()
     procedure(IPD_func1d_proc), pointer :: Func1d => NULL()
     !
@@ -289,13 +290,12 @@ subroutine update_atmos_radiation_physics (Atmos)
 !--- execute the IPD atmospheric setup step
       call mpp_clock_begin(setupClock)
 #ifdef CCPP
-      call CCPP_step (step="time_vary", nblks=Atm_block%nblks, ierr=ierr)
-      if (ierr/=0)  call mpp_error(FATAL, 'Call to CCPP time_vary step failed')
+      call CCPP_step (step="timestep_init", nblks=Atm_block%nblks, ierr=ierr)
+      if (ierr/=0)  call mpp_error(FATAL, 'Call to CCPP timestep_init step failed')
 
 !--- call stochastic physics pattern generation / cellular automata
       call stochastic_physics_wrapper(IPD_Control, IPD_Data, Atm_block, ierr)
       if (ierr/=0)  call mpp_error(FATAL, 'Call to stochastic_physics_wrapper failed')
-
 #else
       Func1d => time_vary_step
       call IPD_step (IPD_Control, IPD_Data(:), IPD_Diag, IPD_Restart, IPD_func1d=Func1d)
@@ -419,6 +419,15 @@ subroutine update_atmos_radiation_physics (Atmos)
       endif
       call getiauforcing(IPD_Control,IAU_data)
       if (mpp_pe() == mpp_root_pe() .and. debug) write(6,*) "end of radiation and physics step"
+
+#ifdef CCPP
+!--- execute the IPD atmospheric timestep finalize step
+      call mpp_clock_begin(setupClock)
+      call CCPP_step (step="timestep_finalize", nblks=Atm_block%nblks, ierr=ierr)
+      if (ierr/=0)  call mpp_error(FATAL, 'Call to CCPP timestep_finalize step failed')
+      call mpp_clock_end(setupClock)
+#endif
+
     endif
 
 #ifdef CCPP

ccpp/config/ccpp_prebuild_config.py

@@ -17,6 +17,8 @@
     # actual variable definition files
     'ccpp/physics/physics/machine.F',
     'ccpp/physics/physics/radsw_param.f',
+    'ccpp/physics/physics/h2o_def.f',
+    'ccpp/physics/physics/ozne_def.f',
     'ccpp/physics/physics/radlw_param.f',
     'gfsphysics/CCPP_layer/CCPP_typedefs.F90',
     'gfsphysics/GFS_layer/GFS_typedefs.F90',

ccpp/driver/CCPP_driver.F90

@@ -3,14 +3,17 @@ module CCPP_driver
   use ccpp_api,           only: ccpp_t
 
   use ccpp_static_api,    only: ccpp_physics_init,                   &
+                                ccpp_physics_timestep_init,          &
                                 ccpp_physics_run,                    &
+                                ccpp_physics_timestep_finalize,      &
                                 ccpp_physics_finalize
 
   use CCPP_data,          only: cdata_tile,                          &
                                 cdata_domain,                        &
                                 cdata_block,                         &
                                 ccpp_suite,                          &
-                                GFS_control
+                                GFS_control,                         &
+                                GFS_data
 
   implicit none
 
@@ -55,6 +58,8 @@ subroutine CCPP_step (step, nblks, ierr)
     ! Local variables
     integer :: nb, nt, ntX
     integer :: ierr2
+    ! DH* 20210104 - remove kdt_rad when code to clear diagnostic buckets is removed
+    integer :: kdt_rad
 
     ierr = 0
 
@@ -95,9 +100,9 @@ subroutine CCPP_step (step, nblks, ierr)
 
     else if (trim(step)=="physics_init") then
 
-      ! Since the physics init steps are independent of the blocking structure,
+      ! Since the physics init step is independent of the blocking structure,
       ! we can use cdata_domain here. Since we don't use threading on the outside,
-      ! we can allow threading inside the time_vary routines.
+      ! we can allow threading inside the physics init routines.
       GFS_control%nthreads = nthrds
 
       call ccpp_physics_init(cdata_domain, suite_name=trim(ccpp_suite), ierr=ierr)
@@ -107,22 +112,53 @@ subroutine CCPP_step (step, nblks, ierr)
         return
       end if
 
-    else if (trim(step)=="time_vary") then
+    ! Timestep init = time_vary
+    else if (trim(step)=="timestep_init") then
 
-      ! Since the time_vary steps only use data structures for all blocks (except the
-      ! CCPP-internal variables ccpp_error_flag and ccpp_error_message, which are defined
-      ! for all cdata structures independently), we can use cdata_domain here.
-      ! Since we don't use threading on the outside, we can allow threading
-      ! inside the time_vary routines.
+      ! Since the physics timestep init step is independent of the blocking structure,
+      ! we can use cdata_domain here. Since we don't use threading on the outside,
+      ! we can allow threading inside the timestep init (time_vary) routines.
       GFS_control%nthreads = nthrds
 
-      call ccpp_physics_run(cdata_domain, suite_name=trim(ccpp_suite), group_name="time_vary", ierr=ierr)
+      call ccpp_physics_timestep_init(cdata_domain, suite_name=trim(ccpp_suite), group_name="time_vary", ierr=ierr)
       if (ierr/=0) then
-        write(0,'(a)') "An error occurred in ccpp_physics_run for group time_vary"
+        write(0,'(a)') "An error occurred in ccpp_physics_timestep_init for group time_vary"
         write(0,'(a)') trim(cdata_domain%errmsg)
         return
       end if
 
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      ! DH* 20210104 - this block of code will be removed once the CCPP framework    !
+      ! fully supports handling diagnostics through its metadata, work in progress   !
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+      !--- determine if radiation diagnostics buckets need to be cleared
+      if (nint(GFS_control%fhzero*3600) >= nint(max(GFS_control%fhswr,GFS_control%fhlwr))) then
+        if (mod(GFS_control%kdt,GFS_control%nszero) == 1) then
+          do nb = 1,nblks
+            call GFS_data(nb)%Intdiag%rad_zero(GFS_control)
+          end do
+        endif
+      else
+        kdt_rad = nint(min(GFS_control%fhswr,GFS_control%fhlwr)/GFS_control%dtp)
+        if (mod(GFS_control%kdt,kdt_rad) == 1) then
+          do nb = 1,nblks
+            call GFS_data(nb)%Intdiag%rad_zero(GFS_control)
+          enddo
+        endif
+      endif
+
+      !--- determine if physics diagnostics buckets need to be cleared
+      if (mod(GFS_control%kdt,GFS_control%nszero) == 1) then
+        do nb = 1,nblks
+          call GFS_data(nb)%Intdiag%phys_zero(GFS_control)
+        end do
+      endif
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      ! *DH 20210104                                                                 !
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
     ! Radiation and stochastic physics
     else if (trim(step)=="radiation" .or. trim(step)=="physics" .or. trim(step)=="stochastics") then
 
@@ -162,6 +198,21 @@ subroutine CCPP_step (step, nblks, ierr)
 !$OMP end parallel
       if (ierr/=0) return
 
+    ! Timestep finalize = time_vary
+    else if (trim(step)=="timestep_finalize") then
+
+      ! Since the physics timestep finalize step is independent of the blocking structure,
+      ! we can use cdata_domain here. Since we don't use threading on the outside,
+      ! we can allow threading inside the timestep finalize (time_vary) routines.
+      GFS_control%nthreads = nthrds
+
+      call ccpp_physics_timestep_finalize(cdata_domain, suite_name=trim(ccpp_suite), group_name="time_vary", ierr=ierr)
+      if (ierr/=0) then
+        write(0,'(a)') "An error occurred in ccpp_physics_timestep_finalize for group time_vary"
+        write(0,'(a)') trim(cdata_domain%errmsg)
+        return
+      end if
+
     ! Finalize
     else if (trim(step)=="finalize") then
 

gfsphysics/GFS_layer/GFS_typedefs.F90

@@ -14,6 +14,8 @@ module GFS_typedefs
 
        use module_radsw_parameters,  only: topfsw_type, sfcfsw_type, profsw_type, cmpfsw_type, NBDSW
        use module_radlw_parameters,  only: topflw_type, sfcflw_type, proflw_type, NBDLW
+       use ozne_def,                 only: levozp, oz_coeff
+       use h2o_def,                  only: levh2o, h2o_coeff
        use mo_gas_optics_rrtmgp,     only: ty_gas_optics_rrtmgp
        use mo_optical_props,         only: ty_optical_props_1scl,ty_optical_props_2str
        use mo_cloud_optics,          only: ty_cloud_optics
@@ -48,10 +50,10 @@ module GFS_typedefs
       ! from aerclm_def
       integer, parameter :: ntrcaerm = 15
 
-      ! These will be set later in GFS_Control%initialize,
-      ! since they depend on the runtime config (e.g. Model%ntoz, Model%h2o_phys, Model%aero_in)
-      private :: levozp, oz_coeff, levh2o, h2o_coeff, ntrcaer
-      integer :: levozp, oz_coeff, levh2o, h2o_coeff, ntrcaer
+      ! This will be set later in GFS_Control%initialize, since
+      ! it depends on the runtime config (Model%aero_in)
+      private :: ntrcaer
+      integer :: ntrcaer
 #endif
       ! If these are changed to >99, need to adjust formatting string in GFS_diagnostics.F90 (and names in diag_tables)
       integer, parameter :: naux2dmax = 20 !< maximum number of auxiliary 2d arrays in output (for debugging)
@@ -1171,10 +1173,11 @@ module GFS_typedefs
     integer              :: nqvdelt         !< the index of specific humidity at the previous timestep for Z-C MP in phy_f3d
     integer              :: nps2delt        !< the index of surface air pressure 2 timesteps back for Z-C MP in phy_f2d
     integer              :: npsdelt         !< the index of surface air pressure at the previous timestep for Z-C MP in phy_f2d
+    integer              :: ncnvwind        !< the index of surface wind enhancement due to convection for MYNN SFC and RAS CNV in phy f2d
 #endif
 
 !--- debug flag
-    logical              :: debug         
+    logical              :: debug
     logical              :: pre_rad         !< flag for testing purpose
 
 !--- variables modified at each time step
@@ -1863,8 +1866,7 @@ module GFS_typedefs
     real (kind=kind_phys), pointer      :: graupelmp(:)       => null()  !<
     real (kind=kind_phys), pointer      :: gwdcu(:,:)         => null()  !<
     real (kind=kind_phys), pointer      :: gwdcv(:,:)         => null()  !<
-    integer                             :: h2o_coeff                     !<
-    real (kind=kind_phys), pointer      :: h2o_pres(:)        => null()  !<
+
     real (kind=kind_phys), pointer      :: hefac(:)           => null()  !<
     real (kind=kind_phys), pointer      :: hffac(:)           => null()  !<
     real (kind=kind_phys), pointer      :: hflxq(:)           => null()  !<
@@ -1898,8 +1900,6 @@ module GFS_typedefs
     integer,               pointer      :: ktop(:)            => null()  !<
     integer                             :: latidxprnt                    !<
     integer                             :: levi                          !<
-    integer                             :: levh2o                        !<
-    integer                             :: levozp                        !<
     integer                             :: lmk                           !<
     integer                             :: lmp                           !<
     integer,               pointer      :: mbota(:,:)         => null()  !<
@@ -1929,9 +1929,7 @@ module GFS_typedefs
     real (kind=kind_phys), pointer      :: oc(:)              => null()  !<
     real (kind=kind_phys), pointer      :: olyr(:,:)          => null()  !<
     logical              , pointer      :: otspt(:,:)         => null()  !<
-    integer                             :: oz_coeff                      !<
     integer                             :: oz_coeffp5                    !<
-    real (kind=kind_phys), pointer      :: oz_pres(:)         => null()  !<
     logical                             :: phys_hydrostatic              !<
     real (kind=kind_phys), pointer      :: plvl(:,:)          => null()  !<
     real (kind=kind_phys), pointer      :: plyr(:,:)          => null()  !<
@@ -4358,7 +4356,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
           stop
        else
           levozp   = 1
-          oz_coeff = 0
+          oz_coeff = 1
        end if
     end if
 #endif
@@ -4876,6 +4874,9 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
 #endif
     if(Model%cnvcld) Model%ncnvcld3d = 1
 
+!--- get cnvwind index in phy_f2d; last entry in phy_f2d array
+    Model%ncnvwind = Model%num_p2d
+
 !--- get cnvw and cnvc indices in phy_f3d
     Model%ncnvw = -999
     Model%ncnvc = -999
@@ -5212,6 +5213,8 @@ subroutine control_print(Model)
       print *, ' shocaftcnv        : ', Model%shocaftcnv
       print *, ' shoc_cld          : ', Model%shoc_cld
       print *, ' uni_cld           : ', Model%uni_cld
+      print *, ' oz_phys           : ', Model%oz_phys
+      print *, ' oz_phys_2015      : ', Model%oz_phys_2015
       print *, ' h2o_phys          : ', Model%h2o_phys
       print *, ' pdfcld            : ', Model%pdfcld
       print *, ' shcnvcw           : ', Model%shcnvcw
@@ -5509,8 +5512,6 @@ subroutine tbd_create (Tbd, IM, Model)
     endif
 
 !--- ozone and stratosphere h2o needs
-    ! DH* oz_coeff is set to zero if both ozphys options are false,
-    ! better to use conditional allocations here for ozpl (and h2opl)? *DH
     allocate (Tbd%ozpl  (IM,levozp,oz_coeff))
     allocate (Tbd%h2opl (IM,levh2o,h2o_coeff))
     Tbd%ozpl  = clear_val
@@ -5529,8 +5530,6 @@ subroutine tbd_create (Tbd, IM, Model)
     Tbd%aer_nm = clear_val
 
 #ifdef CCPP
-! DH* TODO - MOVE THIS TO a block-vector dependent structure in GFS_control?
-! e.g. GFS_Control%imap(blk), GFS_Control%jmap(blk), or ii instead if imap etc? *DH
 !--- maps of local index ix to global indices i and j for this block
     allocate (Tbd%imap (IM))
     allocate (Tbd%jmap (IM))
@@ -6552,7 +6551,6 @@ subroutine interstitial_create (Interstitial, IM, Model)
     allocate (Interstitial%gflx_ocean      (IM))
     allocate (Interstitial%gwdcu           (IM,Model%levs))
     allocate (Interstitial%gwdcv           (IM,Model%levs))
-    allocate (Interstitial%h2o_pres        (levh2o))
     allocate (Interstitial%hefac           (IM))
     allocate (Interstitial%hffac           (IM))
     allocate (Interstitial%hflxq           (IM))
@@ -6581,7 +6579,6 @@ subroutine interstitial_create (Interstitial, IM, Model)
     allocate (Interstitial%oa4             (IM,4))
     allocate (Interstitial%oc              (IM))
     allocate (Interstitial%olyr            (IM,Model%levr+LTP))
-    allocate (Interstitial%oz_pres         (levozp))
     allocate (Interstitial%plvl            (IM,Model%levr+1+LTP))
     allocate (Interstitial%plyr            (IM,Model%levr+LTP))
     allocate (Interstitial%prnum           (IM,Model%levs))
@@ -6827,23 +6824,18 @@ subroutine interstitial_create (Interstitial, IM, Model)
     Interstitial%ipr              = min(IM,10)
     Interstitial%latidxprnt       = 1
     Interstitial%levi             = Model%levs+1
-    Interstitial%levh2o           = levh2o
-    Interstitial%levozp           = levozp
     Interstitial%lmk              = Model%levr+LTP
     Interstitial%lmp              = Model%levr+1+LTP
-    Interstitial%h2o_coeff        = h2o_coeff
     Interstitial%nbdlw            = NBDLW
     Interstitial%nbdsw            = NBDSW
     Interstitial%nf_aelw          = NF_AELW
     Interstitial%nf_aesw          = NF_AESW
     Interstitial%nspc1            = NSPC1
-    Interstitial%oz_coeff         = oz_coeff
-    Interstitial%oz_coeffp5       = oz_coeff+5
-    ! h2o_pres and oz_pres do not change during the run, but
-    ! need to be set later in GFS_phys_time_vary_init (after
-    ! h2o_pres/oz_pres are read in read_h2odata/read_o3data)
-    Interstitial%h2o_pres         = clear_val
-    Interstitial%oz_pres          = clear_val
+    if (Model%oz_phys .or. Model%oz_phys_2015) then
+      Interstitial%oz_coeffp5     = oz_coeff+5
+    else
+      Interstitial%oz_coeffp5     = 5
+    endif
     !
     Interstitial%skip_macro       = .false.
     ! The value phys_hydrostatic from dynamics does not match the
@@ -7449,14 +7441,10 @@ subroutine interstitial_print(Interstitial, Model, mpirank, omprank, blkno)
     ! Print static variables
     write (0,'(a,3i6)') 'Interstitial_print for mpirank, omprank, blkno: ', mpirank, omprank, blkno
     write (0,*) 'Interstitial_print: values that do not change'
-    write (0,*) 'Interstitial%h2o_coeff         = ', Interstitial%h2o_coeff
-    write (0,*) 'sum(Interstitial%h2o_pres)     = ', sum(Interstitial%h2o_pres)
     write (0,*) 'Interstitial%ipr               = ', Interstitial%ipr
     write (0,*) 'Interstitial%itc               = ', Interstitial%itc
     write (0,*) 'Interstitial%latidxprnt        = ', Interstitial%latidxprnt
     write (0,*) 'Interstitial%levi              = ', Interstitial%levi
-    write (0,*) 'Interstitial%levh2o            = ', Interstitial%levh2o
-    write (0,*) 'Interstitial%levozp            = ', Interstitial%levozp
     write (0,*) 'Interstitial%lmk               = ', Interstitial%lmk
     write (0,*) 'Interstitial%lmp               = ', Interstitial%lmp
     write (0,*) 'Interstitial%nbdlw             = ', Interstitial%nbdlw
@@ -7468,8 +7456,6 @@ subroutine interstitial_print(Interstitial, Model, mpirank, omprank, blkno)
     write (0,*) 'Interstitial%nspc1             = ', Interstitial%nspc1
     write (0,*) 'Interstitial%ntiwx             = ', Interstitial%ntiwx
     write (0,*) 'Interstitial%nvdiff            = ', Interstitial%nvdiff
-    write (0,*) 'Interstitial%oz_coeff          = ', Interstitial%oz_coeff
-    write (0,*) 'sum(Interstitial%oz_pres)      = ', sum(Interstitial%oz_pres)
     write (0,*) 'Interstitial%phys_hydrostatic  = ', Interstitial%phys_hydrostatic
     write (0,*) 'Interstitial%skip_macro        = ', Interstitial%skip_macro
     write (0,*) 'Interstitial%trans_aero        = ', Interstitial%trans_aero