search

E3SM-Project / scream

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

exner is still incorrectly defined in p3 #1020

Closed PeterCaldwell closed 3 years ago

PeterCaldwell commented 3 years ago

Matt Norman pointed out that https://github.com/E3SM-Project/scream/blob/90c85ecac8f6e1c482f06fb281613b41f7f85124/components/eam/src/physics/cam/micro_p3.F90#L452 and https://github.com/E3SM-Project/scream/blob/90c85ecac8f6e1c482f06fb281613b41f7f85124/components/eam/src/physics/cam/micro_p3.F90#L452 of micro_p3.F90 use inv_exner where the standard definition of exner=T/theta (e.g. https://en.wikipedia.org/wiki/Exner_function) suggests NOT the inverse of exner is needed. I don't think there's a numerical error here because exner is also defined wrong in https://github.com/E3SM-Project/scream/blob/888c8fb92c4a5ed2301c03d7fc798b8789ced47e/components/eam/src/physics/cam/micro_p3_interface.F90#L1038 and everyone knows 2 wrongs make a right. We should fix this. Note that @AaronDonahue and I struggled with EAM misdefining exner in the past... It looks like this screw up managed to snare us yet again.

AaronDonahue commented 3 years ago

On the plus side this error is not replicated in SCREAMv1 - see #1019 definition of exner for example.

AaronDonahue commented 3 years ago

And just for reference, I remember Peter and I talking about this a lot in the early porting of P3. And based on how exner is defined in the dynamics->physics coupling we must have adopted that approach, even though it is also wrong, https://github.com/E3SM-Project/scream/blob/888c8fb92c4a5ed2301c03d7fc798b8789ced47e/components/eam/src/dynamics/se/dp_coupling.F90#L555-L556 I'm surprised that we went with this definition instead of the correct definition. Oh well. Good time to fix it though.

PeterCaldwell commented 3 years ago

I recall @bartgol convincing us both that the dynamics definition was correct... though upon further review it is clearly the reciprocal of what it is supposed to be (noting in particular that cappa=+R/cp on https://github.com/E3SM-Project/scream/blob/888c8fb92c4a5ed2301c03d7fc798b8789ced47e/components/eam/src/control/physconst.F90#L89 . It's also weird that dynamics uses phys_state(lchnk)%pint(i,pver+1) as p0. I think this means that the actual surface pressure rather than 1000 mb is used in dynamics exner. That's definitely non-standard.

AaronDonahue commented 3 years ago

Oh yeah, that is for sure incorrect. We have a git issue that is rather old still open on E3SM about the incorrect usage of pint vs p0. What we found is that the exner definition in this context is only used in the surface coupling and it isn't really exner that is used, instead it is this weird formulation. So essentially, it is named exner but really isn't, it isn't used anywhere so it has never been a problem, and it only became a problem when we started porting P3 and discovered that the state%exner was incorrect and we couldn't use it. see: https://github.com/E3SM-Project/E3SM/issues/2680

AaronDonahue commented 3 years ago

It looks like SHOC and CLUBB use a definition similar to P3. See: https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/shoc_intr.F90#L711 and https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/clubb_intr.F90#L1665

AaronDonahue commented 3 years ago

but CLUBB converts it back to the true exner later: https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/clubb_intr.F90#L1806

AaronDonahue commented 3 years ago

Keeping on the macrophysics assessment of usage of exner. It also looks like CLUBB uses the state%exner variable which is noted as being incorrect, at least from the perspective of being a pure definition of exner. See https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/clubb_intr.F90#L2502-L2509, https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/clubb_intr.F90#L2502-L2509, and https://github.com/E3SM-Project/scream/blob/3bd49c26f3471b036aa211f73df56c15d4c5133f/components/eam/src/physics/cam/clubb_intr.F90#L2502-L2509 for example.

bartgol commented 3 years ago

I recall @bartgol convincing us both that the dynamics definition was correct... though upon further review it is clearly the reciprocal of what it is supposed to be (noting in particular that cappa=+R/cp on

https://github.com/E3SM-Project/scream/blob/888c8fb92c4a5ed2301c03d7fc798b8789ced47e/components/eam/src/control/physconst.F90#L89 . It's also weird that dynamics uses phys_state(lchnk)%pint(i,pver+1) as p0. I think this means that the actual surface pressure rather than 1000 mb is used in dynamics exner. That's definitely non-standard.

My first 4 google hits for "exner function" are: first, second, third, fourth. All of those claim exner=T/Theta (modulo some constants). That's what Homme uses too.

I apoligize if I pushed you to the dark side. I'll blame it on Page and Brin though...

PeterCaldwell commented 3 years ago

@bartgol - maybe you were telling us that the dycore did things the correct way (like the 4 citations you mention do) and we agreed with you (because you were right) without noticing that the physics was doing things wrong. I can't imagine that you would have said something and Aaron and I would have checked it and still got things backwards. Sorry for besmirching your good name!

PeterCaldwell commented 3 years ago

Closed via #1028