package

Author: teseoch

.github/workflows/continuous.yml

@@ -0,0 +1,44 @@
+name: CI
+
+on:
+  push:
+    branches: [main]
+  pull_request:
+
+permissions:
+  contents: read
+
+jobs:
+  test:
+    runs-on: ${{ matrix.os }}
+    name: ${{ matrix.os }} (py${{ matrix.python-version }})
+
+    strategy:
+      fail-fast: false
+      matrix:
+        os: [ubuntu-latest, macos-latest, windows-latest]
+        python-version: ["3.9", "3.11", "3.14"]
+
+    steps:
+      - uses: actions/checkout@v6
+
+      - uses: actions/setup-python@v6
+        with:
+          python-version: ${{ matrix.python-version }}
+
+      - name: Install build tooling
+        run: |
+          python -m pip install --upgrade pip
+          python -m pip install build
+
+      - name: Build wheel to prove binary build works
+        run: |
+          python -m build --wheel --outdir dist
+
+      - name: Install wheel with test extras
+        run: |
+          python -c "import pathlib, subprocess, sys; wheel = next(pathlib.Path('dist').glob('*.whl')); subprocess.check_call([sys.executable, '-m', 'pip', 'install', f'{wheel}[test]'])"
+
+      - name: Run test suite
+        run: |
+          python -m pytest -q

.github/workflows/release.yml

@@ -0,0 +1,77 @@
+name: Build & publish
+
+on:
+  push:
+    tags: ["v[0-9]+.[0-9]+.[0-9]+*"]
+
+permissions:
+  contents: read
+  id-token: write
+
+jobs:
+  build-wheels:
+    runs-on: ${{ matrix.os }}
+    strategy:
+      fail-fast: false
+      matrix:
+        os: [ ubuntu-22.04, macos-14, windows-2022 ]
+
+    steps:
+      - uses: actions/checkout@v6
+        with:
+          fetch-depth: 0
+
+      - name: Set up QEMU for ARM64 emulation
+        if: runner.os == 'Linux'
+        uses: docker/setup-qemu-action@v4
+        with:
+          platforms: linux/arm64
+
+      - name: Build wheels with cibuildwheel
+        uses: pypa/cibuildwheel@v3.4.1
+
+      - name: Upload wheel artifacts
+        uses: actions/upload-artifact@v7
+        with:
+          name: wheels-${{ matrix.os }}
+          path: wheelhouse/*
+          if-no-files-found: error
+
+  sdist:
+    needs: build-wheels
+    runs-on: ubuntu-latest
+
+    steps:
+      - uses: actions/checkout@v6
+      - uses: actions/setup-python@v6
+        with:
+          python-version: '3.11'
+
+      - run: python -m pip install --upgrade pip build
+      - run: python -m build --sdist
+
+      - uses: actions/upload-artifact@v7
+        with:
+          name: sdist
+          path: dist/*.tar.gz
+
+  publish:
+    needs: [ build-wheels, sdist ]
+    runs-on: ubuntu-latest
+
+    steps:
+      - uses: actions/download-artifact@v8
+        with:
+          pattern: wheels-*
+          path: dist
+          merge-multiple: true
+
+      - uses: actions/download-artifact@v8
+        with:
+          name: sdist
+          path: dist
+
+      - name: Publish to PyPI
+        uses: pypa/gh-action-pypi-publish@release/v1
+        with:
+          packages_dir: dist

.gitignore

@@ -1 +1,4 @@
 build
+_skbuild
+__pycache__/
+.pytest_cache/

pyproject.toml

@@ -1,38 +1,68 @@
 [build-system]
 requires = [
-  "scikit-build-core>=0.8.0"
+  "scikit-build-core>=0.8.0",
 ]
 build-backend = "scikit_build_core.build"
 
 [project]
 name = "polysolve"
-version = "v0.0.1"
-description = "A Python binding for PolySolve"
+version = "0.0.1"
+description = "Python bindings for PolySolve nonlinear optimization and linear solvers"
 authors = [
   { name = "Teseo Schneider", email = "teseo@uvic.ca" }
 ]
 license = { text = "MIT" }
 readme = { file = "README.md", content-type = "text/markdown" }
-requires-python = ">=3.8"
-keywords = ["Optimization", "Newton", "Python Bindings"]
+requires-python = ">=3.9"
+keywords = ["Optimization", "Newton", "Linear Solvers", "Python Bindings"]
 classifiers = [
+  "Development Status :: 3 - Alpha",
+  "Intended Audience :: Science/Research",
   "License :: OSI Approved :: MIT License",
   "License :: OSI Approved :: Mozilla Public License 2.0 (MPL 2.0)",
   "Programming Language :: C++",
   "Programming Language :: Python :: 3",
-  "Programming Language :: Python :: 3.8",
+  "Programming Language :: Python :: 3.9",
+  "Programming Language :: Python :: 3.10",
+  "Programming Language :: Python :: 3.11",
+  "Programming Language :: Python :: 3.12",
+  "Programming Language :: Python :: 3.13",
+  "Programming Language :: Python :: 3.14",
   "Operating System :: OS Independent"
 ]
 dependencies = [
   "numpy",
   "scipy",
 ]
 
+[project.optional-dependencies]
+test = [
+  "pytest>=8",
+]
+
+[project.urls]
+Homepage = "https://github.com/polyfem/polysolve-python"
+Repository = "https://github.com/polyfem/polysolve-python"
+Issues = "https://github.com/polyfem/polysolve-python/issues"
+
 # ---------- scikit-build-core ----------
 [tool.scikit-build]
 cmake.version = ">=3.18"
-build-dir = "build/{wheel_tag}"
+build-dir = "_skbuild/{wheel_tag}"
 cmake.args = ["-DPYBIND11_FINDPYTHON=ON", "-DCMAKE_BUILD_TYPE=Release", "-DPOLYSOLVE_WITH_AMGCL=OFF"]
 wheel.exclude = [
   "**/.git/**/*"                       # no embedded Git repos
 ]
