Merge pull request #229 from olyson/nextsw_cday

Change nextsw_cday calculation for cam_dev (CAM7) in cplhist mode

### Description of changes
The radiation re-ordering in cam_dev (CAM7) requires a change to the calculation of nextsw_cday in cplhist mode. This PR introduces a new namelist variable, nextsw_cday_calc, which controls this behavior. The default value is appropriate for cplhist files generated with CAM7. Note that, by default, this namelist flag has no effect on non-cplhist cases (because of the default value of iradsw in non-cplhist cases).

### Specific notes

Contributors other than yourself, if any:  @billsacks @adamrher

CDEPS Issues Fixed (include github issue #):

Are there dependencies on other component PRs (if so list):

Are changes expected to change answers (bfb, different to roundoff, more substantial): Greater than roundoff-level changes in cases using datm in CPLHIST mode.

Any User Interface Changes (namelist or namelist defaults changes): New namelist variable, nextsw_cday_calc.

Testing performed (e.g. aux_cdeps, CESM prealpha, etc):
I have plugged this branch into the latest version of CTSM and ran the following tests:
- SMS_D_Ld3.f10_f10_mg37.I1850Clm50BgcCrop.derecho_intel.clm-default - passes and bit-for-bit
- SMS_D_Ld1.ne30pg3_t232.I1850Clm50BgcSpinup.derecho_intel.clm-cplhist - passes but changes answers as expected
- SMS_D_Ld1.ne30pg3_t232.I1850Clm50BgcSpinup.derecho_intel.clm-cplhist, with the new namelist flag changed to "cam6" - passes and bit-for-bit

Hashes used for testing:
ctsm5.3.009

datm/atm_comp_nuopc.F90

@@ -111,6 +111,7 @@ module cdeps_datm_comp
   character(CL)                :: model_meshfile = nullstr            ! full pathname to model meshfile
   character(CL)                :: model_maskfile = nullstr            ! full pathname to obtain mask from
   integer                      :: iradsw = 0                          ! radiation interval (input namelist)
+  logical                      :: nextsw_cday_calc_cam7               ! true => use logic appropriate to cam7 (and later) for calculating nextsw_cday
   character(CL)                :: factorFn_mesh = 'null'              ! file containing correction factors mesh
   character(CL)                :: factorFn_data = 'null'              ! file containing correction factors data
   logical                      :: flds_presaero = .false.             ! true => send valid prescribed aero fields to mediator
@@ -213,6 +214,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
     integer           :: nu         ! unit number
     integer           :: ierr       ! error code
     integer           :: bcasttmp(10)
+    character(CL)     :: nextsw_cday_calc
     type(ESMF_VM)     :: vm
     character(len=*),parameter :: subname=trim(modName) // ':(InitializeAdvertise) '
     character(*)    ,parameter :: F00 = "('(" // trim(modName) // ") ',8a)"
@@ -228,6 +230,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
          ny_global, &
          restfilm, &
          iradsw, &
+         nextsw_cday_calc, &
          factorFn_data, &
          factorFn_mesh, &
          flds_presaero, &
@@ -242,6 +245,9 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
 
     rc = ESMF_SUCCESS
 
+    ! Initialize locally-declared namelist items to default values
+    nextsw_cday_calc = nullstr
+
     call NUOPC_CompAttributeGet(gcomp, name='case_name', value=case_name, rc=rc)
     if (ChkErr(rc,__LINE__,u_FILE_u)) return
 
@@ -301,6 +307,8 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
     if (ChkErr(rc,__LINE__,u_FILE_u)) return
     call ESMF_VMBroadcast(vm, restfilm, CL, main_task, rc=rc)
     if (ChkErr(rc,__LINE__,u_FILE_u)) return
+    call ESMF_VMBroadcast(vm, nextsw_cday_calc, CL, main_task, rc=rc)
+    if (ChkErr(rc,__LINE__,u_FILE_u)) return
     call ESMF_VMBroadcast(vm, bcasttmp, 10, main_task, rc=rc)
     if (ChkErr(rc,__LINE__,u_FILE_u)) return
     nx_global         = bcasttmp(1)
@@ -314,6 +322,14 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
     skip_restart_read = (bcasttmp(9) == 1)
     export_all        = (bcasttmp(10) == 1)
 
+    if (nextsw_cday_calc == 'cam7') then
+       nextsw_cday_calc_cam7 = .true.
+    else if (nextsw_cday_calc == 'cam6') then
+       nextsw_cday_calc_cam7 = .false.
+    else
+       call shr_sys_abort(' ERROR illegal datm nextsw_cday_calc = '//trim(nextsw_cday_calc))
+    end if
+
     ! write namelist input to standard out
     if (my_task == main_task) then
        write(logunit,F00)' case_name      = ',trim(case_name)
@@ -324,6 +340,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
        write(logunit,F01)' ny_global      = ',ny_global
        write(logunit,F00)' restfilm       = ',trim(restfilm)
        write(logunit,F01)' iradsw         = ',iradsw
+       write(logunit,F00)' nextsw_cday_calc = ', trim(nextsw_cday_calc)
        write(logunit,F00)' factorFn_data  = ',trim(factorFn_data)
        write(logunit,F00)' factorFn_mesh  = ',trim(factorFn_mesh)
        write(logunit,F02)' flds_presaero  = ',flds_presaero
@@ -825,15 +842,14 @@ end subroutine datm_comp_run
   real(R8) function getNextRadCDay( julday, tod, stepno, dtime, iradsw )
 
     ! Return the calendar day of the next radiation time-step.
-    ! General Usage: nextswday = getNextRadCDay(curr_date) iradsw is
-    ! the frequency to update the next shortwave.  in number of steps
-    ! (or hours if negative) Julian date.
-    ! -- values greater than 1 set
-    !    the next radiation to the present time plus 2 timesteps every iradsw
-    ! -- values less than 0 turn set the next radiation to the  present time
-    !    plus two timesteps every -iradsw hours.
-    ! -- if iradsw is zero, the next radiation time is the
-    !    present time plus 1 timestep.
+    ! General Usage: nextswday = getNextRadCDay(curr_date). iradsw is
+    ! the frequency to update the next shortwave in number of steps
+    ! (or hours if negative).
+    ! -- values greater than 1 set the next radiation to the present time plus either 1 or
+    !    2 timesteps (depending on the value of nextsw_cday_calc_cam7) every iradsw timesteps.
+    ! -- values less than 0 set the next radiation to the present time plus either 1 or 2
+    !    timesteps (depending on the value of nextsw_cday_calc_cam7) every -iradsw hours.
+    ! -- if iradsw is either 0 or 1, the next radiation time is the present time plus 1 timestep.
 
     ! input/output variables
     real(r8)    , intent(in) :: julday
@@ -862,10 +878,21 @@ real(R8) function getNextRadCDay( julday, tod, stepno, dtime, iradsw )
 
     if (liradsw > 1) then
        delta_radsw = liradsw * dtime
-       if (mod(tod+dtime,delta_radsw) == 0 .and. stepno > 0) then
-          nextsw_cday = julday + 2*dtime/shr_const_cday
+       if (nextsw_cday_calc_cam7) then
+          ! The logic in this block is consistent with the driver ordering in CAM7 and
+          ! later. So this is appropriate when using cplhist forcings generated from CAM7
+          ! or later.
+          if (mod(tod,delta_radsw) == 0 .and. stepno > 0) then
+             nextsw_cday = julday + 1*dtime/shr_const_cday
+          else
+             nextsw_cday = -1._r8
+          end if
        else
-          nextsw_cday = -1._r8
+          if (mod(tod+dtime,delta_radsw) == 0 .and. stepno > 0) then
+             nextsw_cday = julday + 2*dtime/shr_const_cday
+          else
+             nextsw_cday = -1._r8
+          end if
        end if
     else
        nextsw_cday = julday + dtime/shr_const_cday

datm/cime_config/namelist_definition_datm.xml

@@ -335,20 +335,47 @@
     <category>datm</category>
     <group>datm_nml</group>
     <desc>
-      Frequency to update radiation in number of steps (or hours if negative)
-      irdasw is the radiation setting used to compute the next shortwave
-      Julian date.  values greater than 1 set the next radiation to the
-      present time plus 2 timesteps every iradsw.  values less than 0 turn
-      set the next radiation to the present time plus two timesteps every
-      -iradsw hours.  if iradsw is zero, the next radiation time is the
-      present time plus 1 timestep.  (default=0.)
+      Frequency to update radiation in number of steps (or hours if negative).
+
+      irdasw is the radiation setting used to compute the next shortwave Julian date.
+      Values greater than 1 set the next radiation to the present time plus either 1 or 2
+      timesteps (depending on the value of nextsw_cday_calc) every iradsw timesteps.
+      Values less than 0 set the next radiation to the present time plus either 1 or 2
+      timesteps (depending on the value of nextsw_cday_calc) every -iradsw hours. If
+      iradsw is either 0 or 1, the next radiation time is the present time plus 1
+      timestep. (default=1 for non-CPLHIST cases, and -1 for CPLHIST cases.)
     </desc>
     <values>
       <value>1</value>
       <value datm_mode="CPLHIST">-1</value>
     </values>
   </entry>
 
+  <entry id="nextsw_cday_calc">
+    <type>char</type>
+    <category>datm</category>
+    <group>datm_nml</group>
+    <valid_values>cam6,cam7</valid_values>
+    <desc>
+      Logic to use for calculating nextsw_cday. This variable has no effect when iradsw is
+      0 or 1, so by default this only applies in CPLHIST cases.
+
+      For CPLHIST cases, this should agree with the version of CAM (or other atmosphere
+      model) used to generate the CPLHIST forcings; the valid values for this variable are
+      based on this: 'cam6' is appropriate for cases generated with the driver ordering in
+      CAM6 and earlier, and 'cam7' is appropriate for cases generated with the driver
+      ordering in CAM7 and later.
+
+      For 'cam6', the next radiation timestep is set to the present time plus 2 timesteps
+      when mod(tod+dtime,delta_radsw)==0. For 'cam7', the next radiation timestep is set
+      to the present time plus 1 timestep when mod(tod,delta_radsw)==0.
+    </desc>
+    <values>
+      <value>cam6</value>
+      <value datm_mode="CPLHIST">cam7</value>
+    </values>
+  </entry>
+
   <entry id="restfilm">
     <type>char</type>
     <category>datm</category>