Skip to content

Add the option for inactive top cells #825

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.

Sign up for GitHub

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

Merged
merged 51 commits into from
Mar 25, 2021
Merged

Add the option for inactive top cells #825

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
51 commits
Select commit Hold shift + click to select a range
f9260b2
Add minLevelCell variable
xylar Feb 14, 2021
15a899e
Add minLevelCell and z-level coord to init mode
xylar Feb 14, 2021
b0ea910
Compute minLevelCell and maxLevelCell for z-level
xylar Feb 15, 2021
74157ae
Add fill values above ice shelf in z-level test
xylar Feb 16, 2021
fbefa9e
Fix full partial top cells
xylar Feb 18, 2021
18c917a
Fix top density in landIcePressure computation
xylar Mar 3, 2021
ce27e0d
Unstretch z-level layerThickness to get restingThickness
xylar Mar 6, 2021
478338c
Add minLevelEdge,minLevelVertex variables to the Registry
cbegeman Feb 15, 2021
c19c05b
Add minLevelCell,Edge,Vertex variables to ocn_init_routines
cbegeman Feb 15, 2021
f23e315
Add minLevel* variables to ocn_mesh
cbegeman Feb 16, 2021
6d88e0f
Test k-loop range minLevelCell,maxLevelCell
cbegeman Feb 16, 2021
a9aca5b
Add minLevel* to ocn_compute_Haney_number
xylar Feb 17, 2021
13e3c7c
add minLevelCell loops to ocn_forcing_build_fraction_absorbed_array
xylar Feb 18, 2021
0a7bb0b
MinLevel* changes to ocn_validate_state:
cbegeman Feb 16, 2021
714698e
Set inactive top cells to -1e34, ocn_init_routines_block
cbegeman Feb 17, 2021
61dd613
MinLevel* changes to ocn_time_integration_split
cbegeman Feb 17, 2021
165db9d
Remove comments flagging do minLevel* loops
cbegeman Feb 19, 2021
f1774c1
Add minLevel* to ocn_thick_hadv_tend
xylar Feb 19, 2021
56b1a47
Add minLevel* to ocn_thick_vadv_tend
xylar Feb 19, 2021
811fed1
Add minLevel* to ocn_thick_surface_flux_tend
xylar Feb 19, 2021
95162b1
Add minLevel* to ocn_vel_forcing_surface_stress_tend
xylar Feb 20, 2021
9f185e1
add minLevel* to ocn_vel_hadv_coriolis_tend
xylar Feb 20, 2021
09d8b77
Add maxLevel* to ocn_vel_vadv_tend
xylar Feb 20, 2021
2697c71
add minLevel* to ocn_vel_pressure_grad_tend
xylar Feb 20, 2021
8ddb271
MinLevel* changes to ocn_diagnostic_solve (pick tracer,velocity surfa…
cbegeman Feb 17, 2021
df75462
MinLevel* changes to ocn_relativeVorticity_Circulation
cbegeman Feb 17, 2021
80e16de
MinLevel* changes to ocn_diagnostic_solve_layerThicknessEdge
cbegeman Feb 17, 2021
7022da4
MinLevel* changes to ocn_diagnostic_solve_richardson
cbegeman Feb 17, 2021
bfc9fce
MinLevel* changes to ocn_diagnostic_solve_surfaceLayer
cbegeman Feb 18, 2021
e622456
MinLevel* changes to ocn_diagnostic_solve_vortVel
cbegeman Feb 18, 2021
22dd934
MinLevel* changes to ocn_diagnostic_solve_z_coordinates
cbegeman Feb 18, 2021
9fbb93f
MinLevel* changes to ocn_diagnostic_solve_pressure
cbegeman Feb 18, 2021
acaa3fc
MinLevel* changes to ocn_fuperp
cbegeman Feb 18, 2021
dd5903f
MinLevel* changes to ocn_ALE_thickness
cbegeman Feb 18, 2021
46d9047
MinLevel* changes to ocn_compute_land_ice_flux_input_fields
cbegeman Feb 18, 2021
4712d6c
MinLevel* changes to ocn_vert_transport_velocity_top
cbegeman Feb 18, 2021
4c1d6f8
MinLevel* changes to ocn_diagnostic_solve_vorticity:
cbegeman Feb 17, 2021
5f518fe
Add maxLevel* to ocn_tracer_advection_mono_tend
xylar Feb 20, 2021
2b1f4ca
MinLevel* changes to ocn_vmix_implicit
cbegeman Feb 19, 2021
624f23b
MinLevel* changes to ocn_vel_hmix_del2_tend
cbegeman Feb 19, 2021
78a5aef
Minlevel* changes to ocn_vel_vmix_tend_implicit
cbegeman Feb 20, 2021
c8ce359
MinLevel* changes to ocn_tracer_vmix_tend_implicit:
cbegeman Feb 19, 2021
906c019
Minlevel* changes to ocn_vmix_cvmix
cbegeman Feb 22, 2021
33f1927
Check for layerThicknessVertex=0 in vorticity solve
cbegeman Mar 1, 2021
8eb7296
Do not use dry cells to determine minLevelEdge*
cbegeman Feb 26, 2021
99c24d7
Do not use dry cells to determine minLevelVertex*
cbegeman Feb 26, 2021
14feaaf
Set vertical velocity top to 0 from k=1:minLevelCell-1
cbegeman Mar 10, 2021
d127ca0
Change k-limits of tridiagonal_solve_mult for tracers
cbegeman Mar 12, 2021
5f696ea
Fixup ocn_time_integration_split
cbegeman Mar 15, 2021
d9fcf4b
Fixup land ice drag in ocn_diagnostics
cbegeman Mar 15, 2021
8329b9c
Fixup ocn_time_integration_split_init
cbegeman Mar 19, 2021
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
22 changes: 22 additions & 0 deletions src/core_ocean/Registry.xml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1412,6 +1412,7 @@
<var name="refBottomDepth"/>
<var name="bottomDepth"/>
<var name="bottomDepthObserved"/>
<var name="minLevelCell"/>
<var name="maxLevelCell"/>
<var name="vertCoordMovementWeights"/>
<var name="edgeMask"/>
Expand Down Expand Up @@ -1506,6 +1507,7 @@
<var name="fVertex"/>
<var name="fCell"/>
<var name="bottomDepth"/>
<var name="minLevelCell"/>
<var name="maxLevelCell"/>
<var name="refBottomDepth" packages="forwardMode;analysisMode"/>
<var name="restingThickness" packages="forwardMode;analysisMode"/>
Expand Down Expand Up @@ -1646,6 +1648,7 @@
<var name="layerThickness"/>
<var name="ssh"/>
<var name="maxLevelEdgeTop"/>
<var name="minLevelEdgeBot"/>
<var name="vertCoordMovementWeights"/>
<var name="edgeMask"/>
<var name="vertexMask"/>
Expand Down Expand Up @@ -2163,21 +2166,40 @@
<var name="coeffs_reconstruct" type="real" dimensions="R3 maxEdges nCells" units="unitless"
description="Coefficients to reconstruct velocity vectors at cells centers."
/>
<var name="minLevelCell" type="integer" dimensions="nCells" units="unitless" default_value="1"
description="Index to the first active ocean cell in each column."
/>
<var name="maxLevelCell" type="integer" dimensions="nCells" units="unitless"
description="Index to the last active ocean cell in each column."
/>
<var name="minLevelEdgeTop" type="integer" dimensions="nEdges" units="unitless"
description="Index to the first edge in a column with at least one active ocean cell on either side of it."
packages="forwardMode;analysisMode"
/>
<var name="maxLevelEdgeTop" type="integer" dimensions="nEdges" units="unitless"
description="Index to the last edge in a column with active ocean cells on both sides of it."
packages="forwardMode;analysisMode"
/>
<var name="minLevelEdgeBot" type="integer" dimensions="nEdges" units="unitless"
description="Index to the first edge in a column with active ocean cells on both sides of it."
packages="forwardMode;analysisMode"
/>
<var name="maxLevelEdgeBot" type="integer" dimensions="nEdges" units="unitless"
description="Index to the last edge in a column with at least one active ocean cell on either side of it."
packages="forwardMode;analysisMode"
/>
<var name="minLevelVertexTop" type="integer" dimensions="nVertices" units="unitless"
description="Index to the first vertex in a column with at least one active ocean cell around it."
packages="forwardMode;analysisMode"
/>
<var name="maxLevelVertexTop" type="integer" dimensions="nVertices" units="unitless"
description="Index to the last vertex in a column with all active cells around it."
packages="forwardMode;analysisMode"
/>
<var name="minLevelVertexBot" type="integer" dimensions="nVertices" units="unitless"
description="Index to the first vertex in a column with all active cells around it."
packages="forwardMode;analysisMode"
/>
<var name="maxLevelVertexBot" type="integer" dimensions="nVertices" units="unitless"
description="Index to the last vertex in a column with at least one active ocean cell around it."
packages="forwardMode;analysisMode"
Expand Down
2 changes: 1 addition & 1 deletion src/core_ocean/mode_analysis/mpas_ocn_analysis_mode.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -112,7 +112,7 @@ function ocn_analysis_mode_init(domain, startTimeStamp) result(ierr)!{{{

call ocn_init_routines_vert_coord(domain)

call ocn_init_routines_compute_max_level(domain)
call ocn_init_routines_compute_min_max_level(domain)

timeStep = mpas_get_clock_timestep(domain % clock, ierr=err_tmp)
call mpas_get_timeInterval(timeStep, dt=dt)
Expand Down
2 changes: 1 addition & 1 deletion src/core_ocean/mode_forward/mpas_ocn_forward_mode.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -253,7 +253,7 @@ function ocn_forward_mode_init(domain, startTimeStamp) result(ierr)!{{{

call ocn_init_routines_vert_coord(domain)

call ocn_init_routines_compute_max_level(domain)
call ocn_init_routines_compute_min_max_level(domain)

call mpas_pool_get_config(domain % configs, 'config_time_integrator', config_time_integrator)

Expand Down
72 changes: 37 additions & 35 deletions src/core_ocean/mode_forward/mpas_ocn_time_integration_split.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -387,7 +387,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
sshNew(iCell) = sshCur(iCell)
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
layerThicknessNew(k,iCell) = layerThicknessCur(k,iCell)
end do
end do
Expand All @@ -408,7 +408,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
tracersGroupNew(:,k,iCell) = &
tracersGroupCur(:,k,iCell)
end do
Expand Down Expand Up @@ -506,7 +506,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
! this is h^{hf}_{n+1}
highFreqThicknessNew(k,iCell) = &
highFreqThicknessCur(k,iCell) + dt* &
Expand Down Expand Up @@ -573,7 +573,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
cell2 = cellsOnEdge(2,iEdge)
! could put this after with uTemp(maxleveledgetop+1:nvertlevels)=0
uTemp = 0.0_RKIND
do k = 1, maxLevelEdgeTop(iEdge)
do k = minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)

! normalBaroclinicVelocityNew =
! normalBaroclinicVelocityOld +
Expand All @@ -593,11 +593,11 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
! on land boundaries maxLevelEdgeTop=0, but we want to
! initialize thicknessSum with a nonzero value to avoid
! a NaN.
normalThicknessFluxSum = layerThickEdge(1,iEdge)* &
uTemp(1)
thicknessSum = layerThickEdge(1,iEdge)
normalThicknessFluxSum = layerThickEdge(minLevelEdgeBot(iEdge),iEdge)* &
uTemp(minLevelEdgeBot(iEdge))
thicknessSum = layerThickEdge(minLevelEdgeBot(iEdge),iEdge)

do k = 2, maxLevelEdgeTop(iEdge)
do k = minLevelEdgeBot(iEdge)+1, maxLevelEdgeTop(iEdge)
normalThicknessFluxSum = normalThicknessFluxSum + &
layerThickEdge(k,iEdge)*uTemp(k)
thicknessSum = thicknessSum + &
Expand All @@ -606,7 +606,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
barotropicForcing(iEdge) = splitFact* &
normalThicknessFluxSum/thicknessSum/dt

do k = 1, maxLevelEdgeTop(iEdge)
do k = minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)
! These two steps are together here:
!{\bf u}'_{k,n+1} =
! {\bf u}'_{k,n} - \Delta t {\overline {\bf G}}
Expand Down Expand Up @@ -1395,11 +1395,11 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
! on land boundaries maxLevelEdgeTop=0, but I want to
! initialize thicknessSum with a nonzero value to avoid
! a NaN.
normalThicknessFluxSum = layerThickEdge(1,iEdge)* &
uTemp(1)
thicknessSum = layerThickEdge(1,iEdge)
normalThicknessFluxSum = layerThickEdge(minLevelEdgeBot(iEdge),iEdge)* &
uTemp(minLevelEdgeBot(iEdge))
thicknessSum = layerThickEdge(minLevelEdgeBot(iEdge),iEdge)

do k = 2, maxLevelEdgeTop(iEdge)
do k = minLevelEdgeBot(iEdge)+1, maxLevelEdgeTop(iEdge)
normalThicknessFluxSum = normalThicknessFluxSum + &
layerThickEdge(k,iEdge)*&
uTemp(k)
Expand Down Expand Up @@ -1542,7 +1542,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(i, k, temp_h, temp)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)

! this is h_{n+1}
temp_h = layerThicknessCur(k,iCell) + dt* &
Expand Down Expand Up @@ -1572,7 +1572,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k, temp)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)

! h^{hf}_{n+1} was computed in Stage 1

Expand Down Expand Up @@ -1628,7 +1628,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
! this is h_{n+1}
layerThicknessNew(k,iCell) = &
layerThicknessCur(k,iCell) + &
Expand All @@ -1643,7 +1643,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iEdge = 1, nEdgesAll
do k= 1, maxLevelEdgeTop(iEdge)
do k= minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)
activeTracerHorizontalAdvectionEdgeFlux(:,k,iEdge) = &
activeTracerHorizontalAdvectionEdgeFlux(:,k,iEdge) / &
layerThickEdge(k,iEdge)
Expand All @@ -1653,7 +1653,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{

!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k= 1, maxLevelCell(iCell)
do k= minLevelCell(iCell), maxLevelCell(iCell)
activeTracerHorizontalAdvectionTendency(:,k,iCell) = &
activeTracerHorizontalAdvectionTendency(:,k,iCell) / &
layerThicknessNew(k,iCell)
Expand Down Expand Up @@ -1708,7 +1708,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
tracersGroupNew(:,k,iCell) = &
(tracersGroupCur(:,k,iCell) * &
layerThicknessCur(k,iCell) + dt* &
Expand All @@ -1725,7 +1725,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)
tracersGroupNew(indexSalinity,k,iCell) = &
max(0.001_RKIND, &
tracersGroupNew(indexSalinity,k,iCell))
Expand All @@ -1745,7 +1745,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
do iCell = 1, nCellsAll

! Reset tracer1 to 2 in top n layers
do k = 1, config_reset_debugTracers_top_nLayers
do k = minLevelCell(iCell), minLevelCell(iCell)+config_reset_debugTracers_top_nLayers-1
tracersGroupNew(1,k,iCell) = 2.0_RKIND
end do

Expand All @@ -1757,11 +1757,11 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
.or.lat>- 2.5.and.lat< 2.5 &
.or.lat> 35.0.and.lat< 40.0 &
.or.lat> 55.0.and.lat< 60.0 ) then
do k = 1, config_reset_debugTracers_top_nLayers
do k = minLevelCell(iCell), minLevelCell(iCell)+config_reset_debugTracers_top_nLayers-1
tracersGroupNew(2,k,iCell) = 2.0_RKIND
end do
else
do k = 1, config_reset_debugTracers_top_nLayers
do k = minLevelCell(iCell), minLevelCell(iCell)+config_reset_debugTracers_top_nLayers-1
tracersGroupNew(2,k,iCell) = 1.0_RKIND
end do
end if
Expand All @@ -1775,11 +1775,11 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
.or.lat> 10.0.and.lat< 15.0 &
.or.lat> 30.0.and.lat< 35.0 &
.or.lat> 50.0.and.lat< 55.0 ) then
do k = 1, config_reset_debugTracers_top_nLayers
do k = minLevelCell(iCell), minLevelCell(iCell)+config_reset_debugTracers_top_nLayers-1
tracersGroupNew(3,k,iCell) = 2.0_RKIND
end do
else
do k = 1, config_reset_debugTracers_top_nLayers
do k = minLevelCell(iCell), minLevelCell(iCell)+config_reset_debugTracers_top_nLayers-1
tracersGroupNew(3,k,iCell) = 1.0_RKIND
end do
end if
Expand All @@ -1795,7 +1795,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCellsAll
do k = 1, maxLevelCell(iCell)
do k = minLevelCell(iCell), maxLevelCell(iCell)

