NOAA-GFDL / MOM6

Modular Ocean Model
Other
27 stars 59 forks source link

SET_DTBT using wrong face area #660

Open herrwang0 opened 5 months ago

herrwang0 commented 5 months ago

Background

Subroutine set_dtbt() is used to estimate barotropic time step size. The subroutine is called in

In the subroutine, there are currently three options in calculating/estimating face area

  1. With optional input argument BT_cont
  2. With optional input argument eta, passing to subroutine find_face_areas()
  3. With optional input argument add_SSH, passing to subroutine find_face_areas()

And all of them seem to have an issue.

Issues

1. BT_cont

For reasons unknown, BT_cont option is never used in this subroutine.

2. eta

Instead, eta is used in recalculating dtbt. The problem is that eta is almost never used because of this buggy block (NONLINEAR_BT_CONTINUITY is by default False and not used with USE_BT_CONT_TYPE): https://github.com/NOAA-GFDL/MOM6/blob/7305528fd99bc09bda648693cd7130e4627de25c/src/core/MOM_barotropic.F90#L2869-L2875

3. add_SSH

The bug in the if-block above effectively let the model use zero sea surface height (add_SSH is zero by default) to calculate face area, which may underestimate gravity wave speed.

Moreover, without eta, subroutine find_face_areas() will use (Z_ref-bathymetry) to estimate thickness. This does not apply to cases like the Great Lakes, where 1) Within the domain, the reference height is spatially varying 2) portions of the water bodies are above the sea level reference height . (and runtime parameter SSH_EXTRA, which is currently only used to estimate an upper bound of eta during initialization would not help.)

Solutions

  1. Use BT_cont when possible
  2. Fix the if-block bug related to eta
  3. *Make Z_ref a 2D field, rather than a scalar. I put an asterisk because the scalar should be fine for most cases, even with spatially varying reference heights. But I wonder if something weird may come up with OBCs...

The fixes may make the model slower with shorter dtbt.

Hallberg-NOAA commented 2 months ago

After looking through the logic of when dtbt is being used and reset, I think that we can safely add the option to use BT_cont in the call to set_dtbt() inside of step_MOM_dyn_split_RK2, because the inconsistently value that is set in barotropic_init() is never actually used when set_dtbt() is called from step_MOM_dyn_split_RK2. Doing so should go a long way toward addressing this issue.