Fix energy leak in GWD #3957

oksanaguba · 2020-11-24T03:41:57Z

Add an energy brute force fix to the GWD. Off by default.

The fix works only if orographic GW are on.

[bfb]

golaz · 2020-11-24T23:45:28Z

@oksanaguba : thanks for starting to work on this.

@rljacob : this is something we would like to consider for v2 when it is ready.

oksanaguba · 2020-11-25T03:42:54Z

Preliminary: Here are 10 years of the run with this PR

I am not sure such variability of restom-ressurf is expected. Green curve is this PR and so, master from 2 days ago + GW fix #1. Blue curve is master from September and the same fix #1.

oksanaguba · 2020-11-25T17:46:01Z

This now is only waiting for more years to complete, along with a baseline run.

I will also put in a small test for this fix.

oksanaguba · 2020-11-29T06:27:35Z

Results for 20 years with the PR and 18 years with the baseline are at https://acme-climate.atlassian.net/wiki/spaces/NGDNA/pages/2064941070/GWD+PR+1 .

singhbalwinder · 2020-12-03T02:18:24Z

Thanks @oksanaguba. The changes look good to me. I have assigned reviewers and I will merge it after it is reviewed.

wlin7 · 2020-12-03T05:17:17Z

components/eam/src/physics/cam/gw_drag.F90

+           ttgw(:,k) = ttgw(:,k) / cpairv(:ncol, k, lchnk)
+        end do
+     else
+        !new code, brute force to conserve energy in each level locally


Hi @oksanaguba , compared to original, with fix1, it seems to amount to throwing away temperature tendency ttgw (which is of diffusion of dry static energy and conversion of tke) from all gw sources. Is this interpretation correct?

i am not sure what you mean. the change in this pr overwrites ttgw at the end of gw_tend call but at that point ttgw is not actually used anywhere below (as it is not used in the original). even for diagnostics, ptend%s is used instead:
call outfld ('TTGW', ptend%s/cpairv(:,:,lchnk), pcols, lchnk)
also, from what i saw ttgw is always ptend%s/c_p .

assuming ptend%u and ptend%v are given, to conserve energy one needs to modify ptend%s (which modifies temperature later, in the interface). there are a few ways to conserve energy, but if ttgw=ptend%s carries some other meaning/quantity that is needed, please, let me know. thanks.

Sorry for not making it clearer. I was focusing on the resulting ptend%s that will be returned to the caller.

When I mentioned ttgw, I didn't mean the final one. I meant the ttgw returned from each call to gw_drag_pro before this conditional block. In the original codes, ptend%s accumulates these ttgw, plus -(udu+vdv). With fix1, the accumulation in ptend%s up to this block is not retained.

ttgw calculated inside gw_drag_pro is the sum of diffusion of dry static energy and conversion of kinetic energy, as the comments and code flow in gw_drag_pro suggest. I have to admit I am not an expert on GWD; just following the codes

yes, the original values in ttgw and ptend%s are lost.

components/eam/src/physics/cam/check_energy.F90

oksanaguba · 2020-12-03T05:50:21Z

@singhbalwinder and @wlin7 just now i pushed a minor fix about namelists logic. before the PR was inserting the new parameter into each namelist, so, the PR would have NLFAILs for each test. Now there is no isse when i run atm_integration suite:

[ac.onguba@blueslogin4 as]$ ./cs.status.20201202_223931_9g7sxu | grep FAIL
    FAIL SMS.ne4_ne4.FC5AV1C-L.anvil_intel.eam-preqx_ftype4 TPUTCOMP Error: Computation time increase > 10 pct from baseline

oksanaguba · 2020-12-03T05:53:06Z

Another detail is that this fix works only if orographic GW are on (use_gw_oro == true). I assume this is the default, but I am not sure about all cases. I will add this to description.

oksanaguba · 2020-12-03T05:58:37Z

Seems that use_gw_oro is always true:

