Merge pull request #4 from PMEAL/dev

Merging dev into main and setting version to 0.1.0, #feature

Author: jgostick

CLAUDE.md .claude/CLAUDE.mdCLAUDE.md renamed to .claude/CLAUDE.md

@@ -6,7 +6,7 @@
 
 ### Source layout
 - `kabs/_single_phase_solver.py` — `SinglePhaseSolver` class: D3Q19 MRT-LBM Taichi kernels (collision, streaming, BCs)
-- `kabs/_solve_flow.py` — `solve_flow()` entry point and `FlowResult` container class; sets BCs, runs loop, handles convergence/VTK export; returns a `FlowResult`
+- `kabs/_solve_flow.py` — `solve_flow()` entry point and `FlowResult` container class; sets BCs, runs loop, handles convergence; returns a `FlowResult`
 - `kabs/_compute_permeability.py` — `compute_permeability()`: accepts a `FlowResult`, applies Darcy's law
 - `kabs/_compute_hydraulic_conductance.py` — `compute_hydraulic_conductance()`: accepts a `FlowResult`, computes Q/ΔP conductance
 - `kabs/utils.py` — `read_flow_vtr()` / `write_flow_vtr()`: VTR file I/O (pure `numpy` + `struct`, no pyvista dependency); embeds `direction` and `nu` as an XML comment for round-trip fidelity
@@ -31,8 +31,8 @@
 - For the cylinders test image (`r=10`, `ν=1/6`): converges around 500–1000 steps
 
 ### VTR output
-- Written via `result.export_to_vtk(prefix)` or automatically when `export_vtk=True` (default)
-- Filename: `{prefix}-{final_step}-{axis}.vtr`  (step count reflects actual converged step, not `n_steps`)
+- Written via `result.export_to_vtk(prefix)` on the `FlowResult`, or via `write_flow_vtr(prefix, result)` from `kabs.utils`
+- Filename: `{prefix}.vtr`  (caller controls the full prefix, including any step/axis suffix if desired)
 - Flow VTR contains: `Solid` (int8), `rho` (float32), `velocity` (3-component float32)
 - `direction` and `nu` embedded as an XML comment (`<!-- kabs-meta ... -->`) for round-trip fidelity
 - Read back with `read_flow_vtr(path)` from `kabs.utils`; uses only `numpy` + `struct` (no pyvista dependency)

.codecov.yml

@@ -0,0 +1,27 @@
+codecov:
+  branch: dev
+  max_report_age: off
+
+coverage:
+  precision: 1
+  round: down
+  range: "50...100"
+
+  status:
+    project:
+      default:
+        target: auto
+        threshold: 2%
+        branches: null
+        carryforward: true
+
+    patch:
+      default:
+        target: auto
+        threshold: 2%
+        branches: null
+
+comment:
+  layout: "header, diff, changes, sunburst, uncovered"
+  branches: null
+  behavior: default

.github/workflows/bump-version.yml

@@ -10,12 +10,12 @@ jobs:
     runs-on: ubuntu-latest
 
     steps:
-      - uses: actions/checkout@v6
+      - uses: actions/checkout@v4
         with:
           token: ${{ secrets.PAT }}
           fetch-depth: 0
 
-      - uses: astral-sh/setup-uv@v7
+      - uses: astral-sh/setup-uv@v5
         with:
           python-version: "3.13"
           enable-cache: false
@@ -29,6 +29,9 @@ jobs:
         if: github.ref == 'refs/heads/dev'
         run: |
           uvx bump-my-version bump pre_n --no-tag
+          uv lock
+          git add uv.lock
+          git commit --amend --no-edit
           git push
 
       - name: Bump version (release)
@@ -41,6 +44,9 @@ jobs:
           elif [[ "${{ github.event.head_commit.message }}" == *"#major"* ]]; then
             uvx bump-my-version bump major
           fi
+          uv lock
+          git add uv.lock
+          git commit --amend --no-edit
           git push --follow-tags
 
       - name: Merge release back into dev

