Skip to content

Commit

Permalink
Merge pull request #330 from rust-ndarray/katexit
Browse files Browse the repository at this point in the history
Revise lax document with katexit
  • Loading branch information
termoshtt committed Sep 8, 2022
2 parents a4b3118 + c25cc23 commit 991a744
Show file tree
Hide file tree
Showing 16 changed files with 257 additions and 121 deletions.
2 changes: 0 additions & 2 deletions .github/workflows/gh-pages.yml
Expand Up @@ -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
Expand Down
4 changes: 1 addition & 3 deletions lax/Cargo.toml
Expand Up @@ -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"
Expand All @@ -50,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"]
16 changes: 0 additions & 16 deletions lax/katex-header.html

This file was deleted.

46 changes: 39 additions & 7 deletions lax/src/cholesky.rs
@@ -1,21 +1,53 @@
//! 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 |
///
/// **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 |
///
/// **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 |
///
fn solve_cholesky(l: MatrixLayout, uplo: UPLO, a: &[Self], b: &mut [Self]) -> Result<()>;
}

Expand Down
15 changes: 11 additions & 4 deletions 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,
Expand Down
24 changes: 20 additions & 4 deletions lax/src/eigh.rs
@@ -1,20 +1,36 @@
//! 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,
uplo: UPLO,
a: &mut [Self],
) -> Result<Vec<Self::Real>>;

/// 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,
Expand Down
12 changes: 7 additions & 5 deletions lax/src/flags.rs
Expand Up @@ -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',
}

Expand Down
6 changes: 5 additions & 1 deletion lax/src/least_squares.rs
Expand Up @@ -12,14 +12,18 @@ pub struct LeastSquaresOutput<A: Scalar> {
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<LeastSquaresOutput<Self>>;

/// Solve least square problems $\argmin_X \| AX - B\|$
fn least_squares_nrhs(
a_layout: MatrixLayout,
a: &mut [Self],
Expand Down
98 changes: 55 additions & 43 deletions lax/src/lib.rs
@@ -1,63 +1,75 @@
//! Linear Algebra eXtension (LAX)
//! ===============================
//!
//! ndarray-free safe Rust wrapper for LAPACK FFI
//!
//! Linear equation, Inverse matrix, Condition number
//! --------------------------------------------------
//! `Lapack` trait and sub-traits
//! -------------------------------
//!
//! As the property of $A$, several types of triangular factorization are used:
//! This crates provides LAPACK wrapper as `impl` of traits to base scalar types.
//! For example, LU decomposition to double-precision matrix is provided like:
//!
//! - 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.
//! ```ignore
//! impl Solve_ for f64 {
//! fn lu(l: MatrixLayout, a: &mut [Self]) -> Result<Pivot> { ... }
//! }
//! ```
//!
//! | 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] | - |
//! see [Solve_] for detail. You can use it like `f64::lu`:
//!
//! [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
//! ```
//! use lax::{Solve_, layout::MatrixLayout, Transpose};
//!
//! [bk]: solveh/trait.Solveh_.html#tymethod.bk
//! [solveh]: solveh/trait.Solveh_.html#tymethod.solveh
//! [invh]: solveh/trait.Solveh_.html#tymethod.invh
//! 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();
//! ```
//!
//! Eigenvalue Problem
//! -------------------
//! 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};
//!
//! Solve eigenvalue problem for a matrix $A$
//! fn solve_at_once<T: Solve_>(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(())
//! }
//! ```
//!
//! $$ Av_i = \lambda_i v_i $$
//! There are several similar traits as described below to keep development easy.
//! They are merged into a single trait, [Lapack].
//!
//! or generalized eigenvalue problem
//! Linear equation, Inverse matrix, Condition number
//! --------------------------------------------------
//!
//! According to the property input metrix, several types of triangular decomposition are used:
//!
//! $$ Av_i = \lambda_i B v_i $$
//! - [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.
//! - [Cholesky_] triat provides methods for Cholesky decomposition for symmetric/hermite positive dinite matrix.
//!
//! Eigenvalue Problem
//! -------------------
//!
//! | matrix type | Eigenvalue (EV) | Generalized Eigenvalue Problem (EG) |
//! |:--------------------------------|:----------------|:------------------------------------|
//! | General (GE) |[eig] | - |
//! | Symmetric (SY) / Hermitian (HE) |[eigh] |[eigh_generalized] |
//! According to the property input metrix,
//! there are several types of eigenvalue problem API
//!
//! [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
//! ----------------------------------------------------------
//! 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;
Expand Down

0 comments on commit 991a744

Please sign in to comment.