As @olyson pointed out, Sa_dens is currently passed from atm to the mediator in CESM, but is only passed from the mediator to ice, not lnd. Instead, CTSM recomputes atmospheric density:
This leads to some differences in the atmospheric density used in CTSM vs. the one sent from CAM – at least when running CTSM in a datm-forced CPLHIST mode (which is what @olyson was investigating). (We have not checked if the same discrepancy exists when running coupled to CAM. It could be that the discrepancy in CPLHIST mode arises from the fact that we're using 3-hour-average fields and so are calculating density from these 3-hour averages rather than instantaneous values; but it could also be that there are small differences in CTSM's calculation of density relative to CAM's.)
It seems like it would be best to ensure consistency between atmospheric density in atm and lnd by passing the atm values to lnd. I thought we might need to make some changes in CDEPS's datm to calculate atmospheric density for modes that don't provide it, but it looks like this is already done for the clmncep mode that is used to force land cases:
We cannot totally remove the above calculation of forc_rho_not_downscaled_grc because we will still need it in the LILAC coupling and (until it's removed) the MCT coupling. So I suggest moving that calculation to its own subroutine that can be called from the LILAC and MCT caps after the call to derive_quantities.
There could be small differences from making this change, even in I compsets where the calculation of density is effectively the same as is currently done inside CTSM (due to changes in order of operations). Unless we reconcile the differences in order of operations, I would suggest the following path forward:
Add this coupling field in CMEPS (ESCOMP/CMEPS#354)
Add some temporary code in CTSM to verify that the received Sa_dens is the same as the internally-calculated value within some tolerance (calling endrun if not).
With the above temporary code in place, run a few short cases:
An F compset
I compsets with different atmospheric forcing datasets (probably doing it with one forcing dataset is sufficient, but it feels safest to test this for our different forcing datasets in case there's something I'm missing that would cause behavior to differ for these different forcing datasets).
Change CTSM code to get atmospheric density from atm when running with the nuopc mediator, by:
Changing the nuopc lnd_import to get forc_rho_not_downscaled_grc from Sa_dens
Moving the calculation of forc_rho_not_downscaled_grc out of derive_quantities and into a new subroutine
Call that new subroutine from the lilac and mct caps after the call to derive_quantities.
As @olyson pointed out, Sa_dens is currently passed from atm to the mediator in CESM, but is only passed from the mediator to ice, not lnd. Instead, CTSM recomputes atmospheric density:
https://github.com/ESCOMP/CTSM/blob/4e126c75a2eaebc6c1c32da4f728bd0891599ff5/src/cpl/utils/lnd_import_export_utils.F90#L72-L73
This leads to some differences in the atmospheric density used in CTSM vs. the one sent from CAM – at least when running CTSM in a datm-forced CPLHIST mode (which is what @olyson was investigating). (We have not checked if the same discrepancy exists when running coupled to CAM. It could be that the discrepancy in CPLHIST mode arises from the fact that we're using 3-hour-average fields and so are calculating density from these 3-hour averages rather than instantaneous values; but it could also be that there are small differences in CTSM's calculation of density relative to CAM's.)
It seems like it would be best to ensure consistency between atmospheric density in atm and lnd by passing the atm values to lnd. I thought we might need to make some changes in CDEPS's datm to calculate atmospheric density for modes that don't provide it, but it looks like this is already done for the clmncep mode that is used to force land cases:
https://github.com/ESCOMP/CDEPS/blob/b5d4929ae2e6bdb8fa3f459ee3adc1a357e1b19f/datm/datm_datamode_clmncep_mod.F90#L467
We cannot totally remove the above calculation of
forc_rho_not_downscaled_grc
because we will still need it in the LILAC coupling and (until it's removed) the MCT coupling. So I suggest moving that calculation to its own subroutine that can be called from the LILAC and MCT caps after the call toderive_quantities
.There could be small differences from making this change, even in I compsets where the calculation of density is effectively the same as is currently done inside CTSM (due to changes in order of operations). Unless we reconcile the differences in order of operations, I would suggest the following path forward:
forc_rho_not_downscaled_grc
fromSa_dens
forc_rho_not_downscaled_grc
out ofderive_quantities
and into a new subroutinederive_quantities
.