.github/workflows/generate-docs.yml

@@ -29,7 +29,7 @@ jobs:
 
       - name: Install dependencies
         run: |
-          uv sync
+          uv sync --group dev
           uv pip install "quartodoc" "griffe<1.0"
 
       - name: Build the documentation

.github/workflows/tests.yml

@@ -10,20 +10,28 @@ jobs:
   test:
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v4
+      - uses: actions/checkout@v6
 
       - name: Install uv
-        uses: astral-sh/setup-uv@v5
+        uses: astral-sh/setup-uv@v7
         with:
           enable-cache: true
 
       - name: Set up Python
-        uses: actions/setup-python@v5
+        uses: actions/setup-python@v6
         with:
           python-version: "3.12"
 
       - name: Install dependencies
         run: uv sync --group dev
 
       - name: Run tests
-        run: uv run pytest tests/ -v
+        run: uv run pytest tests/ -v --cov=kabs --cov-report=xml
+
+      - name: Upload coverage to Codecov
+        uses: codecov/codecov-action@v5
+        with:
+          token: ${{ secrets.CODECOV_TOKEN }}
+          files: ./coverage.xml
+          flags: unittests
+          name: codecov-umbrella

.gitignore

@@ -1,3 +1,5 @@
 __pycache__
 *.vtr
 .DS_Store
+.coverage
+coverage.xml

README.md

@@ -1,5 +1,8 @@
 # `kabs`
 
