From 2372f7fc2fc1c71900f57ce1131c2b323f8803c3 Mon Sep 17 00:00:00 2001 From: Haipeng Lin Date: Wed, 2 Oct 2024 17:52:57 -0400 Subject: [PATCH] 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 | 15 ++++++++------- src/data/cam_thermo.F90 | 4 ++++ src/dynamics/se/dp_coupling.F90 | 16 ++++++++++++++-- src/dynamics/se/dycore/prim_advance_mod.F90 | 6 +++++- 4 files changed, 31 insertions(+), 10 deletions(-) diff --git a/src/data/air_composition.F90 b/src/data/air_composition.F90 index 43ce33d0..443565eb 100644 --- a/src/data/air_composition.F90 +++ b/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 ! diff --git a/src/data/cam_thermo.F90 b/src/data/cam_thermo.F90 index 2c7f62a5..a6de3d2b 100644 --- a/src/data/cam_thermo.F90 +++ b/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)) diff --git a/src/dynamics/se/dp_coupling.F90 b/src/dynamics/se/dp_coupling.F90 index ea251163..1eca0982 100644 --- a/src/dynamics/se/dp_coupling.F90 +++ b/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, & diff --git a/src/dynamics/se/dycore/prim_advance_mod.F90 b/src/dynamics/se/dycore/prim_advance_mod.F90 index 1341b9f4..cf88d7f9 100644 --- a/src/dynamics/se/dycore/prim_advance_mod.F90 +++ b/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