Initialise with true gravity#511
Conversation
Thomas Bendall (tommbendall)
left a comment
There was a problem hiding this comment.
Thanks Chris, sorry if I've jumped the gun here with a review since the PR is still in "draft" state. I'm really excited about this potentially facilitating the use deep gravity.
Along with some minor in-line comments, I have two major suggestions:
-
Do we need to still include all of the different
isotherm*options? I don't think any of theisotherm*options are tested, and since thegeopotoption seems to perform best, should we remove all of theisotherm*options? Or possibly keep only one that preserves the old behaviour? -
Do we need unit-tests for the new kernels?
| =Used when the shallow switch is false and reading a start dump | ||
| =that has been created by a constant gravity model, such as the UM. | ||
| =All methods calculate an Exner field in hydrostatic balance with | ||
| =height-varying gravity. This is done by upward integr |
There was a problem hiding this comment.
there seems to be some of the help section missing here!
| @@ -3247,23 +3221,39 @@ help=Horizontal winds are read in on a W2h space from a start dump. If | |||
| sort-key=02b | |||
| type=logical | |||
|
|
|||
| [namelist:initialization=regravitate] | |||
| compulsory=true | |||
| description=Method for switching gravity | |||
There was a problem hiding this comment.
Could you make this slightly clearer?
| description=Method for switching gravity | |
| description=Method for setting up initial hydrostatic balance when using deep gravity but a start dump that assumed shallow gravity |
| =that has been created by a constant gravity model, such as the UM. | ||
| =All methods calculate an Exner field in hydrostatic balance with | ||
| =height-varying gravity. This is done by upward integr | ||
| =Isothermodynamic: Preserve absolute temperature and pressure on |
There was a problem hiding this comment.
If we're to keep the different isothermal methods (I wonder if we don't!), should we explain the differences between them here?
| ! under which the code may be used. | ||
| !----------------------------------------------------------------------------- | ||
|
|
||
| !> @brief Need to say something here |
There was a problem hiding this comment.
I think you might have forgotten to fill these in!
| end do | ||
|
|
||
| ! Map temperature from newly computed grid back onto current model grid | ||
| call interp( nlayers+1, ht_wt, th, nlayers+1, height_wth(map_wt(1)), & |
There was a problem hiding this comment.
I'm a bit confused about this call, based on the interp function, shouldn't this be height_wth(map_wt(1):map_wt(1)+nlayers)?
It looks like we are passing in an element of an array instead of a section of an array
|
|
||
| ! Map temperature from newly computed grid back onto current model grid | ||
| call interp( nlayers+1, ht_wt, th, nlayers+1, height_wth(map_wt(1)), & | ||
| temperature(map_wt(1)) ) |
There was a problem hiding this comment.
Similarly, I'm a bit confused about this call, based on the interp function, shouldn't this be temperature(map_wt(1):map_wt(1)+nlayers)?
It looks like we are passing in an element of an array instead of a section of an array
|
|
||
| end subroutine interp | ||
|
|
||
| ! Does the input UM data, T(z), pi(z), have the correct z information? |
There was a problem hiding this comment.
Just checking whether you intended to leave this comment in there?
| ! Geopotential height | ||
| do k = 0, nlayers | ||
| ht_wt(k) = planet_radius * height_wth(map_wt(1)+k) / & | ||
| ( planet_radius - height_wth(map_wt(1)+k) ) |
There was a problem hiding this comment.
I'm wondering if we should enforce that this calculation is done at double precision, since planet_radius could be very large and height_wth could be very small
| weight1 = height_w3( map_w3(1) ) - height_wth( map_wt(1) ) | ||
| exner_surf = exner( map_w3(1) ) * ( cp * temp_virt ) / & | ||
| ( cp * temp_virt - ( gravity - coriolis_term( map_wt(1) ) ) * weight1 ) | ||
| theta( map_wt(1) ) = temperature( map_wt(1) ) / exner_surf |
There was a problem hiding this comment.
I'm a bit confused at the surface theta being updated but no other values. Should theta actually be passed into this kernel?
| call invoke( & | ||
|
|
||
| select case( regravitate ) | ||
| case( regravitate_isotherm1 ) |
There was a problem hiding this comment.
I wonder if we need all of these isotherm methods? Or if you want to preserve the current behaviour, we could just keep your isotherm1 method?
PR Summary
Sci/Tech Reviewer: Thomas Bendall (@tommbendall)
Code Reviewer:
Although the shallow atmosphere approximation is not used by UM(ENDgame) , it still retains the constant gravity field of earlier models. If LFRic is to run with true height-varying gravity, initialisation from a UM start-dump must take account of the change in gravity between the two models to avoid a significant change to the meteorology of the data.
The existing code to do this has a known defect, which is to be corrected here. An additional variation of this method is added, along with the option to use a method that takes another approach to the problem.
== Code Changes ==
geo_on_pres_kernel_mod.F90 and pres_lev_diags_alg_mod.x90: Diagnostic of geopotential height on a constant pressure level changed to use the model geopotential, rather than height_w3 field, so that it is valid for shallow=.false. as well as shallow=.true..
rose-meta.conf: New setting in the initialization namelist, regravitate, which selects the method to be used for recalculating theta when switching from constant to height-varying gravity.
map_fd_to_prognostics_alg_mod.x90: Place where regravitate is used.
== New Code ==
initialisation/regrav_isotherm_kernel_mod.F90:
Code Quality Checklist
Testing
trac.log
Test Suite Results - lfric_apps - test_reset_gravity/run2
Suite Information
Task Information
✅ succeeded tasks - 1164
Security Considerations
Performance Impact
AI Assistance and Attribution
Documentation
PSyclone Approval
Sci/Tech Review
(Please alert the code reviewer via a tag when you have approved the SR)
Code Review