Skip to content

Commit

Permalink
Use LrEdge in meke tracer tend
Browse files Browse the repository at this point in the history
  • Loading branch information
cbegeman committed Sep 5, 2022
1 parent c548a2a commit 0445aeb
Show file tree
Hide file tree
Showing 2 changed files with 36 additions and 20 deletions.
15 changes: 15 additions & 0 deletions components/mpas-ocean/src/shared/mpas_ocn_tendency.F
Original file line number Diff line number Diff line change
Expand Up @@ -1371,6 +1371,7 @@ subroutine ocn_tend_meke(tendPool, statePool, &
real (kind=RKIND), dimension(:,:), allocatable :: &
normalThicknessFlux ! Flux of thickness through edge (m^2/s)

character(len=100) :: log_string
! these variables input from shared diagnostics module:
! normalTransportVelocity - transport velocity across edge
! layerThickEdgeFlux - flux-related layer thickness on edge
Expand All @@ -1379,6 +1380,7 @@ subroutine ocn_tend_meke(tendPool, statePool, &
!-----------------------------------------------------------------
! Begin code

call mpas_log_write("ocn_tend_meke")
call mpas_timer_start("ocn_tend_meke")

! set time level to default 1 or override with input arg
Expand All @@ -1395,6 +1397,9 @@ subroutine ocn_tend_meke(tendPool, statePool, &
mekeTend)
call mpas_pool_get_array(statePool, 'layerThickness', &
layerThickness, timeLevel)
iCell = 10
write(log_string,*) 'meke, mekeTend =',meke(1,iCell),mekeTend(1,iCell)
call mpas_log_write(log_string)

allocate(meke3d (1, nVertLevels, nCellsAll+1), &
mekeTend3d(1, nVertLevels, nCellsAll+1))
Expand Down Expand Up @@ -1472,13 +1477,23 @@ subroutine ocn_tend_meke(tendPool, statePool, &
end do
!$omp end do
!$omp end parallel
iCell = 10
write(log_string,*) 'After adv mekeTend =',mekeTend(1,iCell)
call mpas_log_write(log_string)
call mpas_log_write('MEKE source')
call ocn_tracer_meke_source_tend(layerThickness, dt, mekeTend, err)
call mpas_log_write('... completed')
iCell = 10
write(log_string,*) 'After source mekeTend =',mekeTend(1,iCell)
call mpas_log_write(log_string)
call mpas_log_write('MEKE sink')
call ocn_tracer_meke_sink_tend(meke, layerThickness, dt, mekeTend, err)
call mpas_log_write('... completed')
iCell = 10
write(log_string,*) 'After sink meke, mekeTend =',meke(1,iCell),mekeTend(1,iCell)
call mpas_log_write(log_string)
call mpas_timer_stop("ocn_tend_meke")
Expand Down
41 changes: 21 additions & 20 deletions components/mpas-ocean/src/shared/mpas_ocn_tracer_meke_tend.F
Original file line number Diff line number Diff line change
Expand Up @@ -126,7 +126,6 @@ subroutine ocn_tracer_meke_source_tend(layerThickness, dt, tend, err)!{{{

! mekeSourceTopOfEdge comes from diagnostics_variables
! first, interpolate mekeSourceTopOfEdge to mekeSourceEdge
call mpas_log_write('Interpolate meke source to mid-edge')
!$omp parallel
!$omp do schedule(runtime) private(k)
do iEdge = 1, nEdges
Expand All @@ -139,10 +138,8 @@ subroutine ocn_tracer_meke_source_tend(layerThickness, dt, tend, err)!{{{
end do
!$omp end do
!$omp end parallel
call mpas_log_write('... completed')

! second, interpolate mekeSourceEdge to cell centers
call mpas_log_write('Interpolate meke source to cell center')
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCells
Expand All @@ -156,9 +153,7 @@ subroutine ocn_tracer_meke_source_tend(layerThickness, dt, tend, err)!{{{
end do
!$omp end do
!$omp end parallel
call mpas_log_write('... completed')

call mpas_log_write('Compute meke source tend')
!$omp parallel
!$omp do schedule(runtime) private(k)
do iCell = 1, nCells
Expand All @@ -170,7 +165,6 @@ subroutine ocn_tracer_meke_source_tend(layerThickness, dt, tend, err)!{{{
end do
!$omp end do
!$omp end parallel
call mpas_log_write('... completed')

call mpas_timer_stop('MEKE source')

Expand Down Expand Up @@ -227,11 +221,12 @@ subroutine ocn_tracer_meke_sink_tend(meke, layerThickness, dt, tend, err)!{{{
integer :: cell1, cell2, iCell, iEdge, j, k, nCells, nEdges
real(kind=RKIND), parameter :: c_diss = 0.1_RKIND
real(kind=RKIND) :: Lr, Velocity, Length, betaCell
real (kind=RKIND), dimension(:), allocatable :: &
LrEdge
character(len=100) :: log_string

call mpas_log_write('MEKE sink start')
iCell = 10
write(log_string,*) 'meke =',meke(1,iCell)
call mpas_log_write(log_string)
err = 0

if (.not. config_prognostic_meke) return
Expand All @@ -241,28 +236,31 @@ subroutine ocn_tracer_meke_sink_tend(meke, layerThickness, dt, tend, err)!{{{
nCells = nCellsOwned
nEdges = nEdgesOwned

allocate(LrEdge(nEdges))

!$omp parallel
!$omp do schedule(runtime) private(iEdge, j, k, Lr, betaCell, Velocity, Length)
do iCell=1,nCells
Lr = 0.0_RKIND
do j = 1,nEdgesOnCell(iCell)
iEdge = edgesOnCell(j,iCell)
Lr = Lr + LrEdge(iEdge)*edgeAreaFractionOfCell(j, iCell)
enddo
if (sphereRadius == 0.0_RKIND) then
betaCell = 1e-20_RKIND
else
betaCell = 2.0_RKIND*omega*cos(latCell(iCell)) / sphereRadius
endif
do k = minLevelCell(iCell),maxLevelCell(iCell)
Lr = 0.0_RKIND
do j = 1,nEdgesOnCell(iCell)
iEdge = edgesOnCell(j,iCell)
Lr = Lr + LrEdge(iEdge)*edgeAreaFractionOfCell(j, iCell)
enddo
Velocity = meke(k,iCell)**0.5_RKIND
if (sphereRadius == 0.0_RKIND) then
betaCell = 1e-20_RKIND
Length = Lr
else
betaCell = 2.0_RKIND*omega*cos(latCell(iCell)) / sphereRadius
Length = min((Velocity/betaCell)**0.5_RKIND, Lr)
endif
Velocity = meke(k,iCell)**0.5_RKIND
Length = max(1e-10_RKIND, min((Velocity/betaCell)**0.5_RKIND, Lr))
Length = max(1e-10_RKIND, Length)
mekeSink(k,iCell) = c_diss * &
(meke(k,iCell)/layerThickness(k,iCell))**1.5_RKIND / &
Length
if (mekeSink(k,iCell) > 1.0_RKIND) then
if (mekeSink(k,iCell) > 1.0_RKIND .or. (iCell == 10 .and. k==1)) then
write(log_string,*) 'Lr', Lr
call mpas_log_write(log_string)
write(log_string,*) 'omega, sphereRadius', omega, sphereRadius
Expand All @@ -280,6 +278,9 @@ subroutine ocn_tracer_meke_sink_tend(meke, layerThickness, dt, tend, err)!{{{
enddo
!$omp end do
!$omp end parallel
iCell = 10
write(log_string,*) 'mekeTend =',tend(1,iCell)
call mpas_log_write(log_string)
call mpas_timer_stop('MEKE sink')
call mpas_log_write('MEKE sink end')

Expand Down

0 comments on commit 0445aeb

Please sign in to comment.