From 244583587cc2d2fcb375831dc721e87d63387e22 Mon Sep 17 00:00:00 2001 From: wpbonelli Date: Wed, 8 Jan 2025 18:46:05 -0500 Subject: [PATCH] fix(prt): don't reset transform in subcell rect method (#2131) #2066 introduced a bug where particle coordinates were prematurely transformed back from cell-local coords to model coords. This produced incorrect particle positions in/near quad-refined cells. --- src/Solution/ParticleTracker/CellDefn.f90 | 2 -- src/Solution/ParticleTracker/Method.f90 | 3 +-- src/Solution/ParticleTracker/MethodDisv.f90 | 11 ++++++++++- src/Solution/ParticleTracker/MethodSubcellPollock.f90 | 1 - 4 files changed, 11 insertions(+), 6 deletions(-) diff --git a/src/Solution/ParticleTracker/CellDefn.f90 b/src/Solution/ParticleTracker/CellDefn.f90 index 836f0421633..a2f0fd5c5f4 100644 --- a/src/Solution/ParticleTracker/CellDefn.f90 +++ b/src/Solution/ParticleTracker/CellDefn.f90 @@ -1,7 +1,5 @@ module CellDefnModule use KindModule, only: DP, I4B, LGP - use ConstantsModule, only: DZERO - use MathUtilModule, only: is_close implicit none private diff --git a/src/Solution/ParticleTracker/Method.f90 b/src/Solution/ParticleTracker/Method.f90 index 061a627e1ad..4035a0ad8a1 100644 --- a/src/Solution/ParticleTracker/Method.f90 +++ b/src/Solution/ParticleTracker/Method.f90 @@ -134,9 +134,8 @@ subroutine try_pass(this, particle, nextlevel, advancing) ! otherwise pass the particle to the next subdomain. ! if that leaves it on a boundary, stop advancing. call this%pass(particle) - if (particle%iboundary(nextlevel - 1) .ne. 0) then + if (particle%iboundary(nextlevel - 1) .ne. 0) & advancing = .false. - end if end if end subroutine try_pass diff --git a/src/Solution/ParticleTracker/MethodDisv.f90 b/src/Solution/ParticleTracker/MethodDisv.f90 index de65bbbf134..d6fefc4118a 100644 --- a/src/Solution/ParticleTracker/MethodDisv.f90 +++ b/src/Solution/ParticleTracker/MethodDisv.f90 @@ -88,6 +88,8 @@ subroutine load_disv(this, particle, next_level, submethod) ic = particle%idomain(next_level) call this%load_cell_defn(ic, cell%defn) if (this%fmi%ibdgwfsat0(ic) == 0) then + ! Cell is active but dry, so select and initialize pass-to-bottom + ! cell method and set cell method pointer call method_cell_ptb%init( & fmi=this%fmi, & cell=this%cell, & @@ -95,6 +97,7 @@ subroutine load_disv(this, particle, next_level, submethod) tracktimes=this%tracktimes) submethod => method_cell_ptb else if (particle%ifrctrn > 0) then + ! Force the ternary method call method_cell_tern%init( & fmi=this%fmi, & cell=this%cell, & @@ -102,6 +105,8 @@ subroutine load_disv(this, particle, next_level, submethod) tracktimes=this%tracktimes) submethod => method_cell_tern else if (cell%defn%can_be_rect) then + ! Cell is a rectangle, convert it to a rectangular cell type and + ! initialize Pollock's method call cell_poly_to_rect(cell, rect) base => rect call method_cell_plck%init( & @@ -111,6 +116,8 @@ subroutine load_disv(this, particle, next_level, submethod) tracktimes=this%tracktimes) submethod => method_cell_plck else if (cell%defn%can_be_quad) then + ! Cell is quad-refined, convert to a quad rect cell type and + ! initialize the corresponding method call cell_poly_to_quad(cell, quad) base => quad call method_cell_quad%init( & @@ -120,6 +127,7 @@ subroutine load_disv(this, particle, next_level, submethod) tracktimes=this%tracktimes) submethod => method_cell_quad else + ! Default to the ternary method call method_cell_tern%init( & fmi=this%fmi, & cell=this%cell, & @@ -232,7 +240,8 @@ subroutine pass_disv(this, particle) particle%advancing = .false. call this%save(particle, reason=3) else - ! Update old to new cell properties + ! Otherwise, load cell properties into the + ! particle. It may be marked to terminate. call this%load_particle(cell, particle) if (.not. particle%advancing) return diff --git a/src/Solution/ParticleTracker/MethodSubcellPollock.f90 b/src/Solution/ParticleTracker/MethodSubcellPollock.f90 index 1cbd0cce45b..55e8989a925 100644 --- a/src/Solution/ParticleTracker/MethodSubcellPollock.f90 +++ b/src/Solution/ParticleTracker/MethodSubcellPollock.f90 @@ -73,7 +73,6 @@ subroutine apply_msp(this, particle, tmax) call particle%transform(xOrigin, yOrigin) call this%track_subcell(subcell, particle, tmax) call particle%transform(xOrigin, yOrigin, invert=.true.) - call particle%reset_transform() end select end subroutine apply_msp