[ac.onguba@blueslogin4 acme-b4]$ grep -r use_gw_oro *
components/eam/doc/ChangeLog: - Three new namelist variables, use_gw_oro, use_gw_front, and
components/eam/bld/namelist_files/namelist_definition.xml:<entry id="use_gw_oro" type="logical" category="gw_drag"
components/eam/bld/build-namelist:    add_default($nl, 'use_gw_oro'    , 'val'=>'.true.');
components/eam/bld/build-namelist:    add_default($nl, 'use_gw_oro'    , 'val'=>'.true.');
components/eam/bld/build-namelist:if ($nl->get_value('use_gw_oro') =~ /$TRUE/io) {
components/eam/src/physics/cam/phys_control.F90:logical, public, protected :: use_gw_oro = .true.
components/eam/src/physics/cam/phys_control.F90:      use_hetfrz_classnuc, use_gw_oro, use_gw_front, use_gw_convect, &
components/eam/src/physics/cam/phys_control.F90:   call mpibcast(use_gw_oro,                      1 , mpilog,  0, mpicom)
components/eam/src/physics/cam/gw_drag.F90:  use phys_control,  only: use_gw_oro, use_gw_front, use_gw_convect, use_gw_energy_fix1
components/eam/src/physics/cam/gw_drag.F90:  orographic_only = (use_gw_oro .and. .not. do_spectral_waves)
components/eam/src/physics/cam/gw_drag.F90:  if (use_gw_oro) then
components/eam/src/physics/cam/gw_drag.F90:  if (use_gw_oro) then

This fix should be revisited in the future at least because of this issue.

oksanaguba · 2020-12-03T23:56:42Z

@singhbalwinder and @wlin7 -- are you ok to merge this? WC wanted it soon.

wlin7 · 2020-12-04T03:30:20Z

@singhbalwinder and @wlin7 -- are you ok to merge this? WC wanted it soon.

Sorry, @oksanaguba , I still have some reservation because the fix appears like discounting part of the gwd processes (in a systematic way, I would say). It would not only remove contribution to ptend%s from ttgw due to gw_oro, but also gw_front and gw_convect. All three sources are active in our typical compsets. Your results no doubt shows that part of calculation is responsible for the energy conservation issue.

Have you checked with Yaga if she is comfortable with the fix1 approach?

oksanaguba · 2020-12-04T03:53:12Z

We are talking to Yaga and her co-workers (i do not see her on github, but i will ask her about this). It will take me 1-2 weeks to understand what would be a good energy fix here. Please, also note that the fix is off by default and that this PR is here because WC group requested it.

Regarding overwriting ptend%s, in convective and frontal GW code there is a call to momentum_energy_conservation(). In there, the second part of the routine does essentially the same thing as the fix introduced here (up to minor details). I do not understand why the routine wasn't called in the oro GW part (aside from tend_levels assigned probably in a wrong way for oro part of the code). That is, in convect and front GW ptend%s was computed for energy conservation similarly to this fix.

Essentially, there is no way to conserve energy without recomputing ptend%s the way the fix here is implemented. There are minor details whether to peanutbutter energy and whether to do it below or above some regions. These details are what I am going to ask Yaga.

Thanks for your thoughts, tagging @golaz in case there is any interest in this discussion.

wlin7 · 2020-12-04T04:55:56Z

Hi @oksanaguba , thanks for pointing out the call to momentum_energy_conservation for gw_convect and gw_front; and it is true doing the conservation constraint in similar spirit. There is still a difference there, though. By calling momentum_energy_conservation, it just subtracts a net energy change below source level, from both ttgw (for output) and ptend%s (which will then be carried forward).

That said, I agree it is ok to merge this PR for now, because it is not on by default, and allow for tests to see the impact on other aspects of the simulation in addition to energy conservation. The impact on simulation results is of my major concern because the fix seems to drop some part of the gwd processes. Maybe the net impact on simulation results are small and all will be good.

wlin7 · 2020-12-04T05:43:16Z

@oksanaguba , as for your question above why there wasn't a call to momentum_energy_conservation for gw_oro, I think it is because the source for gw_oro is at bottom level.

There is probably something we can learn from the treatment in that subroutine: it first computes the net energy gain or loss due to the gwd (convection or frontal), then subtracting that amount of energy over the levels below the source level to maintain energy conservation. That means if we add ptend%s accumulated from gw_convect and gw_front to your fix1 assignment for ptend%s, it would still preserve the conservation property as with the current fix1.

And if that is true, it would mean the source of energy leak is solely in the ttgw from gw_oro in the original update of ptend%s. And if we want to account for gw_oro's ttgw in ptend%s, while applying conservation constraint, we could still call momentum_energy_conservation but change to subtract mass weighted net column energy gain/loss from all levels.

Just some more speculation here, I should say.

oksanaguba · 2020-12-04T06:55:09Z

My understanding is that the fix in momentum_.. routine is peanutbuttered where the GW are active (at least, convect GW are active in region 0<z<h according to their papers), not below. For oro waves, i am not sure what's going on, iirc tend_level was set to 72 everywhere in the run i tried. The way i see it, there is a confusion about state%s and ptend%s. The former is DSE, the latter is the change in enthalpy (roughly, h=c_pT, not s=c_pT+gz). When in EAM ptend%s is computed and applied, it is the same as computing and applying dT, not dDSE. Thus i am not sure about comments for DSE in GW code and whether it is important to have DSE tendency (or, maybe, it should be converted to temperature tendency).

That aside, the difference between fix1 and the original code is that the original code peanutbutters the leak, and fix1 does it locally. I agree that the leak is probably in the oro part only (though i did not check numerically), but momentum_... routine also applies taucd to the momentum flux. There was probably a reason not to apply it in the oro part. That is another question that i will try to figure out.

wlin7 · 2020-12-04T13:09:08Z

My understanding is that the fix in momentum_.. routine is peanutbuttered where the GW are active (at least, convect GW are active in region 0<z<h according to their papers), not below.

You are right the GWD effect more active above the source level. So my understanding of the implementation is that the artificial conservation fix, which is expected to be small, is applied to the levels below (from tend_level to pver). I think this is purposely done to preserve the parameterized effect in the active region. For oro_gw, it make senses tend_level=72 because. that is the source level, thanks for the confirmation.

As for enthalpy h vs DSE, as ptend is for local time changes, dh=dDSE. Therefore, dT=ds/cpair.

oksanaguba · 2020-12-06T06:09:13Z

Do not merge yet.

oksanaguba · 2020-12-07T22:41:13Z

Closing this, will do another PR soon.

oksanaguba added 2 commits November 23, 2020 20:44

start on fix1

73c6970

fixing build

51e814d

rljacob added this to the v2.0beta milestone Nov 25, 2020

oksanaguba self-assigned this Nov 25, 2020

rljacob added the Atmosphere label Nov 25, 2020

oksanaguba added the BFB PR leaves answers BFB label Nov 29, 2020

oksanaguba added 2 commits November 30, 2020 21:16

a comment

1e1e7d8

add a comment

c7d1c7e

rljacob assigned singhbalwinder and unassigned oksanaguba Dec 2, 2020

singhbalwinder requested review from kaizhangpnl and wlin7 December 3, 2020 02:16

singhbalwinder approved these changes Dec 3, 2020

View reviewed changes

wlin7 reviewed Dec 3, 2020

View reviewed changes

components/eam/src/physics/cam/check_energy.F90 Show resolved Hide resolved

oksanaguba changed the title ~~[wip] Fix energy leak in GWD~~ Fix energy leak in GWD Dec 3, 2020

removing the parameter from default namelist

ade0350

oksanaguba added the Do not merge label Dec 6, 2020

oksanaguba closed this Dec 7, 2020

singhbalwinder mentioned this pull request Dec 15, 2020

Fix energy leak from GW, fix number 2 #3975

Merged

rljacob deleted the oksanaguba/eam/gw-energy-fix1 branch January 11, 2021 03:50

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Fix energy leak in GWD #3957

Fix energy leak in GWD #3957

oksanaguba commented Nov 24, 2020 •

edited

Loading

golaz commented Nov 24, 2020

oksanaguba commented Nov 25, 2020 •

edited

Loading

oksanaguba commented Nov 25, 2020

oksanaguba commented Nov 29, 2020

singhbalwinder commented Dec 3, 2020

wlin7 Dec 3, 2020

oksanaguba Dec 3, 2020

wlin7 Dec 3, 2020 •

edited

Loading

oksanaguba Dec 3, 2020

oksanaguba commented Dec 3, 2020

oksanaguba commented Dec 3, 2020 •

edited

Loading

oksanaguba commented Dec 3, 2020

oksanaguba commented Dec 3, 2020

wlin7 commented Dec 4, 2020

oksanaguba commented Dec 4, 2020

wlin7 commented Dec 4, 2020 •

edited

Loading

wlin7 commented Dec 4, 2020

oksanaguba commented Dec 4, 2020

wlin7 commented Dec 4, 2020 •

edited

Loading

oksanaguba commented Dec 6, 2020

oksanaguba commented Dec 7, 2020

Fix energy leak in GWD #3957

Fix energy leak in GWD #3957

Conversation

oksanaguba commented Nov 24, 2020 • edited Loading

golaz commented Nov 24, 2020

oksanaguba commented Nov 25, 2020 • edited Loading

oksanaguba commented Nov 25, 2020

oksanaguba commented Nov 29, 2020

singhbalwinder commented Dec 3, 2020

wlin7 Dec 3, 2020

Choose a reason for hiding this comment

oksanaguba Dec 3, 2020

Choose a reason for hiding this comment

wlin7 Dec 3, 2020 • edited Loading

Choose a reason for hiding this comment

oksanaguba Dec 3, 2020

Choose a reason for hiding this comment

oksanaguba commented Dec 3, 2020

oksanaguba commented Dec 3, 2020 • edited Loading

oksanaguba commented Dec 3, 2020

oksanaguba commented Dec 3, 2020

wlin7 commented Dec 4, 2020

oksanaguba commented Dec 4, 2020

wlin7 commented Dec 4, 2020 • edited Loading

wlin7 commented Dec 4, 2020

oksanaguba commented Dec 4, 2020

wlin7 commented Dec 4, 2020 • edited Loading

oksanaguba commented Dec 6, 2020

oksanaguba commented Dec 7, 2020

oksanaguba commented Nov 24, 2020 •

edited

Loading

oksanaguba commented Nov 25, 2020 •

edited

Loading

wlin7 Dec 3, 2020 •

edited

Loading

oksanaguba commented Dec 3, 2020 •

edited

Loading

wlin7 commented Dec 4, 2020 •

edited

Loading

wlin7 commented Dec 4, 2020 •

edited

Loading