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.

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

New ocean variables for floating land ice #5464

Merged
merged 16 commits into from
Apr 5, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
2 changes: 1 addition & 1 deletion components/mpas-ocean/bld/build-namelist
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
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
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
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
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
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
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
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