Skip to content

OTHER: Vendor pentapy's solvers#50

Merged
derb12 merged 11 commits intodevelopmentfrom
vendor_pentapy
Sep 30, 2025
Merged

OTHER: Vendor pentapy's solvers#50
derb12 merged 11 commits intodevelopmentfrom
vendor_pentapy

Conversation

@derb12
Copy link
Copy Markdown
Owner

@derb12 derb12 commented Sep 30, 2025

Description

Vendors the PTRANS-I and PTRANS-II solvers from pentapy as private functions. The solvers were internally changed following the guidelines presented in pentapy issue 11 and pull request 27, as well as to fit some personal preferences. Changes outlined as follows:

  1. Works natively on LAPACK banded format (column-wise) rather than row-wise banded format, which makes it easier to use alongside SciPy's banded solvers, especially when doing intermediate calculations with the banded matrices. Theoretically the extra indexing should require ~2N extra addtitions compared to the 19N-29 total operations in row-wise, but I didn't notice any timing increase when changing the internal logic to work column-wise rather than row-wise.
  2. Supports one or multiple right hand sides. This could be done using broadcasting without having to factorize internally, so it was a nice feature to add without any major changes to the actual solver logic.
  3. Added separate functions for taking the factorization and reusing for subsequent solves.
  4. Actually raises an exception when the solve fails, rather than emitting a warning or silently returning nans.

Points (2) and (3) aren't of use in pybaselines yet, but the optimize_lam branch will benefit from them.

The using_pentapy attribute of _banded_utils.PenalizedSystem was renamed to using_penta. It didn't need to be done, but a lot of my personal code uses the PenalizedSystem object and would otherwise potentially be using banded matrices in an unexpected band arrangement without any warning (why point (1) was done), so I thought it'd be better to cause any code using that attribute to fail rather than silently be in an unexpected band format.

Concerning performance, the direct solver is ~20-30% faster than the current main branch of pentapy thanks to improved indexing and solving without tracking the full factorization, both explained in pentapy issue 11. Both the direct solver and factorized solver produce identical outputs to pentapy's solvers, at least on all of the different Whittaker smoothing systems I've tested. Timings and correctness can be compared with the following gist: https://gist.github.com/derb12/1eff076884a5512fb5e833e47cea91bc

Type of Pull Request

  • Bug Fix
  • New Feature
  • Miscellaneous Changes (refactor, code improvements, etc.)
  • Documentation or Example Programs

Pull Request Checklist

  • New code and/or documentation is valid for use with the BSD 3-clause license.
  • New code is fully documented with docstrings that follow Numpy style,
    if applicable.
  • New code follows PEP 8 standards as closely as possible, if applicable.
  • Added/updated tests and ensured they pass locally, if applicable.
  • Verified that documentation builds locally, if applicable.

derb12 added 11 commits August 10, 2025 17:13
Used the modifications in pentapy issue 11 and PR 27 to improve indexing and error handling of the solvers. These versions are ~20-35% faster than pentapy's current versions.
Much easier to use it as a drop-in replacement for solve_banded and does not incur a significant time penalty.
Makes it easier to follow along when cross-checking with the reference paper.
Also added preliminary tests for the pentadiagonal solvers.
Since using_pentapy is no longer a valid way to tell if matrices are row-wise banded, need to cause a breaking change to any code that uses it. Likely not used by outside users, but this ensures the new behavior is tracked.
Uses 2N less additions within the inner loop (time saving negligible), but more importantly unifies the code for the direct solvers and factorizations.
@derb12
Copy link
Copy Markdown
Owner Author

derb12 commented Sep 30, 2025

The failing doctest was introduced in the development branch, not this one, so I'll fix it after merging this.

@derb12 derb12 merged commit a700a58 into development Sep 30, 2025
9 of 11 checks passed
@derb12 derb12 deleted the vendor_pentapy branch September 30, 2025 01:39
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant