Modify quadrature used for ice shelf viscosity #468
Conversation
Codecov Report
@@ Coverage Diff @@
## dev/gfdl #468 +/- ##
============================================
- Coverage 37.83% 37.81% -0.02%
============================================
Files 270 270
Lines 78340 78381 +41
Branches 14503 14515 +12
============================================
+ Hits 29639 29640 +1
- Misses 43297 43338 +41
+ Partials 5404 5403 -1
... and 1 file with indirect coverage changes 📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more |
marshallward
force-pushed
the
dev/gfdl
branch
from
8a66adc to
9de6ce7
Compare
Previously, when ice_viscosity_compute == MODEL, ice shelf viscosity was calculated using 1 quadrature point per cell; however, this quadrature point was not centered in the cell. -Added a routine bilinear_shape_fn_grid_1qp, which is used to instead calculate viscosity for a single cell-centered quadrature point -Added an option to define parameter ice_viscosity_compute = MODEL_QUADRATURE, where viscosity is calculated at the same four quadrature points used during the SSA solution (the typical approach in finite element codes). Note that when using one quadrature point (ice_viscosity_compute == MODEL), array ice_visc(:,:) is the cell-centered ice viscosity. When using four quadrature points (ice_viscosity_compute == MODEL_QUADRATURE), ice_visc(:,:) is the cell-centered ice viscosity divided by the effective stress, Ee, which varies between each quadrature point and is saved in a separate array CS%Ee(:,:,:). In the SSA, ice viscosity is calculated for a cell with indices i,j as ice_visc(i,j) * Ee, where Ee=1 if using one quadrature point for viscosity and Ee=CS%Ee(i,j,k) if using four quad points (where k the current quad point). Also note the post_data call for ice_visc when using four quad points, where Ice_visc is outputted for visualization from a cell as ice_visc(i,j)*Ee_av(i,j), where Ee_av(i,j) is the average Ee in the cell.
alex-huth
force-pushed
the
IS_viscosity_quad_edits
branch
from
c95b33f to
395fdb6
Compare
There was a problem hiding this comment.
Overall these changes seem sensible to me.
I do have one coding-style suggestion to make, however, @alex-huth. This new code repeatedly compares a string with a specified value within the time stepping of the model. This can be inefficient. Rather than repeatedly comparing the string CS%ice_viscosity_compute with "MODEL_QUADRATURE" or "MODEL", etc., it would be better to follow the practice elsewhere in the MOM6 code and doing this string comparison once when the string is first read in, and then using it to set one or more logical variables in the control structure. These logical variables can then be used in place of all the other string comparisons. This is not urgent, so I am approving the code as written, but it would be I good idea to make this change at some point.
|
This has passed pipeline testing at https://gitlab.gfdl.noaa.gov/ogrp/MOM6/-/pipelines/20896. |
Previously, when ice_viscosity_compute == MODEL, ice shelf viscosity was calculated using 1 quadrature point per cell; however, this quadrature point was not centered in the cell.
-Added a routine bilinear_shape_fn_grid_1qp, which is used to instead calculate viscosity for a single, cell-centered quadrature point
-Also added an option to define parameter ice_viscosity_compute = MODEL_QUADRATURE, where viscosity is instead calculated at the same four quadrature points used during the SSA solution (the typical approach in finite element codes, which should be more accurate).
The effective stress, Ee, varies between each quadrature point. If using one quad point, Ee is incorporated directly into ice_visc(:,:), but if using 4 quad points, it is saved in a separate array CS%Ee(:,:,:). In the latter case, ice viscosity is then incorporated into the SSA for a cell with indices i,j as ice_visc(i,j) * Ee(i,j,k) (where k the current quad point). The post_data call for ice_visc when using 4 quad points is then ice_visc(i,j)*Ee_av(i,j), where Ee_av(i,j) is the average Ee in the cell.