MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
methoddisvmodule Module Reference

Data Types

type  methoddisvtype
 

Functions/Subroutines

subroutine, public create_method_disv (method)
 Create a new vertex grid (DISV) tracking method. More...
 
subroutine deallocate (this)
 Destroy the tracking method. More...
 
subroutine load_disv (this, particle, next_level, submethod)
 Load the cell and the tracking method. More...
 
subroutine load_particle (this, cell, particle)
 
subroutine update_flowja (this, cell, particle)
 
subroutine pass_disv (this, particle)
 Pass a particle to the next cell, if there is one. More...
 
subroutine map_neighbor (this, defn, inface, z)
 Map location on cell face to shared cell face of neighbor. More...
 
subroutine apply_disv (this, particle, tmax)
 Apply the DISV tracking method to a particle. More...
 
subroutine load_cell_defn (this, ic, defn)
 Load cell definition from the grid. More...
 
subroutine load_cell_properties (this, ic, defn)
 Loads cell properties to cell definition from the grid. More...
 
subroutine load_cell_polygon (this, defn)
 
subroutine load_cell_neighbors (this, defn)
 Loads face neighbors to cell definition from the grid Assumes cell index and number of vertices are already loaded. More...
 
subroutine load_cell_flows (this, defn)
 Load flows into the cell definition. These include face, boundary and net distributed flows. Assumes cell index and number of vertices are already loaded. More...
 
subroutine load_cell_face_flows (this, defn)
 
subroutine load_cell_boundary_flows (this, defn)
 Load boundary flows from the grid into a polygonal cell. Assumes cell index and number of vertices are already loaded. More...
 
subroutine load_cell_flags (this, defn)
 Load 180-degree vertex indicator array and set flags indicating how cell can be represented. Assumes cell index and number of vertices are already loaded. More...
 

Function/Subroutine Documentation

◆ apply_disv()

