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.

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

Jump to bottom

Generalize weighted z-star coordinate #100

Draft
wants to merge 6 commits into
base: master
Choose a base branch
from
+31 −19
Draft

Generalize weighted z-star coordinate #100

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
2289074
Re-apply PR#5254 Modify ocean z-star ALE coordinate for inactive top …
cbegeman Oct 27, 2022
5bdb41a
Change sign of restingThickness ALE contribution
cbegeman May 28, 2024
43fe745
Set minimum vertCoordMovementWeight
cbegeman May 28, 2024
5fbc260
fixup openacc
cbegeman Jun 4, 2024
33fe309
Use pressureAdjustedSSH with ssh_gradient pressure tendency
cbegeman Jun 12, 2024
8176242
Make init mode vertical coord consistent with allowing differences be…
cbegeman Aug 8, 2024
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
7 changes: 5 additions & 2 deletions components/mpas-ocean/src/mode_init/mpas_ocn_init_vertical_grids.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -591,15 +591,18 @@ subroutine ocn_compute_z_star_layerThickness(layerThickness, restingThickness, s
integer, intent(out) :: iErr

integer :: k
real (kind=RKIND) :: layerStretch
real (kind=RKIND) :: layerStretch, restingThicknessSum

iErr = 0

layerThickness(:) = 0.0_RKIND ! = MPAS_REAL_FILLVAL

if (maxLevelCell<=0) return

layerStretch = (ssh + bottomDepth)/bottomDepth
do k = minLevelCell, maxLevelCell
restingThicknessSum = restingThicknessSum + restingThickness(k)
end do
layerStretch = (ssh + bottomDepth)/restingThicknessSum
do k = minLevelCell, maxLevelCell
layerThickness(k) = layerStretch * restingThickness(k)
end do
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_mesh.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1568,7 +1568,7 @@ subroutine ocn_meshVertCoordInit(domain)!{{{
select case (trim(config_vert_coord_movement))
case ('fixed')

vertCoordMovementWeights = 0.0_RKIND
vertCoordMovementWeights = 1e-14_RKIND
vertCoordMovementWeights(1) = 1.0_RKIND

case ('uniform_stretching')
Expand All @@ -1587,7 +1587,7 @@ subroutine ocn_meshVertCoordInit(domain)!{{{
(config_vert_taper_weight_depth_2 - &
config_vert_taper_weight_depth_1)
! limit range to (0,1)
vertCoordMovementWeights(k) = max( 0.0_RKIND, &
vertCoordMovementWeights(k) = max( 1e-14_RKIND, &
min( 1.0_RKIND, vertCoordMovementWeights(k) ) )
end do

Expand Down
2 changes: 1 addition & 1 deletion components/mpas-ocean/src/shared/mpas_ocn_tendency.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -419,7 +419,7 @@ subroutine ocn_tend_vel(domain, tendPool, statePool, forcingPool, &
vertAleTransportTop, tendVel, err)

! Add pressure gradient
call ocn_vel_pressure_grad_tend(ssh, pressure, surfacePressure, &
call ocn_vel_pressure_grad_tend(ssh, gradSSH, pressure, surfacePressure, &
montgomeryPotential, zMid, &
density, potentialDensity, &
indxTemp, indxSalt, activeTracers, &
Expand Down
32 changes: 20 additions & 12 deletions components/mpas-ocean/src/shared/mpas_ocn_thick_ale.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -125,10 +125,12 @@ subroutine ocn_ALE_thickness(verticalMeshPool, SSH, ALE_thickness, &
nCells ! number of cells

real (kind=RKIND) :: &
weightSum, &! sum of weights in vertical
thicknessSum, &! total thickness
remainder, &! track remainder in mix/max alteration
newThickness ! temp used during min/max adjustment
weightSum, &! sum of weights in vertical
restingThicknessSum, &! total resting thickness
restingThicknessDiff, &! difference from total resting thickness
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could we call restingThicknessDiff something else (like restingSSH)? It isn't clear what it's the difference with respect to.

weightedThicknessSum, &! total resting thickness, weighted
remainder, &! track remainder in mix/max alteration
newThickness ! temp used during min/max adjustment

real (kind=RKIND), dimension(:), allocatable :: &
prelim_ALE_thickness, & ! ALE thickness at new time
Expand Down Expand Up @@ -170,33 +172,39 @@ subroutine ocn_ALE_thickness(verticalMeshPool, SSH, ALE_thickness, &
#ifdef MPAS_OPENACC
!xacc parallel loop &
!xacc present(ALE_thickness, SSH, restingThickness, &
!xacc minLevelCell, maxLevelCell, &
!xacc minLevelCell, maxLevelCell, bottomDepth, &
!xacc vertCoordMovementWeights) &
!xacc private(k, kMin, kMax, thicknessSum)
!xacc private(k, kMin, kMax, restingThicknessDiff, &
!xacc restingThicknessSum, weightedThicknessSum)
#else
!$omp parallel
!$omp do schedule(runtime) &
!$omp private(k, kMin, kMax, thicknessSum)
!$omp private(k, kMin, kMax, restingThicknessDiff, &
!$omp restingThicknessSum, weightedThicknessSum)
#endif
do iCell = 1, nCells
kMax = maxLevelCell(iCell)
kMin = minLevelCell(iCell)

thicknessSum = 1e-14_RKIND
restingThicknessSum = 0.0_RKIND
weightedThicknessSum = 0.0_RKIND
do k = kMin, kMax
thicknessSum = thicknessSum &
restingThicknessSum = restingThicknessSum &
+ restingThickness(k,iCell)
weightedThicknessSum = weightedThicknessSum &
+ vertCoordMovementWeights(k) &
* restingThickness(k,iCell)
end do

restingThicknessDiff = restingThicknessSum - bottomDepth(iCell)
! Note that restingThickness is nonzero, and remaining
! terms are perturbations about zero.
! This is equation 4 and 6 in Petersen et al 2015,
! but with eqn 6
do k = kMin, kMax
ALE_thickness(k,iCell) = restingThickness(k,iCell) &
+ (SSH(iCell)*vertCoordMovementWeights(k)* &
restingThickness(k,iCell) )/thicknessSum
ALE_thickness(k, iCell) = restingThickness(k, iCell) &
+ ( (SSH(iCell) - restingThicknessDiff) * vertCoordMovementWeights(k) * &
restingThickness(k, iCell) ) / weightedThicknessSum
end do
enddo
#ifndef MPAS_OPENACC
Expand Down
5 changes: 3 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_vel_pressure_grad.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -94,7 +94,7 @@ module ocn_vel_pressure_grad
!
!-----------------------------------------------------------------------

subroutine ocn_vel_pressure_grad_tend(ssh, pressure, surfacePressure, &
subroutine ocn_vel_pressure_grad_tend(ssh, gradSSH, pressure, surfacePressure, &
montgomeryPotential, zMid, &
density, potentialDensity, &
indxT, indxS, tracers, &
Expand All @@ -113,6 +113,7 @@ subroutine ocn_vel_pressure_grad_tend(ssh, pressure, surfacePressure, &

real (kind=RKIND), dimension(:), intent(in) :: &
ssh, &!< [in] sea surface height
gradSSH, &!< [in] sea surface height gradient on edges
surfacePressure !< [in] surface pressure

real (kind=RKIND), dimension(:,:), intent(in) :: &
Expand Down Expand Up @@ -241,7 +242,7 @@ subroutine ocn_vel_pressure_grad_tend(ssh, pressure, surfacePressure, &
do k=kMin,kMax
tend(k,iEdge) = tend(k,iEdge) - &
edgeMask(k,iEdge)*invdcEdge* &
( gravity*(ssh(cell2) - ssh(cell1)) &
( gravity*gradSSH(iEdge) &
+ density0Inv*( surfacePressure(cell2) &
- surfacePressure(cell1)) )
end do
Expand Down
Loading