diff --git a/src/Solution/ParticleTracker/MethodCellPassToBot.f90 b/src/Solution/ParticleTracker/MethodCellPassToBot.f90 index eb417457377..433b3ca9650 100644 --- a/src/Solution/ParticleTracker/MethodCellPassToBot.f90 +++ b/src/Solution/ParticleTracker/MethodCellPassToBot.f90 @@ -9,8 +9,6 @@ module MethodCellPassToBotModule use CellModule, only: CellType use SubcellModule, only: SubcellType use TrackControlModule, only: TrackControlType - use DisModule, only: DisType - use DisvModule, only: DisvType implicit none private @@ -47,8 +45,6 @@ subroutine apply_ptb(this, particle, tmax) class(MethodCellPassToBotType), intent(inout) :: this type(ParticleType), pointer, intent(inout) :: particle real(DP), intent(in) :: tmax - ! local - integer(I4B) :: nlay, reason ! Check termination/reporting conditions call this%check(particle, this%cell%defn) @@ -58,20 +54,6 @@ subroutine apply_ptb(this, particle, tmax) particle%z = this%cell%defn%bot particle%iboundary(2) = this%cell%defn%npolyverts + 2 - ! Terminate if bottom layer - select type (dis => this%fmi%dis) - type is (DisType) - nlay = dis%nlay - type is (DisvType) - nlay = dis%nlay - end select - if (this%cell%defn%ilay == nlay) then - particle%istatus = 5 - reason = 3 - else - reason = 1 - end if - ! Save datum call this%save(particle, reason=1) end subroutine apply_ptb diff --git a/src/Solution/ParticleTracker/MethodDis.f90 b/src/Solution/ParticleTracker/MethodDis.f90 index ace302949e0..61c67430f14 100644 --- a/src/Solution/ParticleTracker/MethodDis.f90 +++ b/src/Solution/ParticleTracker/MethodDis.f90 @@ -195,12 +195,19 @@ subroutine load_particle(this, cell, particle) ic = dis%con%ja(ipos) if (ic == particle%icp) then - particle%istatus = 2 + ! terminate in the previous cell. + ! got here by falling through the + ! bottom of a well.. TODO: check + ! if entry/exit face is top/bot? particle%advancing = .false. + particle%idomain(2) = particle%icp + particle%istatus = 2 + particle%izone = particle%izp call this%save(particle, reason=3) return else particle%icp = particle%idomain(2) + particle%izp = particle%izone end if icu = dis%get_nodeuser(ic) diff --git a/src/Solution/ParticleTracker/MethodDisv.f90 b/src/Solution/ParticleTracker/MethodDisv.f90 index 479b0f355bc..d33a60e1ba0 100644 --- a/src/Solution/ParticleTracker/MethodDisv.f90 +++ b/src/Solution/ParticleTracker/MethodDisv.f90 @@ -157,12 +157,19 @@ subroutine load_particle(this, cell, particle) ic = dis%con%ja(ipos) if (ic == particle%icp) then - particle%istatus = 2 + ! terminate in the previous cell. + ! got here by falling through the + ! bottom of a well.. TODO: check + ! if entry/exit face is top/bot? particle%advancing = .false. + particle%idomain(2) = particle%icp + particle%istatus = 2 + particle%izone = particle%izp call this%save(particle, reason=3) return else particle%icp = particle%idomain(2) + particle%izp = particle%izone end if icu = dis%get_nodeuser(ic) diff --git a/src/Solution/ParticleTracker/Particle.f90 b/src/Solution/ParticleTracker/Particle.f90 index a3e7e73d9ce..add4dbac914 100644 --- a/src/Solution/ParticleTracker/Particle.f90 +++ b/src/Solution/ParticleTracker/Particle.f90 @@ -47,7 +47,8 @@ module ParticleModule integer(I4B), public :: icp !< previous cell number (reduced) integer(I4B), public :: icu !< user cell number integer(I4B), public :: ilay !< grid layer - integer(I4B), public :: izone !< zone number + integer(I4B), public :: izone !< current zone number + integer(I4B), public :: izp !< previous zone number integer(I4B), public :: istatus !< tracking status real(DP), public :: x !< x coordinate real(DP), public :: y !< y coordinate @@ -92,6 +93,7 @@ module ParticleModule integer(I4B), dimension(:), pointer, public, contiguous :: icu !< cell number (user) integer(I4B), dimension(:), pointer, public, contiguous :: ilay !< layer integer(I4B), dimension(:), pointer, public, contiguous :: izone !< current zone number + integer(I4B), dimension(:), pointer, public, contiguous :: izp !< previous zone number integer(I4B), dimension(:), pointer, public, contiguous :: istatus !< particle status real(DP), dimension(:), pointer, public, contiguous :: x !< model x coord of particle real(DP), dimension(:), pointer, public, contiguous :: y !< model y coord of particle @@ -135,6 +137,7 @@ subroutine allocate_particle_store(this, np, mempath) call mem_allocate(this%icu, np, 'PLICU', mempath) call mem_allocate(this%ilay, np, 'PLILAY', mempath) call mem_allocate(this%izone, np, 'PLIZONE', mempath) + call mem_allocate(this%izp, np, 'PLIZP', mempath) call mem_allocate(this%istatus, np, 'PLISTATUS', mempath) call mem_allocate(this%x, np, 'PLX', mempath) call mem_allocate(this%y, np, 'PLY', mempath) @@ -166,6 +169,7 @@ subroutine deallocate (this, mempath) call mem_deallocate(this%icu, 'PLICU', mempath) call mem_deallocate(this%ilay, 'PLILAY', mempath) call mem_deallocate(this%izone, 'PLIZONE', mempath) + call mem_deallocate(this%izp, 'PLIZP', mempath) call mem_deallocate(this%istatus, 'PLISTATUS', mempath) call mem_deallocate(this%x, 'PLX', mempath) call mem_deallocate(this%y, 'PLY', mempath) @@ -200,6 +204,7 @@ subroutine resize(this, np, mempath) call mem_reallocate(this%icu, np, 'PLICU', mempath) call mem_reallocate(this%ilay, np, 'PLILAY', mempath) call mem_reallocate(this%izone, np, 'PLIZONE', mempath) + call mem_reallocate(this%izp, np, 'PLIZP', mempath) call mem_reallocate(this%istatus, np, 'PLISTATUS', mempath) call mem_reallocate(this%x, np, 'PLX', mempath) call mem_reallocate(this%y, np, 'PLY', mempath) @@ -243,6 +248,7 @@ subroutine load_particle(this, store, imdl, iprp, ip) this%icu = store%icu(ip) this%ilay = store%ilay(ip) this%izone = store%izone(ip) + this%izp = store%izp(ip) this%istatus = store%istatus(ip) this%x = store%x(ip) this%y = store%y(ip) @@ -279,6 +285,7 @@ subroutine save_particle(this, particle, ip) this%icu(ip) = particle%icu this%ilay(ip) = particle%ilay this%izone(ip) = particle%izone + this%izp(ip) = particle%izp this%istatus(ip) = particle%istatus this%x(ip) = particle%x this%y(ip) = particle%y