From 7aab3ca738b80f12a92c220ba1b2f8cfed23a7e0 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Fri, 2 Sep 2022 14:48:22 +0900 Subject: [PATCH 01/13] Add katexit in dependencies --- lax/Cargo.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/lax/Cargo.toml b/lax/Cargo.toml index a36c34ab..26a6f35e 100644 --- a/lax/Cargo.toml +++ b/lax/Cargo.toml @@ -33,6 +33,7 @@ thiserror = "1.0.24" cauchy = "0.4.0" num-traits = "0.2.14" lapack-sys = "0.14.0" +katexit = "0.1.2" [dependencies.intel-mkl-src] version = "0.7.0" From 014b5eee3df309e8e40f0857aead02524ca0a8b8 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Sat, 3 Sep 2022 16:09:43 +0900 Subject: [PATCH 02/13] Use katex in JobSvd document --- lax/src/flags.rs | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/lax/src/flags.rs b/lax/src/flags.rs index 37a11b3c..dd9dde3d 100644 --- a/lax/src/flags.rs +++ b/lax/src/flags.rs @@ -92,17 +92,19 @@ impl JobEv { } } -/// Specifies how many of the columns of *U* and rows of *V*ᵀ are computed and returned. +/// Specifies how many singular vectors are computed /// -/// For an input array of shape *m*×*n*, the following are computed: +/// For an input matrix $A$ of shape $m \times n$, +/// the following are computed on the singular value decomposition $A = U\Sigma V^T$: +#[cfg_attr(doc, katexit::katexit)] #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash)] #[repr(u8)] pub enum JobSvd { - /// All *m* columns of *U* and all *n* rows of *V*ᵀ. + /// All $m$ columns of $U$, and/or all $n$ rows of $V^T$. All = b'A', - /// The first min(*m*,*n*) columns of *U* and the first min(*m*,*n*) rows of *V*ᵀ. + /// The first $\min(m, n)$ columns of $U$ and/or the first $\min(m, n)$ rows of $V^T$. Some = b'S', - /// No columns of *U* or rows of *V*ᵀ. + /// No columns of $U$ and/or rows of $V^T$. None = b'N', } From fe5deea72524c09b9b20667fdea31202376e7109 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Sat, 3 Sep 2022 17:28:31 +0900 Subject: [PATCH 03/13] Update document in solve.rs Drop module-level doc in solve.rs --- lax/src/solve.rs | 35 +++++++++++++++++++++++++++-------- 1 file changed, 27 insertions(+), 8 deletions(-) diff --git a/lax/src/solve.rs b/lax/src/solve.rs index ae76f190..ab601031 100644 --- a/lax/src/solve.rs +++ b/lax/src/solve.rs @@ -1,24 +1,43 @@ -//! Solve linear problem using LU decomposition - use crate::{error::*, layout::MatrixLayout, *}; use cauchy::*; use num_traits::{ToPrimitive, Zero}; +#[cfg_attr(doc, katexit::katexit)] +/// Solve linear equations using LU-decomposition +/// +/// For a given matrix $A$, LU decomposition is described as $PA = LU$ where +/// +/// - $L$ is lower matrix +/// - $U$ is upper matrix +/// - $P$ is permutation matrix represented by [Pivot] +/// +/// This is designed as two step computation according to LAPACK API: +/// +/// 1. Factorize input matrix $A$ into $L$, $U$, and $P$. +/// 2. Solve linear equation $Ax = b$ or compute inverse matrix $A^{-1}$ +/// using the output of LU factorization. +/// pub trait Solve_: Scalar + Sized { - /// Computes the LU factorization of a general `m x n` matrix `a` using - /// partial pivoting with row interchanges. + /// Computes the LU factorization of a general $m \times n$ matrix + /// with partial pivoting with row interchanges. /// - /// $ PA = LU $ + /// Output + /// ------- + /// - $U$ and $L$ are stored in `a` after LU factorization has succeeded. + /// - $P$ is returned as [Pivot] /// /// Error /// ------ - /// - `LapackComputationalFailure { return_code }` when the matrix is singular - /// - Division by zero will occur if it is used to solve a system of equations - /// because `U[(return_code-1, return_code-1)]` is exactly zero. + /// - if the matrix is singular + /// - On this case, `return_code` in [Error::LapackComputationalFailure] means + /// `return_code`-th diagonal element of $U$ becomes zero. + /// fn lu(l: MatrixLayout, a: &mut [Self]) -> Result; + /// Compute inverse matrix $A^{-1}$ from the output of LU-factorization fn inv(l: MatrixLayout, a: &mut [Self], p: &Pivot) -> Result<()>; + /// Solve linear equations $Ax = b$ using the output of LU-factorization fn solve(l: MatrixLayout, t: Transpose, a: &[Self], p: &Pivot, b: &mut [Self]) -> Result<()>; } From 1ea3a767d322ef46e09ac30f87036fc0a5ea5639 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Sat, 3 Sep 2022 17:39:59 +0900 Subject: [PATCH 04/13] Update document of solveh.rs --- lax/src/solveh.rs | 67 ++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 60 insertions(+), 7 deletions(-) diff --git a/lax/src/solveh.rs b/lax/src/solveh.rs index 9f65978d..76829982 100644 --- a/lax/src/solveh.rs +++ b/lax/src/solveh.rs @@ -1,17 +1,70 @@ -//! Solve symmetric linear problem using the Bunch-Kaufman diagonal pivoting method. -//! -//! See also [the manual of dsytrf](http://www.netlib.org/lapack/lapack-3.1.1/html/dsytrf.f.html) - use crate::{error::*, layout::MatrixLayout, *}; use cauchy::*; use num_traits::{ToPrimitive, Zero}; +#[cfg_attr(doc, katexit::katexit)] +/// Solve symmetric/hermite indefinite linear problem using the [Bunch-Kaufman diagonal pivoting method][BK]. +/// +/// For a given symmetric matrix $A$, +/// this method factorizes $A = U^T D U$ or $A = L D L^T$ where +/// +/// - $U$ (or $L$) are is a product of permutation and unit upper (lower) triangular matrices +/// - $D$ is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. +/// +/// This takes two-step approach based in LAPACK: +/// +/// 1. Factorize given matrix $A$ into upper ($U$) or lower ($L$) form with diagonal matrix $D$ +/// 2. Then solve linear equation $Ax = b$, and/or calculate inverse matrix $A^{-1}$ +/// +/// [BK]: https://doi.org/10.2307/2005787 +/// pub trait Solveh_: Sized { - /// Bunch-Kaufman: wrapper of `*sytrf` and `*hetrf` + /// Factorize input matrix using Bunch-Kaufman diagonal pivoting method + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [ssytrf] | [dsytrf] | [chetrf] | [zhetrf] | + /// + /// [ssytrf]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_ga12d2e56511cf7df066712c61d9acec45.html + /// [dsytrf]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_gad91bde1212277b3e909eb6af7f64858a.html + /// [chetrf]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_ga081dd1908e46d064c2bf0a1f6b664b86.html + /// [zhetrf]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_gadc84a5c9818ee12ea19944623131bd52.html + /// fn bk(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result; - /// Wrapper of `*sytri` and `*hetri` + + /// Compute inverse matrix $A^{-1}$ from factroized result + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [ssytri] | [dsytri] | [chetri] | [zhetri] | + /// + /// [ssytri]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_gaef378ec0761234aac417f487b43b7a8b.html + /// [dsytri]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_ga75e09b4299b7955044a3bbf84c46b593.html + /// [chetri]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_gad87a6a1ac131c5635d47ac440e672bcf.html + /// [zhetri]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_ga4d9cfa0653de400029b8051996ce1e96.html + /// fn invh(l: MatrixLayout, uplo: UPLO, a: &mut [Self], ipiv: &Pivot) -> Result<()>; - /// Wrapper of `*sytrs` and `*hetrs` + + /// Solve linear equation $Ax = b$ using factroized result + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [ssytrs] | [dsytrs] | [chetrs] | [zhetrs] | + /// + /// [ssytrs]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_gae20133a1119b69a7319783ff982c8c62.html + /// [dsytrs]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_ga6a223e61effac7232e98b422f147f032.html + /// [chetrs]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_ga6f9d8da222ffaa7b7535efc922faa1dc.html + /// [zhetrs]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_gacf697e3bb72c5fd88cd90972999401dd.html + /// fn solveh(l: MatrixLayout, uplo: UPLO, a: &[Self], ipiv: &Pivot, b: &mut [Self]) -> Result<()>; } From 26251681026bd09719a55fffd1bb605f83007a39 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Wed, 7 Sep 2022 14:05:47 +0900 Subject: [PATCH 05/13] About trait design in crate level document --- lax/src/lib.rs | 67 +++++++++++++++++++++++++++++++++----------------- 1 file changed, 45 insertions(+), 22 deletions(-) diff --git a/lax/src/lib.rs b/lax/src/lib.rs index 83cf3658..2c78cf82 100644 --- a/lax/src/lib.rs +++ b/lax/src/lib.rs @@ -1,32 +1,55 @@ -//! Linear Algebra eXtension (LAX) -//! =============================== -//! //! ndarray-free safe Rust wrapper for LAPACK FFI //! -//! Linear equation, Inverse matrix, Condition number -//! -------------------------------------------------- +//! `Lapack` trait and sub-traits +//! ------------------------------- +//! +//! This crates provides LAPACK wrapper as `impl` of traits to base scalar types. +//! For example, LU factorization to double-precision matrix is provided like: +//! +//! ```ignore +//! impl Solve_ for f64 { +//! fn lu(l: MatrixLayout, a: &mut [Self]) -> Result { ... } +//! } +//! ``` +//! +//! see [Solve_] for detail. You can use it like `f64::lu`: //! -//! As the property of $A$, several types of triangular factorization are used: +//! ``` +//! use lax::{Solve_, layout::MatrixLayout, Transpose}; //! -//! - LU-decomposition for general matrix -//! - $PA = LU$, where $L$ is lower matrix, $U$ is upper matrix, and $P$ is permutation matrix -//! - Bunch-Kaufman diagonal pivoting method for nonpositive-definite Hermitian matrix -//! - $A = U D U^\dagger$, where $U$ is upper matrix, -//! $D$ is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. +//! let mut a = vec![ +//! 1.0, 2.0, +//! 3.0, 4.0 +//! ]; +//! let mut b = vec![1.0, 2.0]; +//! let layout = MatrixLayout::C { row: 2, lda: 2 }; +//! let pivot = f64::lu(layout, &mut a).unwrap(); +//! f64::solve(layout, Transpose::No, &a, &pivot, &mut b).unwrap(); +//! ``` //! -//! | matrix type | Triangler factorization (TRF) | Solve (TRS) | Inverse matrix (TRI) | Reciprocal condition number (CON) | -//! |:--------------------------------|:------------------------------|:------------|:---------------------|:----------------------------------| -//! | General (GE) | [lu] | [solve] | [inv] | [rcond] | -//! | Symmetric (SY) / Hermitian (HE) | [bk] | [solveh] | [invh] | - | +//! When you want to write generic algorithm for real and complex matrices, +//! this trait can be used as a trait bound: +//! +//! ``` +//! use lax::{Solve_, layout::MatrixLayout, Transpose}; +//! +//! fn solve_at_once(layout: MatrixLayout, a: &mut [T], b: &mut [T]) -> Result<(), lax::error::Error> { +//! let pivot = T::lu(layout, a)?; +//! T::solve(layout, Transpose::No, a, &pivot, b)?; +//! Ok(()) +//! } +//! ``` +//! +//! There are several similar traits as described below to keep development easy. +//! They are merged into a single trait, [Lapack]. +//! +//! Linear equation, Inverse matrix, Condition number +//! -------------------------------------------------- //! -//! [lu]: solve/trait.Solve_.html#tymethod.lu -//! [solve]: solve/trait.Solve_.html#tymethod.solve -//! [inv]: solve/trait.Solve_.html#tymethod.inv -//! [rcond]: solve/trait.Solve_.html#tymethod.rcond +//! According to the property input metrix, several types of triangular factorization are used: //! -//! [bk]: solveh/trait.Solveh_.html#tymethod.bk -//! [solveh]: solveh/trait.Solveh_.html#tymethod.solveh -//! [invh]: solveh/trait.Solveh_.html#tymethod.invh +//! - [Solve_] trait provides methods for LU-decomposition for general matrix. +//! - [Solveh_] triat provides methods for Bunch-Kaufman diagonal pivoting method for symmetric/hermite indefinite matrix //! //! Eigenvalue Problem //! ------------------- From e94e03d7467ccb3431804318b9bc546a8d61e1c6 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 12:21:13 +0900 Subject: [PATCH 06/13] Document for Cholesky decomposition --- lax/src/cholesky.rs | 61 +++++++++++++++++++++++++++++++++++++++------ lax/src/lib.rs | 3 ++- 2 files changed, 56 insertions(+), 8 deletions(-) diff --git a/lax/src/cholesky.rs b/lax/src/cholesky.rs index 94e683ff..04875b64 100644 --- a/lax/src/cholesky.rs +++ b/lax/src/cholesky.rs @@ -1,21 +1,68 @@ -//! Cholesky decomposition - use super::*; use crate::{error::*, layout::*}; use cauchy::*; +#[cfg_attr(doc, katexit::katexit)] +/// Solve symmetric/hermite positive-definite linear equations using Cholesky decomposition +/// +/// For a given positive definite matrix $A$, +/// Cholesky decomposition is described as $A = U^T U$ or $A = LL^T$ where +/// +/// - $L$ is lower matrix +/// - $U$ is upper matrix +/// +/// This is designed as two step computation according to LAPACK API +/// +/// 1. Factorize input matrix $A$ into $L$ or $U$ +/// 2. Solve linear equation $Ax = b$ or compute inverse matrix $A^{-1}$ +/// using $U$ or $L$. pub trait Cholesky_: Sized { - /// Cholesky: wrapper of `*potrf` + /// Compute Cholesky decomposition $A = U^T U$ or $A = L L^T$ according to [UPLO] + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [spotrf] | [dpotrf] | [cpotrf] | [zpotrf] | + /// + /// [spotrf]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_gaaf31db7ab15b4f4ba527a3d31a15a58e.html + /// [dpotrf]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga2f55f604a6003d03b5cd4a0adcfb74d6.html + /// [cpotrf]: https://netlib.org/lapack/explore-html/d8/d6c/group__variants_p_ocomputational_ga4e85f48dbd837ccbbf76aa077f33de19.html + /// [zpotrf]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_ga93e22b682170873efb50df5a79c5e4eb.html /// - /// **Warning: Only the portion of `a` corresponding to `UPLO` is written.** fn cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()>; - /// Wrapper of `*potri` + /// Compute inverse matrix $A^{-1}$ using $U$ or $L$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [spotri] | [dpotri] | [cpotri] | [zpotri] | + /// + /// [spotri]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_ga4c381894bb34b1583fcc0dceafc5bea1.html + /// [dpotri]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga9dfc04beae56a3b1c1f75eebc838c14c.html + /// [cpotri]: https://netlib.org/lapack/explore-html/d6/df6/group__complex_p_ocomputational_ga52b8da4d314abefaee93dd5c1ed7739e.html + /// [zpotri]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_gaf37e3b8bbacd3332e83ffb3f1018bcf1.html /// - /// **Warning: Only the portion of `a` corresponding to `UPLO` is written.** fn inv_cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()>; - /// Wrapper of `*potrs` + /// Solve linear equation $Ax = b$ using $U$ or $L$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:---------|:---------|:---------|:---------| + /// | [spotrs] | [dpotrs] | [cpotrs] | [zpotrs] | + /// + /// [spotrs]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_gaf5cc1531aa5ffe706533fbca343d55dd.html + /// [dpotrs]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga167aa0166c4ce726385f65e4ab05e7c1.html + /// [cpotrs]: https://netlib.org/lapack/explore-html/d6/df6/group__complex_p_ocomputational_gad9052b4b70569dfd6e8943971c9b38b2.html + /// [zpotrs]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_gaa2116ea574b01efda584dff0b74c9fcd.html + /// fn solve_cholesky(l: MatrixLayout, uplo: UPLO, a: &[Self], b: &mut [Self]) -> Result<()>; } diff --git a/lax/src/lib.rs b/lax/src/lib.rs index 2c78cf82..d8d2d831 100644 --- a/lax/src/lib.rs +++ b/lax/src/lib.rs @@ -49,7 +49,8 @@ //! According to the property input metrix, several types of triangular factorization are used: //! //! - [Solve_] trait provides methods for LU-decomposition for general matrix. -//! - [Solveh_] triat provides methods for Bunch-Kaufman diagonal pivoting method for symmetric/hermite indefinite matrix +//! - [Solveh_] triat provides methods for Bunch-Kaufman diagonal pivoting method for symmetric/hermite indefinite matrix. +//! - [Cholesky_] triat provides methods for Cholesky decomposition for symmetric/hermite positive dinite matrix. //! //! Eigenvalue Problem //! ------------------- From 7432f0237569327cec5561a463d2cfcb6b7ad099 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 14:31:13 +0900 Subject: [PATCH 07/13] Unify LAPACK correspondance table --- lax/src/cholesky.rs | 33 +++++++++------------------------ lax/src/solve.rs | 27 +++++++++++++++++++++++++-- lax/src/solveh.rs | 33 +++++++++------------------------ 3 files changed, 43 insertions(+), 50 deletions(-) diff --git a/lax/src/cholesky.rs b/lax/src/cholesky.rs index 04875b64..9b213246 100644 --- a/lax/src/cholesky.rs +++ b/lax/src/cholesky.rs @@ -22,14 +22,9 @@ pub trait Cholesky_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [spotrf] | [dpotrf] | [cpotrf] | [zpotrf] | - /// - /// [spotrf]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_gaaf31db7ab15b4f4ba527a3d31a15a58e.html - /// [dpotrf]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga2f55f604a6003d03b5cd4a0adcfb74d6.html - /// [cpotrf]: https://netlib.org/lapack/explore-html/d8/d6c/group__variants_p_ocomputational_ga4e85f48dbd837ccbbf76aa077f33de19.html - /// [zpotrf]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_ga93e22b682170873efb50df5a79c5e4eb.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | spotrf | dpotrf | cpotrf | zpotrf | /// fn cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()>; @@ -38,14 +33,9 @@ pub trait Cholesky_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [spotri] | [dpotri] | [cpotri] | [zpotri] | - /// - /// [spotri]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_ga4c381894bb34b1583fcc0dceafc5bea1.html - /// [dpotri]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga9dfc04beae56a3b1c1f75eebc838c14c.html - /// [cpotri]: https://netlib.org/lapack/explore-html/d6/df6/group__complex_p_ocomputational_ga52b8da4d314abefaee93dd5c1ed7739e.html - /// [zpotri]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_gaf37e3b8bbacd3332e83ffb3f1018bcf1.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | spotri | dpotri | cpotri | zpotri | /// fn inv_cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()>; @@ -54,14 +44,9 @@ pub trait Cholesky_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [spotrs] | [dpotrs] | [cpotrs] | [zpotrs] | - /// - /// [spotrs]: https://netlib.org/lapack/explore-html/d8/db2/group__real_p_ocomputational_gaf5cc1531aa5ffe706533fbca343d55dd.html - /// [dpotrs]: https://netlib.org/lapack/explore-html/d1/d7a/group__double_p_ocomputational_ga167aa0166c4ce726385f65e4ab05e7c1.html - /// [cpotrs]: https://netlib.org/lapack/explore-html/d6/df6/group__complex_p_ocomputational_gad9052b4b70569dfd6e8943971c9b38b2.html - /// [zpotrs]: https://netlib.org/lapack/explore-html/d3/d8d/group__complex16_p_ocomputational_gaa2116ea574b01efda584dff0b74c9fcd.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | spotrs | dpotrs | cpotrs | zpotrs | /// fn solve_cholesky(l: MatrixLayout, uplo: UPLO, a: &[Self], b: &mut [Self]) -> Result<()>; } diff --git a/lax/src/solve.rs b/lax/src/solve.rs index ab601031..565aab59 100644 --- a/lax/src/solve.rs +++ b/lax/src/solve.rs @@ -32,12 +32,35 @@ pub trait Solve_: Scalar + Sized { /// - On this case, `return_code` in [Error::LapackComputationalFailure] means /// `return_code`-th diagonal element of $U$ becomes zero. /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | sgetrf | dgetrf | cgetrf | zgetrf | + /// fn lu(l: MatrixLayout, a: &mut [Self]) -> Result; - /// Compute inverse matrix $A^{-1}$ from the output of LU-factorization + /// Compute inverse matrix $A^{-1}$ from the output of LU-decomposition + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | sgetri | dgetri | cgetri | zgetri | + /// fn inv(l: MatrixLayout, a: &mut [Self], p: &Pivot) -> Result<()>; - /// Solve linear equations $Ax = b$ using the output of LU-factorization + /// Solve linear equations $Ax = b$ using the output of LU-decomposition + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | sgetrs | dgetrs | cgetrs | zgetrs | + /// fn solve(l: MatrixLayout, t: Transpose, a: &[Self], p: &Pivot, b: &mut [Self]) -> Result<()>; } diff --git a/lax/src/solveh.rs b/lax/src/solveh.rs index 76829982..5d1a3d02 100644 --- a/lax/src/solveh.rs +++ b/lax/src/solveh.rs @@ -24,14 +24,9 @@ pub trait Solveh_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [ssytrf] | [dsytrf] | [chetrf] | [zhetrf] | - /// - /// [ssytrf]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_ga12d2e56511cf7df066712c61d9acec45.html - /// [dsytrf]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_gad91bde1212277b3e909eb6af7f64858a.html - /// [chetrf]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_ga081dd1908e46d064c2bf0a1f6b664b86.html - /// [zhetrf]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_gadc84a5c9818ee12ea19944623131bd52.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | ssytrf | dsytrf | chetrf | zhetrf | /// fn bk(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result; @@ -40,14 +35,9 @@ pub trait Solveh_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [ssytri] | [dsytri] | [chetri] | [zhetri] | - /// - /// [ssytri]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_gaef378ec0761234aac417f487b43b7a8b.html - /// [dsytri]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_ga75e09b4299b7955044a3bbf84c46b593.html - /// [chetri]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_gad87a6a1ac131c5635d47ac440e672bcf.html - /// [zhetri]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_ga4d9cfa0653de400029b8051996ce1e96.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | ssytri | dsytri | chetri | zhetri | /// fn invh(l: MatrixLayout, uplo: UPLO, a: &mut [Self], ipiv: &Pivot) -> Result<()>; @@ -56,14 +46,9 @@ pub trait Solveh_: Sized { /// LAPACK correspondance /// ---------------------- /// - /// | f32 | f64 | c32 | c64 | - /// |:---------|:---------|:---------|:---------| - /// | [ssytrs] | [dsytrs] | [chetrs] | [zhetrs] | - /// - /// [ssytrs]: https://netlib.org/lapack/explore-html/d0/d14/group__real_s_ycomputational_gae20133a1119b69a7319783ff982c8c62.html - /// [dsytrs]: https://netlib.org/lapack/explore-html/d3/db6/group__double_s_ycomputational_ga6a223e61effac7232e98b422f147f032.html - /// [chetrs]: https://netlib.org/lapack/explore-html/d4/d74/group__complex_h_ecomputational_ga6f9d8da222ffaa7b7535efc922faa1dc.html - /// [zhetrs]: https://netlib.org/lapack/explore-html/d3/d80/group__complex16_h_ecomputational_gacf697e3bb72c5fd88cd90972999401dd.html + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | ssytrs | dsytrs | chetrs | zhetrs | /// fn solveh(l: MatrixLayout, uplo: UPLO, a: &[Self], ipiv: &Pivot, b: &mut [Self]) -> Result<()>; } From 4fe45f991ff83e107bb95307c735279d49dfd678 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 14:31:28 +0900 Subject: [PATCH 08/13] s/factorization/decomposition/g --- lax/src/lib.rs | 4 ++-- lax/src/solve.rs | 8 ++++---- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/lax/src/lib.rs b/lax/src/lib.rs index d8d2d831..d7b1be4a 100644 --- a/lax/src/lib.rs +++ b/lax/src/lib.rs @@ -4,7 +4,7 @@ //! ------------------------------- //! //! This crates provides LAPACK wrapper as `impl` of traits to base scalar types. -//! For example, LU factorization to double-precision matrix is provided like: +//! For example, LU decomposition to double-precision matrix is provided like: //! //! ```ignore //! impl Solve_ for f64 { @@ -46,7 +46,7 @@ //! Linear equation, Inverse matrix, Condition number //! -------------------------------------------------- //! -//! According to the property input metrix, several types of triangular factorization are used: +//! According to the property input metrix, several types of triangular decomposition are used: //! //! - [Solve_] trait provides methods for LU-decomposition for general matrix. //! - [Solveh_] triat provides methods for Bunch-Kaufman diagonal pivoting method for symmetric/hermite indefinite matrix. diff --git a/lax/src/solve.rs b/lax/src/solve.rs index 565aab59..9f25c9bf 100644 --- a/lax/src/solve.rs +++ b/lax/src/solve.rs @@ -5,7 +5,7 @@ use num_traits::{ToPrimitive, Zero}; #[cfg_attr(doc, katexit::katexit)] /// Solve linear equations using LU-decomposition /// -/// For a given matrix $A$, LU decomposition is described as $PA = LU$ where +/// For a given matrix $A$, LU decomposition is described as $A = PLU$ where: /// /// - $L$ is lower matrix /// - $U$ is upper matrix @@ -15,15 +15,15 @@ use num_traits::{ToPrimitive, Zero}; /// /// 1. Factorize input matrix $A$ into $L$, $U$, and $P$. /// 2. Solve linear equation $Ax = b$ or compute inverse matrix $A^{-1}$ -/// using the output of LU factorization. +/// using the output of LU decomposition. /// pub trait Solve_: Scalar + Sized { - /// Computes the LU factorization of a general $m \times n$ matrix + /// Computes the LU decomposition of a general $m \times n$ matrix /// with partial pivoting with row interchanges. /// /// Output /// ------- - /// - $U$ and $L$ are stored in `a` after LU factorization has succeeded. + /// - $U$ and $L$ are stored in `a` after LU decomposition has succeeded. /// - $P$ is returned as [Pivot] /// /// Error From 0d89696569e247e9fc2e5e11445bafb2f4563d7f Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 14:43:50 +0900 Subject: [PATCH 09/13] Update document of `Eig_` and `Eigh_` traits --- lax/src/eig.rs | 15 +++++++++++---- lax/src/eigh.rs | 24 ++++++++++++++++++++---- lax/src/lib.rs | 19 ++++--------------- 3 files changed, 35 insertions(+), 23 deletions(-) diff --git a/lax/src/eig.rs b/lax/src/eig.rs index d3c07222..deafba2f 100644 --- a/lax/src/eig.rs +++ b/lax/src/eig.rs @@ -1,12 +1,19 @@ -//! Eigenvalue decomposition for general matrices - use crate::{error::*, layout::MatrixLayout, *}; use cauchy::*; use num_traits::{ToPrimitive, Zero}; -/// Wraps `*geev` for general matrices +#[cfg_attr(doc, katexit::katexit)] +/// Eigenvalue problem for general matrix pub trait Eig_: Scalar { - /// Calculate Right eigenvalue + /// Compute right eigenvalue and eigenvectors $Ax = \lambda x$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:------|:------|:------|:------| + /// | sgeev | dgeev | cgeev | zgeev | + /// fn eig( calc_v: bool, l: MatrixLayout, diff --git a/lax/src/eigh.rs b/lax/src/eigh.rs index a9406ee6..055b2f47 100644 --- a/lax/src/eigh.rs +++ b/lax/src/eigh.rs @@ -1,12 +1,20 @@ -//! Eigenvalue decomposition for Symmetric/Hermite matrices - use super::*; use crate::{error::*, layout::MatrixLayout}; use cauchy::*; use num_traits::{ToPrimitive, Zero}; +#[cfg_attr(doc, katexit::katexit)] +/// Eigenvalue problem for symmetric/hermite matrix pub trait Eigh_: Scalar { - /// Wraps `*syev` for real and `*heev` for complex + /// Compute right eigenvalue and eigenvectors $Ax = \lambda x$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:------|:------|:------|:------| + /// | ssyev | dsyev | cheev | zheev | + /// fn eigh( calc_eigenvec: bool, layout: MatrixLayout, @@ -14,7 +22,15 @@ pub trait Eigh_: Scalar { a: &mut [Self], ) -> Result>; - /// Wraps `*sygv` for real and `*hegv` for complex + /// Compute generalized right eigenvalue and eigenvectors $Ax = \lambda B x$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:------|:------|:------|:------| + /// | ssygv | dsygv | chegv | zhegv | + /// fn eigh_generalized( calc_eigenvec: bool, layout: MatrixLayout, diff --git a/lax/src/lib.rs b/lax/src/lib.rs index d7b1be4a..9a8bf1cf 100644 --- a/lax/src/lib.rs +++ b/lax/src/lib.rs @@ -55,22 +55,11 @@ //! Eigenvalue Problem //! ------------------- //! -//! Solve eigenvalue problem for a matrix $A$ +//! According to the property input metrix, +//! there are several types of eigenvalue problem API //! -//! $$ Av_i = \lambda_i v_i $$ -//! -//! or generalized eigenvalue problem -//! -//! $$ Av_i = \lambda_i B v_i $$ -//! -//! | matrix type | Eigenvalue (EV) | Generalized Eigenvalue Problem (EG) | -//! |:--------------------------------|:----------------|:------------------------------------| -//! | General (GE) |[eig] | - | -//! | Symmetric (SY) / Hermitian (HE) |[eigh] |[eigh_generalized] | -//! -//! [eig]: eig/trait.Eig_.html#tymethod.eig -//! [eigh]: eigh/trait.Eigh_.html#tymethod.eigh -//! [eigh_generalized]: eigh/trait.Eigh_.html#tymethod.eigh_generalized +//! - [Eig_] trait provides methods for eigenvalue problem for general matrix. +//! - [Eigh_] trait provides methods for eigenvalue problem for symmetric/hermite matrix. //! //! Singular Value Decomposition (SVD), Least square problem //! ---------------------------------------------------------- From 94d884346cbca4baed01f5800adbd61d43b33581 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 15:58:11 +0900 Subject: [PATCH 10/13] Update document about SVD --- lax/src/least_squares.rs | 6 +++++- lax/src/lib.rs | 15 +++++++-------- lax/src/svd.rs | 13 +++++++++++-- lax/src/svddc.rs | 11 +++++++++++ 4 files changed, 34 insertions(+), 11 deletions(-) diff --git a/lax/src/least_squares.rs b/lax/src/least_squares.rs index 6be44f33..5699257d 100644 --- a/lax/src/least_squares.rs +++ b/lax/src/least_squares.rs @@ -12,14 +12,18 @@ pub struct LeastSquaresOutput { pub rank: i32, } -/// Wraps `*gelsd` +#[cfg_attr(doc, katexit::katexit)] +/// Solve least square problem pub trait LeastSquaresSvdDivideConquer_: Scalar { + /// Compute a vector $x$ which minimizes Euclidian norm $\| Ax - b\|$ + /// for a given matrix $A$ and a vector $b$. fn least_squares( a_layout: MatrixLayout, a: &mut [Self], b: &mut [Self], ) -> Result>; + /// Solve least square problems $\argmin_X \| AX - B\|$ fn least_squares_nrhs( a_layout: MatrixLayout, a: &mut [Self], diff --git a/lax/src/lib.rs b/lax/src/lib.rs index 9a8bf1cf..47f3eedc 100644 --- a/lax/src/lib.rs +++ b/lax/src/lib.rs @@ -61,16 +61,15 @@ //! - [Eig_] trait provides methods for eigenvalue problem for general matrix. //! - [Eigh_] trait provides methods for eigenvalue problem for symmetric/hermite matrix. //! -//! Singular Value Decomposition (SVD), Least square problem -//! ---------------------------------------------------------- +//! Singular Value Decomposition +//! ----------------------------- //! -//! | matrix type | Singular Value Decomposition (SVD) | SVD with divided-and-conquer (SDD) | Least square problem (LSD) | -//! |:-------------|:-----------------------------------|:-----------------------------------|:---------------------------| -//! | General (GE) | [svd] | [svddc] | [least_squares] | +//! - [SVD_] trait provides methods for singular value decomposition for general matrix +//! - [SVDDC_] trait provides methods for singular value decomposition for general matrix +//! with divided-and-conquer algorithm +//! - [LeastSquaresSvdDivideConquer_] trait provides methods +//! for solving least square problem by SVD //! -//! [svd]: svd/trait.SVD_.html#tymethod.svd -//! [svddc]: svddck/trait.SVDDC_.html#tymethod.svddc -//! [least_squares]: least_squares/trait.LeastSquaresSvdDivideConquer_.html#tymethod.least_squares #[cfg(any(feature = "intel-mkl-system", feature = "intel-mkl-static"))] extern crate intel_mkl_src as _src; diff --git a/lax/src/svd.rs b/lax/src/svd.rs index 0a509a0e..0f0f81ab 100644 --- a/lax/src/svd.rs +++ b/lax/src/svd.rs @@ -14,9 +14,18 @@ pub struct SVDOutput { pub vt: Option>, } -/// Wraps `*gesvd` +#[cfg_attr(doc, katexit::katexit)] +/// Singular value decomposition pub trait SVD_: Scalar { - /// Calculate singular value decomposition $ A = U \Sigma V^T $ + /// Compute singular value decomposition $A = U \Sigma V^T$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | sgesvd | dgesvd | cgesvd | zgesvd | + /// fn svd(l: MatrixLayout, calc_u: bool, calc_vt: bool, a: &mut [Self]) -> Result>; } diff --git a/lax/src/svddc.rs b/lax/src/svddc.rs index bb59348f..55ddfd46 100644 --- a/lax/src/svddc.rs +++ b/lax/src/svddc.rs @@ -2,7 +2,18 @@ use crate::{error::*, layout::MatrixLayout, *}; use cauchy::*; use num_traits::{ToPrimitive, Zero}; +#[cfg_attr(doc, katexit::katexit)] +/// Singular value decomposition with divide-and-conquer method pub trait SVDDC_: Scalar { + /// Compute singular value decomposition $A = U \Sigma V^T$ + /// + /// LAPACK correspondance + /// ---------------------- + /// + /// | f32 | f64 | c32 | c64 | + /// |:-------|:-------|:-------|:-------| + /// | sgesdd | dgesdd | cgesdd | zgesdd | + /// fn svddc(l: MatrixLayout, jobz: JobSvd, a: &mut [Self]) -> Result>; } From eaec8e9dd02af794fefa8a90960b7682846365e8 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 17:09:55 +0900 Subject: [PATCH 11/13] Remove KaTeX auto render header for lax --- lax/Cargo.toml | 3 --- lax/katex-header.html | 16 ---------------- ndarray-linalg/Cargo.toml | 3 --- ndarray-linalg/katex-header.html | 16 ---------------- 4 files changed, 38 deletions(-) delete mode 100644 lax/katex-header.html delete mode 100644 ndarray-linalg/katex-header.html diff --git a/lax/Cargo.toml b/lax/Cargo.toml index 26a6f35e..f46bdbb8 100644 --- a/lax/Cargo.toml +++ b/lax/Cargo.toml @@ -51,6 +51,3 @@ version = "0.10.4" optional = true default-features = false features = ["cblas"] - -[package.metadata.docs.rs] -rustdoc-args = ["--html-in-header", "katex-header.html"] diff --git a/lax/katex-header.html b/lax/katex-header.html deleted file mode 100644 index 6e10c052..00000000 --- a/lax/katex-header.html +++ /dev/null @@ -1,16 +0,0 @@ - - - - - diff --git a/ndarray-linalg/Cargo.toml b/ndarray-linalg/Cargo.toml index 6147a2eb..4ef8c17d 100644 --- a/ndarray-linalg/Cargo.toml +++ b/ndarray-linalg/Cargo.toml @@ -78,6 +78,3 @@ harness = false [[bench]] name = "solveh" harness = false - -[package.metadata.docs.rs] -rustdoc-args = ["--html-in-header", "katex-header.html"] diff --git a/ndarray-linalg/katex-header.html b/ndarray-linalg/katex-header.html deleted file mode 100644 index 6e10c052..00000000 --- a/ndarray-linalg/katex-header.html +++ /dev/null @@ -1,16 +0,0 @@ - - - - - From 25703b83a1ebb6a831a4e727658550a4881e0caf Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 17:20:11 +0900 Subject: [PATCH 12/13] Add katexit for ndarray-linalg --- ndarray-linalg/Cargo.toml | 1 + ndarray-linalg/src/eig.rs | 1 + 2 files changed, 2 insertions(+) diff --git a/ndarray-linalg/Cargo.toml b/ndarray-linalg/Cargo.toml index 4ef8c17d..bfab605f 100644 --- a/ndarray-linalg/Cargo.toml +++ b/ndarray-linalg/Cargo.toml @@ -30,6 +30,7 @@ intel-mkl-system = ["lax/intel-mkl-system"] [dependencies] cauchy = "0.4.0" +katexit = "0.1.2" num-complex = "0.4.0" num-traits = "0.2.14" rand = "0.8.3" diff --git a/ndarray-linalg/src/eig.rs b/ndarray-linalg/src/eig.rs index 833f6d49..bde5bd7d 100644 --- a/ndarray-linalg/src/eig.rs +++ b/ndarray-linalg/src/eig.rs @@ -5,6 +5,7 @@ use crate::layout::*; use crate::types::*; use ndarray::*; +#[cfg_attr(doc, katexit::katexit)] /// Eigenvalue decomposition of general matrix reference pub trait Eig { /// EigVec is the right eivenvector From c25cc23da52b146e8579eed6fee3f3866b3a1ca4 Mon Sep 17 00:00:00 2001 From: Toshiki Teramura Date: Thu, 8 Sep 2022 17:20:54 +0900 Subject: [PATCH 13/13] Drop RUSTDOCFLAGS from GitHub Pages --- .github/workflows/gh-pages.yml | 2 -- 1 file changed, 2 deletions(-) diff --git a/.github/workflows/gh-pages.yml b/.github/workflows/gh-pages.yml index 90ad7b90..f8509459 100644 --- a/.github/workflows/gh-pages.yml +++ b/.github/workflows/gh-pages.yml @@ -40,8 +40,6 @@ jobs: with: command: doc args: --no-deps - env: - RUSTDOCFLAGS: "--html-in-header ndarray-linalg/katex-header.html" - name: Setup Pages uses: actions/configure-pages@v2