Ufs/dev ugwp improvements

[mdtoyNOAA]

This pull request has three main components:

1. To fix some issues with the orographic drag schemes of the "drag_suite" module which have caused model crashes to occur in GFSv17 p8 testing. These issues are described in Issue #20.
2. To add optional diagnostic outputs for each of the components of the drag suite for tuning purposes. This is also described in Issue #20. The optional diagnostic outputs can be selected at run time by the already-existing ldiag_ugwp namelist flag.
3. As part of the 1st component above, at the suggestion of Philip Pegion, we have added the capability to update the u and v winds in each model column based on the tendencies calculated by the mesoscale gravity wave drag and blocking schemes, and using these for the calculation of the tendencies by the small-scale gravity wave drag and turbulent orographic form drag schemes. This sequential updating of the wind profile has been shown to improve the numerical stability. The option to update the winds in this manner is controlled by a new logic flag ugwp_seq_update, which for now, is not a namelist option, but is hardwired by default to .false. in GFS_typedefs.F90. The user can change this to .true. at compile time.

The above modifications involve 14 modules in this repo in order to maintain "backwards compatibility" amongst the following CCPP schemes: unified_ugwp, drag_suite, ugwpv1_gsldrag and cires_ugwp.

RT's and ORT's were conducted on Hera with passing results. Log files are at:
/scratch1/BMC/wrfruc/mtoy/git_local/ufs-weather-model/tests/RegressionTests_hera.intel.log
/scratch1/BMC/wrfruc/mtoy/git_local/ufs-weather-model/tests/OpnReqTests_control_p8_hera.intel.log

Issue: ufs-community/ccpp-physics/#20

The added optional diagnostic outputs are:
dws3dt_ogw -- time-averaged wind speed tendency due to mesoscale gravity wave drag
dws3dt_obl -- time-averaged wind speed tendency due to blocking drag
dws3dt_oss -- time-averaged wind speed tendency due to small-scale gravity wave drag
dws3dt_ofd -- time-averaged wind speed tendency due to turbulent orographic form drag
ldu3dt_ogw -- time-averaged x wind tendency due to mesoscale orographic gravity wave drag
ldu3dt_obl -- time-averaged x wind tendency due to blocking drag
ldu3dt_oss -- time-averaged x wind tendency due to small scale gravity wave drag
ldu3dt_ofd -- time-averaged x wind tendency due to form drag
ldu3dt_ngw -- time-averaged u momentum tendency due to non-stationary gravity wave drag
ldv3dt_ngw -- time-averaged v momentum tendency due to non-stationary gravity wave drag
ldt3dt_ngw -- time-averaged temperature tendency due to non-stationary gravity wave drag
dudt_ogw -- instantaneous x wind tendency from mesoscale orographic gravity wave drag
dvdt_ogw -- instantaneous y wind tendency from mesoscale orographic gravity wave drag
dudt_obl -- instantaneous x wind tendency from blocking drag
dvdt_obl -- instantaneous y wind tendency from blocking drag
dudt_oss -- instantaneous x wind tendency from small scale GWD
dvdt_oss -- instantaneous y wind tendency from small scale GWD
dudt_ofd -- instantaneous x wind tendency from form drag
dvdt_ofd -- instantaneous y wind tendency from form drag
du_ogwcol -- instantaneous integrated x momentum flux from mesoscale orographic gravity wave drag
dv_ogwcol -- instantaneous integrated y momentum flux from mesoscale orographic gravity wave drag
du_oblcol -- instantaneous integrated x momentum flux from blocking drag
dv_oblcol -- instantaneous integrated y momentum flux from blocking drag
du_osscol -- instantaneous integrated x momentum flux from small scale gwd
dv_osscol -- instantaneous integrated y momentum flux from small scale gwd
du_ofdcol -- instantaneous integrated x momentum flux from form drag
dv_ofdcol -- instantaneous integrated y momentum flux from form drag
du3_ogwcol -- time-averaged surface x momentum flux from mesoscale orographic gravity wave drag
dv3_ogwcol -- time-averaged surface y momentum flux from mesoscale orographic gravity wave drag
du3_oblcol -- time-averaged surface x momentum flux from blocking drag
dv3_oblcol -- time-averaged surface y momentum flux from blocking drag
du3_osscol -- time-averaged surface x momentum flux from small scale gravity wave drag
dv3_osscol -- time-averaged surface y momentum flux from small scale gravity wave drag
du3_ofdcol -- time-averaged surface x momentum flux from form drag
dv3_ofdcol -- time-averaged surface y momentum flux from form drag

Note that the relevant diag_table entries for these variables are included in:
ufs-weather-model/tests/parm/diag_table/diag_table_rap

[dustinswales]

Looks good. Only a few minor things.

[mdtoyNOAA]

I have added two plots to show examples of the optional diagnostic outputs available:
Figure 1 shows the variable du_ogwcol
Figure 2 shows the variable ldu3dt_ngw

