Skip to content

Commit

Permalink
Merge branch 'trhille/mali/develop-231220' (PR #6126)
Browse files Browse the repository at this point in the history 
Update MALI version to include higher-order advection and time integration

This merge includes numerous updates to MALI
The major updates include:
* Higher-order, flux-corrected thickness and tracer transport.
* Three strong stability-preserving Runge-Kutta time integrations
  schemes (a two-stage, second-order scheme; three- and four-stage
  third-order schemes)
Minor updates include:
* Earlier update of clock within time step
* Prevent unrealistic calving into holes within an ice shelf
* Change default start and stop times to avoid unpredictable I/O
  behavior
* Higher-order grounding-line flux calculation
* Several small bug fixes

[NML] - only for cases with MALI
[BFB]
  • Loading branch information
@jonbob
jonbob committed Jan 25, 2024
2 parents fd65dbe + 9af3b1e commit f6f79a2
Show file tree
Hide file tree
Showing 19 changed files with 3,103 additions and 189 deletions.
6 changes: 6 additions & 0 deletions components/mpas-albany-landice/bld/build-namelist
View file
Original file line number Diff line number Diff line change
Expand Up @@ -456,6 +456,8 @@ if ($MALI_DYNAMIC eq 'TRUE') {
} else {
add_default($nl, 'config_tracer_advection', 'val'=>"none");
}
add_default($nl, 'config_horiz_tracer_adv_order');
add_default($nl, 'config_advection_coef_3rd_order');
add_default($nl, 'config_restore_thickness_after_advection');
add_default($nl, 'config_zero_sfcMassBalApplied_over_bare_land');

Expand Down Expand Up @@ -509,6 +511,7 @@ add_default($nl, 'config_damagecalvingParameter');
add_default($nl, 'config_ismip6_retreat_k');
add_default($nl, 'config_calving_error_threshold');
add_default($nl, 'config_distribute_unablatedVolumeDynCell');
add_default($nl, 'config_update_velocity_before_calving');

##################################
# Namelist group: thermal_solver #
Expand Down Expand Up @@ -574,6 +577,8 @@ add_default($nl, 'config_dynamic_thickness');

add_default($nl, 'config_dt');
add_default($nl, 'config_time_integration');
add_default($nl, 'config_rk_order');
add_default($nl, 'config_rk3_stages');
add_default($nl, 'config_adaptive_timestep');
add_default($nl, 'config_min_adaptive_timestep');
add_default($nl, 'config_max_adaptive_timestep');
Expand Down Expand Up @@ -635,6 +640,7 @@ add_default($nl, 'config_print_calving_info');
add_default($nl, 'config_print_thermal_info');
add_default($nl, 'config_always_compute_fem_grid');
add_default($nl, 'config_print_velocity_cleanup_details');
add_default($nl, 'config_check_tracer_monotonicity');

####################################
# Namelist group: subglacial_hydro #
Expand Down
6 changes: 6 additions & 0 deletions components/mpas-albany-landice/bld/build-namelist-section
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,8 @@ add_default($nl, 'config_effective_pressure_max');

add_default($nl, 'config_thickness_advection');
add_default($nl, 'config_tracer_advection');
add_default($nl, 'config_horiz_tracer_adv_order');
add_default($nl, 'config_advection_coef_3rd_order');
add_default($nl, 'config_restore_thickness_after_advection');
add_default($nl, 'config_zero_sfcMassBalApplied_over_bare_land');

Expand Down Expand Up @@ -71,6 +73,7 @@ add_default($nl, 'config_damagecalvingParameter');
add_default($nl, 'config_ismip6_retreat_k');
add_default($nl, 'config_calving_error_threshold');
add_default($nl, 'config_distribute_unablatedVolumeDynCell');
add_default($nl, 'config_update_velocity_before_calving');

##################################
# Namelist group: thermal_solver #
Expand Down Expand Up @@ -136,6 +139,8 @@ add_default($nl, 'config_dynamic_thickness');

add_default($nl, 'config_dt');
add_default($nl, 'config_time_integration');
add_default($nl, 'config_rk_order');
add_default($nl, 'config_rk3_stages');
add_default($nl, 'config_adaptive_timestep');
add_default($nl, 'config_min_adaptive_timestep');
add_default($nl, 'config_max_adaptive_timestep');
Expand Down Expand Up @@ -197,6 +202,7 @@ add_default($nl, 'config_print_calving_info');
add_default($nl, 'config_print_thermal_info');
add_default($nl, 'config_always_compute_fem_grid');
add_default($nl, 'config_print_velocity_cleanup_details');
add_default($nl, 'config_check_tracer_monotonicity');

####################################
# Namelist group: subglacial_hydro #
Expand Down
10 changes: 8 additions & 2 deletions components/mpas-albany-landice/bld/namelist_files/namelist_defaults_mali.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,8 @@
<!-- advection -->
<config_thickness_advection>'fo'</config_thickness_advection>
<config_tracer_advection>'none'</config_tracer_advection>
<config_horiz_tracer_adv_order>3</config_horiz_tracer_adv_order>
<config_advection_coef_3rd_order>0.25</config_advection_coef_3rd_order>
<config_restore_thickness_after_advection>.false.</config_restore_thickness_after_advection>
<config_zero_sfcMassBalApplied_over_bare_land>.true.</config_zero_sfcMassBalApplied_over_bare_land>

Expand Down Expand Up @@ -62,6 +64,7 @@
<config_ismip6_retreat_k>-170.0</config_ismip6_retreat_k>
<config_calving_error_threshold>0.1</config_calving_error_threshold>
<config_distribute_unablatedVolumeDynCell>.false.</config_distribute_unablatedVolumeDynCell>
<config_update_velocity_before_calving>.false.</config_update_velocity_before_calving>

<!-- thermal_solver -->
<config_thermal_solver>'none'</config_thermal_solver>
Expand Down Expand Up @@ -115,6 +118,8 @@
<!-- time_integration -->
<config_dt>'0000-00-01_00:00:00'</config_dt>
<config_time_integration>'forward_euler'</config_time_integration>
<config_rk_order>2</config_rk_order>
<config_rk3_stages>3</config_rk3_stages>
<config_adaptive_timestep>.false.</config_adaptive_timestep>
<config_min_adaptive_timestep>3600.0</config_min_adaptive_timestep>
<config_max_adaptive_timestep>3.15e9</config_max_adaptive_timestep>
Expand All @@ -129,8 +134,8 @@
<!-- time_management -->
<config_do_restart>.false.</config_do_restart>
<config_restart_timestamp_name>'rpointer.glc'</config_restart_timestamp_name>
<config_start_time>'0000-01-01_00:00:00'</config_start_time>
<config_stop_time>'0000-01-01_00:00:00'</config_stop_time>
<config_start_time>'0001-01-01_00:00:00'</config_start_time>
<config_stop_time>'0001-01-01_00:00:00'</config_stop_time>
<config_run_duration>'none'</config_run_duration>
<config_calendar_type CALENDAR="NO_LEAP">'noleap'</config_calendar_type>
<config_calendar_type>'gregorian'</config_calendar_type>
Expand Down Expand Up @@ -159,6 +164,7 @@
<config_print_thermal_info>.false.</config_print_thermal_info>
<config_always_compute_fem_grid>.true.</config_always_compute_fem_grid>
<config_print_velocity_cleanup_details>.false.</config_print_velocity_cleanup_details>
<config_check_tracer_monotonicity>.false.</config_check_tracer_monotonicity>

<!-- subglacial_hydro -->
<config_SGH>.false.</config_SGH>
Expand Down
56 changes: 52 additions & 4 deletions components/mpas-albany-landice/bld/namelist_files/namelist_definition_mali.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -139,15 +139,31 @@ Default: Defined in namelist_defaults.xml
category="advection" group="advection">
Selection of the method for advecting thickness ('fo' = first-order upwinding).

Valid values: 'fo', 'none'
Valid values: 'fo', 'fct', 'none'
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_tracer_advection" type="char*1024"
category="advection" group="advection">
Selection of the method for advecting tracers.

Valid values: 'fo', 'none'
Valid values: 'fo', 'fct', 'none'
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_horiz_tracer_adv_order" type="integer"
category="advection" group="advection">
Order of polynomial used for tracer reconstruction at cell edges

Valid values: 2, 3 and 4
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_advection_coef_3rd_order" type="real"
category="advection" group="advection">
Reconstruction of 3rd-order reconstruction to blend with 4th-order reconstuction. Equivalent to beta in Skamarock and Gassmann (2011) eq. 7. 0 is fully 4th order, 1 is fully 3rd order.

Valid values: any real between 0 and 1
Default: Defined in namelist_defaults.xml
</entry>

Expand Down Expand Up @@ -469,6 +485,14 @@ Valid values: .true. or .false.
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_update_velocity_before_calving" type="logical"
category="calving" group="calving">
If true, add an additional velocity solve between advection and calving. If false, use velocity field from beginning of time step to calculate calving rate. For certain calving laws, like damage threshold calving, it is not necessary to update the velocity before calving, while for von Mises stress and eigencalving, it is more accurate to have an updated velocity state before solving for calvingThickness.

Valid values: .true. or .false.
Default: Defined in namelist_defaults.xml
</entry>


<!-- thermal_solver -->

Expand Down Expand Up @@ -835,9 +859,25 @@ Default: Defined in namelist_defaults.xml

<entry id="config_time_integration" type="char*1024"
category="time_integration" group="time_integration">
Time integration method (currently, only forward Euler is supported).
Time integration method.

Valid values: 'forward_euler'
Valid values: 'forward_euler' or 'runge_kutta'
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_rk_order" type="integer"
category="time_integration" group="time_integration">
Order of Runge-Kutta time integration to use. A value of 1 would be equivalent to forward euler, but will cause an error to avoid unnecessary redundancy. Values of 2 and 3 indicate strong-stability preserving RK2 and RK3. There is currently no support for classical RK2 or RK4 methods.

Valid values: 2 or 3
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_rk3_stages" type="integer"
category="time_integration" group="time_integration">
Number of stages for strong stability preserving RK3 time integration. If set to 3, this involves 3 velocity solves and a maximum CFL fraction of 1. If set to 4, this involves 4 velocity solves, but the maximum CFL fraction is 2.

Valid values: 3 or 4
Default: Defined in namelist_defaults.xml
</entry>

Expand Down Expand Up @@ -1133,6 +1173,14 @@ Valid values: .true. or .false.
Default: Defined in namelist_defaults.xml
</entry>

<entry id="config_check_tracer_monotonicity" type="logical"
category="debug" group="debug">
Check tracer monotonicity at the end of the monotonic advection routine and write warnings to log file if not monotonic.

Valid values: .true. or .false.
Default: Defined in namelist_defaults.xml
</entry>


<!-- subglacial_hydro -->

Expand Down
53 changes: 44 additions & 9 deletions components/mpas-albany-landice/src/Registry.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@
<dim name="nISMIP6OceanLayers" units="unitless"
description="The number of layers in the ISMIP6 ocean temperature dataset."
/>
<dim name="maxTracersAdvect" definition="2" units="unitless"
<dim name="maxTracersAdvect" definition="4" units="unitless"
description="The maximum number of tracers to be advected."
/>
<dim name="nRegions" definition="1" units="unitless"
Expand Down Expand Up @@ -110,12 +110,20 @@
<nml_record name="advection" in_defaults="true">
<nml_option name="config_thickness_advection" type="character" default_value="fo" units="unitless"
description="Selection of the method for advecting thickness ('fo' = first-order upwinding)."
possible_values="'fo', 'none'"
possible_values="'fo', 'fct', 'none'"
/>
<nml_option name="config_tracer_advection" type="character" default_value="none" units="unitless"
description="Selection of the method for advecting tracers."
possible_values="'fo', 'none'"
/>
possible_values="'fo', 'fct', 'none'"
/>
<nml_option name="config_horiz_tracer_adv_order" type="integer" default_value="3" units="unitless"
description="Order of polynomial used for tracer reconstruction at cell edges"
possible_values="2, 3 and 4"
/>
<nml_option name="config_advection_coef_3rd_order" type="real" default_value="0.25" units="unitless"
description="Reconstruction of 3rd-order reconstruction to blend with 4th-order reconstuction. Equivalent to beta in Skamarock and Gassmann (2011) eq. 7. 0 is fully 4th order, 1 is fully 3rd order."
possible_values="any real between 0 and 1"
/>
<nml_option name="config_restore_thickness_after_advection" type="logical" default_value=".false." units="unitless"
description="If true, reset thickness to values at previous timestep after advection occurs. This is used for spinning up tracer fields such as damage. When this is true, geometry changes from surface and basal mass balance (grounded or floating) and facemelting are not retained, but changes from calving are."
possible_values=".true. or .false."
Expand Down Expand Up @@ -284,6 +292,10 @@
description="If true, then distribute unablatedVolumeDynCell among dynamic neighbors when converting ablation velocity to ablation thickness. This should only be used as a clean-up measure, while limiting the timestep based on ablation velocity should be used as the primary method of getting accurate ablation thickness from ablation velocity. If you choose to set config_adaptive_timestep_calvingCFL_fraction much larger than 1.0 (which is not recommended), setting this option to true usually results in more accurate calving behavior. "
possible_values=".true. or .false."
/>
<nml_option name="config_update_velocity_before_calving" type="logical" default_value=".false." units="unitless"
description="If true, add an additional velocity solve between advection and calving. If false, use velocity field from beginning of time step to calculate calving rate. For certain calving laws, like damage threshold calving, it is not necessary to update the velocity before calving, while for von Mises stress and eigencalving, it is more accurate to have an updated velocity state before solving for calvingThickness."
possible_values=".true. or .false."
/>
</nml_record>

<nml_record name="thermal_solver" in_defaults="true">
Expand Down Expand Up @@ -477,8 +489,16 @@
possible_values="Any time interval of the format 'YYYY-MM-DD_HH:MM:SS', but limited by CFL condition."
/>
<nml_option name="config_time_integration" type="character" default_value="forward_euler" units="unitless"
description="Time integration method (currently, only forward Euler is supported)."
possible_values="'forward_euler'"
description="Time integration method."
possible_values="'forward_euler' or 'runge_kutta'"
/>
<nml_option name="config_rk_order" type="integer" default_value="2" units="unitless"
description="Order of Runge-Kutta time integration to use. A value of 1 would be equivalent to forward euler, but will cause an error to avoid unnecessary redundancy. Values of 2 and 3 indicate strong-stability preserving RK2 and RK3. There is currently no support for classical RK2 or RK4 methods."
possible_values="2 or 3"
/>
<nml_option name="config_rk3_stages" type="integer" default_value="3" units="stages"
description="Number of stages for strong stability preserving RK3 time integration. If set to 3, this involves 3 velocity solves and a maximum CFL fraction of 1. If set to 4, this involves 4 velocity solves, but the maximum CFL fraction is 2."
possible_values="3 or 4"
/>
<nml_option name="config_adaptive_timestep" type="logical" default_value=".false." units="unitless"
description="Determines if the time step should be adjusted based on the CFL condition or should be steady in time. If true, the config_dt_* options are ignored."
Expand Down Expand Up @@ -533,11 +553,11 @@
description="Path to the filename for restart timestamps to be read and written from."
possible_values="Path to a file."
/>
<nml_option name="config_start_time" type="character" default_value="0000-01-01_00:00:00" units="unitless"
<nml_option name="config_start_time" type="character" default_value="0001-01-01_00:00:00" units="unitless"
description="Timestamp describing the initial time of the simulation. If it is set to 'file', the initial time is read from the filename specified by config_restart_timestamp_name (defaults to 'restart_timestamp')."
possible_values="'YYYY-MM-DD_HH:MM:SS' (items in the format string may be dropped from the left if not needed, and the components on either side of the underscore may be replaced with a single integer representing the rightmost unit)"
/>
<nml_option name="config_stop_time" type="character" default_value="0000-01-01_00:00:00" units="unitless"
<nml_option name="config_stop_time" type="character" default_value="0001-01-01_00:00:00" units="unitless"
description="Timestamp describing the final time of the simulation. If it is set to 'none' the final time is determined from config_start_time and config_run_duration. If config_run_duration is also specified, it takes precedence over config_stop_time. Set config_stop_time to be equal to config_start_time (and config_run_duration to 'none') to perform a diagnostic solve only."
possible_values="'YYYY-MM-DD_HH:MM:SS' or 'none' (items in the format string may be dropped from the left if not needed, and the components on either side of the underscore may be replaced with a single integer representing the rightmost unit)"
/>
Expand Down Expand Up @@ -638,6 +658,10 @@
description="After velocity is calculated there are a few checks for appropriate values in certain geometric configurations. Setting this option to .true. will cause detailed information about those adjustments to be printed."
possible_values=".true. or .false."
/>
<nml_option name="config_check_tracer_monotonicity" type="logical" default_value=".false."
description="Check tracer monotonicity at the end of the monotonic advection routine and write warnings to log file if not monotonic."
possible_values=".true. or .false."
/>
</nml_record>


Expand Down Expand Up @@ -741,6 +765,7 @@
<var name="uReconstructX" packages="higherOrderVelocity"/>
<var name="uReconstructY" packages="higherOrderVelocity"/>
<var name="basalFrictionFlux"/>
<var name="passiveTracer2d"/>
<var name="damage"/>
<var name="calvingVelocityData"/>
<var name="basalMeltInput" packages="hydro"/>
Expand Down Expand Up @@ -850,6 +875,7 @@
<var name="waterPressure" packages="hydro"/>
<var name="channelArea" packages="hydro"/> <!-- this is only needed if running with channels - could be package-controlled -->
<var name="waterFluxMask" packages="hydro"/>
<var name="passiveTracer2d"/>
<var name="damage"/>
<var name="calvingVelocityData"/>
<var name="damageMax"/>
Expand Down Expand Up @@ -1207,9 +1233,12 @@ is the value of that variable from the *previous* time level!
<var name="floatingBasalMassBalApplied" type="real" dimensions="nCells Time" units="kg m^{-2} s^{-1}"
description="Applied basal mass balance on floating regions"
/>
<var name="calvingMask" type="integer" dimensions="nCells Time" units="m" time_levs="1"
<var name="calvingMask" type="integer" dimensions="nCells Time" units="" time_levs="1"
description="mask of grid cells that should be calved. 0=no calving, 1=should be calved"
/>
<var name="openOceanMask" type="integer" dimensions="nCells Time" units="" time_levs="1"
description="mask of open ocean cells, including floating non-dynamic cells. Used to eliminate holes in ice shelf from calving"
/>
<var name="calvingFrontMask" type="integer" dimensions="nCells Time" units="none" time_levs="1"
description="mask of the calving front position on cells"
/>
Expand Down Expand Up @@ -1270,6 +1299,9 @@ is the value of that variable from the *previous* time level!
<var name="calvingVelocityData" type="real" dimensions="nCells Time" units="m s^{-1}" time_levs="1"
description="rate of calving front retreat due to calving, represented as a velocity normal to the calving front (in the x-y plane), given by the input netCDF file."
/>
<var name="passiveTracer2d" type="real" dimensions="nCells Time" units="none" time_levs="1"
description="passive tracer used to verify advection routines"
/>
<var name="damage" type="real" dimensions="nCells Time" units="none" time_levs="1"
description="Damage is parameterized as the local ice thickness divided by the fracture depth, such that unfractured ice has damage=0 and ice fractured through its full thickness has damage=1"
/>
Expand Down Expand Up @@ -1413,6 +1445,9 @@ is the value of that variable from the *previous* time level!
<var name="layerNormalVelocity" type="real" dimensions="nVertLevels nEdges Time" units="m s^{-1}"
description="horizonal velocity, normal component to an edge, layer midpoint"
/>
<var name="layerThicknessEdgeFlux" type="real" dimensions="nVertLevels nEdges Time" units="m^2 s^{-1}"
description="layer-normal thickness flux on edges"
/>
<var name="normalVelocityInitial" type="real" dimensions="nVertInterfaces nEdges Time" units="m s^{-1}"
description="horizonal velocity, normal component to an edge, computed at initialization"
/>
Expand Down
Loading

0 comments on commit f6f79a2

Please sign in to comment.