Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update depth-dependent ice-shelf basal melt parameterization and add variability input field #103

Open
wants to merge 19 commits into
base: develop
Choose a base branch
from
Open

Update depth-dependent ice-shelf basal melt parameterization and add variability input field #103

Show file tree
Hide file tree
Changes from 15 commits
Commits
Show all changes
19 commits
Select commit Hold shift + click to select a range
07cac84
Add basal melt draft dependence subroutine
Feb 6, 2024
cfe4264
add driver changes for calculate_aislens_melt_variability_adjustment
matthewhoffman Feb 7, 2024
8611707
Define draft dependence function for basal melt adjustment
mshiv Feb 15, 2024
b64b61b
Update draft dependence function for basalMassBalAdjustment
mshiv Feb 15, 2024
8f9f6cb
Add draft dependent co-efficient variables for the aislens ice shelf …
mshiv Feb 16, 2024
9f124fc
Update calculate_aislens_melt_variability_adjustment subroutine
mshiv Feb 16, 2024
2dc1067
Add aislens namelist and streams variables and package
mshiv Feb 16, 2024
bb469c4
Add draft dependence subroutine
mshiv Feb 16, 2024
b9c50e3
Update description of draft dependent melt rate subroutine
mshiv Feb 26, 2024
ec75f7a
Update comments for draft dependent parameterization
mshiv Mar 1, 2024
dce872d
Update config option description for basal mass balance methods
mshiv Mar 1, 2024
9b3cdd0
Cleanup comments
mshiv Mar 8, 2024
23440c8
Rename draft dependent melt parameterization variables
mshiv Mar 8, 2024
c884a25
Update draft dependent melt parameterization variables and definitions
mshiv Mar 8, 2024
0e3a4fb
Update basal melt config options description
mshiv Mar 8, 2024
f50f131
Remove package assignment for floatingBasalMassBalAdjustment
mshiv Mar 8, 2024
b3c1eb0
Clean up draft_dependence ice shelf basal melt method
mshiv Mar 18, 2024
50abf8d
Clean up draft_dependence ice shelf basal melt method
mshiv Mar 18, 2024
a84cd6a
Update floatingBasalMassBal calculation for consistent SI units
mshiv Jul 16, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
37 changes: 21 additions & 16 deletions components/mpas-albany-landice/src/Registry.xml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -348,9 +348,9 @@

<nml_record name="iceshelf_melt" in_defaults="true">
<nml_option name="config_basal_mass_bal_float" type="character" default_value="none" units="unitless"
description="Selection of the method for computing the basal mass balance of floating ice. 'none' sets the basalMassBal field to 0 everywhere. 'file' uses without modification whatever value was read in through an input or forcing file or the value set by an ESM coupler. 'constant', 'mismip', 'seroussi' use hardcoded fields defined in the code applicable to Thwaites Glacier. 'temperature_profile' generates a depth-melt relation based on an ocean temperature profile and sill depth. ISMIP6 is the method prescribed by ISMIP6."
possible_values="'none', 'file', 'constant', 'mismip', 'seroussi', 'temperature_profile', 'ismip6'"
/>
description="Selection of the method for computing the basal mass balance of floating ice. 'none' sets the basalMassBal field to 0 everywhere. 'file' uses without modification whatever value was read in through an input or forcing file or the value set by an ESM coupler. 'constant' uses the constant value defined by the heat flux config_bmlt_float_flux and restricted to the region of the domain defined by config_bmlt_float_xlimit. 'mismip' uses the method prescribed for MISMIP+. 'draft_dependence' calculates draft dependent floatingBasalMassBal using draftDepenBasalMeltAlpha1 and draftDepenBasalMeltAlpha0. 'temperature_profile' generates a depth-melt relation based on an ocean temperature profile and sill depth. ISMIP6 is the method prescribed by ISMIP6."
possible_values="'none', 'file', 'constant', 'mismip', 'draft_dependence', 'temperature_profile', 'ismip6'"
/>
<nml_option name="config_bmlt_float_flux" type="real" default_value="0.0" units="W m^{-2}"
description="Value of the constant heat flux applied to the base of floating ice (positive upward)."
possible_values="Any positive real value"
Expand All @@ -359,18 +359,6 @@
description="x value defining region where bmlt_float_flux is applied; melt only where abs(x) is greater than xlimit."
possible_values="Any positive real value"
/>
<nml_option name="config_basal_mass_bal_seroussi_amplitude" type="real" default_value="0.0" units="m"
description="amplitude on the depth adjustment applied to the Seroussi subglacial melt parameterization"
possible_values="any positive real value"
/>
<nml_option name="config_basal_mass_bal_seroussi_period" type="real" default_value="1.0" units="a"
description="period of the periodic depth adjustment applied to the Seroussi subglacial melt parameterization"
possible_values="any positive real value"
/>
<nml_option name="config_basal_mass_bal_seroussi_phase" type="real" default_value="0.0" units="cycles"
description="phase of the periodic depth adjustment applied to the Seroussi subglacial melt parameterization. Units are cycles, i.e., 0-1"
possible_values="any positive real value"
/>
<!-- Options related to temperature profile ice shelf basal melt param. -->
<nml_option name="config_temperature_profile_melt_scale_factor" type="real" default_value="6.0" units="m yr^-1 (deg C)^-2"
description="The scale factor in the 'temperature_profile' melt parameterization that converts a product of two ocean temperatures to a melt rate. Called kappa in code."
Expand Down Expand Up @@ -681,6 +669,7 @@

