This adds the 2D Leith viscosity. #507
Conversation
Viscosity is proportional to the magnitude of the vertical vorticity gradient. This commit adds this viscosity to MOM_hor_visc.F90 by mimicking the Smagorinsky code.
No details.
There was a problem hiding this comment.
This passes are internal tests (adds new inactive code).
There are two minor updates I'd like: remove commented code and correct two log messages. The log messages propagate to all the experiments in MOM6-examples so it's worth getting them right before we merge in the code.
For GFDL, this updates MOM_parameter_doc files so needs to be handled on the command line.
| Shear_mag = sqrt(sh_xx(i,j)*sh_xx(i,j) + & | ||
| 0.25*((sh_xy(I-1,J-1)*sh_xy(I-1,J-1) + sh_xy(I,J)*sh_xy(I,J)) + & | ||
| (sh_xy(I-1,J)*sh_xy(I-1,J) + sh_xy(I,J-1)*sh_xy(I,J-1)))) | ||
| ! Vort_mag = sqrt( & |
There was a problem hiding this comment.
I assume these commented out lines are not needed?
| if (CS%Leith_Ah .or. get_all) then | ||
| call get_param(param_file, mod, "LEITH_BI_CONST",Leith_bi_const, & | ||
| "The nondimensional biharmonic Leith constant, \n"//& | ||
| "typically ??", units="nondim", default=0.0, & |
| if (CS%Leith_Kh .or. get_all) & | ||
| call get_param(param_file, mod, "LEITH_LAP_CONST", Leith_Lap_const, & | ||
| "The nondimensional Laplacian Leith constant, \n"//& | ||
| "often ??", units="nondim", default=0.0, & |
There was a problem hiding this comment.
?? (referring to symbols in log message). If there isn't a typical value to recommend then let's remove the sentence from the documentation.
Deleted some commented lines, per request.
|
Thanks for this contribution, Scott, it looks great. I have three suggestions, though, before we do this pull request. First, although F90 is case insensitive, in MOM6 we use lower and upper case index arguments in a very effective soft-convention to indicate whether a variable is at a tracer/thickness point (e.g., h(i,j)) or a vorticity point (e.g., q(I,J)) or a velocity point (e.g., u(I,j) and v(i,J)). Please consider following this convention for all of the new expressions for the Leith code - it really makes it much easier to understand once you get used to it. Secondly, the size declarations for vort_xy, vort_xy_dx and vort_xy_dy should not be the same, because they are staggered at different points and should have different sizes. Specifically, they Thirdly, consider using |
Changed all indices pertaining to Leith to adhere to capitalization soft convention. Changed the array sizes for vort_xy_** arrays. Changed **1.5 exponents to sqrt()*().
|
|
||
| ! Vorticity gradient | ||
| do J=js-1,Jeq+1 ; do I=is-1,Ieq+1 | ||
| vort_xy_dx(i,J) = CS%DY_dxBu(I,J)*(vort_xy(I,J)*G%IdyCu(I,j) - vort_xy(I-1,J)*G%IdyCu(I-1,j)) |
There was a problem hiding this comment.
Because vort_xy_dx and vort_xy_dy have different valid sizes and are used below with different stencils, they can not be calculated in the same loop-blocks, as is occurring in line 513. Please compare these loop extents with the extents used below in lines 552 and 679 (keeping in mind that in symmetric mode Jsq=js-1) to make sure there are no uninitialized values that are being used. (It looks to me like there will be a problem in symmetric memory mode.)
Fixed loop ranges for vort_xy_dx and vort_xy_dy.
Fixed an inconsistency with loop range when calculating vorticity gradients in non-symmetric memory.
Merged in the code that adds the option to use a Leith viscosity. These changes do not change answers by default, but they do add new run-time options LEITH_KH and LEITH_AH and associated sub-options that could, and the comments describing some other existing options also change, so the MOM_parameter_doc.all and MOM_parameter_doc.short files change.
* +Fix for issue mom-ocean#506, tracer OBC bug - it only happens in the advection for certain flow directions, inflow from OBC plus along-boundary flow. * Tracer OBCs need more of an h halo update. - This one should finally fix the processor count issues with OBCs. * Correct the "if" statement. * +Adding more halo points to an exchange - This will change answers if you start with a non-zero velocity. You need three halo points (or maybe cont_stencil) for the continuity solver. - Also trying to put in some initial DEBUG_OBC code. * Fix some DEBUG_OBC logic * Writing to temporary arrays for tres_xy - Way to trick some compiler. * Another fix for DEBUG_OBC * Fixing the whalo troubles for grids that are 2 wide/long. * Exchange all the h_new points. - without this, because we're remapping all the tres points, it dies in debug mode on bad h_new values. * Trying a different exchange - as it was, it was passing my tests, failing the auto-tests. * Fixing the DEBUG_OBC logging * Putting the logging statement back. - Making an error more verbose too.
Viscosity is proportional to the magnitude of the vertical vorticity
gradient. This commit adds this viscosity to MOM_hor_visc.F90 by
mimicking the Smagorinsky code.