E3SM-Project / scream

Fork of E3SM used to develop exascale global atmosphere model written in C++
https://e3sm-project.github.io/scream/
Other
76 stars 54 forks source link

Inconsistent use of dry/wet variables in SHOC #2168

Closed quantheory closed 1 year ago

quantheory commented 1 year ago

Figured I'd post the other issue related to the conservation errors @oksanaguba has been tracking down.

For SHOC to conserve water mass, it really needs to conserve $\rho r_t$, where $r_t$ is total water mixing ratio and $\rho$ is dry air density if $r_t$ is a dry mixing ratio, or total density if $r_t$ is a wet mixing ratio. In particular, the net flux of water vapor across an interface due to sub-grid-scale mixing is

$$ \rho \overline{w'r'_t} $$

where $\rho$ again has to be dry mixing ratio if $r_t$ is a dry mixing ratio. However, the surface value of $\overline{w'r'_t}$ passed to SHOC is currently calculated using

      wprtp_sfc(i)  = surf_evap(i)/rrho_i(i,nlevi_v)[nlevi_p];

(and something equivalent in the Fortran code), where rrho_i is total ("wet") density at the interface, and thus inconsistent with the dry mixing ratios currently being sent to SHOC.

This flux miscalculation is only an issue at the surface, but there is another problem, which is with the turbulent advection term affecting water mass, defined as:

$$ \frac{1}{\rho} \frac{\partial \rho \overline{w'r'_t}}{\partial z} $$

Using a "wet" $\rho$ here could cause conservation violations internal to the column, wherever there is a vertical gradient of humidity, not just at the interface.

It's not really clear to me what the correct way to deal with this is. A plausible solution to this would be to simply use a "dry" rho everywhere in SHOC, but there might be side effects of this because of other uses of rho in the code. For instance, I think that changing the rho used in the momentum budget would amount to neglecting momentum and kinetic energy associated with water mass (would this be more or less consistent with what the dycore does?). I don't know what the effect is of changing the rho used in equations related to potential temperature or buoyancy, or anywhere else it might show up. It seems like it would affect energy conservation.

Another plausible solution would be to just switch everything back to wet mixing ratios, although if that's not consistent with the energy we want the model to conserve, that feels like just kicking the can down the road.

oksanaguba commented 1 year ago

Just to reiterate -- switching to "wet" fixes most errors in shoc.

Energy can be formulated with wet or dry rho, but it is the same energy, so, it is independent of whether shoc uses wet or dry mmrs. Clubb (and shoc, as far as i know) uses a brute force fix for energy for the turbulence part anyway.

quantheory commented 1 year ago

It's true that if you swap out wet vs. dry rho, you can reformulate the model to conserve "the same" energy, but if you just change rho from wet to dry (or vice versa) without making any other changes to the model, then the set of conserved quantities is not going to be the same, and that includes the energy. So that's all I'm saying. (That said, if SHOC is relying on a fixer to conserve energy anyway, maybe it doesn't matter.)

PeterCaldwell commented 1 year ago

I'm pretty sure that SHOC doesn't really care whether it uses dry or wet mixing ratio (so long as the user is consistent about this choice). My recollection is that we chose to use dry mixing ratio for SHOC just for consistency with P3, which needs dry mixing ratio to be strictly correct. If reverting SHOC to use wet mixing ratio fixes our conservation issues and/or simplifies the code, we should just do that (and add a code comment explaining that everything is formulated in terms of wet mixing ratio and total pressure). For provenance, here's the git issue where we hashed all this out previously: https://github.com/E3SM-Project/scream/issues/1066.

PeterCaldwell commented 1 year ago

So I think the plan here is to just revert the wet/dry SHOC conversion, right? In other words, SHOC will just do everything in wet coordinates (which is what it gets by default) and will use full pressure and density instead of dry versions... I think this should be easy to git revert. Will someone volunteer?

bartgol commented 1 year ago

Given that it happened long ago, and there are prob conflicts doing a revert, it might be easier to just go in and modify the code. But I guess one can always try first...

PeterCaldwell commented 1 year ago

Good point. Though now that I think of it, @oksanaguba was somehow able to do runs without wet/dry conversion, so maybe she has a version that works?

oksanaguba commented 1 year ago

should i make a PR on this? my plan would be to remove dry conversion in shoc (Peter said before it is ok to keep it wet, correct?) and use pdel/pdeldry for P3 conversion instead of qv. I think using qv means that someone always need to keep track of qv_init (at the beginning of time step), it seems like too much work. I will see if there are other potential issues in p3 (like dz), and then people could review it. i would make sure diagnostics for the new version looks ok, too.

PeterCaldwell commented 1 year ago

That would be great. I was thinking of making a PR for each of the 3 issues, but one PR probably makes sense. I agree with your list of tasks: 1). remove dry conversion in SHOC, 2). pdel/pdeldry for P3 conversion, and dz using dry mass and density instead of virtual temperature. There may be more (like that additional conservation error @bogensch was going to talk with you about), but that's a good start.