Skip to content

Commit

Permalink
Merge branch 'cbegeman/ocn/new-floating-land-ice-vars' (PR #5464)
Browse files Browse the repository at this point in the history 
New ocean variables for floating land ice

This PR adds new variables to represent floating land ice. This is a
necessary step toward coupling with an ice sheet model component (MALI).
* landIceMask and landIceFraction (existing variables) now represent
  both grounded and floating land ice.
* landIceFloatingMask and landIceFloatingFraction (new variables)
  represent only floating land ice.
The landIceFloating* variables are used to modulate land ice fluxes,
whereas the landIce* variables are used to mask out certain fields, as
before. landIce* variables also modulate top drag, as we want to include
the effect of top drag on thin film regions to slow flow.

[NML]
[non-BFB]
  • Loading branch information
@jonbob
jonbob committed Apr 5, 2023
2 parents a3b334c + aa144a6 commit c292bec
Show file tree
Hide file tree
Showing 15 changed files with 140 additions and 82 deletions.
2 changes: 1 addition & 1 deletion components/mpas-ocean/bld/build-namelist
View file
Original file line number Diff line number Diff line change
Expand Up @@ -782,7 +782,7 @@ if ($OCN_ISMF eq 'coupled') {
} elsif ($OCN_ISMF eq 'internal') {
add_default($nl, 'config_land_ice_flux_mode', 'val'=>"standalone");
} else {
add_default($nl, 'config_land_ice_flux_mode', 'val'=>"pressure_only");
add_default($nl, 'config_land_ice_flux_mode');
}
add_default($nl, 'config_land_ice_flux_formulation');
add_default($nl, 'config_land_ice_flux_useHollandJenkinsAdvDiff');
Expand Down
12 changes: 6 additions & 6 deletions components/mpas-ocean/bld/namelist_files/namelist_defaults_mpaso.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -333,12 +333,12 @@

<!-- land_ice_fluxes -->
<config_land_ice_flux_mode>'off'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oQU240wLI">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oEC60to30v3wLI">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="ECwISC30to60E1r2">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oRRS30to10v3wLI">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="SOwISC12to60E2r4">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="ECwISC30to60E2r1">'standalone'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oQU240wLI">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oEC60to30v3wLI">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="ECwISC30to60E1r2">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="oRRS30to10v3wLI">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="SOwISC12to60E2r4">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_mode ocn_grid="ECwISC30to60E2r1">'pressure_only'</config_land_ice_flux_mode>
<config_land_ice_flux_formulation>'Jenkins'</config_land_ice_flux_formulation>
<config_land_ice_flux_useHollandJenkinsAdvDiff>.false.</config_land_ice_flux_useHollandJenkinsAdvDiff>
<config_land_ice_flux_attenuation_coefficient>10.0</config_land_ice_flux_attenuation_coefficient>
Expand Down
18 changes: 9 additions & 9 deletions components/mpas-ocean/cime_config/buildnml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -90,18 +90,18 @@ def buildnml(case, caseroot, compname):
decomp_prefix = 'mpas-o.graph.info.'
restoring_file = 'sss.monthlyClimatology.PHC2_salx_040803.oEC60to30v3wLI.nc'
analysis_mask_file = 'masks_SingleRegionAtlanticWTransportTransects_EC60to30v3wLI_171116.nc'
ic_date = '171031'
ic_date = '230220'
ic_prefix = 'oEC60to30v3wLI60lev'
if ocn_ic_mode == 'spunup':
ic_date = '171116'
ic_date = '230220'
ic_prefix = 'oEC60to30v3wLI60lev.restart_theta_year26'

elif ocn_grid == 'ECwISC30to60E1r2':
decomp_date = '200408'
decomp_prefix = 'mpas-o.graph.info.'
restoring_file = 'sss.PHC2_monthlyClimatology.ECwISC30to60E1r02.200408.nc'
analysis_mask_file = 'ECwISC30to60E1r02_transportTransects.nc'
ic_date = '200408'
ic_date = '230220'
ic_prefix = 'ocean.ECwISC30to60E1r2'
if ocn_ic_mode == 'spunup':
logger.warning("WARNING: The specified compset is requesting ocean ICs spunup from a G-case")
Expand All @@ -128,7 +128,7 @@ def buildnml(case, caseroot, compname):
elif ocn_grid == 'oQU240wLI':
decomp_date = '160929'
decomp_prefix = 'mpas-o.graph.info.'
ic_date = '160929'
ic_date = '230220'
ic_prefix = 'ocean.QU.240wLI'
if ocn_ic_mode == 'spunup':
logger.warning("WARNING: The specified compset is requesting ocean ICs spunup from a G-case")
Expand Down Expand Up @@ -159,7 +159,7 @@ def buildnml(case, caseroot, compname):
decomp_prefix = 'mpas-o.graph.info.'
restoring_file = 'sss.monthlyClimatology.PHC2_salx.2004_08_03.180503.nc'
analysis_mask_file = 'masks_SingleRegionAtlanticWTransportTransects_RRS30to10v3wLI.nc'
ic_date = '171109'
ic_date = '230220'
ic_prefix = 'oRRS30to10v3wLI'
if ocn_ic_mode == 'spunup':
logger.warning("WARNING: The specified compset is requesting ocean ICs spunup from a G-case")
Expand Down Expand Up @@ -261,21 +261,21 @@ def buildnml(case, caseroot, compname):
decomp_prefix = 'mpas-o.graph.info.'
restoring_file = 'sss.PHC2_monthlyClimatology.SOwISC12to60E2r4_nomask.210120.nc'
analysis_mask_file = 'SOwISC12to60E2r4_mocBasinsAndTransects20210623.nc'
ic_date = '210107'
ic_date = '230220'
ic_prefix = 'ocean.SOwISC12to60E2r4'
if ocn_ic_mode == 'spunup':
ic_date = '210203'
ic_date = '230220'
ic_prefix = 'mpaso.SOwISC12to60E2r4.rstFromG-anvil'

elif ocn_grid == 'ECwISC30to60E2r1':
decomp_date = '200915'
decomp_prefix = 'mpas-o.graph.info.'
restoring_file = 'sss.PHC2_monthlyClimatology.ECwISC30to60E2r1_nomask.201221.nc'
analysis_mask_file = 'ECwISC30to60E2r1_mocBasinsAndTransects20210623.nc'
ic_date = '210413'
ic_date = '230220'
ic_prefix = 'ocean.ECwISC30to60E2r1'
if ocn_ic_mode == 'spunup':
ic_date = '210414'
ic_date = '230220'
ic_prefix = 'mpaso.ECwISC30to60E2r1.rstFromG-anvil'


Expand Down
10 changes: 5 additions & 5 deletions components/mpas-ocean/driver/ocn_comp_mct.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -790,10 +790,9 @@ end subroutine xml_stream_get_attributes
call t_stopf ('mpaso_mct_init')

call mpas_pool_get_config(domain % configs, 'config_land_ice_flux_mode', config_land_ice_flux_mode)
if ( trim(config_land_ice_flux_mode) .eq. 'pressure_only' ) then
call seq_infodata_PutData( infodata, ocn_prognostic=.true., ocnrof_prognostic=.true., &
ocn_c2_glcshelf=.false.)
else if ( trim(config_land_ice_flux_mode) .eq. 'standalone' ) then
if ( trim(config_land_ice_flux_mode) == 'off' .or. &
trim(config_land_ice_flux_mode) == 'pressure_only' .or. &
trim(config_land_ice_flux_mode) == 'standalone' ) then
call seq_infodata_PutData( infodata, ocn_prognostic=.true., ocnrof_prognostic=.true., &
ocn_c2_glcshelf=.false.)
else if ( trim(config_land_ice_flux_mode) .eq. 'coupled' ) then
Expand Down Expand Up @@ -2784,7 +2783,8 @@ subroutine ocn_export_mct(o2x_o, errorCode) !{{{
!JW o2x_o % rAttr(index_o2x_So_htv, n) = landIceHeatTransferVelocity(i)
!JW o2x_o % rAttr(index_o2x_So_stv, n) = landIceSaltTransferVelocity(i)

if ( trim(config_land_ice_flux_mode) .ne. 'pressure_only' ) then
if ( trim(config_land_ice_flux_mode) .eq. 'standalone' .or. &
trim(config_land_ice_flux_mode) .eq. 'coupled' ) then
o2x_o % rAttr(index_o2x_So_blt, n) = landIceBoundaryLayerTracers(indexBLT,i)
o2x_o % rAttr(index_o2x_So_bls, n) = landIceBoundaryLayerTracers(indexBLS,i)
o2x_o % rAttr(index_o2x_So_htv, n) = landIceTracerTransferVelocities(indexHeatTrans,i)
Expand Down
21 changes: 21 additions & 0 deletions components/mpas-ocean/src/Registry.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1640,6 +1640,8 @@
<var name="landIcePressure"/>
<var name="landIceFraction"/>
<var name="landIceDraft"/>
<var name="landIceFloatingMask"/>
<var name="landIceFloatingFraction"/>
</stream>

<stream name="forcing_data_init"
Expand Down Expand Up @@ -1749,6 +1751,8 @@
<var name="landIcePressure"/>
<var name="landIceFraction"/>
<var name="landIceDraft"/>
<var name="landIceFloatingMask"/>
<var name="landIceFloatingFraction"/>
</stream>
<stream name="block_debug_output"
type="output"
Expand Down Expand Up @@ -1848,6 +1852,8 @@
<var name="landIceDraft"/>
<var name="landIceFraction"/>
<var name="landIceMask"/>
<var name="landIceFloatingMask"/>
<var name="landIceFloatingFraction"/>
<var name="seaIcePressure"/>
<var name="atmosphericPressure"/>
<var_array name="tracersSurfaceValue"/>
Expand Down Expand Up @@ -1990,6 +1996,7 @@
immutable="true">
<var name="landIcePressureForcing"/>
<var name="landIceFractionForcing"/>
<var name="landIceFloatingFractionForcing"/>
<var name="landIceDraftForcing"/>
</stream>

Expand Down Expand Up @@ -2106,6 +2113,8 @@
<var name="nAccumulatedCoupled"/>
<var name="landIceFraction"/>
<var name="landIceMask"/>
<var name="landIceFloatingMask"/>
<var name="landIceFloatingFraction"/>
<var_array name="landIceInterfaceTracers"/>
<var_array name="landIceBoundaryLayerTracers"/>
<var name="landIceFrictionVelocity"/>
Expand Down Expand Up @@ -3612,6 +3621,14 @@
description="Mask indicating where land-ice is present (1) or absent (0)"
packages="landIcePressurePKG"
/>
<var name="landIceFloatingFraction" type="real" dimensions="nCells Time" units="1"
description="The fraction of each cell covered by an ice shelf"
packages="landIcePressurePKG"
/>
<var name="landIceFloatingMask" type="integer" dimensions="nCells Time" default_value="-1"
description="Mask indicating where an ice shelf is present (1) or absent (0)"
packages="landIcePressurePKG"
/>
<var name="landIcePressure" type="real" dimensions="nCells Time" units="Pa"
description="Pressure defined at the sea surface due to land ice."
packages="landIcePressurePKG"
Expand Down Expand Up @@ -3809,6 +3826,10 @@
description="The fraction of each cell covered by land ice"
packages="timeVaryingLandIceForcingPKG"
/>
<var name="landIceFloatingFractionForcing" type="real" dimensions="nCells Time" units="1"
description="The fraction of each cell covered by land ice"
packages="timeVaryingLandIceForcingPKG"
/>
<var name="landIcePressureForcing" type="real" dimensions="nCells Time" units="Pa"
description="Pressure defined at the sea surface due to land ice."
packages="timeVaryingLandIceForcingPKG"
Expand Down
8 changes: 4 additions & 4 deletions components/mpas-ocean/src/analysis_members/mpas_ocn_global_stats.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -234,7 +234,7 @@ subroutine ocn_compute_global_stats(domain, timeLevel, err)!{{{

integer, parameter :: kMaxVariables = 1024 ! this must be a little more than double the number of variables to be reduced
integer, dimension(:), pointer :: minLevelCell, minLevelEdgeBot, minLevelVertexTop, maxLevelCell, &
maxLevelEdgeTop, maxLevelVertexBot, landiceMask
maxLevelEdgeTop, maxLevelVertexBot, landIceFloatingMask

type (MPAS_Time_type) :: xtime_timeType, simulationStartTime_timeType
type (MPAS_TimeInterval_type) :: timeStep
Expand Down Expand Up @@ -364,7 +364,7 @@ subroutine ocn_compute_global_stats(domain, timeLevel, err)!{{{
call mpas_pool_get_array(forcingPool, 'rainFlux', rainFlux)
call mpas_pool_get_array(forcingPool, 'frazilLayerThicknessTendency', frazilLayerThicknessTendency)
call mpas_pool_get_array(forcingPool, 'landIceFreshwaterFlux', landIceFreshwaterFlux)
call mpas_pool_get_array(forcingPool, 'landIceMask', landIceMask)
call mpas_pool_get_array(forcingPool, 'landIceFloatingMask', landIceFloatingMask)

call mpas_pool_get_array(tracersSurfaceFluxPool, 'activeTracersSurfaceFlux', activeTracersSurfaceFlux)
call mpas_pool_get_array(tracersSurfaceRestoringFieldsPool, 'activeTracersPistonVelocity', &
Expand All @@ -381,8 +381,8 @@ subroutine ocn_compute_global_stats(domain, timeLevel, err)!{{{
areaEdge = dcEdge(1:nEdgesSolve)*dvEdge(1:nEdgesSolve)

allocate(landIceFloatingArea(1:nCellsSolve))
if ( associated(landIceMask) ) then
landIceFloatingArea = landIceMask(1:nCellsSolve)*areaCell(1:nCellsSolve)
if ( associated(landIceFloatingMask) ) then
landIceFloatingArea = landIceFloatingMask(1:nCellsSolve)*areaCell(1:nCellsSolve)
else
landIceFloatingArea = 0.
end if
Expand Down
14 changes: 6 additions & 8 deletions components/mpas-ocean/src/analysis_members/mpas_ocn_mixed_layer_depths.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -176,7 +176,6 @@ subroutine ocn_compute_mixed_layer_depths(domain, timeLevel, err)!{{{

real (kind=RKIND), dimension(:), pointer :: tThreshMLD, tGradientMLD, landIceDraft
real (kind=RKIND), dimension(:), pointer :: dGradientMLD, landIcePressure
integer, dimension(:), pointer :: landIceMask
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
integer :: interp_local
real (kind=RKIND), allocatable, dimension(:,:) :: densityGradient, temperatureGradient
Expand Down Expand Up @@ -214,7 +213,6 @@ subroutine ocn_compute_mixed_layer_depths(domain, timeLevel, err)!{{{
call mpas_pool_get_array(tracersPool, 'activeTracers', tracers, timeLevel)
call mpas_pool_get_array(forcingPool, 'landIcePressure', landIcePressure)
call mpas_pool_get_array(forcingPool, 'landIceDraft', landIceDraft)
call mpas_pool_get_array(forcingPool, 'landIceMask', landIceMask)
call mpas_pool_get_array(meshPool, 'latCell', latCell)
call mpas_pool_get_array(meshPool, 'lonCell', lonCell)

Expand Down Expand Up @@ -286,11 +284,11 @@ subroutine ocn_compute_mixed_layer_depths(domain, timeLevel, err)!{{{
!$omp end do
!$omp end parallel
if (associated(landIceMask) ) then
if (associated(landIceDraft) ) then
!$omp parallel
!$omp do schedule(runtime)
do iCell = 1, nCells
tThreshMLD(iCell) = tThreshMLD(iCell) - abs(landIceDraft(iCell))*landIceMask(iCell)
tThreshMLD(iCell) = tThreshMLD(iCell) - abs(landIceDraft(iCell))
end do
!$omp end do
!$omp end parallel
Expand Down Expand Up @@ -352,11 +350,11 @@ subroutine ocn_compute_mixed_layer_depths(domain, timeLevel, err)!{{{
!$omp end parallel

!normalize MLD to top of ice cavity
if (associated(landIceMask) ) then
if (associated(landIceDraft) ) then
!$omp parallel
!$omp do schedule(runtime)
do iCell = 1, nCells
tGradientMLD(iCell) = tGradientMLD(iCell) - abs(landIceDraft(iCell))*landIceMask(iCell)
tGradientMLD(iCell) = tGradientMLD(iCell) - abs(landIceDraft(iCell))
end do
!$omp end do
!$omp end parallel
Expand Down Expand Up @@ -417,11 +415,11 @@ subroutine ocn_compute_mixed_layer_depths(domain, timeLevel, err)!{{{
!$omp end do
!$omp end parallel
if (associated(landIceMask) ) then
if (associated(landIceDraft) ) then
!$omp parallel
!$omp do schedule(runtime)
do iCell = 1, nCells
dGradientMLD(iCell) = dGradientMLD(iCell) - abs(landIceDraft(iCell))*landIceMask(iCell)
dGradientMLD(iCell) = dGradientMLD(iCell) - abs(landIceDraft(iCell))
end do
!$omp end do
!$omp end parallel
Expand Down
29 changes: 23 additions & 6 deletions components/mpas-ocean/src/mode_init/mpas_ocn_init_global_ocean.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1073,13 +1073,16 @@ subroutine ocn_init_setup_global_ocean_interpolate_land_ice_topography(domain, i

real (kind=RKIND), dimension(:), pointer :: latCell, lonCell
real (kind=RKIND), dimension(:), pointer :: landIceThkObserved, landIceDraftObserved, &
landIceFracObserved, landIceGroundedFracObserved
landIceFracObserved, &
landIceFloatingFracObserved, &
landIceGroundedFracObserved

real (kind=RKIND), dimension(:), pointer :: landIcePressure, landIceFraction, ssh, &
real (kind=RKIND), dimension(:), pointer :: landIcePressure, landIceFraction, &
landIceFloatingFraction, ssh, &
bottomDepth

integer, pointer :: nCells
integer, dimension(:), pointer :: maxLevelCell, landIceMask
integer, dimension(:), pointer :: maxLevelCell, landIceMask, landIceFloatingMask

integer :: iCell

Expand Down Expand Up @@ -1123,6 +1126,7 @@ subroutine ocn_init_setup_global_ocean_interpolate_land_ice_topography(domain, i
inXPeriod = 2.0_RKIND * pii)

elseif (config_global_ocean_topography_method .eq. "bilinear_interpolation") then

call ocn_init_interpolation_bilinear_horiz(landIceThkLon % array, landIceThkLat % array, &
landIceThkIC % array, nLonLandIceThk, nLatLandIceThk, &
lonCell, latCell, landIceThkObserved, nCells, &
Expand Down Expand Up @@ -1186,28 +1190,42 @@ subroutine ocn_init_setup_global_ocean_interpolate_land_ice_topography(domain, i
call mpas_pool_get_array(landIceInitPool, 'landIceDraftObserved', landIceDraftObserved)
call mpas_pool_get_array(landIceInitPool, 'landIceThkObserved', landIceThkObserved)
call mpas_pool_get_array(landIceInitPool, 'landIceFracObserved', landIceFracObserved)
call mpas_pool_get_array(landIceInitPool, 'landIceFracObserved', landIceFloatingFracObserved)
call mpas_pool_get_array(landIceInitPool, 'landIceGroundedFracObserved', landIceGroundedFracObserved)
call mpas_pool_get_array(statePool, 'ssh', ssh, 1)
call mpas_pool_get_array(forcingPool, 'landIceMask', landIceMask)
call mpas_pool_get_array(forcingPool, 'landIceFloatingMask', landIceFloatingMask)
call mpas_pool_get_array(forcingPool, 'landIcePressure', landIcePressure)
call mpas_pool_get_array(forcingPool, 'landIceFraction', landIceFraction)
call mpas_pool_get_array(forcingPool, 'landIceFloatingFraction', landIceFloatingFraction)

if (landIceFloatingMask(1) == -1) then
call mpas_log_write('landIceFloatingMask contains the default value which likely ' // &
'indicates that this field is missing in the initial condition file (e.g. ' // &
'because it is meant for an older E3SM version).', MPAS_LOG_CRIT)
endif
ssh(:) = 0.0_RKIND
landIceFraction(:) = 0.0_RKIND
modifyLandIcePressureMask(:) = 0
landIcePressure(:) = 0.0_RKIND
do iCell = 1, nCells

if(landIceMask(iCell) == 1) then
if (landIceMask(iCell) == 1) then
landIceFraction(iCell) = landIceFracObserved(iCell)
end if

! this implicitly assumes that when a cell is considered floating, the
! full landIceFracObserved is floating
if (landIceFloatingMask(iCell) == 1) then
landIceFloatingFraction(iCell) = landIceFloatingFracObserved(iCell)
end if

! nothing to do here if the cell is land
if (maxLevelCell(iCell) <= 0) cycle

! we compute the SSH first and find out the land-ice pressure
ssh(iCell) = min(0.0_RKIND, landIceDraftObserved(iCell))
if(ssh(iCell) < 0.0_RKIND) then
if (ssh(iCell) < 0.0_RKIND) then
modifyLandIcePressureMask(iCell) = 1
end if
end do
Expand Down Expand Up @@ -2056,7 +2074,6 @@ subroutine ocn_init_setup_global_ocean_interpolate_tracers(domain, tracerArray,
inXPeriod = 2.0_RKIND * pii)
elseif (config_global_ocean_tracer_method .eq. "bilinear_interpolation") then
call ocn_init_interpolation_bilinear_horiz(tracerLon % array, tracerLat % array, &
tracerIC % array, nLonTracer, nLatTracer, &
lonCell, latCell, interpTracer, nCells, &
Expand Down
Loading

0 comments on commit c292bec

Please sign in to comment.