<package name="thermal" description="This package includes variables required for the thermal solver."/>

<package name="draftDependentMelt" description="This package includes variables required for the draft dependent basal melt parameterization"/>
</packages>


Expand Down Expand Up @@ -793,7 +782,10 @@
<!-- these variables just for ISMIP6 grounded glacier forcing -->
<var name="ismip6Runoff" packages="ismip6GroundedFaceMelt" />
<var name="ismip6_2dThermalForcing" packages="ismip6GroundedFaceMelt" />
</stream>
mshiv marked this conversation as resolved.
Show resolved Hide resolved
<!-- The following variables are defined for the draft dependent ice-shelf basal melt parameterization -->
<var name="floatingBasalMassBalAdjustment" packages="draftDependentMelt"/>
mshiv marked this conversation as resolved.
Show resolved Hide resolved
<var name="draftDepenBasalMeltAlpha1" packages="draftDependentMelt"/>
<var name="draftDepenBasalMeltAlpha0" packages="draftDependentMelt"/>

<!-- An alternate way to allow the HO variables to exist in a separate file.
<stream name="inputHigherOrderVelocity"
Expand Down Expand Up @@ -894,6 +886,10 @@
<var name="ismip6RunoffCurrent" packages="ismip6GroundedFaceMelt" />
<var name="ismip6_2dThermalForcingCurrent" packages="ismip6GroundedFaceMelt" />
<var name="forcingTimeStamp" packages="ismip6GroundedFaceMelt" />
<!-- the following variables just for the draft dependent ice-shelf basal melt parameterization -->
<var name="floatingBasalMassBalAdjustment" packages="draftDependentMelt"/>
mshiv marked this conversation as resolved.
Show resolved Hide resolved
<var name="draftDepenBasalMeltAlpha1" packages="draftDependentMelt"/>
<var name="draftDepenBasalMeltAlpha0" packages="draftDependentMelt"/>
</stream>


Expand Down Expand Up @@ -1224,6 +1220,15 @@ is the value of that variable from the *previous* time level!
<var name="floatingBasalMassBal" type="real" dimensions="nCells Time" units="kg m^{-2} s^{-1}"
description="Potential basal mass balance on floating regions"
/>
<var name="floatingBasalMassBalAdjustment" type="real" dimensions="nCells Time" units="kg m^{-2} s^{-1}" default_value="0.0"
mshiv marked this conversation as resolved.
Show resolved Hide resolved
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's not clear what the difference is between floatingBasalMassBalAdjustment and draftDepenBasalMeltAlpha0. The descriptions need more detail. For which cases is floatingBasalMassBalAdjustment used?

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Technically, they are similar in that they are both added to the product of draft and draftDepenBasalMeltAlpha1. Main difference is that draftDepenBasalMeltAlpha0 (and draftDepenBasalMeltAlpha1) is time independent while floatingBasalMassBalAdjustment can be time dependent, with the latter being read in via a forcing file input stream.

draftDepenBasalMeltAlpha0 and draftDepenBasalMeltAlpha1 are maps derived via a regression of basal melt with draft.

floatingBasalMassBalAdjustment is treated as an optional additional basal melt forcing value - in the context of the AISLENS experiments, it includes time varying basal melt forcing that can have secular/seasonal trends and variability, or a combination of them. Will update the main description above as well.

