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

Avoid copying the tracer registry in advect_tracer #85

Merged
merged 2 commits into from
Mar 16, 2022
Merged
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
68 changes: 32 additions & 36 deletions src/tracer/MOM_tracer_advect.F90
Original file line number Diff line number Diff line change
Expand Up @@ -83,7 +83,6 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
optional, intent(out) :: vhr_out !< Remaining accumulated volume or mass fluxes
!! through the meridional faces [H L2 ~> m3 or kg]

type(tracer_type) :: Tr(MAX_FIELDS_) ! The array of registered tracers
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: &
hprev ! cell volume at the end of previous tracer change [H L2 ~> m3 or kg]
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)) :: &
Expand Down Expand Up @@ -127,7 +126,6 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
if (CS%usePPM .and. .not. CS%useHuynh) stencil = 3

ntr = Reg%ntr
do m=1,ntr ; Tr(m) = Reg%Tr(m) ; enddo
Idt = 1.0 / dt

max_iter = 2*INT(CEILING(dt/CS%dt)) + 1
Expand All @@ -138,7 +136,7 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
call create_group_pass(CS%pass_uhr_vhr_t_hprev, uhr, vhr, G%Domain)
call create_group_pass(CS%pass_uhr_vhr_t_hprev, hprev, G%Domain)
do m=1,ntr
call create_group_pass(CS%pass_uhr_vhr_t_hprev, Tr(m)%t, G%Domain)
call create_group_pass(CS%pass_uhr_vhr_t_hprev, Reg%Tr(m)%t, G%Domain)
enddo
call cpu_clock_end(id_clock_pass)

Expand Down Expand Up @@ -188,27 +186,11 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
! initialize diagnostic fluxes and tendencies
!$OMP do
do m=1,ntr
if (associated(Tr(m)%ad_x)) then
do k=1,nz ; do j=jsd,jed ; do i=isd,ied
Tr(m)%ad_x(I,j,k) = 0.0
enddo ; enddo ; enddo
endif
if (associated(Tr(m)%ad_y)) then
do k=1,nz ; do J=jsd,jed ; do i=isd,ied
Tr(m)%ad_y(i,J,k) = 0.0
enddo ; enddo ; enddo
endif
if (associated(Tr(m)%advection_xy)) then
do k=1,nz ; do j=jsd,jed ; do i=isd,ied
Tr(m)%advection_xy(i,j,k) = 0.0
enddo ; enddo ; enddo
endif
if (associated(Tr(m)%ad2d_x)) then
do j=jsd,jed ; do i=isd,ied ; Tr(m)%ad2d_x(I,j) = 0.0 ; enddo ; enddo
endif
if (associated(Tr(m)%ad2d_y)) then
do J=jsd,jed ; do i=isd,ied ; Tr(m)%ad2d_y(i,J) = 0.0 ; enddo ; enddo
endif
if (associated(Reg%Tr(m)%ad_x)) Reg%Tr(m)%ad_x(:,:,:) = 0.0
if (associated(Reg%Tr(m)%ad_y)) Reg%Tr(m)%ad_y(:,:,:) = 0.0
if (associated(Reg%Tr(m)%advection_xy)) Reg%Tr(m)%advection_xy(:,:,:) = 0.0
if (associated(Reg%Tr(m)%ad2d_x)) Reg%Tr(m)%ad2d_x(:,:) = 0.0
if (associated(Reg%Tr(m)%ad2d_y)) Reg%Tr(m)%ad2d_y(:,:) = 0.0
enddo
!$OMP end parallel

