Skip to content

Commit

Permalink
refactor(prt): introduce cell method type/module (#2203)
Browse files Browse the repository at this point in the history
Introduce CellMethodType and move the check() routine there from the base MethodType. These checks occur immediately before tracking begins across a cell, so this is the proper scope — they didn't belong in the base method type.
  • Loading branch information
wpbonelli authored Feb 6, 2025
1 parent 2d01e20 commit 0bc5c4a
Show file tree
Hide file tree
Showing 9 changed files with 189 additions and 163 deletions.
1 change: 1 addition & 0 deletions make/makefile
Original file line number Diff line number Diff line change
Expand Up @@ -287,6 +287,7 @@ $(OBJDIR)/SfrCrossSectionUtils.o \
$(OBJDIR)/TernarySolveTrack.o \
$(OBJDIR)/SubcellTri.o \
$(OBJDIR)/Method.o \
$(OBJDIR)/MethodCell.o \
$(OBJDIR)/SubcellRect.o \
$(OBJDIR)/VirtualBase.o \
$(OBJDIR)/STLVecInt.o \
Expand Down
1 change: 1 addition & 0 deletions msvs/mf6core.vfproj
Original file line number Diff line number Diff line change
Expand Up @@ -396,6 +396,7 @@
<File RelativePath="..\src\Solution\ParticleTracker\CellRectQuad.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\CellUtil.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\Method.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\MethodCell.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\MethodCellPassToBot.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\MethodCellPollock.f90"/>
<File RelativePath="..\src\Solution\ParticleTracker\MethodCellPollockQuad.f90"/>
Expand Down
157 changes: 0 additions & 157 deletions src/Solution/ParticleTracker/Method.f90
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,6 @@ module MethodModule
procedure :: save
procedure :: track
procedure :: try_pass
procedure :: check
end type MethodType

abstract interface
Expand Down Expand Up @@ -185,160 +184,4 @@ subroutine save(this, particle, reason)
call this%trackctl%save(particle, kper=per, kstp=stp, reason=reason)
end subroutine save

!> @brief Check reporting/terminating conditions before tracking.
!!
!! Check a number of conditions determining whether to continue
!! tracking the particle or terminate it, as well as whether to
!! record any output data as per selected reporting conditions.
!<
subroutine check(this, particle, cell_defn, tmax)
! modules
use TdisModule, only: endofsimulation, totimc, totim
use ParticleModule, only: TERM_WEAKSINK, TERM_NO_EXITS, &
TERM_STOPZONE, TERM_INACTIVE
! dummy
class(MethodType), intent(inout) :: this
type(ParticleType), pointer, intent(inout) :: particle
type(CellDefnType), pointer, intent(inout) :: cell_defn
real(DP), intent(in) :: tmax
! local
logical(LGP) :: dry_cell, dry_particle, no_exit_face, stop_zone, weak_sink
integer(I4B) :: i
real(DP) :: t, ttrackmax

dry_cell = this%fmi%ibdgwfsat0(cell_defn%icell) == 0
dry_particle = particle%z > cell_defn%top
no_exit_face = cell_defn%inoexitface > 0
stop_zone = cell_defn%izone > 0 .and. particle%istopzone == cell_defn%izone
weak_sink = cell_defn%iweaksink > 0

particle%izone = cell_defn%izone
if (stop_zone) then
particle%advancing = .false.
particle%istatus = TERM_STOPZONE
call this%save(particle, reason=3)
return
end if

if (no_exit_face .and. .not. dry_cell) then
particle%advancing = .false.
particle%istatus = TERM_NO_EXITS
call this%save(particle, reason=3)
return
end if

if (weak_sink) then
if (particle%istopweaksink > 0) then
particle%advancing = .false.
particle%istatus = TERM_WEAKSINK
call this%save(particle, reason=3)
return
else
call this%save(particle, reason=4)
end if
end if

if (dry_cell) then
if (particle%idrymeth == 0) then
! drop to cell bottom. handled by pass
! to bottom method, nothing to do here
no_exit_face = .false.
else if (particle%idrymeth == 1) then
! stop
particle%advancing = .false.
particle%istatus = TERM_INACTIVE
call this%save(particle, reason=3)
return
else if (particle%idrymeth == 2) then
! stay
particle%advancing = .false.
no_exit_face = .false.

! we might report tracking times
! out of order here, but we want
! the particle termination event
! (if this is the last time step)
! to have the maximum tracking t,
! so we need to keep tabs on it.
ttrackmax = totim

! update tracking time to time
! step end time and save record
particle%ttrack = totim
call this%save(particle, reason=2)

! record user tracking times
call this%tracktimes%advance()
if (this%tracktimes%any()) then
do i = this%tracktimes%selection(1), this%tracktimes%selection(2)
t = this%tracktimes%times(i)
if (t < totimc) cycle
if (t >= tmax) exit
particle%ttrack = t
call this%save(particle, reason=5)
if (t > ttrackmax) ttrackmax = t
end do
end if

! terminate if last period/step
if (endofsimulation) then
particle%istatus = TERM_NO_EXITS
particle%ttrack = ttrackmax
call this%save(particle, reason=3)
return
end if
end if
else if (dry_particle .and. this%name /= "passtobottom") then
if (particle%idrymeth == 0) then
! drop to water table
particle%z = cell_defn%top
call this%save(particle, reason=1)
else if (particle%idrymeth == 1) then
! stop
particle%advancing = .false.
particle%istatus = TERM_INACTIVE
call this%save(particle, reason=3)
return
else if (particle%idrymeth == 2) then
! stay
particle%advancing = .false.
no_exit_face = .false.

! we might report tracking times
! out of order here, but we want
! the particle termination event
! (if this is the last time step)
! to have the maximum tracking t,
! so we need to keep tabs on it.
ttrackmax = totim

! update tracking time to time
! step end time and save record
particle%ttrack = totim
call this%save(particle, reason=2)

! record user tracking times
call this%tracktimes%advance()
if (this%tracktimes%any()) then
do i = this%tracktimes%selection(1), this%tracktimes%selection(2)
t = this%tracktimes%times(i)
if (t < totimc) cycle
if (t >= tmax) exit
particle%ttrack = t
call this%save(particle, reason=5)
if (t > ttrackmax) ttrackmax = t
end do
end if
end if
end if

if (no_exit_face) then
particle%advancing = .false.
particle%istatus = TERM_NO_EXITS
call this%save(particle, reason=3)
return
end if

end subroutine check

end module MethodModule
177 changes: 177 additions & 0 deletions src/Solution/ParticleTracker/MethodCell.f90
Original file line number Diff line number Diff line change
@@ -0,0 +1,177 @@
module MethodCellModule

use KindModule, only: DP, I4B, LGP
use ConstantsModule, only: DONE, DZERO
use MethodModule, only: MethodType
use ParticleModule, only: ParticleType
use CellDefnModule, only: CellDefnType
implicit none

private
public :: MethodCellType

type, abstract, extends(MethodType) :: MethodCellType
contains
procedure, public :: check
end type MethodCellType

contains

!> @brief Check reporting/terminating conditions before tracking
!! the particle across the cell.
!!
!! Check a number of conditions determining whether to continue
!! tracking the particle or terminate it, as well as whether to
!! record any output data as per selected reporting conditions.
!<
subroutine check(this, particle, cell_defn, tmax)
! modules
use TdisModule, only: endofsimulation, totimc, totim
use ParticleModule, only: TERM_WEAKSINK, TERM_NO_EXITS, &
TERM_STOPZONE, TERM_INACTIVE
! dummy
class(MethodCellType), intent(inout) :: this
type(ParticleType), pointer, intent(inout) :: particle
type(CellDefnType), pointer, intent(inout) :: cell_defn
real(DP), intent(in) :: tmax
! local
logical(LGP) :: dry_cell, dry_particle, no_exit_face, stop_zone, weak_sink
integer(I4B) :: i
real(DP) :: t, ttrackmax

dry_cell = this%fmi%ibdgwfsat0(cell_defn%icell) == 0
dry_particle = particle%z > cell_defn%top
no_exit_face = cell_defn%inoexitface > 0
stop_zone = cell_defn%izone > 0 .and. particle%istopzone == cell_defn%izone
weak_sink = cell_defn%iweaksink > 0

particle%izone = cell_defn%izone
if (stop_zone) then
particle%advancing = .false.
particle%istatus = TERM_STOPZONE
call this%save(particle, reason=3)
return
end if

if (no_exit_face .and. .not. dry_cell) then
particle%advancing = .false.
particle%istatus = TERM_NO_EXITS
call this%save(particle, reason=3)
return
end if

if (weak_sink) then
if (particle%istopweaksink > 0) then
particle%advancing = .false.
particle%istatus = TERM_WEAKSINK
call this%save(particle, reason=3)
return
else
call this%save(particle, reason=4)
end if
end if

if (dry_cell) then
if (particle%idrymeth == 0) then
! drop to cell bottom. handled by pass
! to bottom method, nothing to do here
no_exit_face = .false.
else if (particle%idrymeth == 1) then
! stop
particle%advancing = .false.
particle%istatus = TERM_INACTIVE
call this%save(particle, reason=3)
return
else if (particle%idrymeth == 2) then
! stay
particle%advancing = .false.
no_exit_face = .false.

! we might report tracking times
! out of order here, but we want
! the particle termination event
! (if this is the last time step)
! to have the maximum tracking t,
! so we need to keep tabs on it.
ttrackmax = totim

! update tracking time to time
! step end time and save record
particle%ttrack = totim
call this%save(particle, reason=2)

! record user tracking times
call this%tracktimes%advance()
if (this%tracktimes%any()) then
do i = this%tracktimes%selection(1), this%tracktimes%selection(2)
t = this%tracktimes%times(i)
if (t < totimc) cycle
if (t >= tmax) exit
particle%ttrack = t
call this%save(particle, reason=5)
if (t > ttrackmax) ttrackmax = t
end do
end if

! terminate if last period/step
if (endofsimulation) then
particle%istatus = TERM_NO_EXITS
particle%ttrack = ttrackmax
call this%save(particle, reason=3)
return
end if
end if
else if (dry_particle .and. this%name /= "passtobottom") then
if (particle%idrymeth == 0) then
! drop to water table
particle%z = cell_defn%top
call this%save(particle, reason=1)
else if (particle%idrymeth == 1) then
! stop
particle%advancing = .false.
particle%istatus = TERM_INACTIVE
call this%save(particle, reason=3)
return
else if (particle%idrymeth == 2) then
! stay
particle%advancing = .false.
no_exit_face = .false.

! we might report tracking times
! out of order here, but we want
! the particle termination event
! (if this is the last time step)
! to have the maximum tracking t,
! so we need to keep tabs on it.
ttrackmax = totim

! update tracking time to time
! step end time and save record
particle%ttrack = totim
call this%save(particle, reason=2)

! record user tracking times
call this%tracktimes%advance()
if (this%tracktimes%any()) then
do i = this%tracktimes%selection(1), this%tracktimes%selection(2)
t = this%tracktimes%times(i)
if (t < totimc) cycle
if (t >= tmax) exit
particle%ttrack = t
call this%save(particle, reason=5)
if (t > ttrackmax) ttrackmax = t
end do
end if
end if
end if

if (no_exit_face) then
particle%advancing = .false.
particle%istatus = TERM_NO_EXITS
call this%save(particle, reason=3)
return
end if

end subroutine check

end module MethodCellModule
4 changes: 2 additions & 2 deletions src/Solution/ParticleTracker/MethodCellPassToBot.f90
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
module MethodCellPassToBotModule

use KindModule, only: DP, I4B
use MethodModule, only: MethodType
use MethodCellModule, only: MethodCellType
use CellDefnModule, only: CellDefnType, create_defn
use PrtFmiModule, only: PrtFmiType
use BaseDisModule, only: DisBaseType
Expand All @@ -17,7 +17,7 @@ module MethodCellPassToBotModule
public :: MethodCellPassToBotType
public :: create_method_cell_ptb

type, extends(MethodType) :: MethodCellPassToBotType
type, extends(MethodCellType) :: MethodCellPassToBotType
private
contains
procedure, public :: apply => apply_ptb
Expand Down
3 changes: 2 additions & 1 deletion src/Solution/ParticleTracker/MethodCellPollock.f90
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@ module MethodCellPollockModule
use KindModule, only: DP, I4B
use ConstantsModule, only: DONE, DZERO
use MethodModule, only: MethodType
use MethodCellModule, only: MethodCellType
use MethodSubcellPoolModule, only: method_subcell_plck, &
method_subcell_tern
use CellRectModule, only: CellRectType, create_cell_rect
Expand All @@ -14,7 +15,7 @@ module MethodCellPollockModule
public :: MethodCellPollockType
public :: create_method_cell_pollock

type, extends(MethodType) :: MethodCellPollockType
type, extends(MethodCellType) :: MethodCellPollockType
contains
procedure, public :: apply => apply_mcp
procedure, public :: deallocate => destroy_mcp
Expand Down
Loading

0 comments on commit 0bc5c4a

Please sign in to comment.