+
+[tool.pytest.ini_options]
+testpaths = ["tests"]
+addopts = "-ra"
+
+[tool.cibuildwheel]
+build = "cp39-* cp310-* cp311-* cp312-* cp313-* cp314-*"
+skip = "pp* *-musllinux*"
+test-extras = ["test"]
+test-command = "pytest {project}/tests -q"
+
+[tool.cibuildwheel.linux]
+archs = ["x86_64", "aarch64"]

tests/test_linear.py

@@ -0,0 +1,118 @@
+import numpy as np
+import pytest
+import scipy.sparse
+
+import polysolve
+
+
+SPARSE_SOLVER = {"solver": "Eigen::SimplicialLDLT"}
+DENSE_SOLVER = {"solver": "Eigen::PartialPivLU"}
+
+
+def sparse_system():
+    A = scipy.sparse.csc_matrix(
+        [
+            [4.0, 1.0, 0.0],
+            [1.0, 3.0, 1.0],
+            [0.0, 1.0, 2.0],
+        ]
+    )
+    b = np.array([1.0, 2.0, 3.0])
+    return A, b
+
+
+def assert_solution(A, x, b):
+    np.testing.assert_allclose(A @ x, b, atol=1e-10, rtol=1e-10)
+
+
+def test_solve_sparse_system_one_shot():
+    A, b = sparse_system()
+
+    x = polysolve.solve(
+        A,
+        b,
+        SPARSE_SOLVER,
+        log_level=polysolve.LogLevel.off,
+    )
+
+    assert_solution(A, x, b)
+
+
+def test_linear_solver_reuses_factorization_for_multiple_rhs():
+    A, b = sparse_system()
+    solver = polysolve.LinearSolver(SPARSE_SOLVER, log_level=polysolve.LogLevel.off)
+
+    solver.analyse_pattern(A)
+    with pytest.raises(RuntimeError, match="factorize"):
+        solver.solve(b)
+
+    solver.factorise(A)
+    assert solver.name == "Eigen::SimplicialLDLT"
+    assert not solver.is_dense
+
+    assert_solution(A, solver.solve(b), b)
+
+    b2 = np.array([2.0, -1.0, 0.5])
+    assert_solution(A, solver.solve(b2), b2)
+    assert solver.info()["solver_info"] == "Success"
+
+
+def test_solve_dense_system_one_shot():
+    A = np.array(
+        [
+            [4.0, 1.0, 0.0],
+            [1.0, 3.0, 1.0],
+            [0.0, 1.0, 2.0],
+        ]
+    )
+    b = np.array([1.0, 2.0, 3.0])
+
+    x = polysolve.solve(
+        A,
+        b,
+        DENSE_SOLVER,
+        log_level=polysolve.LogLevel.off,
+    )
+
+    assert_solution(A, x, b)
+
+
+def test_dense_linear_solver_factorize_alias():
+    A = np.array([[4.0, 1.0], [1.0, 3.0]])
+    b = np.array([1.0, 2.0])
+    solver = polysolve.LinearSolver(DENSE_SOLVER, log_level=polysolve.LogLevel.off)
+
+    solver.factorize(A)
+
+    assert solver.is_dense
+    assert_solution(A, solver.solve(b), b)
+
+
+def test_linear_solver_rejects_invalid_shapes():
+    with pytest.raises(RuntimeError, match="square"):
+        polysolve.solve(
+            np.ones((2, 3)),
+            np.ones(2),
+            DENSE_SOLVER,
+            log_level=polysolve.LogLevel.off,
+        )
+
+    solver = polysolve.LinearSolver(DENSE_SOLVER, log_level=polysolve.LogLevel.off)
+    solver.factorize(np.eye(2))
+
+    with pytest.raises(RuntimeError, match="right-hand side"):
+        solver.solve(np.ones(3))
+
+    with pytest.raises(RuntimeError, match="initial guess"):
+        solver.solve(np.ones(2), x0=np.ones(3))
+
+
+def test_linear_solver_introspection_helpers():
+    available = polysolve.available_linear_solvers()
+
+    assert "Eigen::SimplicialLDLT" in available
+    assert polysolve.default_linear_solver() in available
+    assert polysolve.LinearSolver.available_solvers() == available
+    assert polysolve.LinearSolver.default_solver() == polysolve.default_linear_solver()
+    assert isinstance(polysolve.available_linear_preconditioners(), list)
+    assert isinstance(polysolve.default_linear_preconditioner(), str)

