Skip to content

Commit

Permalink
update docstrings, simplify and make private MethodDis cell load rout…
Browse files Browse the repository at this point in the history 
…ines
  • Loading branch information
@wpbonelli
wpbonelli committed Oct 29, 2023
1 parent 96dbcc2 commit 307e661
Show file tree
Hide file tree
Showing 2 changed files with 35 additions and 64 deletions.
3 changes: 2 additions & 1 deletion src/Model/ParticleTracking/prt1prp1.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ module PrtPrpModule
use ArrayHandlersModule, only: expandarray
use GlobalDataModule
use TrackModule, only: TrackControlType
use GeomUtilModule, only: point_in_polygon

implicit none

Expand Down Expand Up @@ -463,7 +464,7 @@ subroutine prp_ad(this)
call get_jk(icu, dis%ncpl, dis%nlay, icpl, ilay)
end select

! -- todo: check that location is within the specified cell
! -- todo: get cell polygon and check that location is within it

! -- Set particle boundname
if (size(this%boundname) /= 0) &
Expand Down
96 changes: 33 additions & 63 deletions src/Solution/ParticleTracker/MethodDis.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,15 +32,14 @@ module MethodDisModule
procedure, public :: apply => apply_mGD ! applies the DIS-grid method
procedure, public :: pass => pass_mGD ! passes the particle to the next cell
procedure, public :: loadsub => loadsub_mGD ! loads the cell method
procedure, private :: get_npolyverts ! returns the number of polygon vertices for a cell in the grid
procedure, private :: get_iatop ! returns index used to locate top elevation of cell in the grid
procedure, private :: get_top ! returns top elevation based on index iatop
procedure, public :: load_cellDefn ! loads cell definition from the grid
procedure, public :: load_cellDefn_facenbr ! loads face neighbors to a cell object
procedure, public :: load_cellDefn_polyverts ! loads polygon vertices to a cell object
procedure, public :: load_cellDefn_flows ! loads flows to a cell object
procedure, public :: load_cellDefn_ispv180 ! loads 180-degree vertex indicator to a cell object
procedure, public :: load_cellDefn_basic ! loads basic components to a cell object from its grid
procedure, private :: load_cellDefn_basic ! loads basic components to a cell object from its grid
procedure, private :: load_cellDefn_polyverts ! loads polygon vertices to a cell object
procedure, private :: load_cellDefn_facenbr ! loads face neighbors to a cell object
procedure, private :: load_cellDefn_ispv180 ! loads 180-degree vertex indicator to a cell object
procedure, private :: load_cellDefn_flows ! loads flows to a cell object
procedure, public :: addBoundaryFlows_cellRect ! adds BoundaryFlows from the grid to the faceflow array of a rectangular cell
end type MethodDisType

Expand Down Expand Up @@ -117,14 +116,15 @@ subroutine loadsub_mGD(this, particle, levelNext, submethod)
ic = particle%iTrackingDomain(levelNext) ! kluge note: is cell number always known coming in?
! -- load cellDefn
call this%load_cellDefn(ic, this%cellRect%cellDefn)

! -- If cell is active but dry, select and initialize
! -- pass-to-bottom method and set cell method pointer
if (this%fmi%ibdgwfsat0(ic) == 0) then ! kluge note: use cellDefn%sat == DZERO here instead?
! -- Cell is active but dry, so select and initialize pass-to-bottom
! -- cell method and set cell method pointer
call methodCellPassToBot%init(particle, this%cellRect%cellDefn, &
this%trackctl)
submethod => methodCellPassToBot
else
! -- load rectangular cell
! -- load rectangular cell (todo: refactor into separate routine)
select type (dis => this%fmi%dis) ! kluge???
type is (GwfDisType)
icu = dis%get_nodeuser(ic)
Expand Down Expand Up @@ -156,8 +156,8 @@ subroutine loadsub_mGD(this, particle, levelNext, submethod)
term = factor / areaz
this%cellRect%vz1 = this%cellRect%cellDefn%faceflow(6) * term
this%cellRect%vz2 = -this%cellRect%cellDefn%faceflow(7) * term
! -- Select and initialize Pollock's cell method and set cell method
! -- pointer