+[![Tests](https://github.com/PMEAL/kabs/actions/workflows/tests.yml/badge.svg)](https://github.com/PMEAL/kabs/actions/workflows/tests.yml)
+[![codecov](https://codecov.io/gh/PMEAL/kabs/branch/dev/graph/badge.svg)](https://codecov.io/gh/PMEAL/kabs)
+
 `kabs` computes the absolute (Darcy) permeability of a porous material from its 3D tomographic image using the Lattice Boltzmann Method (LBM). Given a binary voxel image of the pore space, it solves single-phase incompressible creeping flow, returning results in lattice units or physical units.
 
 ![](taichi_lbm.png)
@@ -36,7 +39,7 @@ ti.init(arch=ti.cpu)  # use ti.gpu for GPU acceleration
 im = ps.generators.cylinders([200, 200, 200], r=10, porosity=0.7).astype(int)
 
 # Run the LBM simulation and get a FlowResult back
-result = solve_flow(im, direction="x", export_vtk=False)
+result = solve_flow(im, direction="x")
 
 # Optionally save to a VTR file for later inspection
 result.export_to_vtk("sample")
@@ -93,14 +96,21 @@ compute_permeability(result)
 The convergence check fires every `log_every` steps (default 500), so the true stopping
 point is rounded to that interval.
 
+To save the converged result to a VTR file, call `export_to_vtk` on the returned object:
+
+```python
+result = solve_flow(im, direction="x", tol=1e-4)
+result.export_to_vtk("sample")  # writes sample-<step>-x.vtr
+```
+
 ### Full permeability tensor
 
 For anisotropic materials, run all three directions:
 
 ```python
 results = {}
 for ax in ("x", "y", "z"):
-    result = solve_flow(im, direction=ax, export_vtk=False)
+    result = solve_flow(im, direction=ax)
     results[ax] = compute_permeability(result, dx_m=2.85e-6)
 
 print(f"Kx={results['x']['k_mD']:.2f}  Ky={results['y']['k_mD']:.2f}  Kz={results['z']['k_mD']:.2f}  mD")
@@ -124,7 +134,7 @@ For images with a high solid fraction, enable sparse storage so only pore voxels
 allocated in GPU memory:
 
 ```python
-solver = solve_flow(im, direction="x", sparse=True)
+result = solve_flow(im, direction="x", sparse=True)
 ```
 
 ## Return values (permeability)

docs/taichi_lbm.png

@@ -19,19 +19,18 @@
     dx_m    = physical voxel size in metres
 """
 
-
 import numpy as np
 
 
 __all__ = ["compute_hydraulic_conductance"]
 
 
-_RHO_IN  = 1.00
+_RHO_IN = 1.00
 _RHO_OUT = 0.99
 
 
 def compute_hydraulic_conductance(
-    source,
+    soln,
     direction=None,
     nu=None,
     dx_m=None,
@@ -47,18 +46,14 @@ def compute_hydraulic_conductance(
 
     Parameters
     ----------
-    source : FlowResult or str/path-like
-        Either a ``FlowResult`` returned by ``solve_flow()``, or a path to a
-        ``.vtr`` file written by ``SinglePhaseSolver.export_VTK()``.
-        When a ``FlowResult`` is given, ``direction`` and ``nu`` default to the
-        values stored in the result.
+    soln : FlowResult
+        A ``FlowResult`` returned by ``solve_flow()``.  ``direction`` and
+        ``nu`` default to the values stored in the result.
     direction : {'x', 'y', 'z'} or None
-        Flow direction.  If *None* and ``source`` is a ``FlowResult``, taken
-        from ``source.direction``; otherwise defaults to ``'x'``.
+        Flow direction.  If *None*, taken from ``soln.direction``.
     nu : float or None
         Kinematic viscosity used in the LBM simulation (lattice units).
-        If *None* and ``source`` is a ``FlowResult``, taken from ``source.nu``;
-        otherwise defaults to 1/6.
+        If *None*, taken from ``soln.nu``.
     dx_m : float or None
         Physical voxel size in metres.  Required for physical-unit output.
     mu_phys : float or None
@@ -76,11 +71,12 @@ def compute_hydraulic_conductance(
         Q_m3s      – volumetric flow rate in m^3/s  (None if dx_m/mu_phys not given)
         dP_Pa      – total pressure drop in Pa       (None if dx_m/mu_phys not given)
         g_SI       – conductance in m^3/(Pa·s)       (None if dx_m/mu_phys not given)
+        summary    – human-readable result summary string (always populated)
     """
-    _dir = direction if direction is not None else source.direction
-    _nu  = nu        if nu        is not None else source.nu
-    solid    = source.solid
-    velocity = source.velocity
+    _dir = direction if direction is not None else soln.direction
+    _nu = nu if nu is not None else soln.nu
+    solid = soln.solid
+    velocity = soln.velocity
 
     _dir = _dir.lower()
     if _dir not in ("x", "y", "z"):
@@ -94,53 +90,68 @@ def compute_hydraulic_conductance(
     vz = velocity[..., 2]
     nx, ny, nz = solid.shape
 
-    v_flow  = {"x": vx,      "y": vy,      "z": vz     }[_dir]
-    L_flow  = {"x": nx,      "y": ny,      "z": nz     }[_dir]
+    v_flow = {"x": vx, "y": vy, "z": vz}[_dir]
+    L_flow = {"x": nx, "y": ny, "z": nz}[_dir]
     A_cross = {"x": ny * nz, "y": nx * nz, "z": nx * ny}[_dir]
 
     u_darcy = float(np.mean(v_flow))
-    Q_lu    = u_darcy * A_cross
-    dP_lu   = (_RHO_IN - _RHO_OUT) * cs2
-    g_lu    = Q_lu / dP_lu
+    Q_lu = u_darcy * A_cross
+    dP_lu = (_RHO_IN - _RHO_OUT) * cs2
+    g_lu = Q_lu / dP_lu
 
     Q_m3s = dP_Pa = g_SI = None
     can_convert = (dx_m is not None) and (mu_phys is not None)
     if can_convert:
         g_SI = g_lu * _nu * dx_m**3 / mu_phys
 
+    n_pore = int(pore_mask.sum())
+    porosity = n_pore / pore_mask.size
+    ax_idx = {"x": 0, "y": 1, "z": 2}[_dir]
+    slices = np.array([np.sum(np.take(v_flow, i, axis=ax_idx)) for i in range(L_flow)])
+    lines = [
+        f"",
+        f"Flow direction           = {_dir}",
+        f"Conduit length           = {L_flow}  [lu]",
+        f"Cross-section area       = {A_cross}  [lu^2]",
+        f"Pore voxels              = {n_pore}  (porosity = {porosity:.4f})",
+        f"",
+        f"Darcy velocity  u_D      = {u_darcy:.6e}  [lu/ts]",
+        f"Volumetric flow Q        = {Q_lu:.6e}  [lu^3/ts]",
+        f"Pressure drop   dP       = {dP_lu:.6f}  [lu pressure]",
+        f"",
+        f"Conductance     g        = {g_lu:.6e}  [lu^3/ts / lu_pressure]",
+    ]
+    if can_convert:
+        lines += [
+            f"",
+            f"With dx = {dx_m:.4e} m  and  mu = {mu_phys:.4e} Pa·s:",
+            f"  g = {g_SI:.4e}  m^3/(Pa·s)",
+        ]
+    elif dx_m is None:
+        lines.append(
+            "To get physical units: pass dx_m (voxel size in metres) and mu_phys (dynamic viscosity in Pa·s)."
+        )
+    else:
+        lines.append(
+            "To get physical units: also pass mu_phys (dynamic viscosity in Pa·s)."
+        )
+    lines += [
+        f"",
+        f"--- Sanity check: per-slice Q (should be constant) ---",
+        f"  Q_slice min={slices.min():.4e}  max={slices.max():.4e}"
+        f"  mean={slices.mean():.4e}  std={slices.std():.4e}",
+    ]
+    summary = "\n".join(lines)
+
     if verbose:
-        n_pore   = int(pore_mask.sum())
-        porosity = n_pore / pore_mask.size
-        print(f"\nFlow direction           = {_dir}")
-        print(f"Conduit length           = {L_flow}  [lu]")
-        print(f"Cross-section area       = {A_cross}  [lu^2]")
-        print(f"Pore voxels              = {n_pore}  (porosity = {porosity:.4f})")
-        print(f"\nDarcy velocity  u_D      = {u_darcy:.6e}  [lu/ts]")
-        print(f"Volumetric flow Q        = {Q_lu:.6e}  [lu^3/ts]")
-        print(f"Pressure drop   dP       = {dP_lu:.6f}  [lu pressure]")
-        print(f"\nConductance     g        = {g_lu:.6e}  [lu^3/ts / lu_pressure]")
-        if can_convert:
-            print(f"\nWith dx = {dx_m:.4e} m  and  mu = {mu_phys:.4e} Pa·s:")
-            print(f"  g = {g_SI:.4e}  m^3/(Pa·s)")
-        elif dx_m is None:
-            print("\nTo get physical units: pass dx_m (voxel size in metres)"
-                  " and mu_phys (dynamic viscosity in Pa·s).")
-        else:
-            print("\nTo get physical units: also pass mu_phys (dynamic viscosity in Pa·s).")
-
-        ax_idx = {"x": 0, "y": 1, "z": 2}[_dir]
-        slices = np.array([np.sum(np.take(v_flow, i, axis=ax_idx)) for i in range(L_flow)])
-        print(f"\n--- Sanity check: per-slice Q (should be constant) ---")
-        print(f"  Q_slice min={slices.min():.4e}  max={slices.max():.4e}"
-              f"  mean={slices.mean():.4e}  std={slices.std():.4e}")
+        print(summary)
 
     return {
-        "Q_lu":  Q_lu,
+        "Q_lu": Q_lu,
         "dP_lu": dP_lu,
-        "g_lu":  g_lu,
+        "g_lu": g_lu,
         "Q_m3s": Q_m3s,
         "dP_Pa": dP_Pa,
-        "g_SI":  g_SI,
+        "g_SI": g_SI,
+        "summary": summary,
     }
-
-