Skip to content

Commit

Permalink
Merge branch 'cbegeman/ocn/new-floating-land-ice-vars' into next (PR #…
Browse files Browse the repository at this point in the history 
…5464)

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 3, 2023
2 parents e3f4936 + aa144a6 commit 57dada5
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 57dada5

Please sign in to comment.