! -- Select and initialize Pollock's method and set method pointer
call methodCellPollock%init(particle, this%cellRect, this%trackctl)
submethod => methodCellPollock
end if
Expand Down Expand Up @@ -264,21 +264,6 @@ subroutine apply_mGD(this, particle, tmax)
!
end subroutine apply_mGD

!> @brief Return the number of polygon vertices for a cell in the grid
function get_npolyverts(this, ic) result(npolyverts)
use GwfDisModule
implicit none
! -- dummy
class(MethodDisType), intent(inout) :: this
integer, intent(in) :: ic
! -- result
integer :: npolyverts
!
npolyverts = 4
!
return
end function get_npolyverts

!> @brief Get the index used to locate top elevation of a cell in the grid
function get_iatop(this, ic) result(iatop)
use GwfDisModule
Expand Down Expand Up @@ -340,15 +325,15 @@ subroutine load_cellDefn(this, ic, cellDefn)
! -- Load 180-degree indicator
call this%load_cellDefn_ispv180(cellDefn)
!
! -- Load flows
! -- (assumes face neighbors already loaded)
! -- Load flows (assumes face neighbors already loaded)
call this%load_cellDefn_flows(cellDefn)
!
return
!
end subroutine load_cellDefn

!> @brief Loads basic cell definition components from the grid
!> @brief Loads basic cell definition components from the grid.
!! These include cell index, vertex count, and (ia)top and botm.
subroutine load_cellDefn_basic(this, ic, cellDefn)
implicit none
! -- dummy
Expand All @@ -359,15 +344,15 @@ subroutine load_cellDefn_basic(this, ic, cellDefn)
integer :: iatop
real(DP) :: top, bot, sat
!
! -- cell index
cellDefn%icell = ic
!
! -- Load npolyverts
cellDefn%npolyverts = this%get_npolyverts(ic)
! -- number of vertices
cellDefn%npolyverts = 4 ! rectangular cell always has 4 vertices
!
! -- Load iatop, top, and bot
! -- iatop, top, and bot
iatop = this%get_iatop(ic)
cellDefn%iatop = iatop
! cellDefn%top = this%get_top(iatop)
top = this%fmi%dis%top(ic)
bot = this%fmi%dis%bot(ic)
sat = this%fmi%gwfsat(ic)
Expand All @@ -376,7 +361,7 @@ subroutine load_cellDefn_basic(this, ic, cellDefn)
cellDefn%bot = bot
cellDefn%sat = sat
!
! -- Load porosity, retfactor, and izone
! -- porosity, retfactor, and izone
cellDefn%porosity = this%porosity(ic)
cellDefn%retfactor = this%retfactor(ic)
cellDefn%izone = this%izone(ic)
Expand All @@ -386,10 +371,7 @@ subroutine load_cellDefn_basic(this, ic, cellDefn)
end subroutine load_cellDefn_basic

!> @brief Load polygon vertices into cell definition from the grid
!!
!! kluge note: are polyverts even needed for MethodDis???
!!
!<
!! Assumes cell index and number of vertices are already loaded.
subroutine load_cellDefn_polyverts(this, cellDefn)
use InputOutputModule ! kluge
use GwfDisModule ! kluge???
Expand All @@ -404,8 +386,6 @@ subroutine load_cellDefn_polyverts(this, cellDefn)
double precision :: cellx, celly, dxhalf, dyhalf
!
ic = cellDefn%icell
!
if (cellDefn%npolyverts .eq. 0) cellDefn%npolyverts = this%get_npolyverts(ic) ! kluge note: just set it to 4 ???
npolyverts = cellDefn%npolyverts
!
! -- Load polygon vertices. Note that the polyverts array
Expand Down Expand Up @@ -440,13 +420,14 @@ subroutine load_cellDefn_polyverts(this, cellDefn)
cellDefn%polyvert(5)%ivert = cellDefn%polyvert(1)%ivert
cellDefn%polyvert(5)%x = cellDefn%polyvert(1)%x
cellDefn%polyvert(5)%y = cellDefn%polyvert(1)%y
end select ! ... kluge???
end select
!
return
!
end subroutine load_cellDefn_polyverts

