feat(PRT): add MF6 PRT utilities

* prt -> mp7 pathline/endpoint data format conversion adapters
* to_coords() method for ParticleData classes

autotest/test_binaryfile.py

@@ -1,4 +1,5 @@
 import os
+from typing import List
 
 import numpy as np
 import pytest

autotest/test_particledata.py

@@ -1,8 +1,17 @@
 import numpy as np
+import pytest
+from autotest.test_grid_cases import GridCases
 
-from flopy.modpath import ParticleData
+from flopy.discretization import StructuredGrid
+from flopy.modpath import (
+    CellDataType,
+    FaceDataType,
+    LRCParticleData,
+    NodeParticleData,
+    ParticleData,
+)
 
-structured_plocs = [(1, 1, 1), (1, 1, 2)]
+structured_cells = [(0, 1, 1), (0, 1, 2)]
 structured_dtype = np.dtype(
     [
         ("k", "<i4"),
@@ -17,35 +26,227 @@
 )
 structured_array = np.core.records.fromrecords(
     [
-        (1, 1, 1, 0.5, 0.5, 0.5, 0.0, 0),
-        (1, 1, 2, 0.5, 0.5, 0.5, 0.0, 0),
+        (0, 1, 1, 0.5, 0.5, 0.5, 0.0, 0),
+        (0, 1, 2, 0.5, 0.5, 0.5, 0.0, 0),
     ],
     dtype=structured_dtype,
 )
 
 
-def test_particledata_structured_partlocs_as_list_of_tuples():
-    locs = structured_plocs
+def test_particledata_structured_ctor_with_partlocs_as_list_of_tuples():
+    locs = structured_cells
     data = ParticleData(partlocs=locs, structured=True)
 
     assert data.particlecount == 2
     assert data.dtype == structured_dtype
     assert np.array_equal(data.particledata, structured_array)
 
 
-def test_particledata_structured_partlocs_as_ndarray():
-    locs = np.array(structured_plocs)
+def test_particledata_structured_ctor_with_partlocs_as_ndarray():
+    locs = np.array(structured_cells)
     data = ParticleData(partlocs=locs, structured=True)
 
     assert data.particlecount == 2
     assert data.dtype == structured_dtype
     assert np.array_equal(data.particledata, structured_array)
 
 
-def test_particledata_structured_partlocs_as_list_of_lists():
-    locs = [list(p) for p in structured_plocs]
+def test_particledata_structured_ctor_with_partlocs_as_list_of_lists():
+    locs = [list(p) for p in structured_cells]
     data = ParticleData(partlocs=locs, structured=True)
 
     assert data.particlecount == 2
     assert data.dtype == structured_dtype
     assert np.array_equal(data.particledata, structured_array)
+
+
+def test_particledata_structured_to_coords_with_ids():
+    data = ParticleData(
+        partlocs=structured_cells,
+        structured=True,
+        particleids=range(len(structured_cells)),
+    )
+    grid = GridCases().structured_small()
+    coords = data.to_coords(grid)
+
+    assert len(coords) == len(structured_cells)
+    assert all(
+        len(c) == 4 for c in coords
+    )  # each coord should be a tuple (id, x, y, z)
+    for ci, c in enumerate(coords):
+        assert ci == c[0]
+
+
+def test_particledata_structured_to_coords_with_no_ids():
+    data = ParticleData(partlocs=structured_cells, structured=True)
+    grid = GridCases().structured_small()
+    coords = data.to_coords(grid)
+
+    assert len(coords) == len(structured_cells)
+    assert all(
+        len(c) == 3 for c in coords
+    )  # each coord should be a tuple (x, y, z)
+
+    for ci, cell in enumerate(structured_cells):
+        # check containing cell is correct
+        coord = coords[ci]
+        assert cell == grid.intersect(*coord)
+
+        # check global coords are equivalent to local
+        k, i, j = cell
+        verts = grid.get_cell_vertices(i, j)
+        xs, ys = list(zip(*verts))
+        minx, maxx = min(xs), max(xs)
+        miny, maxy = min(ys), max(ys)
+        minz, maxz = grid.botm[k, i, j], grid.top[i, j]
+        assert np.isclose(
+            coord[2], minz + (structured_array[ci][3] * (maxz - minz))
+        )
+        assert np.isclose(
+            coord[0], minx + (structured_array[ci][4] * (maxx - minx))
+        )
+        assert np.isclose(
+            coord[1], miny + (structured_array[ci][5] * (maxy - miny))
+        )
+
+
+def test_lrcparticledata_to_coords_with_particle_on_top_and_bottom_faces():
+    rd = 1
+    cd = 1
+    sddata = FaceDataType(
+        horizontaldivisions1=0,
+        verticaldivisions1=0,
+        horizontaldivisions2=0,
+        verticaldivisions2=0,
+        horizontaldivisions3=0,
+        verticaldivisions3=0,
+        horizontaldivisions4=0,
+        verticaldivisions4=0,
+        rowdivisions5=rd,
+        columndivisions5=cd,
+        rowdivisions6=rd,
+        columndivisions6=cd,
+    )
+    lrcregions = [[0, 1, 1, 0, 1, 1]]
+    data = LRCParticleData(subdivisiondata=sddata, lrcregions=lrcregions)
+    grid = GridCases().structured_small()
+    coords = data.to_coords(grid)
+
+    num_cells = len(
+        [
+            (lrc[3] - lrc[0]) * (lrc[4] - lrc[1]) * (lrc[5] - lrc[2])
+            for lrc in lrcregions
+        ]
+    )
+    assert (
+        len(coords) == num_cells * rd * cd * 2
+    )  # 1 particle each on top and bottom faces
+
+    # particle should be centered on each face
+    verts = grid.get_cell_vertices(1, 1)
+    xs, ys = list(zip(*verts))
+    for coord in coords:
+        assert np.isclose(coord[0], np.mean(xs))
+        assert np.isclose(coord[1], np.mean(ys))
+
+    # check elevation
+    assert coords[0][2] == grid.top_botm[1, 1, 1]
+    assert coords[1][2] == grid.top_botm[0, 1, 1]
+
+
+def test_lrcparticledata_to_coords_with_particle_on_all_faces():
+    sddata = FaceDataType(
+        horizontaldivisions1=1,
+        verticaldivisions1=1,
+        horizontaldivisions2=1,
+        verticaldivisions2=1,
+        horizontaldivisions3=1,
+        verticaldivisions3=1,
+        horizontaldivisions4=1,
+        verticaldivisions4=1,
+        rowdivisions5=1,
+        columndivisions5=1,
+        rowdivisions6=1,
+        columndivisions6=1,
+    )
+    lrcregions = [[0, 1, 1, 0, 1, 1]]
+    data = LRCParticleData(subdivisiondata=sddata, lrcregions=lrcregions)
+    grid = GridCases().structured_small()
+    coords = data.to_coords(grid)
+
+    num_cells = len(
+        [
+            (lrc[3] - lrc[0]) * (lrc[4] - lrc[1]) * (lrc[5] - lrc[2])
+            for lrc in lrcregions
+        ]
+    )
+    assert len(coords) == num_cells * 6  # 1 particle on each face
+
+
+def test_lrcparticledata_to_coords_with_particles_within_cell():
+    sddata = CellDataType()
+    lrcregions = [[0, 1, 1, 0, 1, 1]]
+    data = LRCParticleData(subdivisiondata=sddata, lrcregions=lrcregions)
+    grid = GridCases().structured_small()
+    coords = data.to_coords(grid)
+
+    num_cells = len(
+        [
+            (lrc[3] - lrc[0]) * (lrc[4] - lrc[1]) * (lrc[5] - lrc[2])
+            for lrc in lrcregions
+        ]
+    )
+    assert (
+        len(coords)
+        == num_cells
+        * sddata.rowcelldivisions
+        * sddata.columncelldivisions
+        * sddata.layercelldivisions
+    )
+
+
+def get_nn(grid: StructuredGrid, k, i, j):
+    return k * grid.nrow * grid.ncol + i * grid.ncol + j
+
+
+def test_nodeparticledata_to_coords_with_particle_on_all_faces():
+    sddata = FaceDataType(
+        horizontaldivisions1=1,
+        verticaldivisions1=1,
+        horizontaldivisions2=1,
+        verticaldivisions2=1,
+        horizontaldivisions3=1,
+        verticaldivisions3=1,
+        horizontaldivisions4=1,
+        verticaldivisions4=1,
+        rowdivisions5=1,
+        columndivisions5=1,
+        rowdivisions6=1,
+        columndivisions6=1,
+    )
+    grid = GridCases().structured_small()
+    nodes = [get_nn(grid, 0, 1, 1)]
+    data = NodeParticleData(subdivisiondata=sddata, nodes=[nodes])
+
+    coords = data.to_coords(grid)
+
+    num_cells = len(nodes)
+    assert len(coords) == num_cells * 6
+
+
+def test_nodeparticledata_to_coords_with_particles_within_cell():
+    sddata = CellDataType()
+    grid = GridCases().structured_small()
+    nodes = [get_nn(grid, 0, 1, 1)]
+    data = NodeParticleData(subdivisiondata=sddata, nodes=[nodes])
+
+    coords = data.to_coords(grid)
+
+    num_cells = len(nodes)
+    assert (
+        len(coords)
+        == num_cells
+        * sddata.rowcelldivisions
+        * sddata.columncelldivisions
+        * sddata.layercelldivisions
+    )