Feature/qr and colpiv qr improvements and solver (#1561)

* add coupled blas and lapack dependency

* factor out some qr code so I can reuse it from the unpivoted qr implementation

* update internal docs, add safety checks

* factor out qr solver as well

* start reworking the (colpiv) qr traits

* factor out error

* rework qr traits for better error handling

* start unifying colpiv qr and qr trait interfaces to lapack

* unification of qr traits complete

* kill some useless trait bounds

* add solver and q multiplication methods to the non-pivoted qr factorization

* more cleanups and refactorings

* add trmm interface

* start with r multiplication stuff

* start with unifying qr implementation

* unify calculation logic for most of qr logic in one trait

* doc updates

* export decomposition from crate root

* start with left R multiplication functions

* R*B and R^T*B multiplication

* r mul mut test failing

* r multiplication rework

* restrict QR decomposition to overdetermined system m>=n

* docs and changelog

* fix snafu

* more exhaustive tests for QR decomposition

* update docs

* update changelog

* add tests for Q^T, Q multiplication from left and right

* improve trmm impl

* add R, R^T multiplication functionality from left and right

* add proptests for Q multiplication

* update changelog

* [ai assisted, manually reviewed] updates to documentation

* better safeguards for lapack calls

Author: geo-ant

nalgebra-lapack/CHANGELOG.md

@@ -5,31 +5,38 @@ For the **changes between versions 0.4.0 and 0.27.0** see the main
 
 ## Unreleased
 
-* bugfixes in Schur decomposition
-* bugfixes in LU decomposition
-* fix failing tests for Cholesky decomposition
-* fix compilation with `serde-serialize` feature enabled
-* add column-pivoting QR decomposition and solver
-* fix logic error in calculation of complex eigenvalues in eigen-decomposition.
-* change the feature flags for choosing the lapack backend, update docs accordingly
+* Bugfixes in Schur decomposition
+* Bugfixes in LU decomposition
+* Fix failing tests for Cholesky decomposition
+* Fix compilation with `serde-serialize` feature enabled
+* Add column-pivoting QR decomposition and solver
+* Introduce a common, trait-based interface to QR and column-pivoted QR
+* Extend QR decompositions to allow multiplication with Q, Q^T, R, R^T efficiently
+  from both sides
+* Fix logic error in calculation of complex eigenvalues in eigen-decomposition
+* Change the feature flags for choosing the LAPACK backend, update docs accordingly
+* Remove untested complex support in QR decomposition
+* Remove `Qr::unpack`, since it brings no practical benefit over `(qr.q(), qr.r())`,
+  but looks like it would bring an efficiency gain
+* QR decomposition restricted to m >= n
 
 ## [0.4.0] - 2016-09-07
 
 * Made all traits use associated types for their output type parameters. This
   simplifies usage of the traits and is consistent with the concept of
   associated types used as output type parameters (not input type parameters) as
   described in [the associated type
-  RFC](https://github.com/rust-lang/rfcs/blob/master/text/0195-associated-items.md).
-* Implemented `check_info!` macro to check all LAPACK calls.
-* Implemented error handling with [error_chain](https://crates.io/crates/error-chain).
+  RFC](https://github.com/rust-lang/rfcs/blob/master/text/0195-associated-items.md)
+* Implemented `check_info!` macro to check all LAPACK calls
+* Implemented error handling with [error_chain](https://crates.io/crates/error-chain)
 
 ## [0.3.0] - 2016-09-06
 
 * Documentation is hosted at https://docs.rs/nalgebra-lapack/
-* Updated `nalgebra` to 0.10.
-* Rename traits `HasSVD` to `SVD` and `HasEigensystem` to `Eigensystem`.
-* Added `Solve` trait for solving a linear matrix equation.
-* Added `Inverse` for computing the multiplicative inverse of a matrix.
-* Added `Cholesky` for decomposing a positive-definite matrix.
+* Updated `nalgebra` to 0.10
+* Rename traits `HasSVD` to `SVD` and `HasEigensystem` to `Eigensystem`
+* Added `Solve` trait for solving a linear matrix equation
+* Added `Inverse` for computing the multiplicative inverse of a matrix
+* Added `Cholesky` for decomposing a positive-definite matrix
 * The `Eigensystem` and `SVD` traits are now generic over types. The
-  associated types have been removed.
+  associated types have been removed

nalgebra-lapack/Cargo.toml

@@ -6,7 +6,7 @@ authors = ["Sébastien Crozet <developer@crozet.re>", "Andrew Straw <strawman@as
 description = "Matrix decompositions using nalgebra matrices and Lapack bindings."
 homepage = "https://nalgebra.rs"
 repository = "https://github.com/dimforge/nalgebra"
-readme = "../README.md"
+readme = "README.md"
 categories = ["science", "mathematics"]
 keywords = ["linear", "algebra", "matrix", "vector", "lapack"]
 license = "MIT"
@@ -23,14 +23,14 @@ arbitrary = ["nalgebra/arbitrary"]
 slow-tests = []
 
 # For BLAS/LAPACK
-lapack-openblas = ["lapack-src/openblas"]
-lapack-netlib = ["lapack-src/netlib"]
-lapack-accelerate = ["lapack-src/accelerate"]
+lapack-openblas = ["lapack-src/openblas", "blas-src/openblas"]
+lapack-netlib = ["lapack-src/netlib", "blas-src/netlib"]
+lapack-accelerate = ["lapack-src/accelerate", "blas-src/accelerate"]
 lapack-mkl = ["lapack-mkl-static-seq"]
-lapack-mkl-static-seq = ["lapack-src/intel-mkl-static-sequential"]
-lapack-mkl-static-par = ["lapack-src/intel-mkl-static-parallel"]
-lapack-mkl-dynamic-seq = ["lapack-src/intel-mkl-dynamic-sequential"]
-lapack-mkl-dynamic-par = ["lapack-src/intel-mkl-dynamic-parallel"]
+lapack-mkl-static-seq = ["lapack-src/intel-mkl-static-sequential", "blas-src/intel-mkl-static-sequential"]
+lapack-mkl-static-par = ["lapack-src/intel-mkl-static-parallel", "blas-src/intel-mkl-static-parallel"]
+lapack-mkl-dynamic-seq = ["lapack-src/intel-mkl-dynamic-sequential", "blas-src/intel-mkl-dynamic-sequential"]
+lapack-mkl-dynamic-par = ["lapack-src/intel-mkl-dynamic-parallel", "blas-src/intel-mkl-dynamic-parallel"]
 lapack-custom = []
 
 [dependencies]
@@ -40,9 +40,10 @@ num-complex = { version = "0.4", default-features = false }
 simba = "0.9"
 serde = { version = "1.0", features = ["derive"], optional = true }
 lapack = { version = "0.20", default-features = false }
-lapack-src = { version = "0.13", optional = true}
+blas = { version = "0.23", default-features = false }
+lapack-src = { version = "0.13", optional = true }
+blas-src = { version = "0.14", optional = true, default-features = false}
 thiserror = "1.0"
-# clippy = "*"
 
 [dev-dependencies]
 nalgebra = { version = "0.34", features = ["arbitrary", "rand"], path = ".." }

nalgebra-lapack/README.md

@@ -1,4 +1,4 @@
-# nalgebra-lapack [![Version][version-img]][version-url] [![Status][status-img]][status-url] [![Doc][doc-img]][doc-url]
+# nalgebra-lapack [![Version][version-img]][version-url] [![Doc][doc-img]][doc-url]
 
 Rust library for linear algebra using nalgebra and LAPACK.
 
@@ -10,45 +10,43 @@ Documentation is available [here](https://docs.rs/nalgebra-lapack/).
 
 MIT
 
-## Cargo features to select lapack provider
+## Cargo features to select LAPACK provider
 
 Like the [lapack crate](https://crates.io/crates/lapack) from which this
 behavior is inherited, nalgebra-lapack uses [cargo
 features](https://doc.crates.io/manifest.html#the-[features]-section) to select
-which lapack provider (or implementation) is used. Command line arguments to
+which LAPACK provider (or implementation) is used. Command line arguments to
 cargo are the easiest way to do this, and the best provider depends on your
 particular system. In some cases, the providers can be further tuned with
 environment variables.
 
 Below are given examples of how to invoke `cargo build` on two different systems
-using two different providers. The `--no-default-features --features "provider"`
+using two different providers. The `--no-default-features --features "lapack-*"`
 arguments will be consistent for other `cargo` commands.
 
 ### Ubuntu
 
-As tested on Ubuntu 12.04, do this to build the lapack package against
+As tested on Ubuntu 24.04, do this to build the LAPACK package against
 the system installation of netlib without LAPACKE (note the E) or
 CBLAS:
 
-    sudo apt-get install gfortran libblas3gf liblapack3gf
+    sudo apt-get install gfortran libblas-dev liblapack-dev
     export CARGO_FEATURE_SYSTEM_NETLIB=1
     export CARGO_FEATURE_EXCLUDE_LAPACKE=1
     export CARGO_FEATURE_EXCLUDE_CBLAS=1
 
-    export CARGO_FEATURES="--no-default-features --features netlib"
+    export CARGO_FEATURES="--no-default-features --features lapack-netlib"
     cargo build ${CARGO_FEATURES}
 
 ### macOS
 
 On macOS, do this to use Apple's Accelerate framework:
 
-    export CARGO_FEATURES="--no-default-features --features accelerate"
+    export CARGO_FEATURES="--no-default-features --features lapack-accelerate"
     cargo build ${CARGO_FEATURES}
 
 [version-img]: https://img.shields.io/crates/v/nalgebra-lapack.svg
 [version-url]: https://crates.io/crates/nalgebra-lapack
-[status-img]: https://travis-ci.org/strawlab/nalgebra-lapack.svg?branch=master
-[status-url]: https://travis-ci.org/strawlab/nalgebra-lapack
 [doc-img]: https://docs.rs/nalgebra-lapack/badge.svg
 [doc-url]: https://docs.rs/nalgebra-lapack/
 

nalgebra-lapack/src/cholesky.rs

@@ -11,7 +11,7 @@ use na::{DefaultAllocator, Matrix, OMatrix, Scalar};
 
 use lapack;
 
-/// The cholesky decomposition of a symmetric-definite-positive matrix.
+/// The Cholesky decomposition of a symmetric positive-definite matrix.
 #[cfg_attr(feature = "serde-serialize", derive(Serialize, Deserialize))]
 #[cfg_attr(
     feature = "serde-serialize",
@@ -42,7 +42,7 @@ impl<T: CholeskyScalar + Zero, D: Dim> Cholesky<T, D>
 where
     DefaultAllocator: Allocator<D, D>,
 {
-    /// Computes the cholesky decomposition of the given symmetric-definite-positive square
+    /// Computes the Cholesky decomposition of the given symmetric positive-definite square
     /// matrix.
     ///
     /// Only the lower-triangular part of the input matrix is considered.
@@ -51,7 +51,7 @@ where
         // TODO: check symmetry as well?
         assert!(
             m.is_square(),
-            "Unable to compute the cholesky decomposition of a non-square matrix."
+            "Unable to compute the Cholesky decomposition of a non-square matrix."
         );
 
         let uplo = b'L';
@@ -64,13 +64,13 @@ where
         Some(Self { l: m })
     }
 
-    /// Retrieves the lower-triangular factor of the cholesky decomposition.
+    /// Retrieves the lower-triangular factor of the Cholesky decomposition.
     pub fn unpack(mut self) -> OMatrix<T, D, D> {
         self.l.fill_upper_triangle(Zero::zero(), 1);
         self.l
     }
 
-    /// Retrieves the lower-triangular factor of che cholesky decomposition, without zeroing-out
+    /// Retrieves the lower-triangular factor of the Cholesky decomposition, without zeroing-out
     /// its strict upper-triangular part.
     ///
     /// This is an allocation-less version of `self.l()`. The values of the strict upper-triangular
@@ -79,15 +79,15 @@ where
         self.l
     }
 
-    /// Retrieves the lower-triangular factor of the cholesky decomposition.
+    /// Retrieves the lower-triangular factor of the Cholesky decomposition.
     #[must_use]
     pub fn l(&self) -> OMatrix<T, D, D> {
         let mut res = self.l.clone();
         res.fill_upper_triangle(Zero::zero(), 1);
         res
     }
 
-    /// Retrieves the lower-triangular factor of the cholesky decomposition, without zeroing-out
+    /// Retrieves the lower-triangular factor of the Cholesky decomposition, without zeroing-out
     /// its strict upper-triangular part.
     ///
     /// This is an allocation-less version of `self.l()`. The values of the strict upper-triangular
@@ -177,7 +177,7 @@ where
  *
  */
 /// Trait implemented by floats (`f32`, `f64`) and complex floats (`Complex<f32>`, `Complex<f64>`)
-/// supported by the cholesky decomposition.
+/// supported by the Cholesky decomposition.
 pub trait CholeskyScalar: Scalar + Copy {
     #[allow(missing_docs)]
     fn xpotrf(uplo: u8, n: i32, a: &mut [Self], lda: i32, info: &mut i32);