!> @brief Loads face neighbors to cell definition from the grid
!> @brief Loads face neighbors to cell definition from the grid.
!! Assumes cell index and number of vertices are already loaded.
subroutine load_cellDefn_facenbr(this, cellDefn)
use InputOutputModule
use GwfDisModule
Expand All @@ -462,8 +443,6 @@ subroutine load_cellDefn_facenbr(this, cellDefn)
integer :: ncpl
!
ic = cellDefn%icell
!
if (cellDefn%npolyverts .eq. 0) cellDefn%npolyverts = this%get_npolyverts(ic)
npolyverts = cellDefn%npolyverts
!
! -- Load face neighbors. Note that the facenbr array
Expand Down Expand Up @@ -516,12 +495,9 @@ subroutine load_cellDefn_facenbr(this, cellDefn)
!
end subroutine load_cellDefn_facenbr

!> @brief Load flows into the cell definition
!!
!! Loads polygon face, top and bottom, and net distributed
!! flows from the grid into the cell definition.
!!
!>
!> @brief Load flows into the cell definition.
!! These include face flows and net distributed flows.
!! Assumes cell index and number of vertices are already loaded.
subroutine load_cellDefn_flows(this, cellDefn)
implicit none
! -- dummy
Expand All @@ -531,8 +507,6 @@ subroutine load_cellDefn_flows(this, cellDefn)
integer :: ic, npolyverts, m, n
!
ic = cellDefn%icell
!
if (cellDefn%npolyverts .eq. 0) cellDefn%npolyverts = this%get_npolyverts(ic)
npolyverts = cellDefn%npolyverts
!
! -- Load face flows. Note that the faceflow array
Expand Down Expand Up @@ -574,7 +548,8 @@ subroutine load_cellDefn_flows(this, cellDefn)
!
end subroutine load_cellDefn_flows

!> @brief Add BoundaryFlows to the cell definition faceflow array
!> @brief Add boundary flows to the cell definition faceflow array.
!! Assumes cell index and number of vertices are already loaded.
subroutine addBoundaryFlows_cellRect(this, cellDefn)
implicit none
! -- dummy
Expand Down Expand Up @@ -606,13 +581,9 @@ subroutine addBoundaryFlows_cellRect(this, cellDefn)
!
end subroutine addBoundaryFlows_cellRect

!> @brief Load vertex 180-degree indicator array into cell definition
!!
!! Loads 180-degree vertex indicator array to cell
!! definition and sets flags that indicate how cell
!! can be represented. Todo: rename?
!!
!<
!> @brief Load 180-degree vertex indicator array and set flags
!! indicating how cell can be represented (kluge: latter needed?).
!! Assumes cell index and number of vertices are already loaded.
subroutine load_cellDefn_ispv180(this, cellDefn)
implicit none
! -- dummy
Expand All @@ -621,16 +592,15 @@ subroutine load_cellDefn_ispv180(this, cellDefn)
! -- local
integer :: npolyverts
!
! -- Set up polyverts if not done previously
if (cellDefn%npolyverts .eq. 0) call this%load_cellDefn_polyverts(cellDefn)
npolyverts = cellDefn%npolyverts
!
! -- Load 180-degree indicator. Note that the ispv180 array
! -- does not get reallocated if it is already allocated
! -- to a size greater than or equal to npolyverts+1.
! -- Also, set flags that indicate how cell can be represented.
call allocate_as_needed(cellDefn%ispv180, npolyverts + 1)
cellDefn%ispv180(1:npolyverts + 1) = .false.
!
! -- Set flags that indicate how cell can be represented.
cellDefn%canBeCellRect = .true.
cellDefn%canBeCellRectQuad = .false.
!
Expand Down

0 comments on commit 307e661

Please sign in to comment.