! h^{hf}_{n+1} was computed in Stage 1

Expand Down Expand Up @@ -1828,7 +1828,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp parallel
!$omp do schedule(runtime) private(k)
do iEdge = 1, nEdgesAll
do k = 1, maxLevelEdgeTop(iEdge)
do k = minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)
normalVelocityNew(k,iEdge) = &
normalBarotropicVelocityNew(iEdge) + &
2*normalBaroclinicVelocityNew(k,iEdge) - &
Expand Down Expand Up @@ -1956,9 +1956,9 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{
!$omp do schedule(runtime)
do iCell = 1, nCellsAll
surfaceVelocity(indexSurfaceVelocityZonal,iCell) = &
velocityZonal(1,iCell)
velocityZonal(minLevelCell(iCell),iCell)
surfaceVelocity(indexSurfaceVelocityMeridional,iCell) = &
velocityMeridional(1,iCell)
velocityMeridional(minLevelCell(iCell),iCell)

SSHGradient(indexSSHGradientZonal,iCell) = gradSSHZonal(iCell)
SSHGradient(indexSSHGradientMeridional,iCell) = &
Expand Down Expand Up @@ -2152,6 +2152,8 @@ subroutine ocn_time_integration_split_init(domain)!{{{

call mpas_pool_get_array(meshPool, 'refBottomDepth', refBottomDepth)
call mpas_pool_get_array(meshPool, 'cellsOnEdge', cellsOnEdge)
call mpas_pool_get_array(meshPool, 'minLevelCell', minLevelCell)
call mpas_pool_get_array(meshPool, 'minLevelEdgeBot', minLevelEdgeBot)
call mpas_pool_get_array(meshPool, 'maxLevelEdgeTop', maxLevelEdgeTop)

!*** Compute barotropic velocity at first timestep
Expand All @@ -2168,7 +2170,7 @@ subroutine ocn_time_integration_split_init(domain)!{{{

if (config_filter_btr_mode) then
do iCell = 1, nCells
layerThickness(1,iCell) = refBottomDepth(1)
layerThickness(minLevelCell(iCell),iCell) = refBottomDepth(minLevelCell(iCell))
enddo
endif

Expand All @@ -2186,13 +2188,13 @@ subroutine ocn_time_integration_split_init(domain)!{{{
! initialize thicknessSum with a nonzero value to avoid
! a NaN.
layerThicknessEdge1 = 0.5_RKIND* &
(layerThickness(1,cell1) + &
layerThickness(1,cell2) )
(layerThickness(minLevelCell(cell1),cell1) + &
layerThickness(minLevelCell(cell2),cell2) )
normalThicknessFluxSum = layerThicknessEdge1* &
normalVelocity(1,iEdge)
normalVelocity(minLevelEdgeBot(iEdge),iEdge)
layerThicknessSum = layerThicknessEdge1

do k=2, kmax
do k=minLevelEdgeBot(iEdge)+1, kmax
layerThicknessEdge1 = 0.5_RKIND* &
(layerThickness(k,cell1) + &
layerThickness(k,cell2))
Expand All @@ -2209,7 +2211,7 @@ subroutine ocn_time_integration_split_init(domain)!{{{

! normalBaroclinicVelocity = normalVelocity -
! normalBarotropicVelocity
do k = 1, kmax
do k = minLevelEdgeBot(iEdge), kmax
normalBaroclinicVelocity(k,iEdge) = &
normalVelocity(k,iEdge) - &
normalBarotropicVelocity(iEdge)
Expand Down
Loading