diff --git a/components/mpas-albany-landice/src/Registry.xml b/components/mpas-albany-landice/src/Registry.xml index 06bf6c6c9399..18d3ab67c751 100644 --- a/components/mpas-albany-landice/src/Registry.xml +++ b/components/mpas-albany-landice/src/Registry.xml @@ -367,10 +367,10 @@ - + @@ -378,18 +378,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" /> - - - + + @@ -852,8 +842,14 @@ + + + + + + + - + + + + + + @@ -1301,7 +1304,25 @@ is the value of that variable from the *previous* time level! - + + + + + + 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) @@ -586,27 +614,34 @@ 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) +!> The linear relationship is modified by: +!> minDraft: a depth upto which melt is a constant (melt = alpha0 + (minDraft * alpha1)) +!> For select ice shelves, a constant melt rate is set instead of the draft dependence. +!> This constant melt rate is set by adding a selector field for the +!> parameterization type (constant, linear function): paramType +!> and a constantValue field for the constant melt rate values. +!> These last two parameters can be modified to be more universal +!> outside of this subroutine if relevant. !----------------------------------------------------------------------- - 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, draftDepenBasalMeltAlpha1, draftDepenBasalMeltAlpha0, & + draftDepenBasalMelt_minDraft, draftDepenBasalMelt_paramType, draftDepenBasalMelt_constantMeltValue, err) + + use mpas_log, only: mpas_log_write + implicit none !----------------------------------------------------------------- ! input variables @@ -620,7 +655,11 @@ 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 + real (kind=RKIND), dimension(:), pointer :: draftDepenBasalMelt_minDraft ! min draft at which (linear) draft dependent parameterization of basal melt is valid + real (kind=RKIND), dimension(:), pointer :: draftDepenBasalMelt_paramType ! parameterization selector field (linear dependence = 0, constant value = 1) + real (kind=RKIND), dimension(:), pointer :: draftDepenBasalMelt_constantMeltValue ! constant melt rate value for selected ice shelves/regions !----------------------------------------------------------------- ! input/output variables !----------------------------------------------------------------- @@ -635,64 +674,96 @@ 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 - + logical :: print_debug + + ! Ranges for sanity checks (SI units) and debugging, min, max + ! values defined from original parameter fields before + ! interpolation onto MPAS grid. + real(kind=RKIND), parameter :: alpha0Min = -0.0008_RKIND + real(kind=RKIND), parameter :: alpha0Max = 0.002_RKIND + real(kind=RKIND), parameter :: alpha1Min = -5e-06_RKIND + real(kind=RKIND), parameter :: alpha1Max = 10.e-07_RKIND + real(kind=RKIND), parameter :: constantMeltMin = -0.0008_RKIND + real(kind=RKIND), parameter :: constantMeltMax = 0.0005_RKIND + real(kind=RKIND), parameter :: minDraftMin = -2500._RKIND + real(kind=RKIND), parameter :: minDraftMax = -0.0_RKIND + real(kind=RKIND), parameter :: meltMaxAbs = 0.01_RKIND ! kg/m2/s 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. - if ( (li_mask_is_floating_ice(cellMask(iCell))) .and. (lowerSurface(iCell) < -10.0_RKIND) ) then + if ( li_mask_is_floating_ice(cellMask(iCell))) then ! 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 + ! Logs for debugging. + if (nint(draftDepenBasalMelt_paramType(iCell)) /= 0 .and. nint(draftDepenBasalMelt_paramType(iCell)) /= 1) then + call mpas_log_write('Warning: paramType not 0 or 1 at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/draftDepenBasalMelt_paramType(iCell)/)) + end if + + if (draftDepenBasalMeltAlpha0(iCell) < alpha0Min .or. draftDepenBasalMeltAlpha0(iCell) > alpha0Max) then + call mpas_log_write('Warning: alpha0 out of range at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/draftDepenBasalMeltAlpha0(iCell)/)) + end if + + if (draftDepenBasalMeltAlpha1(iCell) < alpha1Min .or. draftDepenBasalMeltAlpha1(iCell) > alpha1Max) then + call mpas_log_write('Warning: alpha1 out of range at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/draftDepenBasalMeltAlpha1(iCell)/)) + end if + + if (draftDepenBasalMelt_constantMeltValue(iCell) < constantMeltMin .or. & + draftDepenBasalMelt_constantMeltValue(iCell) > constantMeltMax) then + call mpas_log_write('Warning: constantMeltValue out of range at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/draftDepenBasalMelt_constantMeltValue(iCell)/)) + end if + + if (draftDepenBasalMelt_minDraft(iCell) > minDraftMax .or. draftDepenBasalMelt_minDraft(iCell) < minDraftMin) then + call mpas_log_write('Warning: minDraft out of expected range at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/draftDepenBasalMelt_minDraft(iCell)/)) + end if + + if (zDraft > 0.0_RKIND) then + call mpas_log_write('Warning: zDraft is positive at iCell=$i value=$r', & + intArgs=(/iCell/), realArgs=(/zDraft/)) + end if + + ! --- Compute melt --- + ! If the draft is deeper than or equal to the minimum threshold, apply linear draft dependence + if (nint(draftDepenBasalMelt_paramType(iCell)) == 1) then + ! No linear draft dependence, set a constant melt rate + ! value melt = constantMeltValue + floatingBasalMassBal(iCell) = draftDepenBasalMelt_constantMeltValue(iCell) + else + ! Linear draft dependence: melt = alpha0 + alpha1 * draft + ! at deeper drafts (zDraft > minDraft) + ! Calculation is performed by comparing absolute + ! magnitudes of zDraft and minDraft, + ! TODO: Clean up for final version after diagnostic checks. + ! For shallower draft, use a melt rate calculated for the minDraft threshold + floatingBasalMassBal(iCell) = draftDepenBasalMeltAlpha0(iCell) + max(abs(zDraft), abs(draftDepenBasalMelt_minDraft(iCell))) * draftDepenBasalMeltAlpha1(iCell) + end if ! linear dependence selector + + ! --- Check for unreasonable melt rates --- + if (abs(floatingBasalMassBal(iCell)) > meltMaxAbs) then + call mpas_log_write('Warning: melt rate out of range at iCell=$i melt=$r', & + intArgs=(/iCell/), realArgs=(/floatingBasalMassBal(iCell)/)) + end if + endif ! ice is floating enddo ! iCell - ! change units from m/yr to kg/m2/s - floatingBasalMassBal(:) = floatingBasalMassBal(:) * config_ice_density / scyr + + ! change units from m/yr to kg/m2/s. Uncomment below if using + ! m/yr. Not required if using SI units. + ! floatingBasalMassBal(:) = floatingBasalMassBal(:) * config_ice_density / scyr - end subroutine basal_melt_thwaites_seroussi + end subroutine basal_melt_draft_dependence !-----------------------------------------------------------------------