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

Enables and modifies the submesoscale eddy parameterization #5172

Merged
Original file line number Diff line number Diff line change
Expand Up @@ -147,9 +147,9 @@
<config_Redi_horizontal_ramp_max ocn_grid="WCAtl12to45E2r4">40e3</config_Redi_horizontal_ramp_max>

<!-- submesoscale_eddy_parameterization -->
<config_submesoscale_enable>.false.</config_submesoscale_enable>
<config_submesoscale_enable>.true.</config_submesoscale_enable>
<config_submesoscale_tau>172800</config_submesoscale_tau>
<config_submesoscale_Ce>0.06</config_submesoscale_Ce>
<config_submesoscale_Ce>0.08</config_submesoscale_Ce>
<config_submesoscale_Lfmin>1000.0</config_submesoscale_Lfmin>
<config_submesoscale_ds_max>100000.0</config_submesoscale_ds_max>

Expand Down
2 changes: 2 additions & 0 deletions components/mpas-ocean/cime_config/buildnml
Original file line number Diff line number Diff line change
Expand Up @@ -1118,6 +1118,8 @@ def buildnml(case, caseroot, compname):
lines.append(' <var name="vertVelocityTop"/>')
lines.append(' <var name="zMid"/>')
lines.append(' <var name="atmosphericPressure"/>')
lines.append(' <var name="normalMLEvelocity"/>')
lines.append(' <var name="vertMLEBolusVelocityTop"/>')
if not ocn_grid.startswith("oRRS1"):
lines.append(' <var name="normalGMBolusVelocity"/>')
lines.append(' <var name="vertGMBolusVelocityTop"/>')
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/Registry.xml
Original file line number Diff line number Diff line change
Expand Up @@ -2875,8 +2875,8 @@
description="count of times redi limiter is invoked on a timestep"
packages="gm"
/>
<var name="gradBuoyEddy" type="real" dimensions="nVertLevelsP1 nEdges Time" units="s^-2"
description="horizontal gradient of buoyancy at cell edge and interfaces in vertical. this is used for the GM parameterization and the submesoscale eddy parameterization"
<var name="gradDensityEddy" type="real" dimensions="nVertLevelsP1 nEdges Time" units="kg m^-4"
description="horizontal gradient of density at cell edge and interfaces in vertical. this is used for the GM parameterization and the submesoscale eddy parameterization"
packages="gm;submeso"
/>
<var name="mleVelocityZonal" type="real" dimensions="nVertLevels nCells Time" units="m s^-1"
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -134,6 +134,8 @@
<var name="atmosphericPressure"/>
<var name="normalGMBolusVelocity"/>
<var name="vertGMBolusVelocityTop"/>
<var name="normalMLEvelocity"/>
<var name="vertMLEBolusVelocityTop"/>
<var name="GMBolusVelocityZonal"/>
<var name="GMBolusVelocityMeridional"/>
<var name="cGMphaseSpeed"/>
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -122,7 +122,7 @@ module ocn_diagnostics_variables
real(kind=RKIND), dimension(:,:), pointer :: RediKappaScaling
real(kind=RKIND), dimension(:,:), pointer :: RediKappaSfcTaper
real(kind=RKIND), dimension(:,:), pointer :: k33
real(kind=RKIND), dimension(:,:), pointer :: gradBuoyEddy
real(kind=RKIND), dimension(:,:), pointer :: gradDensityEddy
real(kind=RKIND), dimension(:), pointer :: gmBolusKappa
real(kind=RKIND), dimension(:), pointer :: cGMphaseSpeed
real(kind=RKIND), dimension(:,:,:), pointer :: slopeTriadUp, slopeTriadDown
Expand Down Expand Up @@ -271,8 +271,8 @@ subroutine ocn_diagnostics_variables_init(domain, jenkinsOn, hollandJenkinsOn, e
end if

if (config_use_GM .or. config_submesoscale_enable) then
call mpas_pool_get_array(diagnosticsPool, 'gradBuoyEddy', &
gradBuoyEddy)
call mpas_pool_get_array(diagnosticsPool, 'gradDensityEddy', &
gradDensityEddy)
end if

if ( config_compute_active_tracer_budgets ) then
Expand Down Expand Up @@ -781,7 +781,7 @@ subroutine ocn_diagnostics_variables_init(domain, jenkinsOn, hollandJenkinsOn, e
!$acc )
end if
if (config_use_GM.or.config_submesoscale_enable) then
!$acc enter data create(gradBuoyEddy)
!$acc enter data create(gradDensityEddy)
end if
if ( config_compute_active_tracer_budgets ) then
!$acc enter data create( &
Expand Down Expand Up @@ -1029,7 +1029,7 @@ subroutine ocn_diagnostics_variables_destroy(err) !{{{
!$acc )
end if
if (config_use_GM.or.config_submesoscale_enable) then
!$acc exit data delete(gradBuoyEddy)
!$acc exit data delete(gradDensityEddy)
end if
if ( config_compute_active_tracer_budgets ) then
!$acc exit data delete( &
Expand Down Expand Up @@ -1230,7 +1230,7 @@ subroutine ocn_diagnostics_variables_destroy(err) !{{{
RossbyRadius)
end if
if (config_use_GM.or.config_submesoscale_enable) then
nullify(gradBuoyEddy)
nullify(gradDensityEddy)
end if
if ( config_compute_active_tracer_budgets ) then
nullify(activeTracerHorMixTendency, &
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -166,7 +166,7 @@ subroutine ocn_eddy_compute_buoyancy_gradient()
do iEdge = 1, nEdges
if (maxLevelEdgeTop(iEdge) .GE. minLevelEdgeBot(iEdge)) then
do k = minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)+1
gradBuoyEddy(k,iEdge) = (gradDensityTopOfEdge(k,iEdge) - dDensityDzTopOfEdge(k,iEdge) &
gradDensityEddy(k,iEdge) = (gradDensityTopOfEdge(k,iEdge) - dDensityDzTopOfEdge(k,iEdge) &
* gradZMidTopOfEdge(k,iEdge))
end do
end if
Expand Down
6 changes: 3 additions & 3 deletions components/mpas-ocean/src/shared/mpas_ocn_gm.F
Original file line number Diff line number Diff line change
Expand Up @@ -716,7 +716,7 @@ subroutine ocn_GM_compute_Bolus_velocity(statePool, &
tridiagC(k - 1) = 2.0_RKIND*cGMphaseSpeed(iEdge)**2/layerThickEdgeMean(k, iEdge) &
/(layerThickEdgeMean(k - 1, iEdge) + layerThickEdgeMean(k, iEdge))
rightHandSide(k - 1) = gmBolusKappa(iEdge)*gmKappaScaling(k,iEdge)*gravity/rho_sw* &
gradBuoyEddy(k,iEdge)
gradDensityEddy(k,iEdge)

! Second to next to the last rows
do k = minLevelEdgeBot(iEdge) + 2, maxLevelEdgeTop(iEdge) - 1
Expand All @@ -729,7 +729,7 @@ subroutine ocn_GM_compute_Bolus_velocity(statePool, &
tridiagC(k - 1) = 2.0_RKIND*cGMphaseSpeed(iEdge)**2/layerThickEdgeMean(k, iEdge) &
/(layerThickEdgeMean(k - 1, iEdge) + layerThickEdgeMean(k, iEdge))
rightHandSide(k - 1) = gmBolusKappa(iEdge)*gmKappaScaling(k,iEdge)*gravity/rho_sw* &
gradBuoyEddy(k,iEdge)
gradDensityEddy(k,iEdge)
end do

! Last row
Expand All @@ -741,7 +741,7 @@ subroutine ocn_GM_compute_Bolus_velocity(statePool, &
tridiagB(k - 1) = -2.0_RKIND*cGMphaseSpeed(iEdge)**2/(layerThickEdgeMean(k - 1, iEdge) &
*layerThickEdgeMean(k, iEdge)) - BruntVaisalaFreqTopEdge
rightHandSide(k - 1) = gmBolusKappa(iEdge)*gmKappaScaling(k,iEdge)*gravity/rho_sw* &
gradBuoyEddy(k,iEdge)
gradDensityEddy(k,iEdge)
! Total number of rows
N = maxLevelEdgeTop(iEdge) - minLevelEdgeBot(iEdge)

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -96,7 +96,7 @@ subroutine ocn_submesoscale_compute_velocity()
!This is thickness weighted to allow for non uniform grid spacing,
!and centered on layer interfaces. The mixed layer goes from the
!ocean top layer to mid-depth of layer(indMLDedge)
!NOTE, gradBuoyEddy and BFV are defined at the top cell interface
!NOTE, gradDensityEddy and BFV are defined at the top cell interface
!and not valid for minLevelCell, so properties for the top layer
!rely on the values at (minLevelCell+1).
!$omp parallel
Expand All @@ -118,7 +118,7 @@ subroutine ocn_submesoscale_compute_velocity()
hML = 0.5_RKIND*layerThickEdgeMean(minLevelEdgeBot(iEdge),iEdge)
bvfML = bvfML + hML*bvfAv
if (minLevelEdgeTop(iEdge) .ge. 1) then
gradBuoyML = hML*gradBuoyEddy(minLevelEdgeTop(iEdge)+1,iEdge)
gradBuoyML = hML*gradDensityEddy(minLevelEdgeTop(iEdge)+1,iEdge)
bvfML = hML*bvfAv
else
cycle
Expand All @@ -131,10 +131,11 @@ subroutine ocn_submesoscale_compute_velocity()
max(BruntVaisalaFreqTop(k,cell2),1.0E-20_RKIND)))
bvfML = bvfML + hAv*bvfAv
hML = hML + hAv
gradBuoyML = gradBuoyML + hAv*gradBuoyEddy(k,iEdge)
gradBuoyML = gradBuoyML + hAv*gradDensityEddy(k,iEdge)
end do
bvfML = bvfML / (1.0E-20_RKIND + hML)
gradBuoyML = gradBuoyML / (1.0E-20_RKIND + hML)
gradBuoyML = gravity*gradBuoyML/rho_sw !convert from density unit to buoyancy

!compute depths and shape function

Expand Down