Expand Down Expand Up @@ -265,14 +247,14 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
!$OMP do ordered
do k=1,nz ; if (domore_k(k) > 0) then
! First, advect zonally.
call advect_x(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &
call advect_x(Reg%Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &
isv, iev, jsv-stencil, jev+stencil, k, G, GV, US, CS%usePPM, CS%useHuynh)
endif ; enddo

!$OMP do ordered
do k=1,nz ; if (domore_k(k) > 0) then
! Next, advect meridionally.
call advect_y(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &
call advect_y(Reg%Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &
isv, iev, jsv, jev, k, G, GV, US, CS%usePPM, CS%useHuynh)

! Update domore_k(k) for the next iteration
Expand All @@ -287,14 +269,14 @@ subroutine advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, x_first
!$OMP do ordered
do k=1,nz ; if (domore_k(k) > 0) then
! First, advect meridionally.
call advect_y(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &
call advect_y(Reg%Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &
isv-stencil, iev+stencil, jsv, jev, k, G, GV, US, CS%usePPM, CS%useHuynh)
endif ; enddo

!$OMP do ordered
do k=1,nz ; if (domore_k(k) > 0) then
! Next, advect zonally.
call advect_x(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &
call advect_x(Reg%Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &
isv, iev, jsv, jev, k, G, GV, US, CS%usePPM, CS%useHuynh)

! Update domore_k(k) for the next iteration
Expand Down Expand Up @@ -390,13 +372,20 @@ subroutine advect_x(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &
real :: tiny_h ! The smallest numerically invertible thickness [H ~> m or kg m-2].
real :: h_neglect ! A thickness that is so small it is usually lost
! in roundoff and can be neglected [H ~> m or kg m-2].
real :: aR, aL ! Reconstructed tracer concentrations at the right and left edges [conc]
real :: dMx ! Difference between the maximum of the surrounding cell concentrations and
! the value in the cell whose reconstruction is being found [conc]
real :: dMn ! Difference between the tracer concentration in the cell whose reconstruction
! is being found and the minimum of the surrounding values [conc]
real :: Tp, Tc, Tm ! Tracer concentrations around the upstream cell [conc]
real :: dA ! Difference between the reconstruction tracer edge values [conc]
real :: mA ! Average of the reconstruction tracer edge values [conc]
real :: a6 ! Curvature of the reconstruction tracer values [conc]
logical :: do_i(SZIB_(G),SZJ_(G)) ! If true, work on given points.
logical :: do_any_i
logical :: usePLMslope
integer :: i, j, m, n, i_up, stencil
real :: aR, aL, dMx, dMn, Tp, Tc, Tm, dA, mA, a6
real :: fac1,u_L_in,u_L_out ! terms used for time-stepping OBC reservoirs
type(OBC_segment_type), pointer :: segment=>NULL()
logical :: usePLMslope
logical, dimension(SZJ_(G),SZK_(GV)) :: domore_u_initial

! keep a local copy of the initial values of domore_u, which is to be used when computing ad2d_x
Expand Down Expand Up @@ -547,7 +536,7 @@ subroutine advect_x(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, &

dA = aR - aL ; mA = 0.5*( aR + aL )
if (G%mask2dCu(I_up,j)*G%mask2dCu(I_up-1,j)*(Tp-Tc)*(Tc-Tm) <= 0.) then
aL = Tc ; aR = Tc ! PCM for local extremum and bounadry cells
aL = Tc ; aR = Tc ! PCM for local extrema and boundary cells
elseif ( dA*(Tc-mA) > (dA*dA)/6. ) then
aL = 3.*Tc - 2.*aR
elseif ( dA*(Tc-mA) < - (dA*dA)/6. ) then
Expand Down Expand Up @@ -754,14 +743,21 @@ subroutine advect_y(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &
real :: tiny_h ! The smallest numerically invertible thickness [H ~> m or kg m-2].
real :: h_neglect ! A thickness that is so small it is usually lost
! in roundoff and can be neglected [H ~> m or kg m-2].
real :: aR, aL ! Reconstructed tracer concentrations at the right and left edges [conc]
real :: dMx ! Difference between the maximum of the surrounding cell concentrations and
! the value in the cell whose reconstruction is being found [conc]
real :: dMn ! Difference between the tracer average in the cell whose reconstruction
! is being found and the minimum of the surrounding values [conc]
real :: Tp, Tc, Tm ! Tracer concentrations around the upstream cell [conc]
real :: dA ! Difference between the reconstruction tracer edge values [conc]
real :: mA ! Average of the reconstruction tracer edge values [conc]
real :: a6 ! Curvature of the reconstruction tracer values [conc]
logical :: do_j_tr(SZJ_(G)) ! If true, calculate the tracer profiles.
logical :: do_i(SZIB_(G), SZJ_(G)) ! If true, work on given points.
logical :: do_any_i
logical :: usePLMslope
integer :: i, j, j2, m, n, j_up, stencil
real :: aR, aL, dMx, dMn, Tp, Tc, Tm, dA, mA, a6
real :: fac1,v_L_in,v_L_out ! terms used for time-stepping OBC reservoirs
type(OBC_segment_type), pointer :: segment=>NULL()
logical :: usePLMslope

usePLMslope = .not. (usePPM .and. useHuynh)
! stencil for calculating slope values
Expand Down Expand Up @@ -919,7 +915,7 @@ subroutine advect_y(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, &

dA = aR - aL ; mA = 0.5*( aR + aL )
if (G%mask2dCv(i,J_up)*G%mask2dCv(i,J_up-1)*(Tp-Tc)*(Tc-Tm) <= 0.) then
aL = Tc ; aR = Tc ! PCM for local extremum and bounadry cells
aL = Tc ; aR = Tc ! PCM for local extrema and boundary cells
elseif ( dA*(Tc-mA) > (dA*dA)/6. ) then
aL = 3.*Tc - 2.*aR
elseif ( dA*(Tc-mA) < - (dA*dA)/6. ) then
Expand Down