tests/test_nonlinear.py

@@ -0,0 +1,73 @@
+import numpy as np
+import pytest
+import scipy.sparse
+
+import polysolve
+
+
+TARGET = np.array([-2.0, 3.0, 1.0])
+NONLINEAR_SOLVER = {
+    "solver": "Newton",
+    "line_search": {"method": "Backtracking"},
+    "max_iterations": 25,
+}
+LINEAR_SOLVER = {"solver": "Eigen::SimplicialLDLT"}
+
+
+class Quadratic(polysolve.Problem):
+    def __init__(self, *, dense_hessian=False):
+        super().__init__()
+        self.dense_hessian = dense_hessian
+        self.steps = []
+
+    def value(self, x):
+        y = x - TARGET
+        return float(y @ y)
+
+    def gradient(self, x):
+        return 2.0 * (x - TARGET)
+
+    def hessian(self, x):
+        if self.dense_hessian:
+            return 2.0 * np.eye(x.size)
+        return 2.0 * scipy.sparse.eye(x.size, format="csc")
+
+    def post_step(self, iter_num, solver_info, x, grad):
+        self.steps.append((iter_num, solver_info, x.copy(), grad.copy()))
+
+
+def minimize(problem):
+    return polysolve.minimize(
+        problem,
+        np.zeros(3),
+        NONLINEAR_SOLVER,
+        LINEAR_SOLVER,
+        log_level=polysolve.LogLevel.off,
+    )
+
+
+def test_minimize_sparse_hessian_quadratic():
+    problem = Quadratic()
+
+    x, info = minimize(problem)
+
+    np.testing.assert_allclose(x, TARGET, atol=1e-10, rtol=1e-10)
+    assert info["energy"] == pytest.approx(0.0, abs=1e-12)
+    assert info["iterations"] <= 2
+    assert problem.steps
+
+
+def test_minimize_accepts_dense_hessian_return_value():
+    x, info = minimize(Quadratic(dense_hessian=True))
+
+    np.testing.assert_allclose(x, TARGET, atol=1e-10, rtol=1e-10)
+    assert info["energy"] == pytest.approx(0.0, abs=1e-12)
+
+
+def test_minimize_validates_gradient_shape():
+    class BadGradient(Quadratic):
+        def gradient(self, x):
+            return np.zeros(x.size - 1)
+
+    with pytest.raises(RuntimeError, match="gradient returned vector"):
+        minimize(BadGradient())