Fig_1_du_ogwcol.pdf
Fig_2_ldu3dt_ngw.pdf

[ChunxiZhang-NOAA]

@mdtoyNOAA May I ask do you have any particular reason that you didn't add ugwp_seq_update as a namelist option. If not, I would suggest to add it in the namelist.

[mdtoyNOAA]

@ChunxiZhang-NOAA I wanted a chance to further test this feature before making it a "public" namelist option. It's possible that it may not be needed for model stability, as the other fixes I made to the code should solve the problem. For now I thought it would be better to leave it as a compile-time option to avoid confusion, and users getting differing results if they set the option to .true. (I plan on leaving it .false. for my next phase of testing.)

[ChunxiZhang-NOAA]

@mdtoyNOAA Since it is not clear if the sequential updating is needed for model stability and you need further tests, I would suggest to keep this feature in your fork and branch but not commit it into the authoritative repo for now. When it is confirmed that the sequential updating is effective for model stability, you can create a PR to merge this feature into ufs/dev branch.

[mdtoyNOAA]

@mdtoyNOAA Since it is not clear if the sequential updating is needed for model stability and you need further tests, I would suggest to keep this feature in your fork and branch but not commit it into the authoritative repo for now. When it is confirmed that the sequential updating is effective for model stability, you can create a PR to merge this feature into ufs/dev branch.

@ChunxiZhang-NOAA I may need to backtrack on my earlier statement. The sequential updating was found by Phil Pegion to have a stabilizing effect for some tests that he ran. There was consensus amongst EMC participants at recent Unified GWD meetings that it should be included as an option. I'd like to reconsider and suggest having it included as a namelist option.

[yangfanglin]

I agree this can be included, and controlled by a namelist option-- *Fanglin Yang, Ph.D.* *Chief, Model Physics Group* *Modeling and Data Assimilation Branch* *NOAA/NWS/NCEP Environmental Modeling Center* *https://www.emc.ncep.noaa.gov/gmb/wx24fy/fyang/ <https://www.emc.ncep.noaa.gov/gmb/wx24fy/fyang/>*

[mdtoyNOAA]

@ChunxiZhang-NOAA @yangfanglin I'll go ahead and update my branch to include ugwp_seq_update in the namelist

[ChunxiZhang-NOAA]

@mdtoyNOAA Ok, thanks.

[grantfirl]

Looks OK from a CCPP point of view. I'm not qualified to comment much on the science.

[ChunxiZhang-NOAA]

@mdtoyNOAA I noticed that in cires_ugwpv1_module.F90, the 'old' ak and bk method is still used to find the launch level of NGWs. Is it possible to change it to the method in the latest drag_suite.F90. So ak and bk calculated reference level can be completely discarded.

[mdtoyNOAA]

@mdtoyNOAA I noticed that in cires_ugwpv1_module.F90, the 'old' ak and bk method is still used to find the launch level of NGWs. Is it possible to change it to the method in the latest drag_suite.F90. So ak and bk calculated reference level can be completely discarded.

No, the method for finding the launch level for NGWs is not related to the method in the latest drag_suite.F90.

[ChunxiZhang-NOAA]

@mdtoyNOAA Just curious, what if use the actual air pressure for the following code:
do k=1, levs
pmb(k) = ak(k) + prefbk(k) ! use actual air pressure here
zkm(k) = -hpskmalog(pmb(k)/pref) ! use alog(prsi(k)/prsi(1))
enddo
And use height (~5000m) for the following code to find the launch level:
do k=levs, 1, -1
if (pmb(k) .gt. knob_ugwp_palaunch ) exit
enddo

[mdtoyNOAA]

@mdtoyNOAA Just curious, what if use the actual air pressure for the following code:
do k=1, levs
pmb(k) = ak(k) + pref_bk(k) ! use actual air pressure here zkm(k) = -hpskm_alog(pmb(k)/pref) ! use alog(prsi(k)/prsi(1))
enddo
And use height (~5000m) for the following code to find the launch level:
do k=levs, 1, -1
if (pmb(k) .gt. knob_ugwp_palaunch ) exit
enddo

I would like for the code modification to produce identical results. The current method finds a unique model level as the "launch level". ak(k) and bk(k) are only functions of model level -- they essentially provide the vertical coordinate. Using the actual pressure would give different model launch levels depending on the horizontal geographic location. I agree we should replace the "ak,bk" method with a universal method. For now, can we create an "Issue" and solve this in a later PR?

[ChunxiZhang-NOAA]

@mdtoyNOAA Sure, please create an issue for it. The commit process for this PR is scheduled to begin tomorrow afternoon.

[mdtoyNOAA]

@mdtoyNOAA Sure, please create an issue for it. The commit process for this PR is scheduled to begin tomorrow afternoon.

Will do. And thank you for the heads-up about tomorrow.