description="adjustment to floatingBasalMassBal."
/>
<var name="draftDepenBasalMeltAlpha0" type="real" dimensions="nCells" units="kg m^{-2} s^{-1}" package="draftDependentMelt"
description="Alpha0 (intercept) parameter in (linear) draft dependent calculation of floatingBasalMassBal."
/>
<var name="draftDepenBasalMeltAlpha1" type="real" dimensions="nCells" units="kg m^{-3} s^{-1}" package="draftDependentMelt"
description="Alpha1 (slope) parameter in (linear) draft dependent calculation of floatingBasalMassBal."
/>
<var name="basalMassBalApplied" type="real" dimensions="nCells Time" units="kg m^{-2} s^{-1}"
description="applied basal mass balance"
/>
Expand Down
8 changes: 8 additions & 0 deletions components/mpas-albany-landice/src/mode_forward/mpas_li_core_interface.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -113,6 +113,7 @@ function li_setup_packages(configPool, packagePool, iocontext) result(ierr)
logical, pointer :: SIAvelocityActive
logical, pointer :: hydroActive
logical, pointer :: observationsActive
logical, pointer :: draftDependentMeltActive
logical, pointer :: ismip6ShelfMeltActive
logical, pointer :: ismip6GroundedFaceMeltActive
logical, pointer :: thermalActive
Expand All @@ -132,6 +133,7 @@ function li_setup_packages(configPool, packagePool, iocontext) result(ierr)
call mpas_pool_get_package(packagePool, 'higherOrderVelocityActive', higherOrderVelocityActive)
call mpas_pool_get_package(packagePool, 'hydroActive', hydroActive)
call mpas_pool_get_package(packagePool, 'observationsActive', observationsActive)
call mpas_pool_get_package(packagePool, 'draftDependentMeltActive', draftDependentMeltActive)
call mpas_pool_get_package(packagePool, 'ismip6ShelfMeltActive', ismip6ShelfMeltActive)
call mpas_pool_get_package(packagePool, 'ismip6GroundedFaceMeltActive', ismip6GroundedFaceMeltActive)
call mpas_pool_get_package(packagePool, 'thermalActive', thermalActive)
Expand All @@ -157,6 +159,12 @@ function li_setup_packages(configPool, packagePool, iocontext) result(ierr)
"'config_write_albany_ascii_mesh' is set to .true.")
endif

if (trim(config_basal_mass_bal_float) == 'draft_dependence') then
draftDependentMeltActive = .true.
call mpas_log_write("The 'draftDependentMelt' package and assocated variables have been enabled because " // &
"'config_basal_mass_bal_float' is set to 'draft_dependence'")
endif

if ( (trim(config_basal_mass_bal_float) == 'ismip6') .or. &
((trim(config_front_mass_bal_grounded) == 'ismip6') .and. &
(config_use_3d_thermal_forcing_for_face_melt)) ) then
Expand Down
102 changes: 43 additions & 59 deletions components/mpas-albany-landice/src/mode_forward/mpas_li_iceshelf_melt.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -272,10 +272,14 @@ subroutine li_basal_melt_floating_ice(domain, err)

real (kind=RKIND), dimension(:), pointer :: &
floatingBasalMassBal, & ! basal mass balance for floating ice
floatingBasalMassBalAdjustment, & ! possible adjustment to basal mass balance for floating ice
thickness, & ! ice thickness (m)
lowerSurface, & ! lower surface elevation (m)
bedTopography ! bed topography (m; negative below sea level)

real (kind=RKIND), dimension(:), pointer :: &
draftDepenBasalMeltAlpha0, draftDepenBasalMeltAlpha1 ! variables used for basal melt draft dependence

real(kind=RKIND), pointer :: daysSinceStart

integer :: iCell, jCell, iEdge, iNeighbor, err_tmp
Expand Down Expand Up @@ -475,10 +479,12 @@ subroutine li_basal_melt_floating_ice(domain, err)
! change units from m/s to kg/m2/s
floatingBasalMassBal(:) = floatingBasalMassBal(:) * config_ice_density

elseif (trim(config_basal_mass_bal_float) == 'seroussi') then
elseif (trim(config_basal_mass_bal_float) =='draft_dependence') then
mshiv marked this conversation as resolved.
Show resolved Hide resolved

