From b74b279432d6b1a35a1b386724a5ff997360ffe9 Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 13:56:58 +0100 Subject: [PATCH 1/6] codex: restructure the Reference BLAS example guide Lead with what the example builds, a two-command quick start with a usage snippet, and a key-files table; explain the single shared-library build with a diagram before the build scripts; collect the PRIK and f2py return-value differences in one table; and move environment setup, versions, and platforms to the end. The quick start and snippet were run against a fresh build (167 tests passed). Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 3 + docs/user/examples/fortran/blas-wrapper.md | 285 +++++++++++---------- 2 files changed, 152 insertions(+), 136 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 752ea5246..a04307d87 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,9 @@ release tags add a leading `v` to the package version. ## Unreleased +- The Reference BLAS example guide now starts with a quick start, its key + files, and a short explanation of the shared-library build, and lists the + return-value differences between the PRIK and f2py wrappers in one table. - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/blas-wrapper.md b/docs/user/examples/fortran/blas-wrapper.md index abdbe6df3..f65610957 100644 --- a/docs/user/examples/fortran/blas-wrapper.md +++ b/docs/user/examples/fortran/blas-wrapper.md @@ -9,94 +9,94 @@ publication: reviewed # Build and Validate the Reference BLAS with PRIK -This example builds the official Reference BLAS sources as two importable -Python extension modules: +This example wraps all 155 routines of the Reference BLAS twice, once with +PRIK and once with NumPy's f2py, and checks both against independent +mathematical results. BLAS is compiled once into a shared library that both +wrappers link, so a difference between them comes from the wrapper, never +from the numerics. -- one generated by **PRIK** -- one generated by **NumPy’s f2py** +### What you get -It compares both wrappers with independent mathematical results across all 155 -callable Reference BLAS routines. - -### What this example shows - -- Build PRIK and f2py wrappers against the same compiled BLAS library. -- Call vector and matrix routines with NumPy arrays. -- Compare numerical results and the Python interfaces produced by each tool. - -You should already be comfortable with NumPy arrays, basic packaging, and building Fortran extensions. +- Two extension modules, `prik_reference_blas` and `f2py_reference_blas`, each + exposing the same 155 routines: Level 1 vector, Level 2 matrix-vector, and + Level 3 matrix-matrix operations, in real and complex precisions, including + packed, banded, symmetric, Hermitian, and triangular storage. +- A named test for every routine that checks PRIK and f2py against an + independent formula and against each other. --- -## Versions used - -| Component | Version / source | -| --- | --- | -| PRIK | current repository checkout | -| Reference BLAS | snapshot shipped in Netlib LAPACK 3.12.1 | -| Python | 3.12 in the dedicated CI job | -| NumPy / f2py | NumPy 2.5.1 | -| Meson | 1.11.2 | -| Ninja | 1.13.0 | -| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | +## Quick start -> **Note:** f2py is part of NumPy. -> On Python 3.12 it uses the Meson backend, which is why Meson and Ninja are required. +From a PRIK checkout with PRIK installed, GNU Fortran on `PATH`, and the pinned +NumPy, Meson, and Ninja (see [Set up a clean environment](#set-up-a-clean-environment)): -## Tested platforms +```bash +source examples/fortran/blas/build_all.sh +python3 -m pytest -q examples/fortran/blas/tests +``` -The Real Libraries Portability workflow builds and runs this example with -Python 3.12 on: +The first command builds both wrappers and puts them on `PYTHONPATH` for this +shell; use `source`, not `bash`, so that setting survives. The second runs the +complete comparison. -| Operating system | Architectures | Native toolchain | -| --- | --- | --- | -| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | -| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | +After this, both modules import in the same shell: -The ordinary numerical suite runs on all four targets. The maintainer -full-surface audit also runs on Linux x86-64. +```python +import numpy as np +import f2py_reference_blas +import prik_reference_blas -For everyday use of this example, prefer the checked-in sources in -`examples/fortran/blas/native/`. -(See the [Source provenance](#source-provenance) section at the end if you want to verify the upstream archive yourself.) +x = np.array([2.0, -4.0, 1.0]) +y = np.array([3.0, 5.0, -2.0]) +prik_reference_blas.daxpy(np.int32(3), np.float64(-1.5), x, np.int32(1), y, np.int32(1)) +print(y) # [ 0. 11. -3.5] +``` --- -## 1. Prepare the repository and toolchain +## Key files -Clone PRIK, create a virtual environment, and install the same Python build -tools used by the dedicated CI job: +Everything lives under +[`examples/fortran/blas/`](../../../../examples/fortran/blas/): -```bash -git clone https://github.com/PyNumLab/prik.git -cd prik -python3 -m venv .venv -. .venv/bin/activate -python3 -m pip install --upgrade pip -python3 -m pip install -e ".[qa]" \ - "numpy==2.5.1" "meson==1.11.2" "ninja==1.13.0" -``` +| File | What it does | +| --- | --- | +| [`native/`](../../../../examples/fortran/blas/native/) | The 155 Reference BLAS sources from Netlib LAPACK 3.12.1; the build downloads nothing. | +| [`build_prik.sh`](../../../../examples/fortran/blas/build_prik.sh) | Compiles BLAS into one shared library and builds the PRIK wrapper against it. | +| [`blas.pyf`](../../../../examples/fortran/blas/blas.pyf) | The reviewed f2py signature file for all 155 routines. | +| [`build_f2py.sh`](../../../../examples/fortran/blas/build_f2py.sh) | Builds the f2py wrapper from `blas.pyf`, linked to the same library. | +| [`build_all.sh`](../../../../examples/fortran/blas/build_all.sh) | Runs both build scripts and adds both modules to `PYTHONPATH`. | +| [`routine_inventory.py`](../../../../examples/fortran/blas/routine_inventory.py) | The authoritative list of routines, grouped by BLAS level and kind. | +| [`tests/`](../../../../examples/fortran/blas/tests/) | One test file per routine family, plus [`test_routine_coverage.py`](../../../../examples/fortran/blas/tests/test_routine_coverage.py), which fails if a routine is missing from the sources, inventory, exports, or tests. | +| [`tests/helpers.py`](../../../../examples/fortran/blas/tests/helpers.py) | The small comparison helpers used by every test. | -Install GNU Fortran separately. On Ubuntu: +--- -```bash -sudo apt-get update -sudo apt-get install --yes gfortran -gfortran --version -``` +## How the build works -All remaining commands run from the repository root with the virtual -environment active. +BLAS is compiled once. Both wrappers link that one shared library: -The runnable material is self-contained in the repository's -[`examples/` directory](../../../../examples/). After PRIK and the listed tools -are installed, you can copy that directory alone. +```text +155 BLAS sources ──compile once──> libprik_full_blas + │ + PRIK API, read from the sources ─────────┼──> prik_reference_blas + f2py API, from blas.pyf ─────────────────┴──> f2py_reference_blas +``` ---- +1. `examples.native_library` compiles the 155 sources into + `libprik_full_blas` and prints its path. +2. PRIK reads the same sources to generate its Python API, compiles no BLAS + source itself (`--no-compile-input-sources`), and links the library. +3. f2py builds its wrapper from the reviewed `blas.pyf` and links the same + library. -## 2. Compile BLAS once and build the PRIK wrapper +Both wrappers are built with `-O0`, so the comparison checks correctness +rather than optimization-dependent results. This example is not a performance +benchmark. `build_all.sh` runs these two scripts; each can also be reused on +its own. -Run the first build script from the repository root: +**The PRIK build**, from `build_prik.sh`: ```bash @@ -122,20 +122,7 @@ python -m prik "$EXAMPLE_WORKSPACE/examples/fortran/blas/native" \ --wrapper-c-flags="-O0 -g0" ``` -`examples.native_library` compiles all 155 implementations and returns the -resulting shared-library path. PRIK reads the same source directory to build -the Python API, skips native implementation compilation, and links that -library. - -`-O0` keeps the PRIK and f2py correctness builds equivalent and avoids making -optimization-dependent claims. This example focuses on correctness, not -performance. - ---- - -## 3. Build f2py against the same native library - -Run the same f2py build script exercised by the test suite: +**The f2py build**, from `build_f2py.sh`: ```bash @@ -161,54 +148,28 @@ python -m numpy.f2py -c \ --opt=-O0 ``` -The committed [`blas.pyf`](../../../../examples/fortran/blas/blas.pyf) defines the f2py -interface. f2py compiles only its wrapper and links it to -`BLAS_SHARED_LIBRARY`, so both wrappers exercise the same compiled BLAS -implementations. - -A few routines expose scalar writebacks differently through the two wrappers; -the comparison below shows those return-value differences explicitly. - -Import both modules: - -```python -import os -import sys - -build_root = os.environ["BLAS_BUILD_ROOT"] -sys.path.insert(0, f"{build_root}/prik") -sys.path.insert(0, f"{build_root}/f2py") - -import f2py_reference_blas -import prik_reference_blas -``` - -### Important wrapper difference - -PRIK deliberately follows the native scalar contract: - -- For a subroutine such as `DAXPY`, arrays are mutated in place and PRIK returns the visible input scalars because they are treated as inout. -- The f2py comparison module also mutates the output array but returns `None`. -- Function routines such as `DDOT` return their numerical result through both wrappers. +Everything is written to the temporary `BLAS_BUILD_ROOT` directory, not to the +repository. --- -## 4. Run the complete test suite +## PRIK and f2py differences -Build both wrappers and run the 155-routine suite: +Both wrappers mutate output arrays in place and compute identical results. +They differ in what a call returns: -```bash -source examples/fortran/blas/build_all.sh -python3 -m pytest -q examples/fortran/blas/tests -``` +| Routine | PRIK returns | f2py returns | +| --- | --- | --- | +| A subroutine such as `daxpy` | The visible scalar arguments, which it treats as in-out: `(n, alpha, incx, incy)` | `None` | +| A function such as `ddot` | The result and the visible scalars: `(result, n, incx, incy)` | The result alone | +| The 6 rotation routines, whose scalars have no Fortran `intent` | The scalar writebacks directly | Typed NumPy 0-D arrays, because `blas.pyf` records those scalars as `intent(inout)` | -The tests cover vector, matrix, packed, banded, symmetric, Hermitian, and -triangular operations. Each routine is called with representative inputs and -checked against an independent mathematical result. +PRIK follows the native argument list and scalar contract; the f2py signature +file is the reviewed comparison interface. --- -## 5. See how results are validated +## How results are validated Each comparison checks three relationships: @@ -223,10 +184,8 @@ leading dimensions, and unused storage where they are part of a routine's contract. The independent formula or residual remains the primary numerical reference. -The two examples below come directly from the runnable suite and use its small -NumPy comparison helpers. - -#### Test helper conventions +The two examples below come directly from the runnable suite and use its +helpers: - `assert_allclose_for_dtype` compares floating-point results with a tolerance matched to their NumPy dtype. Its optional `operation_size` is a @@ -259,8 +218,8 @@ def test_daxpy(prik_blas, f2py_blas): assert_storage_unchanged(f2py_x, x) ``` -Both wrappers must mutate `y` to the expected value. -The input-only array `x` must remain unchanged. +Both wrappers must mutate `y` to the expected value. The input-only array `x` +must remain unchanged. ### DDOT – scalar function result @@ -288,31 +247,85 @@ def test_ddot(prik_blas, f2py_blas): --- -## 6. Run focused examples +## Run the tests -After building the wrappers, run a family or one routine: +Run the complete suite, one family, one routine, or every test that mentions a +routine name: ```bash +python3 -m pytest -q examples/fortran/blas/tests python3 -m pytest -q examples/fortran/blas/tests/test_level1_real.py python3 -m pytest -q examples/fortran/blas/tests/test_level1_real.py::test_daxpy python3 -m pytest -q examples/fortran/blas/tests -k dgemm ``` -- Complete Level-1 examples → [`test_level1_real.py`](../../../../examples/fortran/blas/tests/test_level1_real.py) -- Matrix / packed / banded / symmetric / Hermitian / triangular examples → files under [`examples/fortran/blas/tests/`](../../../../examples/fortran/blas/tests/) -- Public routine list → [`routine_inventory.py`](../../../../examples/fortran/blas/routine_inventory.py) -- Routine coverage check → [`test_routine_coverage.py`](../../../../examples/fortran/blas/tests/test_routine_coverage.py) - -For the copyable build scripts, test commands, and source provenance, see the -[`examples/fortran/blas` project README](../../../../examples/fortran/blas/README.md). +The tests cover vector, matrix, packed, banded, symmetric, Hermitian, and +triangular operations, each called with representative inputs and checked +against an independent result. --- +## Set up a clean environment + +Clone PRIK, create a virtual environment, and install the same Python build +tools used by the dedicated CI job: + +```bash +git clone https://github.com/PyNumLab/prik.git +cd prik +python3 -m venv .venv +. .venv/bin/activate +python3 -m pip install --upgrade pip +python3 -m pip install -e ".[qa]" \ + "numpy==2.5.1" "meson==1.11.2" "ninja==1.13.0" +``` + +Install GNU Fortran separately. On Ubuntu: + +```bash +sudo apt-get update +sudo apt-get install --yes gfortran +gfortran --version +``` + +Run the example's commands from the repository root with the virtual +environment active. + +## Versions used + +| Component | Version / source | +| --- | --- | +| PRIK | current repository checkout | +| Reference BLAS | snapshot shipped in Netlib LAPACK 3.12.1 | +| Python | 3.12 in the dedicated CI job | +| NumPy / f2py | NumPy 2.5.1 | +| Meson | 1.11.2 | +| Ninja | 1.13.0 | +| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | + +f2py is part of NumPy. On Python 3.12 it uses the Meson backend, which is why +Meson and Ninja are required. + +## Tested platforms + +The Real Libraries Portability workflow builds and runs this example with +Python 3.12 on: + +| Operating system | Architectures | Native toolchain | +| --- | --- | --- | +| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | +| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | + +The numerical suite runs on all four targets. A maintainer full-surface audit +also runs on Linux x86-64. + ## Troubleshooting -- Confirm that `gfortran`, `meson` and `ninja` are on your `PATH`. -- On Python 3.12+, do **not** force the old distutils backend of f2py. Use the - pinned Meson and Ninja setup shown above. +- Confirm that `gfortran`, `meson`, and `ninja` are on your `PATH`. +- On Python 3.12 or later, do not force f2py's old distutils backend. Use the + pinned Meson and Ninja shown above. +- Use `source examples/fortran/blas/build_all.sh`; running it with `bash` + starts a child shell, so the exported `PYTHONPATH` is lost. - Run a single failing test with more detail and keep the build directory: ```bash From e3ffd017b24ad4a8591d24e33ac61e31e8600f76 Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 14:06:03 +0100 Subject: [PATCH 2/6] codex: restructure the LAPACK example guide Lead with what the example builds (1,936 wrapped procedures, 127 validated routines), a quick start, and a key-files table; explain the single shared-library build with a diagram before the build scripts; collect the PRIK, f2py, and SciPy calling differences in one table; and move setup, versions, and platforms to the end. Match the one-line changelog style on main for the BLAS and LAPACK guide entries. The LAPACK build and tests were not run locally, per the repository's instruction to leave LAPACK coverage to CI. Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 5 +- docs/user/examples/fortran/lapack-wrapper.md | 284 ++++++++++--------- 2 files changed, 153 insertions(+), 136 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index a04307d87..7701c7d0b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,9 +7,8 @@ release tags add a leading `v` to the package version. ## Unreleased -- The Reference BLAS example guide now starts with a quick start, its key - files, and a short explanation of the shared-library build, and lists the - return-value differences between the PRIK and f2py wrappers in one table. +- Improve the Reference BLAS example guide +- Improve the LAPACK example guide - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/lapack-wrapper.md b/docs/user/examples/fortran/lapack-wrapper.md index a8e9b3007..a5500ac81 100644 --- a/docs/user/examples/fortran/lapack-wrapper.md +++ b/docs/user/examples/fortran/lapack-wrapper.md @@ -9,89 +9,96 @@ publication: reviewed # Build and Validate LAPACK with PRIK -This example builds the complete Reference LAPACK library and wraps it with -PRIK. It validates the 127 double-precision real routines also available -through `scipy.linalg.lapack` in SciPy 1.18.0. +This example wraps the complete Reference LAPACK with PRIK, and a comparison +surface with NumPy's f2py, then checks both against SciPy and independent +mathematical results. LAPACK is compiled once into a shared library that both +wrappers link, so a difference between them comes from the wrapper, never from +the numerics. -### What this example shows +### What you get -- Build PRIK and f2py wrappers against the same compiled LAPACK library. -- Call linear-system, factorization, eigenvalue, and singular-value routines - with NumPy arrays. -- Compare results with SciPy and check solutions, residuals, reconstructions, - and other mathematical properties. +- `prik_reference_lapack_example`: PRIK's wrapper of all 1,936 procedures in + Reference LAPACK's default, non-XBLAS source set. +- `f2py_reference_lapack_example`: the f2py comparison wrapper. +- A named test for each of the 127 double-precision real routines that SciPy + 1.18.0 also exposes: linear systems, least squares, factorizations, + eigenvalue problems, and singular values. Each test checks PRIK, f2py, and + SciPy against an independent mathematical result. -You should already be comfortable with the BLAS wrapper example, NumPy arrays, and basic packaging. +This example builds on the [Reference BLAS example](blas-wrapper.md), which +explains the same shared-library pattern in a smaller setting. --- -## Versions used +## Quick start -| Component | Version / source | -| --- | --- | -| PRIK | current repository checkout | -| Reference LAPACK | Netlib LAPACK 3.12.1 | -| Reference BLAS | BLAS snapshot shipped in LAPACK 3.12.1 | -| Python | 3.12 or newer | -| NumPy / f2py | NumPy 2.5.1 | -| SciPy | exactly 1.18.0 | -| Meson | 1.11.2 | -| Ninja | 1.13.0 | -| Fortran compiler | compatible `gfortran` | +From a PRIK checkout with PRIK installed, GNU Fortran and the LAPACK and BLAS +development libraries on the system, and the pinned NumPy, SciPy, Meson, and +Ninja (see [Set up a clean environment](#set-up-a-clean-environment)): -## Tested platforms - -The Real Libraries Portability workflow builds and runs this example with -Python 3.12 on: +```bash +source examples/fortran/lapack/build_all.sh +python3 -m pytest -q examples/fortran/lapack/tests +``` -| Operating system | Architectures | Native toolchain | -| --- | --- | --- | -| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | -| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | +The first command builds both wrappers and puts them on `PYTHONPATH` for this +shell; use `source`, not `bash`, so that setting survives. The second runs the +127-routine comparison. -The ordinary 127-routine validation suite runs on all four targets. The -maintainer full-surface audit also runs on Linux x86-64. +After this, `prik_reference_lapack_example` and +`f2py_reference_lapack_example` import in the same shell. The +[DGESV example](#dgesv-solve-a-general-linear-system) below shows a complete +call through each wrapper. --- -## 1. Prepare the repository and toolchain +## Key files -Clone PRIK, create a virtual environment, and install the pinned comparison and -build tools: +Everything lives under +[`examples/fortran/lapack/`](../../../../examples/fortran/lapack/): -```bash -git clone https://github.com/PyNumLab/prik.git -cd prik -python3 -m venv .venv -. .venv/bin/activate -python3 -m pip install --upgrade pip -python3 -m pip install -e ".[qa]" \ - "numpy==2.5.1" "scipy==1.18.0" \ - "meson==1.11.2" "ninja==1.13.0" -``` +| File | What it does | +| --- | --- | +| [`native/`](../../../../examples/fortran/lapack/native/) | The Reference LAPACK 3.12.1 source snapshot; the build downloads nothing. | +| [`xblas_sources.txt`](../../../../examples/fortran/lapack/xblas_sources.txt) | The sources left out of the default build because they need the separate XBLAS library. | +| [`support/`](../../../../examples/fortran/lapack/support/) | Two workspace-rounding helpers from upstream `INSTALL/` that the default build needs. | +| [`build_prik.sh`](../../../../examples/fortran/lapack/build_prik.sh) | Compiles LAPACK into one shared library and builds the PRIK wrapper against it. | +| [`lapack.pyf`](../../../../examples/fortran/lapack/lapack.pyf) | The reviewed f2py signature file for the comparison routines and `la_constants`. | +| [`lapack.f2cmap`](../../../../examples/fortran/lapack/lapack.f2cmap) | Tells f2py that `real(wp)` is a C `double`. | +| [`build_f2py.sh`](../../../../examples/fortran/lapack/build_f2py.sh) | Builds the f2py wrapper from `lapack.pyf`, linked to the same library. | +| [`build_all.sh`](../../../../examples/fortran/lapack/build_all.sh) | Runs both build scripts and adds both modules to `PYTHONPATH`. | +| [`routine_inventory.py`](../../../../examples/fortran/lapack/routine_inventory.py) | The reviewed list of the 127 validated routines, grouped by LAPACK family. | +| [`tests/`](../../../../examples/fortran/lapack/tests/) | One test file per routine family, plus [`test_routine_coverage.py`](../../../../examples/fortran/lapack/tests/test_routine_coverage.py), which checks the inventory against the tests. | +| [`tests/helpers.py`](../../../../examples/fortran/lapack/tests/helpers.py) | The comparison helpers used by the tests. | -Install GNU Fortran and the LAPACK and BLAS development packages. On Ubuntu: +--- -```bash -sudo apt-get update -sudo apt-get install --yes gfortran liblapack-dev libblas-dev -gfortran --version -``` +## How the build works -All remaining commands run from the repository root with the virtual -environment active. +LAPACK is compiled once. Both wrappers link that one shared library: -The runnable material is self-contained in the repository's -[`examples/` directory](../../../../examples/). After PRIK and the listed tools -are installed, you can copy that directory alone. +```text +LAPACK and BLAS sources ──compile once──> libprik_full_lapack + │ + PRIK API, read from the sources ─────────────────┼──> prik_reference_lapack_example + f2py API, from lapack.pyf ───────────────────────┴──> f2py_reference_lapack_example +``` ---- +1. `examples.native_library` compiles the bundled LAPACK and BLAS sources into + `libprik_full_lapack`, links the installed LAPACK and BLAS libraries for + support routines outside the bundled set, and keeps the compiler's module + files in `LAPACK_MODULE_DIR`. +2. PRIK reads the same default, non-XBLAS sources to generate its Python API, + compiles no LAPACK source itself (`--no-compile-input-sources`), and links + the library. +3. f2py builds its wrapper from the reviewed `lapack.pyf` and links the same + library. -## 2. Compile LAPACK once and build the PRIK wrapper +Both wrappers are built with `-O0` and use `LAPACK_MODULE_DIR` when they +compile. `build_all.sh` runs these two scripts; each can also be reused on its +own. -Compile the native files once into a shared `.so` file so both wrappers can -reuse it. The native builder links the installed LAPACK and BLAS development -libraries for companion support symbols: +**The PRIK build**, from `build_prik.sh`: ```bash @@ -120,19 +127,7 @@ python -m prik "$LAPACK_SOURCE_ROOT" \ --wrapper-c-flags="-O0 -g0" ``` -PRIK reads the sources to build the Python API, skips native implementation -compilation, and links the shared library. The include path supplies the module -metadata needed by the generated wrapper. - ---- - -## 3. Build the f2py comparison wrapper - -The committed [`lapack.pyf`](../../../../examples/fortran/lapack/lapack.pyf) contains the -125 selected routines and the `la_constants` module signature. f2py compiles -only this wrapper and links `LAPACK_SHARED_LIBRARY`. - -Run the same direct f2py command exercised by the test suite: +**The f2py build**, from `build_f2py.sh`: ```bash @@ -159,66 +154,37 @@ python -m numpy.f2py -c \ --opt=-O0 ``` -`LAPACK_MODULE_DIR` provides the compiler-generated module files needed to -compile each wrapper. Both wrappers link the existing shared library instead of -recompiling LAPACK. - -The comparison excludes `dgees` and `dgges` because f2py 2.5.1 cannot generate -their callback declarations correctly. Those two routines are still checked -through PRIK, SciPy, and their Schur decompositions. - -Import the two built modules and SciPy's LAPACK module from the repository -root: - -```python -import os -import sys - -sys.path.insert(0, f"{os.environ['LAPACK_BUILD_ROOT']}/prik") -sys.path.insert(0, os.environ["LAPACK_F2PY_ROOT"]) - -import f2py_reference_lapack_example -import prik_reference_lapack_example -from scipy.linalg import lapack as scipy_lapack -``` - -### SciPy comparison - -The tests use the 127 double-precision real LAPACK routines available in SciPy -1.18.0 for `np.float64` arrays. Pinning that version keeps the comparison API -and expected results reproducible. +Everything is written to the temporary `LAPACK_BUILD_ROOT` directory, not to +the repository. --- -## 4. Run the complete test suite +## PRIK, f2py, and SciPy differences -Build both wrappers and run all 127 routine tests: +All three compute the same results. They differ in how a call looks: -```bash -source examples/fortran/lapack/build_all.sh -python3 -m pytest -q examples/fortran/lapack/tests -``` +| Case | PRIK | f2py comparison wrapper | SciPy | +| --- | --- | --- | --- | +| A routine such as `dgesv` | Native argument order; updates arrays in place and returns the visible scalars, for example `(2, 1, 2, 2, 0)` | Updates arrays in place and returns `None` | Returns new arrays and `info` | +| Character selectors such as `"L"` in `dpotrf` | Returned with the other scalars, because LAPACK declares no `intent` | Passed as `b"L"` | A keyword such as `lower=1` | +| The 9 routines whose scalar outputs have no Fortran `intent` | Returns the writebacks directly | Typed NumPy 0-D arrays, because `lapack.pyf` records them as `intent(inout)` | — | +| `dgees` and `dgges` | Wrapped and tested | Not in the comparison: f2py 2.5.1 cannot generate their selection callbacks | Used to verify the Schur decompositions | -The suite covers linear systems, least squares, factorizations, eigenvalue -problems, singular values, and related matrix operations. +SciPy reports zero-based pivots, while LAPACK and both wrappers use one-based +pivots. --- -## 5. See how results are validated +## How results are validated LAPACK outputs are not always unique. Eigenvectors and singular vectors may change sign, repeated eigenspaces may use a different orthonormal basis, and -pivot ties may choose another valid permutation. -Therefore byte-for-byte agreement is not the only oracle. +pivot ties may choose another valid permutation, so byte-for-byte agreement is +not the only oracle. Tests use explicit solutions, residuals, factor reconstructions, orthogonality, -eigen equations, and storage checks. The two reusable checks shown below live -in [`tests/helpers.py`](../../../../examples/fortran/lapack/tests/helpers.py). - -#### Test helper conventions - -The snippets use standard NumPy operations whenever the check is local. The -two helpers in the displayed DPOTRF test keep its repeated checks consistent: +eigen equations, and storage checks. The two helpers used below live in +[`tests/helpers.py`](../../../../examples/fortran/lapack/tests/helpers.py): - `assert_allclose_float64` compares values using a tolerance appropriate for float64 arithmetic. Its `operation_size` argument is a rounding-error scale: @@ -227,8 +193,8 @@ two helpers in the displayed DPOTRF test keep its repeated checks consistent: - `assert_storage_unchanged` compares storage exactly, including `NaN` sentinels in parts of an array LAPACK must not read or overwrite. -The examples below show the PRIK, f2py, and SciPy calls together with a direct -mathematical check. They come from the runnable suite. +The examples below come from the runnable suite and show the PRIK, f2py, and +SciPy calls together with a direct mathematical check. ### DGESV – solve a general linear system @@ -307,35 +273,87 @@ The reconstruction `A = L @ L.T` confirms that the factor is correct. --- -## 6. Run focused examples +## Run the tests -After building the wrappers, run a family or one routine: +Run the complete suite, one family, one routine, or every test that mentions a +routine name: ```bash +python3 -m pytest -q examples/fortran/lapack/tests python3 -m pytest -q examples/fortran/lapack/tests/test_linear_general.py python3 -m pytest -q \ examples/fortran/lapack/tests/test_linear_general.py::test_dgesv_solves_general_system python3 -m pytest -q examples/fortran/lapack/tests -k dgesvd ``` -- Full DGESV and related general-system tests → [`test_linear_general.py`](../../../../examples/fortran/lapack/tests/test_linear_general.py) -- Cholesky and other positive-definite examples → [`test_linear_positive_definite.py`](../../../../examples/fortran/lapack/tests/test_linear_positive_definite.py) -- Other families live under [`examples/fortran/lapack/tests/`](../../../../examples/fortran/lapack/tests/) -- Public routine list → [`routine_inventory.py`](../../../../examples/fortran/lapack/routine_inventory.py) -- Routine coverage check → [`test_routine_coverage.py`](../../../../examples/fortran/lapack/tests/test_routine_coverage.py) +--- -For the copyable build scripts, test commands, and source provenance, see the -[`examples/fortran/lapack` project README](../../../../examples/fortran/lapack/README.md). +## Set up a clean environment ---- +Clone PRIK, create a virtual environment, and install the pinned comparison and +build tools: + +```bash +git clone https://github.com/PyNumLab/prik.git +cd prik +python3 -m venv .venv +. .venv/bin/activate +python3 -m pip install --upgrade pip +python3 -m pip install -e ".[qa]" \ + "numpy==2.5.1" "scipy==1.18.0" \ + "meson==1.11.2" "ninja==1.13.0" +``` + +Install GNU Fortran and the LAPACK and BLAS development packages. On Ubuntu: + +```bash +sudo apt-get update +sudo apt-get install --yes gfortran liblapack-dev libblas-dev +gfortran --version +``` + +Run the example's commands from the repository root with the virtual +environment active. + +## Versions used + +| Component | Version / source | +| --- | --- | +| PRIK | current repository checkout | +| Reference LAPACK | Netlib LAPACK 3.12.1 | +| Reference BLAS | BLAS snapshot shipped in LAPACK 3.12.1 | +| Python | 3.12 or newer | +| NumPy / f2py | NumPy 2.5.1 | +| SciPy | exactly 1.18.0 | +| Meson | 1.11.2 | +| Ninja | 1.13.0 | +| Fortran compiler | compatible `gfortran` | + +SciPy is pinned to exactly 1.18.0 so its low-level comparison API and expected +results stay reproducible. + +## Tested platforms + +The Real Libraries Portability workflow builds and runs this example with +Python 3.12 on: + +| Operating system | Architectures | Native toolchain | +| --- | --- | --- | +| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | +| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | + +The 127-routine validation suite runs on all four targets. A maintainer +full-surface audit also runs on Linux x86-64. ## Troubleshooting -- Confirm that `gfortran`, `ar`, `meson` and `ninja` are on `PATH`. +- Confirm that `gfortran`, `ar`, `meson`, and `ninja` are on `PATH`. - Keep SciPy at **exactly 1.18.0** so its low-level comparison API matches this example. - On Python 3.12 or newer, let f2py use Meson; do not force the removed distutils backend. +- Use `source examples/fortran/lapack/build_all.sh`; running it with `bash` + starts a child shell, so the exported `PYTHONPATH` is lost. - Rerun one named test with more detail and keep the build directory: ```bash From 3326edca649bdb5914893d78615799d3fac6aa91 Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 14:14:03 +0100 Subject: [PATCH 3/6] codex: restructure the MINPACK example guide Lead with what the example builds and the 22 procedures by family, a quick start that shows the minpack_module import, and a key-files table; explain the single-command build before the script; add verified hybrd1 and lmdif1 examples with work-array sizes and info meanings; and correct the claim that results are checked against SciPy, which the suite does not use. Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 1 + docs/user/examples/fortran/minpack-wrapper.md | 263 ++++++++++-------- 2 files changed, 153 insertions(+), 111 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 7701c7d0b..cb7c913ce 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,6 +9,7 @@ release tags add a leading `v` to the package version. - Improve the Reference BLAS example guide - Improve the LAPACK example guide +- Improve the MINPACK example guide - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/minpack-wrapper.md b/docs/user/examples/fortran/minpack-wrapper.md index 3a4d98bc5..0bd6ca4e9 100644 --- a/docs/user/examples/fortran/minpack-wrapper.md +++ b/docs/user/examples/fortran/minpack-wrapper.md @@ -9,83 +9,78 @@ publication: reviewed # Build and Validate MINPACK with PRIK -This example takes the checked-in -[fortran-lang/minpack](https://github.com/fortran-lang/minpack) source and -builds an importable Python extension containing all 22 public MINPACK -procedures. +This example turns [fortran-lang/minpack](https://github.com/fortran-lang/minpack), +the modern Fortran MINPACK, into one Python extension with all 22 public +procedures. MINPACK's solvers call back into Python for each residual, so you +write your problem as an ordinary Python function. The suite checks every +procedure against exact solutions and direct linear-algebra identities. -The example solves known nonlinear and least-squares problems and checks their -results with exact solutions and direct linear-algebra identities. +### What you get -### What this example shows +One extension, `prik_reference_minpack`, whose Fortran module `minpack_module` +holds all 22 procedures: -- Wrap a complete numerical solver library as one Python extension. -- Pass NumPy arrays and ordinary Python functions to MINPACK routines. -- Check root-finding, least-squares, Jacobian, and factorization results. - -You should already be comfortable with NumPy arrays, Python callables, and -building a local Fortran extension. +| Family | Procedures | +| --- | --- | +| Hybrid nonlinear solvers (root finding) | `hybrd`, `hybrd1`, `hybrj`, `hybrj1` | +| Levenberg-Marquardt solvers (least squares) | `lmder`, `lmder1`, `lmdif`, `lmdif1`, `lmstr`, `lmstr1` | +| Diagnostics and finite differences | `chkder`, `enorm`, `fdjac1`, `fdjac2` | +| Factorization and update helpers | `dogleg`, `lmpar`, `qform`, `qrfac`, `qrsolv`, `r1mpyq`, `r1updt`, `rwupdt` | --- -## Versions used - -| Component | Version / source | -| --- | --- | -| PRIK | current repository checkout | -| MINPACK | [fortran-lang/minpack commit `c0b5aea`](https://github.com/fortran-lang/minpack/tree/c0b5aea9fcd2b83865af921a7a7e881904f8d3c2) | -| Python | 3.12 in the dedicated CI job | -| NumPy | 2.5.1 | -| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | +## Quick start -The repository owns the checked-in source snapshot under -`examples/fortran/minpack/native/`, so the example does not download code during its -build. +From a PRIK checkout with PRIK installed and GNU Fortran on `PATH` (see +[Set up a clean environment](#set-up-a-clean-environment)): -## Tested platforms +```bash +source examples/fortran/minpack/build_all.sh +python3 -m pytest -q examples/fortran/minpack/tests +``` -The Real Libraries Portability workflow builds and runs the complete numerical -suite with Python 3.12 on: +The first command builds the extension and puts it on `PYTHONPATH` for this +shell; use `source`, not `bash`, so that setting survives. The second runs the +tests. -| Operating system | Architectures | Native toolchain | -| --- | --- | --- | -| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | -| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | +After this, the solvers import in the same shell: ---- +```python +from prik_reference_minpack import minpack_module as minpack +``` -## 1. Prepare the repository and toolchain +[Use the generated API](#use-the-generated-api) shows complete calls. -Clone PRIK, create a virtual environment, and install the Python tools used by -the dedicated CI job: - -```bash -git clone https://github.com/PyNumLab/prik.git -cd prik -python3 -m venv .venv -. .venv/bin/activate -python3 -m pip install --upgrade pip -python3 -m pip install -e ".[qa]" "numpy==2.5.1" -``` +--- -Install GNU Fortran separately. On Ubuntu: +## Key files -```bash -sudo apt-get update -sudo apt-get install --yes gfortran -gfortran --version -``` +Everything lives under [`examples/fortran/minpack/`](../../../../examples/fortran/minpack/): -All remaining commands run from the repository root with the virtual -environment active. The complete runnable project lives under -[`examples/fortran/minpack/`](../../../../examples/fortran/minpack/). +| File | What it does | +| --- | --- | +| [`native/minpack.f90`](../../../../examples/fortran/minpack/native/minpack.f90) | The MINPACK source: one file holding the public module and its implementation. | +| [`build_prik.sh`](../../../../examples/fortran/minpack/build_prik.sh) | Builds the extension with one PRIK command. | +| [`build_all.sh`](../../../../examples/fortran/minpack/build_all.sh) | Runs `build_prik.sh` and adds the extension to `PYTHONPATH`. | +| [`routine_inventory.py`](../../../../examples/fortran/minpack/routine_inventory.py) | The list of the 22 procedures, grouped by family. | +| [`tests/test_solvers.py`](../../../../examples/fortran/minpack/tests/test_solvers.py) | The root-finding and least-squares solvers, with Python callbacks. | +| [`tests/test_diagnostics.py`](../../../../examples/fortran/minpack/tests/test_diagnostics.py) | The diagnostics and finite-difference helpers. | +| [`tests/test_linear_algebra.py`](../../../../examples/fortran/minpack/tests/test_linear_algebra.py) | The factorization and update helpers. | +| [`tests/test_routine_coverage.py`](../../../../examples/fortran/minpack/tests/test_routine_coverage.py) | Checks that the inventory, the generated exports, and the tests stay in sync. | --- -## 2. Build the PRIK wrapper +## How the build works MINPACK keeps its public declarations and implementations in one source file, -so one command can generate the wrapper and compile the library: +so one PRIK command generates the wrapper and compiles MINPACK into the same +extension. There is no separate native library: + +```text +native/minpack.f90 ──prik──> wrapper + MINPACK, compiled together ──> prik_reference_minpack +``` + +`build_prik.sh` runs that command: ```bash @@ -104,58 +99,68 @@ python3 -m prik "$EXAMPLE_WORKSPACE/examples/fortran/minpack/native/minpack.f90" --wrapper-c-flags="-O0 -g0" ``` -The example uses `-O0` so the tests focus on correct results. PRIK compiles the -native source and generated bridge into one extension. +The example uses `-O0` so the tests focus on correct results rather than +optimization-dependent ones. Everything is written to the temporary +`MINPACK_BUILD_ROOT` directory, not to the repository. -For normal use, source the convenience entrypoint: +--- -```bash -source examples/fortran/minpack/build_all.sh -``` +## Use the generated API -It builds the extension and exports its directory on `PYTHONPATH` for the -current shell. +MINPACK routines keep their documented Fortran argument order, including work +arrays and their lengths. The callback receives the problem sizes, the current +point, and an output array it fills with the residuals. Each solver returns +MINPACK's `info` status and updates `x` in place. ---- +**Root finding.** `hybrd1` finds where the circle x² + y² = 4 meets the curve +y = x³: -## 3. Use the generated Python API +```python +import numpy as np +from prik_reference_minpack import minpack_module as minpack -MINPACK routines keep their documented argument order, including work arrays -and status values. Pass NumPy arrays with the generated dtype, shape, and -layout. Solver callbacks are ordinary Python functions with the generated -callback signature. +def equations(n, x, fvec, iflag): + fvec[0] = x[0] ** 2 + x[1] ** 2 - 4.0 + fvec[1] = x[1] - x[0] ** 3 ---- +x = np.array([1.0, 1.0]) +fvec = np.empty(2) +info = minpack.hybrd1(equations, np.int32(2), x, fvec, np.float64(1e-10), np.empty(19), np.int32(19)) +print(info, x.round(6)) # 1 [1.174222 1.619013] +``` -## 4. Run the complete test suite +The work array needs at least n(3n + 13)/2 entries, 19 for two unknowns. +`info == 1` means MINPACK estimates the relative error in `x` is within the +tolerance. -After the build finishes, run: +**Least squares.** `lmdif1` fits y = a·exp(b·t) to five points: -```bash -python3 -m pytest -q examples/fortran/minpack/tests -``` +```python +t = np.array([0.0, 1.0, 2.0, 3.0, 4.0]) +y = 2.0 * np.exp(-0.5 * t) -The tests cover all 22 public procedures: +def residuals(m, n, p, fvec, iflag): + fvec[:] = p[0] * np.exp(p[1] * t) - y -| Family | Procedures | -| --- | ---: | -| Diagnostics and finite differences | 4 | -| Hybrid nonlinear solvers | 4 | -| Levenberg-Marquardt solvers | 6 | -| Factorization and update helpers | 8 | -| **Total** | **22** | +p = np.array([1.0, 0.0]) +fvec = np.empty(5) +info = minpack.lmdif1(residuals, np.int32(5), np.int32(2), p, fvec, np.float64(1e-10), + np.empty(2, dtype=np.int32), np.empty(25), np.int32(25)) +print(info, p.round(6)) # 2 [ 2. -0.5] +``` -Each procedure is called with representative data and checked against SciPy, a -known solution, or a direct linear-algebra result. +Here `m = 5` residuals fit `n = 2` parameters; the work array needs at least +m·n + 5n + m entries, 25 in this case. `lmdif1` estimates the Jacobian by +finite differences; `lmder1` and `hybrj1` take a callback that also supplies +it. --- -## 5. See how results are validated +## How results are validated -For example, `hybrd1` can solve the two-variable equation -`x - [1, -2] = 0`. MINPACK calls the Python function whenever it needs the -current residual. The example below is the runnable `hybrd1` test; its -`minpack` fixture supplies the generated module: +Each procedure is called with representative data and checked against a known +solution or a direct linear-algebra result. The runnable `hybrd1` test solves +`x - [1, -2] = 0`; its `minpack` fixture supplies `minpack_module`: ```python @@ -186,42 +191,78 @@ def test_hybrd1(minpack): np.testing.assert_allclose(fvec, 0.0, atol=1.0e-10) ``` -The complete suite applies the same pattern to root-finding and least-squares -solvers, then checks their solutions and final residuals against the declared -problem. +The complete suite applies the same pattern to the other root-finding and +least-squares solvers, and checks the helpers with algebraic invariants. It +also verifies callback counts, caller-array writebacks, and Fortran-order +matrices. --- -## 6. Run focused examples +## Run the tests -After building the extension, run a family or one routine: +Run the complete suite, one family, or one routine: ```bash +python3 -m pytest -q examples/fortran/minpack/tests python3 -m pytest -q examples/fortran/minpack/tests/test_solvers.py python3 -m pytest -q \ examples/fortran/minpack/tests/test_solvers.py::test_hybrd1 ``` -- Callback-driven nonlinear solvers → - [`test_solvers.py`](../../../../examples/fortran/minpack/tests/test_solvers.py) -- Diagnostics and finite-difference helpers → - [`test_diagnostics.py`](../../../../examples/fortran/minpack/tests/test_diagnostics.py) -- Factorization and update helpers → - [`test_linear_algebra.py`](../../../../examples/fortran/minpack/tests/test_linear_algebra.py) -- Public routine list → - [`routine_inventory.py`](../../../../examples/fortran/minpack/routine_inventory.py) -- Routine coverage check → - [`test_routine_coverage.py`](../../../../examples/fortran/minpack/tests/test_routine_coverage.py) -- Copyable project instructions → - [`examples/fortran/minpack/README.md`](../../../../examples/fortran/minpack/README.md) - --- +## Set up a clean environment + +Clone PRIK, create a virtual environment, and install the Python tools used by +the dedicated CI job: + +```bash +git clone https://github.com/PyNumLab/prik.git +cd prik +python3 -m venv .venv +. .venv/bin/activate +python3 -m pip install --upgrade pip +python3 -m pip install -e ".[qa]" "numpy==2.5.1" +``` + +Install GNU Fortran separately. On Ubuntu: + +```bash +sudo apt-get update +sudo apt-get install --yes gfortran +gfortran --version +``` + +Run the example's commands from the repository root with the virtual +environment active. + +## Versions used + +| Component | Version / source | +| --- | --- | +| PRIK | current repository checkout | +| MINPACK | [fortran-lang/minpack commit `c0b5aea`](https://github.com/fortran-lang/minpack/tree/c0b5aea9fcd2b83865af921a7a7e881904f8d3c2) | +| Python | 3.12 in the dedicated CI job | +| NumPy | 2.5.1 | +| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | + +## Tested platforms + +The Real Libraries Portability workflow builds and runs the complete numerical +suite with Python 3.12 on: + +| Operating system | Architectures | Native toolchain | +| --- | --- | --- | +| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | +| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | + ## Troubleshooting - Confirm that `gfortran` is available on `PATH`. -- Use `source examples/fortran/minpack/build_all.sh`; executing it in a child shell - does not preserve the exported `PYTHONPATH`. +- Use `source examples/fortran/minpack/build_all.sh`; running it with `bash` + starts a child shell, so the exported `PYTHONPATH` is lost. +- `AttributeError` on a routine such as `hybrd1`: import it from + `prik_reference_minpack.minpack_module`, not from the extension's top level. - Start with one helper or solver test and add `-vv -s` when diagnosing a callback or generated-wrapper failure. From b9bcc998c344bf14dbfd1105ce13c9e18b025e1d Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 14:25:54 +0100 Subject: [PATCH 4/6] codex: compare MINPACK's SciPy-exposed solvers with SciPy SciPy's root(method="hybr") and least_squares(method="lm") are built on MINPACK, so hybrd, hybrd1, hybrj, hybrj1, lmdif, lmdif1, lmder, and lmder1 can be cross-checked through SciPy's public API; the other 14 procedures have no public SciPy counterpart. One parametrized test solves a nonlinear root-finding system and a nonlinear curve fit with each solver, checks PRIK's answer against the known solution first, then checks agreement with SciPy. It skips without SciPy; CI already installs scipy==1.18.0 in the MINPACK job. A perturbed-data control confirmed the agreement check fails when the two differ. Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 1 + docs/user/examples/fortran/minpack-wrapper.md | 15 + examples/fortran/minpack/README.md | 6 + .../minpack/tests/test_scipy_comparison.py | 280 ++++++++++++++++++ 4 files changed, 302 insertions(+) create mode 100644 examples/fortran/minpack/tests/test_scipy_comparison.py diff --git a/CHANGELOG.md b/CHANGELOG.md index cb7c913ce..c0d556e39 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -10,6 +10,7 @@ release tags add a leading `v` to the package version. - Improve the Reference BLAS example guide - Improve the LAPACK example guide - Improve the MINPACK example guide +- Compare the MINPACK example's eight SciPy-exposed solvers with SciPy - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/minpack-wrapper.md b/docs/user/examples/fortran/minpack-wrapper.md index 0bd6ca4e9..d9e55b597 100644 --- a/docs/user/examples/fortran/minpack-wrapper.md +++ b/docs/user/examples/fortran/minpack-wrapper.md @@ -66,6 +66,7 @@ Everything lives under [`examples/fortran/minpack/`](../../../../examples/fortra | [`tests/test_solvers.py`](../../../../examples/fortran/minpack/tests/test_solvers.py) | The root-finding and least-squares solvers, with Python callbacks. | | [`tests/test_diagnostics.py`](../../../../examples/fortran/minpack/tests/test_diagnostics.py) | The diagnostics and finite-difference helpers. | | [`tests/test_linear_algebra.py`](../../../../examples/fortran/minpack/tests/test_linear_algebra.py) | The factorization and update helpers. | +| [`tests/test_scipy_comparison.py`](../../../../examples/fortran/minpack/tests/test_scipy_comparison.py) | Compares the eight solvers SciPy exposes with SciPy's MINPACK-based solvers. | | [`tests/test_routine_coverage.py`](../../../../examples/fortran/minpack/tests/test_routine_coverage.py) | Checks that the inventory, the generated exports, and the tests stay in sync. | --- @@ -196,6 +197,20 @@ least-squares solvers, and checks the helpers with algebraic invariants. It also verifies callback counts, caller-array writebacks, and Fortran-order matrices. +**SciPy cross-check.** SciPy's `root(method="hybr")` and +`least_squares(method="lm")` are built on MINPACK, so the eight solvers they +cover (`hybrd`, `hybrd1`, `hybrj`, `hybrj1`, `lmdif`, `lmdif1`, `lmder`, +`lmder1`) are also compared with SciPy on the two nonlinear problems from +[Use the generated API](#use-the-generated-api). Each case first checks PRIK's +answer independently, then that PRIK and SciPy agree. The other 14 procedures +have no public SciPy counterpart. The comparison skips if SciPy is not +installed: + +```bash +python3 -m pip install "scipy==1.18.0" +python3 -m pytest -q examples/fortran/minpack/tests/test_scipy_comparison.py +``` + --- ## Run the tests diff --git a/examples/fortran/minpack/README.md b/examples/fortran/minpack/README.md index 416fc0bea..4d2d0c2bc 100644 --- a/examples/fortran/minpack/README.md +++ b/examples/fortran/minpack/README.md @@ -83,6 +83,12 @@ caller-array writebacks, and Fortran-order matrices. The public routine list stays in sync with the generated exports, and every public procedure is exercised. +The eight solvers SciPy exposes (`hybrd`, `hybrd1`, `hybrj`, `hybrj1`, +`lmdif`, `lmdif1`, `lmder`, `lmder1`) are also compared with SciPy's +MINPACK-based `root(method="hybr")` and `least_squares(method="lm")` on +nonlinear problems in [`tests/test_scipy_comparison.py`](tests/test_scipy_comparison.py). +That comparison skips if SciPy is not installed. + ## Sources and license [`native/minpack.f90`](native/minpack.f90) matches upstream `src/minpack.f90` diff --git a/examples/fortran/minpack/tests/test_scipy_comparison.py b/examples/fortran/minpack/tests/test_scipy_comparison.py new file mode 100644 index 000000000..960e5479c --- /dev/null +++ b/examples/fortran/minpack/tests/test_scipy_comparison.py @@ -0,0 +1,280 @@ +"""MINPACK solvers agree with SciPy's MINPACK-based solvers on nonlinear problems. + +SciPy exposes eight of the 22 procedures through public solvers: the hybrid +root finders through ``root(method="hybr")`` and the Levenberg-Marquardt +solvers through ``least_squares(method="lm")``. Each case solves a genuinely +nonlinear problem with both, checks PRIK's answer independently, then checks +that the two interfaces agree. The known answer stays the primary check. +""" + +from __future__ import annotations + +import numpy as np +import pytest + +optimize = pytest.importorskip("scipy.optimize") + +pytestmark = [pytest.mark.fortran_end_to_end, pytest.mark.real_library] + +TOLERANCE = np.float64(1.0e-12) +MAX_EVALUATIONS = np.int32(1000) +FACTOR = np.float64(100.0) +ZERO = np.float64(0.0) +N = np.int32(2) + +# Root finding: where the circle x^2 + y^2 = 4 meets the curve y = x^3. +ROOT_START = np.array([1.0, 1.0]) + + +def _system(x): + return np.array([x[0] ** 2 + x[1] ** 2 - 4.0, x[1] - x[0] ** 3]) + + +def _system_jacobian(x): + return np.array([[2.0 * x[0], 2.0 * x[1]], [-3.0 * x[0] ** 2, 1.0]]) + + +# Least squares: recover a = 2 and b = -0.5 from y = a * exp(b * t). +T = np.array([0.0, 1.0, 2.0, 3.0, 4.0]) +Y = 2.0 * np.exp(-0.5 * T) +M = np.int32(T.size) +FIT_START = np.array([1.0, 0.0]) +FIT_ANSWER = np.array([2.0, -0.5]) + + +def _fit_residuals(p): + return p[0] * np.exp(p[1] * T) - Y + + +def _fit_jacobian(p): + growth = np.exp(p[1] * T) + return np.column_stack([growth, p[0] * T * growth]) + + +def _hybrd(minpack): + x = ROOT_START.copy() + fvec = np.empty(2) + + def fcn(_n, x, fvec, _iflag): + fvec[:] = _system(x) + + info, _nfev = minpack.hybrd( + fcn, + N, + x, + fvec, + TOLERANCE, + MAX_EVALUATIONS, + np.int32(1), + np.int32(1), + ZERO, + np.ones(2), + np.int32(1), + FACTOR, + np.int32(0), + np.empty((2, 2), order="F"), + N, + np.empty(3), + np.int32(3), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + ) + return info, x + + +def _hybrd1(minpack): + x = ROOT_START.copy() + + def fcn(_n, x, fvec, _iflag): + fvec[:] = _system(x) + + info = minpack.hybrd1(fcn, N, x, np.empty(2), TOLERANCE, np.empty(19), np.int32(19)) + return info, x + + +def _root_with_jacobian(_n, x, fvec, fjac, _ldfjac, iflag): + if iflag == 1: + fvec[:] = _system(x) + elif iflag == 2: + fjac[:, :] = _system_jacobian(x) + + +def _hybrj(minpack): + x = ROOT_START.copy() + info, _nfev, _njev = minpack.hybrj( + _root_with_jacobian, + N, + x, + np.empty(2), + np.empty((2, 2), order="F"), + N, + TOLERANCE, + MAX_EVALUATIONS, + np.ones(2), + np.int32(1), + FACTOR, + np.int32(0), + np.empty(3), + np.int32(3), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + ) + return info, x + + +def _hybrj1(minpack): + x = ROOT_START.copy() + info = minpack.hybrj1( + _root_with_jacobian, + N, + x, + np.empty(2), + np.empty((2, 2), order="F"), + N, + TOLERANCE, + np.empty(15), + np.int32(15), + ) + return info, x + + +def _fit(_m, _n, p, fvec, _iflag): + fvec[:] = _fit_residuals(p) + + +def _fit_with_jacobian(_m, _n, p, fvec, fjac, _ldfjac, iflag): + if iflag == 1: + fvec[:] = _fit_residuals(p) + elif iflag == 2: + fjac[:, :] = _fit_jacobian(p) + + +def _lmdif(minpack): + p = FIT_START.copy() + info, _nfev = minpack.lmdif( + _fit, + M, + N, + p, + np.empty(M), + TOLERANCE, + TOLERANCE, + ZERO, + MAX_EVALUATIONS, + ZERO, + np.ones(2), + np.int32(1), + FACTOR, + np.int32(0), + np.empty((M, 2), order="F"), + M, + np.empty(2, dtype=np.int32), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(M), + ) + return info, p + + +def _lmdif1(minpack): + p = FIT_START.copy() + info = minpack.lmdif1( + _fit, M, N, p, np.empty(M), TOLERANCE, np.empty(2, dtype=np.int32), np.empty(25), np.int32(25) + ) + return info, p + + +def _lmder(minpack): + p = FIT_START.copy() + info, _nfev, _njev = minpack.lmder( + _fit_with_jacobian, + M, + N, + p, + np.empty(M), + np.empty((M, 2), order="F"), + M, + TOLERANCE, + TOLERANCE, + ZERO, + MAX_EVALUATIONS, + np.ones(2), + np.int32(1), + FACTOR, + np.int32(0), + np.empty(2, dtype=np.int32), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(2), + np.empty(M), + ) + return info, p + + +def _lmder1(minpack): + p = FIT_START.copy() + info = minpack.lmder1( + _fit_with_jacobian, + M, + N, + p, + np.empty(M), + np.empty((M, 2), order="F"), + M, + TOLERANCE, + np.empty(2, dtype=np.int32), + np.empty(15), + np.int32(15), + ) + return info, p + + +def _scipy_root(jacobian): + result = optimize.root(_system, ROOT_START, jac=_system_jacobian if jacobian else None, method="hybr") + assert result.success + return result.x + + +def _scipy_fit(jacobian): + result = optimize.least_squares( + _fit_residuals, FIT_START, jac=_fit_jacobian if jacobian else "2-point", method="lm" + ) + assert result.success + return result.x + + +@pytest.mark.parametrize( + ("solve", "scipy_jacobian", "problem"), + [ + pytest.param(_hybrd, False, "root", id="hybrd-vs-root-hybr"), + pytest.param(_hybrd1, False, "root", id="hybrd1-vs-root-hybr"), + pytest.param(_hybrj, True, "root", id="hybrj-vs-root-hybr-jac"), + pytest.param(_hybrj1, True, "root", id="hybrj1-vs-root-hybr-jac"), + pytest.param(_lmdif, False, "fit", id="lmdif-vs-least-squares-lm"), + pytest.param(_lmdif1, False, "fit", id="lmdif1-vs-least-squares-lm"), + pytest.param(_lmder, True, "fit", id="lmder-vs-least-squares-lm-jac"), + pytest.param(_lmder1, True, "fit", id="lmder1-vs-least-squares-lm-jac"), + ], +) +def test_solver_agrees_with_scipy(minpack, solve, scipy_jacobian, problem): + info, x = solve(minpack) + + if problem == "root": + assert info == np.int32(1) + np.testing.assert_allclose(_system(x), 0.0, atol=1.0e-10) + expected = _scipy_root(scipy_jacobian) + else: + assert np.int32(1) <= info <= np.int32(4) + np.testing.assert_allclose(x, FIT_ANSWER, atol=1.0e-8) + expected = _scipy_fit(scipy_jacobian) + + np.testing.assert_allclose(x, expected, atol=1.0e-8) From 2c87e4c82db4cf8fa0693fe8c6345b3bf8bd0262 Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 14:37:24 +0100 Subject: [PATCH 5/6] codex: restructure the FFTPACK example guide Lead with the 31 procedures by family, a quick start with the fftpack import, and a key-files table; explain the three source roles with a diagram before the build script; and add a verified table of where FFTPACK's conventions differ from numpy.fft (integer fftfreq indices, an unnormalized ifft, and packed rfft output). The quick start and every page snippet were run against a fresh build (33 tests passed). Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 1 + docs/user/examples/fortran/fftpack-wrapper.md | 265 +++++++++--------- 2 files changed, 139 insertions(+), 127 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index c0d556e39..830bbc77c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,7 @@ release tags add a leading `v` to the package version. - Improve the LAPACK example guide - Improve the MINPACK example guide - Compare the MINPACK example's eight SciPy-exposed solvers with SciPy +- Improve the FFTPACK example guide - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/fftpack-wrapper.md b/docs/user/examples/fortran/fftpack-wrapper.md index cfa96408f..9834d4b08 100644 --- a/docs/user/examples/fortran/fftpack-wrapper.md +++ b/docs/user/examples/fortran/fftpack-wrapper.md @@ -9,85 +9,89 @@ publication: reviewed # Build and Validate FFTPACK with PRIK -This example takes the checked-in -[fortran-lang/fftpack](https://github.com/fortran-lang/fftpack) sources and -builds an importable Python extension containing all 31 public procedures from -the `fftpack` module. - -The example compares Fourier, cosine, sine, frequency, and spectrum operations -with NumPy, SciPy, or known transform properties. +This example turns [fortran-lang/fftpack](https://github.com/fortran-lang/fftpack), +the modern Fortran FFTPACK, into one Python extension with all 31 public +procedures of its `fftpack` module. It checks every procedure against NumPy, +SciPy, or a known transform property. -### What this example shows +### What you get -- Wrap a complete multi-file Fortran library as one Python extension. -- Call both low-level and high-level transforms with NumPy arrays. -- Check transform values, normalization, frequency ordering, dtype, and shape. +One extension, `prik_reference_fftpack`, whose `fftpack` namespace holds: -You should already be comfortable with NumPy arrays and building a local -Fortran extension. +| Family | Procedures | +| --- | --- | +| High-level Fourier transforms | `fft`, `ifft`, `rfft`, `irfft` | +| High-level cosine transforms | `dct`, `idct`, `dct_t1i`, `dct_t1`, `dct_t23i`, `dct_t2`, `dct_t3` | +| Frequency and spectrum ordering | `fftfreq`, `rfftfreq`, `fftshift`, `ifftshift` | +| Complex work-array transforms | `zffti`, `zfftf`, `zfftb` | +| Real work-array transforms | `dffti`, `dfftf`, `dfftb`, `dzffti`, `dzfftf`, `dzfftb` | +| Cosine and sine work-array transforms | `dcosqi`, `dcosqf`, `dcosqb`, `dcosti`, `dcost`, `dsinti`, `dsint` | --- -## Versions used +## Quick start -| Component | Version / source | -| --- | --- | -| PRIK | current repository checkout | -| FFTPACK | [fortran-lang/fftpack commit `0fffe7c`](https://github.com/fortran-lang/fftpack/tree/0fffe7c05a918363a7cc12ae138a695afd115f36) | -| Python | 3.12 in the dedicated CI job | -| NumPy | 2.5.1 | -| SciPy | 1.18.0 | -| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | +From a PRIK checkout with PRIK installed, GNU Fortran on `PATH`, and the pinned +NumPy and SciPy (see [Set up a clean environment](#set-up-a-clean-environment)): -The repository owns the checked-in source snapshot under -`examples/fortran/fftpack/native/`, so the example does not download code during its -build. +```bash +source examples/fortran/fftpack/build_all.sh +python3 -m pytest -q examples/fortran/fftpack/tests +``` -## Tested platforms +The first command builds the extension and puts it on `PYTHONPATH` for this +shell; use `source`, not `bash`, so that setting survives. The second runs the +tests. -The Real Libraries Portability workflow builds and runs the complete numerical -suite with Python 3.12 on: +After this, the transforms import in the same shell: -| Operating system | Architectures | Native toolchain | -| --- | --- | --- | -| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | -| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | +```python +import prik_reference_fftpack + +fftpack = prik_reference_fftpack.fftpack +``` + +[Use the generated API](#use-the-generated-api) shows complete calls. --- -## 1. Prepare the repository and toolchain +## Key files -Clone PRIK, create a virtual environment, and install the Python tools used by -the dedicated CI job: +Everything lives under [`examples/fortran/fftpack/`](../../../../examples/fortran/fftpack/): -```bash -git clone https://github.com/PyNumLab/prik.git -cd prik -python3 -m venv .venv -. .venv/bin/activate -python3 -m pip install --upgrade pip -python3 -m pip install -e ".[qa]" "numpy==2.5.1" "scipy==1.18.0" -``` +| File | What it does | +| --- | --- | +| [`native/fftpack.f90`](../../../../examples/fortran/fftpack/native/fftpack.f90) | Declares the public `fftpack` module, which defines the Python API. | +| [`native/rk.f90`](../../../../examples/fortran/fftpack/native/rk.f90) | Defines the real kind used by that API. | +| `native/fftpack_*.f90` | Implement the module's procedures as Fortran submodules. | +| The other files in [`native/`](../../../../examples/fortran/fftpack/native/) | The computational kernels, compiled and linked but not exposed to Python. | +| [`build_prik.sh`](../../../../examples/fortran/fftpack/build_prik.sh) | Builds the extension with one PRIK command that gives each source its role. | +| [`build_all.sh`](../../../../examples/fortran/fftpack/build_all.sh) | Runs `build_prik.sh` and adds the extension to `PYTHONPATH`. | +| [`routine_inventory.py`](../../../../examples/fortran/fftpack/routine_inventory.py) | The list of the 31 procedures, grouped by family. | +| [`tests/test_transforms.py`](../../../../examples/fortran/fftpack/tests/test_transforms.py) | One test per procedure, checked against NumPy, SciPy, or a transform property. | +| [`tests/helpers.py`](../../../../examples/fortran/fftpack/tests/helpers.py) | Helpers such as converting between FFTPACK's and NumPy's real-FFT layouts. | +| [`tests/test_routine_coverage.py`](../../../../examples/fortran/fftpack/tests/test_routine_coverage.py) | Checks that the inventory, the generated exports, and the tests stay in sync. | -Install GNU Fortran separately. On Ubuntu: +--- -```bash -sudo apt-get update -sudo apt-get install --yes gfortran -gfortran --version -``` +## How the build works -All remaining commands run from the repository root with the virtual -environment active. The complete runnable project lives under -[`examples/fortran/fftpack/`](../../../../examples/fortran/fftpack/). +FFTPACK splits into public declarations, submodule implementations, and +low-level kernels. One PRIK command compiles all of them once, but only the +first group shapes the Python API: ---- +```text +rk.f90, fftpack.f90, fftpack_*.f90 ── the Python API ──┐ + ├──prik──> prik_reference_fftpack +the other .f90 kernels ── --native-fortran-sources ────┘ +``` -## 2. Build the PRIK wrapper +| Source group | Passed as | Role | +| --- | --- | --- | +| `rk.f90`, `fftpack.f90`, `fftpack_*.f90` | Positional sources | Define the Python-facing API and its implementation. | +| The remaining `.f90` kernels | `--native-fortran-sources` | Satisfy native dependencies without adding their storage-level signatures to the Python API. | -FFTPACK uses public module declarations, submodule implementations, and -link-only computational kernels. The build command gives each source the role -it needs: +`build_prik.sh` runs that command: ```bash @@ -121,37 +125,15 @@ python3 -m prik "${FFTPACK_PUBLIC_SOURCES[@]}" \ --wrapper-c-flags="-O0 -g0" ``` -The example uses `-O0` so the tests focus on correct results. Every source is -compiled once: positional files define the Python-facing API, while -`--native-fortran-sources` adds implementation code without exposing it to -Python. - -For normal use, source the convenience entrypoint: - -```bash -source examples/fortran/fftpack/build_all.sh -``` - -It builds the extension and exports its directory on `PYTHONPATH` for the -current shell. +The example uses `-O0` so the tests focus on correct results. Everything is +written to the temporary `FFTPACK_BUILD_ROOT` directory, not to the +repository. --- -## 3. Understand how sources define the API - -The source groups have different roles: +## Use the generated API -| Source group | Responsibility | -| --- | --- | -| `rk.f90` | Defines the real kind used by the public API. | -| `fftpack.f90` | Declares the public FFTPACK module. | -| `fftpack_*.f90` | Implements its procedures in Fortran submodules. | -| Remaining `.f90` files | Supply linked computational kernels. | - -The public declarations define the Python types. For example, `zfftf` accepts -an ordinary NumPy `complex128` array. - -High-level transform results that are allocatable in Fortran use PRIK's +High-level transforms whose Fortran results are allocatable return PRIK's `AllocatableArray` handle. Read the NumPy view with `to_numpy()` and release the native allocation with `close()`: @@ -167,41 +149,31 @@ finally: result.close() ``` -Fixed-shape frequency and shift results are returned directly as NumPy arrays. - ---- - -## 4. Run the complete test suite +Fixed-shape frequency and shift results, such as `fftfreq`, return NumPy arrays +directly. The work-array routines (`zffti`, `zfftf`, …) keep FFTPACK's +initialize-then-transform pattern and update the caller's array in place; see +the `zfftf` test [below](#how-results-are-validated). -After the build finishes, run: +**FFTPACK's conventions differ from `numpy.fft`.** The same inputs give: -```bash -python3 -m pytest -q examples/fortran/fftpack/tests -``` +| Call | FFTPACK | `numpy.fft` | +| --- | --- | --- | +| `fftfreq(4)` | `[0, 1, -2, -1]`: integer frequency indices | `[0, 0.25, -0.5, -0.25]`: the indices divided by `n` | +| `ifft(fft(x))` | `n * x`: the inverse is not normalized | `x` | +| `rfft([1, 2, 3, 4])` | `[10, -2, 2, -2]`: packed real and imaginary parts | `[10, -2+2j, -2]`: complex values | -The tests cover all 31 public procedures: - -| Family | Procedures | -| --- | ---: | -| Complex work-array transforms | 3 | -| Real work-array transforms | 6 | -| Cosine and sine work-array transforms | 7 | -| High-level Fourier transforms | 4 | -| High-level cosine transforms | 7 | -| Frequency and spectrum ordering | 4 | -| **Total** | **31** | - -Each procedure is called with representative data and checked against NumPy, -SciPy, or a known transform property. +Divide an `ifft` result by `n` to recover the input. The tests convert the +packed real layout with a helper in `tests/helpers.py`. --- -## 5. See how results are validated +## How results are validated -The suite compares transform results with independent NumPy or SciPy results -and also checks in-place mutation, dtype, shape, normalization, and frequency -ordering. For example, this `zfftf` test comes directly from the runnable -suite: +NumPy is the reference for the Fourier transforms, shifts, and frequency +ordering; SciPy is the reference for the cosine and sine families. The suite +also checks normalization, in-place mutation, preservation of high-level +inputs, dtype, shape, frequency ordering, and release of allocatable results. +For example, this `zfftf` test comes directly from the runnable suite: ```python @@ -221,34 +193,73 @@ in place, and compares the result with NumPy's independently implemented FFT. --- -## 6. Run focused examples +## Run the tests -After building the extension, run a family or one procedure: +Run the complete suite, one procedure, or every test that mentions a name: ```bash -python3 -m pytest -q examples/fortran/fftpack/tests/test_transforms.py +python3 -m pytest -q examples/fortran/fftpack/tests python3 -m pytest -q \ examples/fortran/fftpack/tests/test_transforms.py::test_zfftf python3 -m pytest -q examples/fortran/fftpack/tests -k fftshift ``` -- Complete numerical examples → - [`test_transforms.py`](../../../../examples/fortran/fftpack/tests/test_transforms.py) -- Public routine list → - [`routine_inventory.py`](../../../../examples/fortran/fftpack/routine_inventory.py) -- Routine coverage check → - [`test_routine_coverage.py`](../../../../examples/fortran/fftpack/tests/test_routine_coverage.py) -- Copyable project instructions → - [`examples/fortran/fftpack/README.md`](../../../../examples/fortran/fftpack/README.md) - --- +## Set up a clean environment + +Clone PRIK, create a virtual environment, and install the Python tools used by +the dedicated CI job: + +```bash +git clone https://github.com/PyNumLab/prik.git +cd prik +python3 -m venv .venv +. .venv/bin/activate +python3 -m pip install --upgrade pip +python3 -m pip install -e ".[qa]" "numpy==2.5.1" "scipy==1.18.0" +``` + +Install GNU Fortran separately. On Ubuntu: + +```bash +sudo apt-get update +sudo apt-get install --yes gfortran +gfortran --version +``` + +Run the example's commands from the repository root with the virtual +environment active. + +## Versions used + +| Component | Version / source | +| --- | --- | +| PRIK | current repository checkout | +| FFTPACK | [fortran-lang/fftpack commit `0fffe7c`](https://github.com/fortran-lang/fftpack/tree/0fffe7c05a918363a7cc12ae138a695afd115f36) | +| Python | 3.12 in the dedicated CI job | +| NumPy | 2.5.1 | +| SciPy | 1.18.0 | +| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | + +## Tested platforms + +The Real Libraries Portability workflow builds and runs the complete numerical +suite with Python 3.12 on: + +| Operating system | Architectures | Native toolchain | +| --- | --- | --- | +| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | +| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | + ## Troubleshooting - Confirm that `gfortran` is available on `PATH`. -- Use `source examples/fortran/fftpack/build_all.sh`; executing it in a child shell - does not preserve the exported `PYTHONPATH`. -- Run one failing procedure with `-vv -s` to retain its compiler and wrapper +- Use `source examples/fortran/fftpack/build_all.sh`; running it with `bash` + starts a child shell, so the exported `PYTHONPATH` is lost. +- A result off by a factor of `n`, or in an unexpected order: see the + [convention table](#use-the-generated-api) above. +- Run one failing procedure with `-vv -s` to see its compiler and wrapper diagnostics. --- From 1d0ea2f857c559ab82b32c1f4698709343225db4 Mon Sep 17 00:00:00 2001 From: said Date: Mon, 28 Sep 2026 14:47:04 +0100 Subject: [PATCH 6/6] codex: restructure the BSPLINE-FORTRAN example guide Lead with the two namespaces and their public surface, a quick start, and a key-files table; explain the three-source build with a diagram; add a table of what PRIK maps from modern Fortran (abstract base, extensions, deferred and inherited bindings, generic constructors, private members, constants); and extend the API section with printed values, a derivative, a 2-D surface, and status codes. The quick start and every snippet were run against a fresh build (39 tests passed). Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 1 + docs/user/examples/fortran/bspline-wrapper.md | 289 ++++++++++-------- 2 files changed, 165 insertions(+), 125 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 830bbc77c..cb7ced46c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -12,6 +12,7 @@ release tags add a leading `v` to the package version. - Improve the MINPACK example guide - Compare the MINPACK example's eight SciPy-exposed solvers with SciPy - Improve the FFTPACK example guide +- Improve the BSPLINE-FORTRAN example guide - Improve the Open MPI `mpi_f08` tutorial - Improve the PRIMA example guide - The README and the documentation homepage now open with a short terminal demo. diff --git a/docs/user/examples/fortran/bspline-wrapper.md b/docs/user/examples/fortran/bspline-wrapper.md index d9338eb04..52f2a115c 100644 --- a/docs/user/examples/fortran/bspline-wrapper.md +++ b/docs/user/examples/fortran/bspline-wrapper.md @@ -9,86 +9,79 @@ publication: reviewed # Build and Validate BSPLINE-FORTRAN with PRIK -This example takes the checked-in -[BSPLINE-FORTRAN](https://github.com/jacobwilliams/bspline-fortran) source and -builds an importable Python extension with the complete interpolation surface: -15 public procedural routines, eight order constants, and seven public classes. +This example turns [BSPLINE-FORTRAN](https://github.com/jacobwilliams/bspline-fortran), +a modern Fortran 2008 B-spline interpolation library, into one Python +extension. It wraps the upstream source unmodified, including an abstract +derived type, six concrete subclasses, and generic constructors. The suite +checks interpolation from one to six dimensions against analytic functions and +SciPy. -It evaluates B-splines from one to six dimensions. The tests compare results -with analytic functions and SciPy rather than treating the wrapper as its own -reference. +### What you get -### What this example shows +One extension, `prik_bspline`, with two namespaces: -- Wrap a modern multi-file Fortran library as one Python extension. -- Construct and call derived types over an abstract Fortran base. -- Check procedural and object-oriented interpolation with NumPy arrays. - -You should already be comfortable with NumPy arrays, Python classes, and -building a local Fortran extension. +| Namespace | Public surface | +| --- | --- | +| `bspline_oo_module` | The abstract `Bspline_Class` and six concrete classes, `Bspline_1d` through `Bspline_6d` | +| `bspline_sub_module` | 15 procedural routines: `db1ink`-`db6ink` (setup), `db1val`-`db6val` (evaluation), `db1sqad` and `db1fqad` (definite integrals), and `get_status_message`; plus eight order constants such as `bspline_order_cubic` | --- -## Versions used +## Quick start -| Component | Version / source | -| --- | --- | -| PRIK | current repository checkout | -| BSPLINE-FORTRAN | [version 7.4.0, commit `047c7244`](https://github.com/jacobwilliams/bspline-fortran/tree/047c7244) | -| Python | 3.12 in the dedicated CI job | -| NumPy | 2.5.1 | -| SciPy | 1.18.0 | -| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | +From a PRIK checkout with PRIK installed, GNU Fortran on `PATH`, and the pinned +NumPy and SciPy (see [Set up a clean environment](#set-up-a-clean-environment)): -The repository owns the checked-in source snapshot under -`examples/fortran/bspline/native/`, so the example does not download code during its -build. +```bash +source examples/fortran/bspline/build_all.sh +python3 -m pytest -q examples/fortran/bspline/tests +``` -## Tested platforms +The first command builds the extension and puts it on `PYTHONPATH` for this +shell; use `source`, not `bash`, so that setting survives. The second runs the +tests. -The Real Libraries Portability workflow builds and runs the complete numerical -suite with Python 3.12 on: +After this, the classes import in the same shell: -| Operating system | Architectures | Native toolchain | -| --- | --- | --- | -| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | -| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | +```python +import prik_bspline.bspline_oo_module as bspline +``` ---- +[Use the generated API](#use-the-generated-api) shows complete calls. -## 1. Prepare the repository and toolchain +--- -Clone PRIK, create a virtual environment, and install the Python tools used by -the dedicated CI job: +## Key files -```bash -git clone https://github.com/PyNumLab/prik.git -cd prik -python3 -m venv .venv -. .venv/bin/activate -python3 -m pip install --upgrade pip -python3 -m pip install -e ".[qa]" "numpy==2.5.1" "scipy==1.18.0" -``` +Everything lives under [`examples/fortran/bspline/`](../../../../examples/fortran/bspline/): -Install GNU Fortran separately. On Ubuntu: +| File | What it does | +| --- | --- | +| [`native/bspline_kinds_module.F90`](../../../../examples/fortran/bspline/native/bspline_kinds_module.F90) | Defines the kind parameters used by the other modules. | +| [`native/bspline_sub_module.f90`](../../../../examples/fortran/bspline/native/bspline_sub_module.f90) | The procedural interface: setup, evaluation, and integral routines. | +| [`native/bspline_oo_module.f90`](../../../../examples/fortran/bspline/native/bspline_oo_module.f90) | The object-oriented interface: the abstract base and its six extensions. | +| [`build_prik.sh`](../../../../examples/fortran/bspline/build_prik.sh) | Builds the extension with one PRIK command. | +| [`build_all.sh`](../../../../examples/fortran/bspline/build_all.sh) | Runs `build_prik.sh` and adds the extension to `PYTHONPATH`. | +| [`routine_inventory.py`](../../../../examples/fortran/bspline/routine_inventory.py) | The reviewed public surface: classes, bindings, routines, and constants. | +| [`tests/test_object_oriented_api.py`](../../../../examples/fortran/bspline/tests/test_object_oriented_api.py) | Constructs and evaluates every class, and checks the abstract-base and inheritance behavior. | +| [`tests/test_procedural_api.py`](../../../../examples/fortran/bspline/tests/test_procedural_api.py) | One named numerical test per procedural routine. | +| [`tests/test_routine_coverage.py`](../../../../examples/fortran/bspline/tests/test_routine_coverage.py) | Fails if an expected export disappears, an extra one appears, or a routine has no test. | -```bash -sudo apt-get update -sudo apt-get install --yes gfortran -gfortran --version -``` +--- -All remaining commands run from the repository root with the virtual -environment active. The complete runnable project lives under -[`examples/fortran/bspline/`](../../../../examples/fortran/bspline/). +## How the build works ---- +One PRIK command reads the three sources in dependency order and compiles them, +with the generated wrapper, into one extension. The upstream source is not +edited: -## 2. Build the PRIK wrapper +```text +bspline_kinds_module.F90 ─┐ +bspline_sub_module.f90 ───┼──prik──> prik_bspline ┬─ bspline_sub_module (procedural) +bspline_oo_module.f90 ────┘ └─ bspline_oo_module (classes) +``` -BSPLINE-FORTRAN separates its kind definitions, procedural routines, and -object-oriented types into ordered source files. The build command passes those -three public sources in dependency order: +`build_prik.sh` runs that command: ```bash @@ -110,84 +103,94 @@ python3 -m prik \ --wrapper-c-flags="-O0 -g0" ``` -The example uses `-O0` so the tests focus on correct results. PRIK compiles the -native source and generated bridge into one extension. +The example uses `-O0` so the tests focus on correct results. Everything is +written to the temporary `BSPLINE_BUILD_ROOT` directory, not to the repository. -For normal use, source the convenience entrypoint: +--- -```bash -source examples/fortran/bspline/build_all.sh -``` +## What PRIK maps from modern Fortran -It builds the extension and exports its directory on `PYTHONPATH` for the -current shell. +| Fortran construct | In Python | +| --- | --- | +| Abstract type `bspline_class` with deferred bindings | `Bspline_Class`, which raises `TypeError` if constructed directly | +| Six extensions of that type | Subclasses: `issubclass(bspline.Bspline_1d, bspline.Bspline_Class)` is `True` | +| Bindings the base implements once | Inherited methods: `clear_flag`, `status_message`, `status_ok` | +| Deferred bindings | Methods each class answers itself: `destroy`, `size_of` | +| Generic constructor `interface bspline_1d` | `Bspline_1d(...)`, which accepts both the empty and the data-driven form | +| Generic interfaces such as `db1ink` and `db1val` | One Python name for their specific procedures | +| Private components and private bindings | Not exposed | +| Module constants such as `bspline_order_cubic` | Module attributes (`bspline_order_cubic == 4`) | --- -## 3. Use the generated Python API +## Use the generated API -The object-oriented module exposes an abstract `bspline_class` and six concrete -dimension-specific subclasses. The `bspline_1d` generic constructor accepts an -empty form and a data-driven form: +**A one-dimensional spline.** The data-driven constructor builds the spline; +`evaluate` takes the derivative order as its second argument: ```python import numpy as np import prik_bspline.bspline_oo_module as bspline x = np.linspace(0.0, 2.0 * np.pi, 25) -spline = bspline.Bspline_1d(x, np.sin(x), np.int32(4)) +spline = bspline.Bspline_1d(x, np.sin(x), np.int32(4)) # cubic: order 4 value, iflag = spline.evaluate(np.float64(1.234), np.int32(0)) +print(value) # 0.943811 (sin(1.234) = 0.943818) + +slope, iflag = spline.evaluate(np.float64(1.234), np.int32(1)) +print(slope) # 0.330588 (cos(1.234) = 0.330465) + area, iflag = spline.integral(np.float64(0.0), np.float64(np.pi)) +print(area) # 1.999991 (exact: 2) ``` -The abstract base is exported but cannot be constructed. Its concrete -extensions inherit the base bindings and answer its deferred operations: +`iflag == 0` means success. `spline.status_ok()` reports whether the last call +succeeded, and `spline.status_message(iflag)` turns a code into text: evaluating +at `x = 99`, outside the data, gives `iflag == 601`, "Error in db*val: x value +out of bounds". + +**A two-dimensional surface.** Sample on a grid in Fortran order and pass one +order per dimension: ```python -bspline.Bspline_Class() -# TypeError: bspline_class is an abstract native type and cannot be -# instantiated; create one of its concrete extensions instead +x = np.linspace(0.0, np.pi, 20) +y = np.linspace(0.0, np.pi, 20) +samples = np.asfortranarray(np.sin(x)[:, None] * np.cos(y)[None, :]) -issubclass(bspline.Bspline_1d, bspline.Bspline_Class) # True +surface = bspline.Bspline_2d(x, y, samples, np.int32(4), np.int32(4)) +value, iflag = surface.evaluate(np.float64(1.0), np.float64(0.5), np.int32(0), np.int32(0)) +print(value) # 0.738460 (sin(1) * cos(0.5) = 0.738460) ``` -The procedural module exposes the matching `db1ink` through `db6ink` setup -routines and `db1val` through `db6val` evaluators. Pass ordinary NumPy arrays; -PRIK performs the ABI conversion inside the generated wrapper. +`Bspline_3d` through `Bspline_6d` follow the same pattern. ---- - -## 4. Run the complete test suite +**The abstract base** is exported but cannot be constructed: -After the build finishes, run: - -```bash -python3 -m pytest -q examples/fortran/bspline/tests +```python +bspline.Bspline_Class() +# TypeError: bspline_class is an abstract native type and cannot be +# instantiated; create one of its concrete extensions instead ``` -The tests cover every exported routine and class: - -| Family | Public surface | -| --- | ---: | -| Interpolation setup | 6 routines | -| Evaluation | 6 routines | -| Definite integrals | 2 routines | -| Status reporting | 1 routine | -| Order constants | 8 constants | -| Derived types | 1 abstract base + 6 concrete classes | - -The inventory test fails if an expected export disappears, an extra public -export appears, or a procedural routine has no named numerical test. +**The procedural interface** in `bspline_sub_module` mirrors the Fortran +routines: `db1ink` builds the knots and coefficients, and `db1val` evaluates +them, with the arrays passed explicitly. The `db1sqad` test below shows it. --- -## 5. See how results are validated +## How results are validated + +The suite builds every procedural family from one to six dimensions and +evaluates an affine function through every evaluator. It checks +one-dimensional analytic values, derivatives, definite integrals, and +callback-driven integration, plus a comparison with SciPy's +`make_interp_spline`. The object-oriented tests construct and evaluate every +concrete class and check the abstract-base, inheritance, deferred-binding, and +generic-constructor behavior. -The suite checks interpolation against analytic values and SciPy, along with -constructor behavior, inheritance, abstract-base dispatch, generated status, -and Fortran-order array handling. This test comes directly from the runnable -suite and shows the procedural one-dimensional definite integral: +This test builds a cubic spline for `sin(x)` through the procedural interface, +integrates it from zero to π, and checks the known value of two: ```python @@ -201,41 +204,77 @@ def test_db1sqad(bspline_sub): assert value == pytest.approx(2.0, abs=1.0e-6) ``` -It builds a cubic spline for `sin(x)`, integrates it from zero to π, and checks -the known value of two. - --- -## 6. Run focused examples +## Run the tests -After building the extension, run a family or one routine: +Run the complete suite, one interface family, one routine, or every test that +mentions a name: ```bash +python3 -m pytest -q examples/fortran/bspline/tests python3 -m pytest -q examples/fortran/bspline/tests/test_object_oriented_api.py python3 -m pytest -q \ examples/fortran/bspline/tests/test_procedural_api.py::test_db1ink python3 -m pytest -q examples/fortran/bspline/tests -k db6 ``` -- Derived-type examples → - [`test_object_oriented_api.py`](../../../../examples/fortran/bspline/tests/test_object_oriented_api.py) -- Procedural numerical examples → - [`test_procedural_api.py`](../../../../examples/fortran/bspline/tests/test_procedural_api.py) -- Public surface and coverage check → - [`test_routine_coverage.py`](../../../../examples/fortran/bspline/tests/test_routine_coverage.py) -- Reviewed inventory → - [`routine_inventory.py`](../../../../examples/fortran/bspline/routine_inventory.py) -- Copyable project instructions → - [`examples/fortran/bspline/README.md`](../../../../examples/fortran/bspline/README.md) - --- +## Set up a clean environment + +Clone PRIK, create a virtual environment, and install the Python tools used by +the dedicated CI job: + +```bash +git clone https://github.com/PyNumLab/prik.git +cd prik +python3 -m venv .venv +. .venv/bin/activate +python3 -m pip install --upgrade pip +python3 -m pip install -e ".[qa]" "numpy==2.5.1" "scipy==1.18.0" +``` + +Install GNU Fortran separately. On Ubuntu: + +```bash +sudo apt-get update +sudo apt-get install --yes gfortran +gfortran --version +``` + +Run the example's commands from the repository root with the virtual +environment active. + +## Versions used + +| Component | Version / source | +| --- | --- | +| PRIK | current repository checkout | +| BSPLINE-FORTRAN | [version 7.4.0, commit `047c7244`](https://github.com/jacobwilliams/bspline-fortran/tree/047c7244) | +| Python | 3.12 in the dedicated CI job | +| NumPy | 2.5.1 | +| SciPy | 1.18.0 | +| Fortran compiler | GNU Fortran 13 in CI; a compatible `gfortran` works locally | + +## Tested platforms + +The Real Libraries Portability workflow builds and runs the complete numerical +suite with Python 3.12 on: + +| Operating system | Architectures | Native toolchain | +| --- | --- | --- | +| Linux | x86-64, ARM64 | GNU Fortran 13 + GCC 13 | +| macOS | Intel, ARM64 | GNU Fortran 13 + GNU GCC 13 | + ## Troubleshooting - Confirm that `gfortran` is available on `PATH`. -- Use `source examples/fortran/bspline/build_all.sh`; executing it in a child shell does - not preserve the exported `PYTHONPATH`. -- Run one failing procedure with `-vv -s` to retain its compiler and wrapper +- Use `source examples/fortran/bspline/build_all.sh`; running it with `bash` + starts a child shell, so the exported `PYTHONPATH` is lost. +- A nonzero `iflag`: `status_message(iflag)` on the spline, or + `get_status_message(iflag)` in the procedural interface, explains it. +- Run one failing procedure with `-vv -s` to see its compiler and wrapper diagnostics. ---