subroutine methoddisvmodule::apply_disv ( class(methoddisvtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
real(dp), intent(in)  tmax 
)
private

Definition at line 303 of file MethodDisv.f90.

304  class(MethodDisvType), intent(inout) :: this
305  type(ParticleType), pointer, intent(inout) :: particle
306  real(DP), intent(in) :: tmax
307 
308  call this%track(particle, 1, tmax)

◆ create_method_disv()

subroutine, public methoddisvmodule::create_method_disv ( type(methoddisvtype), pointer  method)

Definition at line 46 of file MethodDisv.f90.

47  ! dummy
48  type(MethodDisvType), pointer :: method
49  ! local
50  type(CellPolyType), pointer :: cell
51 
52  allocate (method)
53  allocate (method%name)
54  call create_cell_poly(cell)
55  method%cell => cell
56  method%name = "disv"
57  method%delegates = .true.
58  call create_defn(method%neighbor)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ deallocate()

subroutine methoddisvmodule::deallocate ( class(methoddisvtype), intent(inout)  this)
private

Definition at line 62 of file MethodDisv.f90.

63  class(MethodDisvType), intent(inout) :: this
64  deallocate (this%name)

◆ load_cell_boundary_flows()

subroutine methoddisvmodule::load_cell_boundary_flows ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 503 of file MethodDisv.f90.

504  ! dummy
505  class(MethodDisvType), intent(inout) :: this
506  type(CellDefnType), pointer, intent(inout) :: defn
507 
508  ! local
509  integer(I4B) :: ic, iv, ioffset, npolyverts, max_faces
510 
511  ic = defn%icell
512  npolyverts = defn%npolyverts
513  max_faces = this%fmi%max_faces
514  ioffset = (ic - 1) * max_faces
515  do iv = 1, npolyverts
516  defn%faceflow(iv) = &
517  defn%faceflow(iv) + &
518  this%fmi%BoundaryFlows(ioffset + iv)
519  end do
520  defn%faceflow(npolyverts + 1) = defn%faceflow(1)
521  defn%faceflow(npolyverts + 2) = &
522  defn%faceflow(npolyverts + 2) + &
523  this%fmi%BoundaryFlows(ioffset + max_faces - 1)
524  defn%faceflow(npolyverts + 3) = &
525  defn%faceflow(npolyverts + 3) + &
526  this%fmi%BoundaryFlows(ioffset + max_faces)
527 

◆ load_cell_defn()

subroutine methoddisvmodule::load_cell_defn ( class(methoddisvtype), intent(inout)  this,
integer(i4b), intent(in)  ic,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 312 of file MethodDisv.f90.

313  ! dummy
314  class(MethodDisvType), intent(inout) :: this
315  integer(I4B), intent(in) :: ic
316  type(CellDefnType), pointer, intent(inout) :: defn
317 
318  call this%load_cell_properties(ic, defn)
319  call this%load_cell_polygon(defn)
320  call this%load_cell_neighbors(defn)
321  call this%load_cell_saturation_status(defn)
322  call this%load_cell_flags(defn)
323  call this%load_cell_flows(defn)

◆ load_cell_face_flows()

subroutine methoddisvmodule::load_cell_face_flows ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 483 of file MethodDisv.f90.

484  ! dummy
485  class(MethodDisvType), intent(inout) :: this
486  type(CellDefnType), pointer, intent(inout) :: defn
487  ! local
488  integer(I4B) :: m, n, nfaces
489  real(DP) :: q
490 
491  nfaces = defn%npolyverts + 3
492  do m = 1, nfaces
493  n = defn%facenbr(m)
494  if (n > 0) then
495  q = this%fmi%gwfflowja(this%fmi%dis%con%ia(defn%icell) + n)
496  defn%faceflow(m) = defn%faceflow(m) + q
497  end if
498  end do

◆ load_cell_flags()

subroutine methoddisvmodule::load_cell_flags ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 534 of file MethodDisv.f90.

535  ! dummy
536  class(MethodDisvType), intent(inout) :: this
537  type(CellDefnType), pointer, intent(inout) :: defn
538  ! local
539  integer(I4B) :: npolyverts
540  integer(I4B) :: m
541  integer(I4B) :: m0
542  integer(I4B) :: m1
543  integer(I4B) :: m2
544  integer(I4B) :: ic
545  integer(I4B) :: num90
546  integer(I4B) :: num180
547  real(DP) :: x0
548  real(DP) :: y0
549  real(DP) :: x1
550  real(DP) :: y1
551  real(DP) :: x2
552  real(DP) :: y2
553  real(DP) :: epsang
554  real(DP) :: s0x
555  real(DP) :: s0y
556  real(DP) :: &
557  s0mag, s2x, s2y, s2mag, sinang, cosang, dotprod
558  logical(LGP) last180
559 
560  ic = defn%icell
561  npolyverts = defn%npolyverts
562 
563  if (size(defn%ispv180) < npolyverts + 3) &
564  call expandarray(defn%ispv180, npolyverts + 1)
565 
566  defn%ispv180(1:npolyverts + 1) = .false.
567  defn%can_be_rect = .false.
568  defn%can_be_quad = .false.
569  epsang = 1d-5
570  num90 = 0
571  num180 = 0
572  last180 = .false.
573 
574  ! assumes non-self-intersecting polygon
575  do m = 1, npolyverts
576  m1 = m
577  if (m1 .eq. 1) then
578  m0 = npolyverts
579  m2 = 2
580  else if (m1 .eq. npolyverts) then
581  m0 = npolyverts - 1
582  m2 = 1
583  else
584  m0 = m1 - 1
585  m2 = m1 + 1
586  end if
587  x0 = defn%polyvert(1, m0)
588  y0 = defn%polyvert(2, m0)
589  x1 = defn%polyvert(1, m1)
590  y1 = defn%polyvert(2, m1)
591  x2 = defn%polyvert(1, m2)
592  y2 = defn%polyvert(2, m2)
593  s0x = x0 - x1
594  s0y = y0 - y1
595  s0mag = dsqrt(s0x * s0x + s0y * s0y)
596  s2x = x2 - x1
597  s2y = y2 - y1
598  s2mag = dsqrt(s2x * s2x + s2y * s2y)
599  sinang = (s0x * s2y - s0y * s2x) / (s0mag * s2mag)
600  cosang = dsqrt(dabs(done - (sinang * sinang)))
601  if (dabs(sinang) .lt. epsang) then
602  dotprod = s0x * s2x + s0y * s2y
603  if (dotprod .lt. dzero) then
604  num180 = num180 + 1
605  last180 = .true.
606  defn%ispv180(m) = .true.
607  end if
608  else
609  if (dabs(cosang) .lt. epsang) num90 = num90 + 1
610  last180 = .false.
611  end if
612  end do
613 
614  ! List of 180-degree indicators wraps around for convenience
615  defn%ispv180(npolyverts + 1) = defn%ispv180(1)
616 
617  ! Set rect/quad flags
618  if (num90 .eq. 4) then
619  if (num180 .eq. 0) then
620  defn%can_be_rect = .true.
621  else
622  defn%can_be_quad = .true.
623  end if
624  end if

◆ load_cell_flows()

subroutine methoddisvmodule::load_cell_flows ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 448 of file MethodDisv.f90.

449  ! dummy
450  class(MethodDisvType), intent(inout) :: this
451  type(CellDefnType), pointer, intent(inout) :: defn
452  ! local
453  integer(I4B) :: nfaces, nslots
454 
455  ! expand faceflow array if needed
456  nfaces = defn%npolyverts + 3
457  nslots = size(defn%faceflow)
458  if (nslots < nfaces) call expandarray(defn%faceflow, nfaces - nslots)
459 
460  ! Load face flows, including boundary flows. As with cell verts,
461  ! the face flow array wraps around. Top and bottom flows make up
462  ! the last two elements, respectively, for size npolyverts + 3.
463  ! If there is no flow through any face, set a no-exit-face flag.
464  defn%faceflow = dzero
465  call this%load_cell_boundary_flows(defn)
466  call this%load_cell_face_flows(defn)
467  call this%cap_cell_wt_flow(defn)
468  call this%load_cell_no_exit_face(defn)
469 
470  ! Add up net distributed flow
471  defn%distflow = this%fmi%SourceFlows(defn%icell) + &
472  this%fmi%SinkFlows(defn%icell) + &
473  this%fmi%StorageFlows(defn%icell)
474 
475  ! Set weak sink flag
476  if (this%fmi%SinkFlows(defn%icell) .ne. dzero) then
477  defn%iweaksink = 1
478  else
479  defn%iweaksink = 0
480  end if

◆ load_cell_neighbors()

subroutine methoddisvmodule::load_cell_neighbors ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 374 of file MethodDisv.f90.

375  ! dummy
376  class(MethodDisvType), intent(inout) :: this
377  type(CellDefnType), pointer, intent(inout) :: defn
378  ! local
379  integer(I4B) :: ic1
380  integer(I4B) :: ic2
381  integer(I4B) :: icu1
382  integer(I4B) :: icu2
383  integer(I4B) :: j1
384  integer(I4B) :: j2
385  integer(I4B) :: k1
386  integer(I4B) :: k2
387  integer(I4B) :: iloc
388  integer(I4B) :: ipos
389  integer(I4B) :: istart1
390  integer(I4B) :: istart2
391  integer(I4B) :: istop1
392  integer(I4B) :: istop2
393  integer(I4B) :: isharedface
394  integer(I4B) :: ncpl
395  integer(I4B) :: nfaces
396  integer(I4B) :: nslots
397 
398  ! expand facenbr array if needed
399  nfaces = defn%npolyverts + 3
400  nslots = size(defn%facenbr)
401  if (nslots < nfaces) call expandarray(defn%facenbr, nfaces - nslots)
402 
403  select type (dis => this%fmi%dis)
404  type is (disvtype)
405  ! Load face neighbors.
406  defn%facenbr = 0
407  ic1 = defn%icell
408  icu1 = dis%get_nodeuser(ic1)
409  ncpl = dis%get_ncpl()
410  call get_jk(icu1, ncpl, dis%nlay, j1, k1)
411  istart1 = dis%iavert(j1)
412  istop1 = dis%iavert(j1 + 1) - 1
413  do iloc = 1, dis%con%ia(ic1 + 1) - dis%con%ia(ic1) - 1
414  ipos = dis%con%ia(ic1) + iloc
415  ! mask could become relevant if PRT uses interface model
416  if (dis%con%mask(ipos) == 0) cycle
417  ic2 = dis%con%ja(ipos)
418  icu2 = dis%get_nodeuser(ic2)
419  call get_jk(icu2, ncpl, dis%nlay, j2, k2)
420  istart2 = dis%iavert(j2)
421  istop2 = dis%iavert(j2 + 1) - 1
422  call shared_face(dis%javert(istart1:istop1), &
423  dis%javert(istart2:istop2), &
424  isharedface)
425  if (isharedface /= 0) then
426  ! Edge (polygon) face neighbor
427  defn%facenbr(isharedface) = int(iloc, 1)
428  else
429  if (k2 > k1) then
430  ! Bottom face neighbor
431  defn%facenbr(defn%npolyverts + 2) = int(iloc, 1)
432  else if (k2 < k1) then
433  ! Top face neighbor
434  defn%facenbr(defn%npolyverts + 3) = int(iloc, 1)
435  else
436  call pstop(1, "k2 should be <> k1, since no shared edge face")
437  end if
438  end if
439  end do
440  end select
441  ! List of edge (polygon) faces wraps around
442  defn%facenbr(defn%npolyverts + 1) = defn%facenbr(1)
Here is the call graph for this function:

◆ load_cell_polygon()

subroutine methoddisvmodule::load_cell_polygon ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 360 of file MethodDisv.f90.

361  ! dummy
362  class(MethodDisvType), intent(inout) :: this
363  type(CellDefnType), pointer, intent(inout) :: defn
364 
365  call this%fmi%dis%get_polyverts( &
366  defn%icell, &
367  defn%polyvert, &
368  closed=.true.)
369  defn%npolyverts = size(defn%polyvert, dim=2) - 1

◆ load_cell_properties()

subroutine methoddisvmodule::load_cell_properties ( class(methoddisvtype), intent(inout)  this,
integer(i4b), intent(in)  ic,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 327 of file MethodDisv.f90.

328  ! dummy
329  class(MethodDisvType), intent(inout) :: this
330  integer(I4B), intent(in) :: ic
331  type(CellDefnType), pointer, intent(inout) :: defn
332  ! local
333  real(DP) :: top
334  real(DP) :: bot
335  real(DP) :: sat
336  integer(I4B) :: icu, icpl, ilay
337 
338  defn%icell = ic
339  defn%iatop = get_iatop(this%fmi%dis%get_ncpl(), &
340  this%fmi%dis%get_nodeuser(ic))
341  top = this%fmi%dis%top(ic)
342  bot = this%fmi%dis%bot(ic)
343  sat = this%fmi%gwfsat(ic)
344  top = bot + sat * (top - bot)
345  defn%top = top
346  defn%bot = bot
347  defn%sat = sat
348  defn%porosity = this%porosity(ic)
349  defn%retfactor = this%retfactor(ic)
350  defn%izone = this%izone(ic)
351  select type (dis => this%fmi%dis)
352  type is (disvtype)
353  icu = dis%get_nodeuser(ic)
354  call get_jk(icu, dis%ncpl, dis%nlay, icpl, ilay)
355  defn%ilay = ilay
356  end select
357 
Here is the call graph for this function:

◆ load_disv()

subroutine methoddisvmodule::load_disv ( class(methoddisvtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
integer(i4b), intent(in)  next_level,
class(methodtype), intent(inout), pointer  submethod 
)
private

Definition at line 68 of file MethodDisv.f90.

69  use cellmodule
70  use cellrectmodule
72  use cellutilmodule
73  ! dummy
74  class(MethodDisvType), intent(inout) :: this
75  type(ParticleType), pointer, intent(inout) :: particle
76  integer(I4B), intent(in) :: next_level
77  class(MethodType), pointer, intent(inout) :: submethod
78  ! local
79  integer(I4B) :: ic
80  class(CellType), pointer :: base
81  type(CellRectType), pointer :: rect
82  type(CellRectQuadType), pointer :: quad
83 
84  select type (cell => this%cell)
85  type is (cellpolytype)
86  ic = particle%itrdomain(next_level)
87  call this%load_cell_defn(ic, cell%defn)
88  if (this%fmi%ibdgwfsat0(ic) == 0) then
89  ! Cell is active but dry, so select and initialize pass-to-bottom
90  ! cell method and set cell method pointer
91  call method_cell_ptb%init( &
92  fmi=this%fmi, &
93  cell=this%cell, &
94  events=this%events, &
95  tracktimes=this%tracktimes)
96  submethod => method_cell_ptb
97  else if (particle%frctrn) then
98  ! Force the ternary method
99  call method_cell_tern%init( &
100  fmi=this%fmi, &
101  cell=this%cell, &
102  events=this%events, &
103  tracktimes=this%tracktimes)
104  submethod => method_cell_tern
105  else if (cell%defn%can_be_rect) then
106  ! Cell is a rectangle, convert it to a rectangular cell type and
107  ! initialize Pollock's method
108  call cell_poly_to_rect(cell, rect)
109  base => rect
110  call method_cell_plck%init( &
111  fmi=this%fmi, &
112  cell=base, &
113  events=this%events, &
114  tracktimes=this%tracktimes)
115  submethod => method_cell_plck
116  else if (cell%defn%can_be_quad) then
117  ! Cell is quad-refined, convert to a quad rect cell type and
118  ! initialize the corresponding method
119  call cell_poly_to_quad(cell, quad)
120  base => quad
121  call method_cell_quad%init( &
122  fmi=this%fmi, &
123  cell=base, &
124  events=this%events, &
125  tracktimes=this%tracktimes)
126  submethod => method_cell_quad
127  else
128  ! Default to the ternary method
129  call method_cell_tern%init( &
130  fmi=this%fmi, &
131  cell=this%cell, &
132  events=this%events, &
133  tracktimes=this%tracktimes)
134  submethod => method_cell_tern
135  end if
136  end select
subroutine, public cell_poly_to_rect(poly, rect)
Convert CellPoly representation to CellRect. Assumes the conversion is possible.
Definition: CellUtil.f90:14
subroutine, public cell_poly_to_quad(poly, quad)
Convert CellPoly representation to CellRectQuad. Assumes the conversion is possible.
Definition: CellUtil.f90:115
Here is the call graph for this function:

◆ load_particle()

subroutine methoddisvmodule::load_particle ( class(methoddisvtype), intent(inout)  this,
type(cellpolytype), intent(inout), pointer  cell,
type(particletype), intent(inout), pointer  particle 
)

Definition at line 139 of file MethodDisv.f90.

140  ! modules
141  use disvmodule, only: disvtype
142  use particlemodule, only: term_boundary
143  use particleeventmodule, only: terminate
144  ! dummy
145  class(MethodDisvType), intent(inout) :: this
146  type(CellPolyType), pointer, intent(inout) :: cell
147  type(ParticleType), pointer, intent(inout) :: particle
148  ! local
149  integer(I4B) :: inface
150  integer(I4B) :: ipos
151  integer(I4B) :: ic
152  integer(I4B) :: icu
153  integer(I4B) :: inbr
154  integer(I4B) :: idiag
155  integer(I4B) :: icpl
156  integer(I4B) :: ilay
157  real(DP) :: z
158 
159  select type (dis => this%fmi%dis)
160  type is (disvtype)
161  inface = particle%iboundary(level_feature)
162  idiag = dis%con%ia(cell%defn%icell)
163  inbr = cell%defn%facenbr(inface)
164  ipos = idiag + inbr
165  ic = dis%con%ja(ipos)
166  icu = dis%get_nodeuser(ic)
167  call get_jk(icu, dis%ncpl, dis%nlay, icpl, ilay)
168 
169  particle%itrdomain(level_feature) = ic
170  particle%icu = icu
171  particle%ilay = ilay
172 
173  z = particle%z
174  call this%map_neighbor(cell%defn, inface, z)
175 
176  particle%iboundary(level_feature) = inface
177  particle%itrdomain(level_subfeature:) = 0
178  particle%iboundary(level_subfeature:) = 0
179  particle%z = z
180  end select
181 
@, public terminate
particle terminated
@ term_boundary
terminated at a boundary face
Definition: Particle.f90:30
Vertex grid discretization.
Definition: Disv.f90:24
Here is the call graph for this function:

◆ map_neighbor()

subroutine methoddisvmodule::map_neighbor ( class(methoddisvtype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn,
integer(i4b), intent(inout)  inface,
double precision, intent(inout)  z 
)

Definition at line 237 of file MethodDisv.f90.

238  ! dummy
239  class(MethodDisvType), intent(inout) :: this
240  type(CellDefnType), pointer, intent(inout) :: defn
241  integer(I4B), intent(inout) :: inface
242  double precision, intent(inout) :: z
243  ! local
244  integer(I4B) :: icin
245  integer(I4B) :: npolyvertsin
246  integer(I4B) :: ic
247  integer(I4B) :: npolyverts
248  integer(I4B) :: inbr
249  integer(I4B) :: inbrnbr
250  integer(I4B) :: j
251  integer(I4B) :: m
252  real(DP) :: zrel
253  real(DP) :: topfrom
254  real(DP) :: botfrom
255  real(DP) :: top
256  real(DP) :: bot
257  real(DP) :: sat
258 
259  ! Map to shared cell face of neighbor
260  inbr = defn%facenbr(inface)
261  if (inbr .eq. 0) then
262  ! Exterior face; no neighbor to map to
263  ! todo AMP: reconsider when multiple models allowed
264  inface = -1
265  else
266  ! Load definition for neighbor cell (neighbor with shared face)
267  icin = defn%icell
268  j = this%fmi%dis%con%ia(icin)
269  ic = this%fmi%dis%con%ja(j + inbr)
270  call this%load_cell_defn(ic, this%neighbor)
271 
272  npolyvertsin = defn%npolyverts
273  npolyverts = this%neighbor%npolyverts
274  if (inface .eq. npolyvertsin + 2) then
275  ! Exits through bot, enters through top
276  inface = npolyverts + 3
277  else if (inface .eq. npolyvertsin + 3) then
278  ! Exits through top, enters through bot
279  inface = npolyverts + 2
280  else
281  ! Exits and enters through shared polygon face
282  j = this%fmi%dis%con%ia(ic)
283  do m = 1, npolyverts + 3
284  inbrnbr = this%neighbor%facenbr(m)
285  if (this%fmi%dis%con%ja(j + inbrnbr) .eq. icin) then
286  inface = m
287  exit
288  end if
289  end do
290  ! Map z between cells
291  topfrom = defn%top
292  botfrom = defn%bot
293  zrel = (z - botfrom) / (topfrom - botfrom)
294  top = this%fmi%dis%top(ic)
295  bot = this%fmi%dis%bot(ic)
296  sat = this%fmi%gwfsat(ic)
297  z = bot + zrel * sat * (top - bot)
298  end if
299  end if

◆ pass_disv()

subroutine methoddisvmodule::pass_disv ( class(methoddisvtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle 
)
private

Definition at line 208 of file MethodDisv.f90.

209  use particlemodule, only: term_boundary
210  use particleeventmodule, only: terminate
211  ! dummy
212  class(MethodDisvType), intent(inout) :: this
213  type(ParticleType), pointer, intent(inout) :: particle
214  ! local
215  type(CellPolyType), pointer :: cell
216  integer(I4B) :: iface
217  logical(LGP) :: no_neighbors
218 
219  select type (c => this%cell)
220  type is (cellpolytype)
221  cell => c
222  iface = particle%iboundary(level_feature)
223  no_neighbors = cell%defn%facenbr(iface) == 0
224 
225  ! todo AMP: reconsider when multiple models supported
226  if (no_neighbors) then
227  call this%terminate(particle, status=term_boundary)
228  return
229  end if
230 
231  call this%load_particle(cell, particle)
232  call this%update_flowja(cell, particle)
233  end select

◆ update_flowja()

subroutine methoddisvmodule::update_flowja ( class(methoddisvtype), intent(inout)  this,
type(cellpolytype), intent(inout), pointer  cell,
type(particletype), intent(inout), pointer  particle 
)

Definition at line 184 of file MethodDisv.f90.

185  ! dummy
186  class(MethodDisvType), intent(inout) :: this
187  type(CellPolyType), pointer, intent(inout) :: cell
188  type(ParticleType), pointer, intent(inout) :: particle
189  ! local
190  integer(I4B) :: inbr
191  integer(I4B) :: idiag
192  integer(I4B) :: ipos
193 
194  idiag = this%fmi%dis%con%ia(cell%defn%icell)
195  inbr = cell%defn%facenbr(particle%iboundary(level_feature))
196  ipos = idiag + inbr
197 
198  ! leaving old cell
199  this%flowja(ipos) = this%flowja(ipos) - done
200 
201  ! entering new cell
202  this%flowja(this%fmi%dis%con%isym(ipos)) &
203  = this%flowja(this%fmi%dis%con%isym(ipos)) + done
204