call basal_melt_thwaites_seroussi(floatingBasalMassBal, daysSinceStart, lowerSurface, cellMask, &
config_sea_level, config_ice_density, nCellsSolve, err_tmp)
call mpas_pool_get_array(geometryPool, 'draftDepenBasalMeltAlpha1', draftDepenBasalMeltAlpha1)
call mpas_pool_get_array(geometryPool, 'draftDepenBasalMeltAlpha0', draftDepenBasalMeltAlpha0)
call basal_melt_draft_dependence(floatingBasalMassBal, daysSinceStart, lowerSurface, cellMask, &
config_sea_level, config_ice_density, nCellsSolve, draftDepenBasalMeltAlpha1, draftDepenBasalMeltAlpha0, err_tmp)
err = ior(err, err_tmp)

elseif (trim(config_basal_mass_bal_float) == 'temperature_profile') then
Expand All @@ -505,6 +511,23 @@ subroutine li_basal_melt_floating_ice(domain, err)

endif

! Apply an adjustment to the BMB field. Default value is 0.
block => domain % blocklist
do while (associated(block))

! get pools
call mpas_pool_get_subpool(block % structs, 'geometry', geometryPool)

! get fields from the geometry pool
call mpas_pool_get_array(geometryPool, 'floatingBasalMassBal', floatingBasalMassBal)
call mpas_pool_get_array(geometryPool, 'floatingBasalMassBalAdjustment', floatingBasalMassBalAdjustment)

floatingBasalMassBal(:) = floatingBasalMassBal(:) + floatingBasalMassBalAdjustment(:)

block => block % next
enddo ! associated(block)


! Allocate scratch fields for flood-fill
! Note: This only supports one block per processor
call mpas_pool_get_subpool(domain % blocklist % structs, 'scratch', scratchPool)
Expand Down Expand Up @@ -585,27 +608,22 @@ end subroutine li_basal_melt_floating_ice

!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
!
! ! routine basal_melt_thwaites_seroussi
! ! routine basal_melt_draft_dependence
!
!> \brief Calculate ice shelf melt rate from depth param.
!> \author William Lipscomb
!> \date November 2015
!> \brief Calculate ice shelf melt rate from depth/draft parameters.
!> \author Shivaprakash Muruganandham
!> \date February 2024
!> \details
!> Melt rate parameterization from:
!> Seroussi, H., Y. Nakayama, E. Larour, D. Menemenlis, M. Morlighem, E. Rignot, and A. Khazendar (2017),
!> Continued retreat of Thwaites Glacier, West Antarctica, controlled by bed topography and ocean circulation,
!> Geophys. Res. Lett., 1-9, doi:10.1002/2017GL072910.
!> for Thwaites Glacier.
!> Specifically, this is a linear fit of melt with shelf draft from the Supplemental Information, Figure S1.
!> The linear relation is modified by a:
!> * depth above which there is no melt (Antarctic Surface Water saturation)
!> * a maximum melt rate (Circumpolar Deep Water saturation)
!> * a depth below which melt stops increasing (minimum sill height)
!> Draft dependent melt rate parameterization:
!> It calculates melt rate as a function of shelf draft.
!> In its current form, this is a linear function form, with the draft
!> dependent melt component calculated as:
!> melt = alpha0 + (draft * alpha1)

!-----------------------------------------------------------------------

subroutine basal_melt_thwaites_seroussi(floatingBasalMassBal, daysSinceStart, lowerSurface, cellMask, &
config_sea_level, config_ice_density, nCellsSolve, err)
subroutine basal_melt_draft_dependence(floatingBasalMassBal, daysSinceStart, lowerSurface, cellMask, &
config_sea_level, config_ice_density, nCellsSolve, draftDepenBasalMeltAlpha0, draftDepenBasalMeltAlpha1, err)

