Fix for b4b to CAM-SIMA (w/ History); add notes on modifications and/…

…or divergence from CAM and/or -SIMA

src/data/air_composition.F90

@@ -559,17 +559,15 @@ subroutine water_composition_update(mmr, ncol, vcoord, to_dry_factor)
       else if (vcoord == vc_dry_pressure) then
          ! SE
 
-         ! TODO hplin 9/17/24: for compatibility with CAM-SIMA that allocates thermodynamic_active_species_idx(0:num_advected)
-         ! (whereas CAM only allocates 1-indexed) I subset it here. But from the meaning of the code arguments
-         ! it is unknown where it was meant to be thermodynamic_active_species_idx_dycore.
-         ! This should be verified during code review.
+         ! **TEMP** TODO hplin 9/17/24:
+         ! for compatibility with CAM-SIMA that allocates thermodynamic_active_species_idx(0:num_advected)
+         ! (whereas CAM only allocates 1-indexed) I subset it here. This should be verified during code
+         ! review.
          call get_cp(mmr(:ncol,:,:), .false., cp_or_cv_dycore(:ncol,:), &
                      factor=to_dry_factor, active_species_idx_dycore=thermodynamic_active_species_idx(1:), &
                      cpdry=cpairv(:ncol,:))
       else if (vcoord == vc_height) then
          ! MPAS
-
-         ! TODO hplin 9/17/24 same here.
          call get_R(mmr(:ncol,:,:), thermodynamic_active_species_idx(1:), &
                     cp_or_cv_dycore(:ncol,:), fact=to_dry_factor, Rdry=rairv(:ncol,:))
 
@@ -745,7 +743,10 @@ subroutine get_cp_1hd(tracer, inv_cp, cp, factor, active_species_idx_dycore, cpd
       end do
 
       if (dry_air_species_num == 0) then
-         ! FIXME hplin 9/13/24: this if-clause is always in CAM but was not in CAM-SIMA
+         ! **TMP** TODO CHECK hplin 9/13/24: this if-clause is always in CAM but was not in CAM-SIMA
+         ! and dry_air_species_num is == 0 in CAM-SIMA as of 10/1/24 which means this clause will
+         ! change answers (?) -- previously in CAM-SIMA there was only a call to get_cp_dry
+         ! and not the two IF-clauses here.
          sum_cp = thermodynamic_active_species_cp(0)
       else if (present(cpdry)) then
          !

src/data/cam_thermo.F90

@@ -233,6 +233,10 @@ subroutine cam_thermo_dry_air_update(mmr, T, ncol, update_thermo_variables, to_d
       real(kind_phys) :: sponge_factor(SIZE(mmr, 2))
       character(len=*), parameter :: subname = 'cam_thermo_dry_air_update: '
 
+      if (.not. update_thermo_variables) then
+         return
+      end if
+
       if (present(to_dry_factor)) then
         if (SIZE(to_dry_factor, 1) /= ncol) then
           call endrun(subname//'DIM 1 of to_dry_factor is'//to_str(SIZE(to_dry_factor,1))//'but should be'//to_str(ncol))

src/dynamics/se/dp_coupling.F90

@@ -690,9 +690,16 @@ subroutine derived_phys_dry(cam_runtime_opts, phys_state, phys_tend)
    ! wet pressure variables (should be removed from physics!)
    factor_array(:,:) = 1.0_kind_phys
    !$omp parallel do num_threads(horz_num_threads) private (k, i, m_cnst)
+   ! **TEMP** TODO CHECK hplin: check indices to use here
+   ! in cam6_3_109 after fix in 6_3_127: loop from dry_air_species_num + 1, thermodynamic_active_species_num
+   ! in cam-sima: factor_array only contains m = ix_q
+   ! hplin -- to get same answers as CAM-SIMA ix_q would need to be used,
+   !          but the CAM version appears to be correct. this should be science checked
+   ! // **TEMP**
    do m_cnst = dry_air_species_num + 1, thermodynamic_active_species_num
       m = thermodynamic_active_species_idx(m_cnst)
-      ! for now m is just == ix_q (one entry)
+   !  write(6,*) "hplin dp_coupling m, m_cnst", m, m_cnst
+   ! m = ix_q
       do k = 1, nlev
          do i = 1, pcols
             ! at this point all q's are dry
@@ -784,6 +791,10 @@ subroutine derived_phys_dry(cam_runtime_opts, phys_state, phys_tend)
    ! Compute molecular viscosity(kmvis) and conductivity(kmcnd).
    ! Update zvirv registry variable; calculated for WACCM-X.
    !-----------------------------------------------------------------------------
+   ! **TEMP** TODO CHECK hplin: CAM has this if-clause for dry_air_species_num > 0
+   ! or otherwise uses zvirv = zvir. CAM-SIMA previously did not have this, and
+   ! instead has a switch for update_thermodynamic_variables. Check if we still want
+   ! this if-clause or change it to something else.
    if (dry_air_species_num > 0) then
       call cam_thermo_dry_air_update( &
            mmr                     = const_data_ptr, & ! dry MMR
@@ -808,6 +819,7 @@ subroutine derived_phys_dry(cam_runtime_opts, phys_state, phys_tend)
    !$omp parallel do num_threads(horz_num_threads) private (k, i)
    do k = 1, nlev
       do i = 1, pcols
+         ! **TEMP** TODO CHECK hplin: CAM and CAM-SIMA version differ
          ! CAM version:
          ! phys_state%exner(i,k) = (phys_state%pint(i,pver+1)/phys_state%pmid(i,k))**cappav(i,k)
 
@@ -839,7 +851,7 @@ subroutine derived_phys_dry(cam_runtime_opts, phys_state, phys_tend)
       end do
    end do
 
-   !Call geopotential_temp CCPP scheme:
+   ! Call geopotential_temp CCPP scheme:
    call geopotential_temp_run(pver, lagrangian_vertical, pver, 1,                        &
                               pverp, 1, num_advected, phys_state%lnpint,                 &
                               phys_state%pint, phys_state%pmid, phys_state%pdel,         &

src/dynamics/se/dycore/prim_advance_mod.F90

@@ -1674,8 +1674,12 @@ subroutine calc_tot_energy_dynamics(elem,fvm,nets,nete,tl,tl_qdp,outfld_name_suf
         call get_dp(elem(ie)%state%Qdp(:,:,:,1:qsize,tl_qdp),MASS_MIXING_RATIO,thermodynamic_active_species_idx_dycore,&
              elem(ie)%state%dp3d(:,:,:,tl),pdel,ps=ps,ptop=hyai(1)*ps0)
         call get_cp(elem(ie)%state%Qdp(:,:,:,1:qsize,tl_qdp),&
-             .false.,cp,dp_dry=elem(ie)%state%dp3d(:,:,:,tl),&
+             .false.,cp,factor=1.0_r8/elem(ie)%state%dp3d(:,:,:,tl),&
              active_species_idx_dycore=thermodynamic_active_species_idx_dycore)
+
+        ! TODO: need to port cam6_3_109 changes to total energy using get_hydrostatic_energy
+        ! https://github.com/ESCOMP/CAM/pull/761/files#diff-946bde17289e2f42e43e64413610aa11d102deda8b5199ddaa5b71e67e5d517a
+
         do k = 1, nlev
           do j=1,np
             do i = 1, np