!-----------------------------------------------------------------
! input variables
Expand All @@ -619,7 +637,8 @@ subroutine basal_melt_thwaites_seroussi(floatingBasalMassBal, daysSinceStart, lo
real (kind=RKIND), pointer, intent(in) :: config_ice_density !< ice density
integer, pointer :: &
nCellsSolve !< number of locally owned cells

real (kind=RKIND), dimension(:), pointer :: draftDepenBasalMeltAlpha1 ! slope for (linear) draft dependent parameterization of basal melt
real (kind=RKIND), dimension(:), pointer :: draftDepenBasalMeltAlpha0 ! intercept for (linear) draft dependent parameterization of basal melt
!-----------------------------------------------------------------
! input/output variables
!-----------------------------------------------------------------
Expand All @@ -634,64 +653,29 @@ subroutine basal_melt_thwaites_seroussi(floatingBasalMassBal, daysSinceStart, lo
!-----------------------------------------------------------------
! local variables
!-----------------------------------------------------------------
real (kind=RKIND) :: slopeSer ! slope of relation between depth and melt rate
real (kind=RKIND) :: interceptSer ! depth at which melting goes to 0
real (kind=RKIND) :: maxMeltSer ! maximum allowable melt rate
real (kind=RKIND) :: sillDepth ! depth below which melt rate no longer increases
real (kind=RKIND), pointer :: config_basal_mass_bal_seroussi_amplitude
real (kind=RKIND), pointer :: config_basal_mass_bal_seroussi_period
real (kind=RKIND), pointer :: config_basal_mass_bal_seroussi_phase
real(kind=RKIND) :: hCavity ! depth of ice cavity beneath floating ice (m)
real(kind=RKIND) :: zDraft ! draft of floating ice (m below sea level)
integer :: iCell


err = 0

call mpas_pool_get_config(liConfigs, 'config_basal_mass_bal_seroussi_amplitude', &
config_basal_mass_bal_seroussi_amplitude) ! meters
call mpas_pool_get_config(liConfigs, 'config_basal_mass_bal_seroussi_period', &
config_basal_mass_bal_seroussi_period) ! years
call mpas_pool_get_config(liConfigs, 'config_basal_mass_bal_seroussi_phase', &
config_basal_mass_bal_seroussi_phase) ! cycles

slopeSer = 0.088_RKIND ! slope of relation between depth and melt rate (melt (m/yr) per depth (m))
interceptSer = -100.0_RKIND ! depth (m) at which melting goes to 0 (negative meaning below sea level)
maxMeltSer = 50.0_RKIND ! maximum allowable melt rate (m/yr) (positive meaning melting)
sillDepth = -650.0_RKIND ! depth below which melt stops increasing (m) (negative meaning below sea level)

if (config_basal_mass_bal_seroussi_period <= 0.0_RKIND) then
call mpas_log_write("Value for config_basal_mass_bal_seroussi_period has to be a positive real value.", MPAS_LOG_ERR)
err = ior(err, 1)
endif

! Modify intercept height for variability parameters
interceptSer = interceptSer + config_basal_mass_bal_seroussi_amplitude * &
sin( (2.0_RKIND * pii / config_basal_mass_bal_seroussi_period) * (daysSinceStart/365.0_RKIND) &
+ 2.0_RKIND * pii * config_basal_mass_bal_seroussi_phase)

! Initialize before computing
floatingBasalMassBal(:) = 0.0_RKIND

do iCell = 1, nCellsSolve

! Shut off melt at an arbitrary shallow depth to discourage ice from disappearing.
mshiv marked this conversation as resolved.
Show resolved Hide resolved
if ( (li_mask_is_floating_ice(cellMask(iCell))) .and. (lowerSurface(iCell) < -10.0_RKIND) ) then
if ( li_mask_is_floating_ice(cellMask(iCell))) then
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

To be completely accurate, this should probably treat dynamic and non-dynamic cells differently. Melt rates for non-dynamic cells should probably be calculated using the draft values from neighboring dynamic cells. @matthewhoffman, what do you think?

Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we should probably be doing that, but if we include that, it should probably be done for all melt options. Looks like we sort of say we are doing that at line 586 but it's not actually doing it. Or maybe it should happen where BMB is actually applied, using the dynamic mask that exists at that point. I think that should be a separate PR, given it would change the application of all melt rate methods.

! ice is present and floating

zDraft = lowerSurface(iCell) - config_sea_level
! Coefficients for m/yr melt rate (in units of Seroussi figure but without negative meaning melting)
floatingBasalMassBal(iCell) = max(-1.0_RKIND * maxMeltSer, min(0.0_RKIND, slopeSer * &
(max(zDraft, sillDepth) - interceptSer)))

endif ! ice is present
floatingBasalMassBal(iCell) = draftDepenBasalMeltAlpha0(iCell) + zDraft * draftDepenBasalMeltAlpha1(iCell)
endif ! ice is floating
enddo ! iCell

! change units from m/yr to kg/m2/s
floatingBasalMassBal(:) = floatingBasalMassBal(:) * config_ice_density / scyr
mshiv marked this conversation as resolved.
Show resolved Hide resolved


end subroutine basal_melt_thwaites_seroussi
end subroutine basal_melt_draft_dependence
!-----------------------------------------------------------------------


Expand Down