Skip to content

Split LAPACK error into Computational failure and invalid value #210

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 1 commit into from
Jul 1, 2020
Merged

Split LAPACK error into Computational failure and invalid value #210

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 11 additions & 3 deletions ndarray-linalg/src/error.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,9 +12,17 @@ pub enum LinalgError {
#[error("Not square: rows({}) != cols({})", rows, cols)]
NotSquare { rows: i32, cols: i32 },

/// LAPACK subroutine returns non-zero code
#[error("LAPACK: return_code = {}", return_code)]
Lapack { return_code: i32 },
#[error(
"Invalid value for LAPACK subroutine {}-th argument",
-return_code
)]
LapackInvalidValue { return_code: i32 },

#[error(
"Comutational failure in LAPACK subroutine: return_code = {}",
return_code
)]
LapackComputationalFailure { return_code: i32 },

/// Strides of the array is not supported
#[error("invalid stride: s0={}, s1={}", s0, s1)]
Expand Down
20 changes: 9 additions & 11 deletions ndarray-linalg/src/lapack/cholesky.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,10 +1,7 @@
//! Cholesky decomposition

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::{into_result, UPLO};
use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};

pub trait Cholesky_: Sized {
/// Cholesky: wrapper of `*potrf`
Expand All @@ -25,14 +22,14 @@ macro_rules! impl_cholesky {
impl Cholesky_ for $scalar {
unsafe fn cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()> {
let (n, _) = l.size();
let info = $trf(l.lapacke_layout(), uplo as u8, n, a, n);
into_result(info, ())
$trf(l.lapacke_layout(), uplo as u8, n, a, n).as_lapack_result()?;
Ok(())
}

unsafe fn inv_cholesky(l: MatrixLayout, uplo: UPLO, a: &mut [Self]) -> Result<()> {
let (n, _) = l.size();
let info = $tri(l.lapacke_layout(), uplo as u8, n, a, l.lda());
into_result(info, ())
$tri(l.lapacke_layout(), uplo as u8, n, a, l.lda()).as_lapack_result()?;
Ok(())
}

unsafe fn solve_cholesky(
Expand All @@ -44,8 +41,9 @@ macro_rules! impl_cholesky {
let (n, _) = l.size();
let nrhs = 1;
let ldb = 1;
let info = $trs(l.lapacke_layout(), uplo as u8, n, nrhs, a, l.lda(), b, ldb);
into_result(info, ())
$trs(l.lapacke_layout(), uplo as u8, n, nrhs, a, l.lda(), b, ldb)
.as_lapack_result()?;
Ok(())
}
}
};
Expand Down
19 changes: 9 additions & 10 deletions ndarray-linalg/src/lapack/eig.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,13 +1,9 @@
//! Eigenvalue decomposition for general matrices

use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};
use num_traits::Zero;

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::into_result;

/// Wraps `*geev` for real/complex
pub trait Eig_: Scalar {
unsafe fn eig(
Expand All @@ -30,7 +26,7 @@ macro_rules! impl_eig_complex {
let mut w = vec![Self::Complex::zero(); n as usize];
let mut vl = Vec::new();
let mut vr = vec![Self::Complex::zero(); (n * n) as usize];
let info = $ev(
$ev(
l.lapacke_layout(),
b'N',
jobvr,
Expand All @@ -42,8 +38,9 @@ macro_rules! impl_eig_complex {
n,
&mut vr,
n,
);
into_result(info, (w, vr))
)
.as_lapack_result()?;
Ok((w, vr))
}
}
};
Expand Down Expand Up @@ -82,6 +79,7 @@ macro_rules! impl_eig_real {
.zip(wi.iter())
.map(|(&r, &i)| Self::Complex::new(r, i))
.collect();

// If the j-th eigenvalue is real, then
// eigenvector = [ vr[j], vr[j+n], vr[j+2*n], ... ].
//
Expand Down Expand Up @@ -121,7 +119,8 @@ macro_rules! impl_eig_real {
})
.collect();

into_result(info, (w, v))
info.as_lapack_result()?;
Ok((w, v))
}
}
};
Expand Down
20 changes: 9 additions & 11 deletions ndarray-linalg/src/lapack/eigh.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,13 +1,9 @@
//! Eigenvalue decomposition for Hermite matrices

use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};
use num_traits::Zero;

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::{into_result, UPLO};

/// Wraps `*syev` for real and `*heev` for complex
pub trait Eigh_: Scalar {
unsafe fn eigh(
Expand Down Expand Up @@ -37,8 +33,9 @@ macro_rules! impl_eigh {
let (n, _) = l.size();
let jobz = if calc_v { b'V' } else { b'N' };
let mut w = vec![Self::Real::zero(); n as usize];
let info = $ev(l.lapacke_layout(), jobz, uplo as u8, n, &mut a, n, &mut w);
into_result(info, w)
$ev(l.lapacke_layout(), jobz, uplo as u8, n, &mut a, n, &mut w)
.as_lapack_result()?;
Ok(w)
}

unsafe fn eigh_generalized(
Expand All @@ -51,7 +48,7 @@ macro_rules! impl_eigh {
let (n, _) = l.size();
let jobz = if calc_v { b'V' } else { b'N' };
let mut w = vec![Self::Real::zero(); n as usize];
let info = $evg(
$evg(
l.lapacke_layout(),
1,
jobz,
Expand All @@ -62,8 +59,9 @@ macro_rules! impl_eigh {
&mut b,
n,
&mut w,
);
into_result(info, w)
)
.as_lapack_result()?;
Ok(w)
}
}
};
Expand Down
39 changes: 15 additions & 24 deletions ndarray-linalg/src/lapack/least_squares.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,14 +1,10 @@
//! Least squares

use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};
use ndarray::{ErrorKind, ShapeError};
use num_traits::Zero;

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::into_result;

/// Result of LeastSquares
pub struct LeastSquaresOutput<A: Scalar> {
/// singular values
Expand Down Expand Up @@ -53,7 +49,7 @@ macro_rules! impl_least_squares {
let mut singular_values: Vec<Self::Real> = vec![Self::Real::zero(); k as usize];
let mut rank: i32 = 0;

let status = $gelsd(
$gelsd(
a_layout.lapacke_layout(),
m,
n,
Expand All @@ -67,15 +63,13 @@ macro_rules! impl_least_squares {
&mut singular_values,
rcond,
&mut rank,
);

into_result(
status,
LeastSquaresOutput {
singular_values,
rank,
},
)
.as_lapack_result()?;

Ok(LeastSquaresOutput {
singular_values,
rank,
})
}

unsafe fn least_squares_nrhs(
Expand All @@ -99,7 +93,7 @@ macro_rules! impl_least_squares {
let mut singular_values: Vec<Self::Real> = vec![Self::Real::zero(); k as usize];
let mut rank: i32 = 0;

let status = $gelsd(
$gelsd(
a_layout.lapacke_layout(),
m,
n,
Expand All @@ -111,15 +105,12 @@ macro_rules! impl_least_squares {
&mut singular_values,
rcond,
&mut rank,
);

into_result(
status,
LeastSquaresOutput {
singular_values,
rank,
},
)
.as_lapack_result()?;
Ok(LeastSquaresOutput {
singular_values,
rank,
})
}
}
};
Expand Down
18 changes: 13 additions & 5 deletions ndarray-linalg/src/lapack/mod.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -54,11 +54,19 @@ impl Lapack for f64 {}
impl Lapack for c32 {}
impl Lapack for c64 {}

pub fn into_result<T>(return_code: i32, val: T) -> Result<T> {
if return_code == 0 {
Ok(val)
} else {
Err(LinalgError::Lapack { return_code })
trait AsLapackResult {
fn as_lapack_result(self) -> Result<()>;
}

impl AsLapackResult for i32 {
fn as_lapack_result(self) -> Result<()> {
if self > 0 {
return Err(LinalgError::LapackComputationalFailure { return_code: self });
}
if self < 0 {
return Err(LinalgError::LapackInvalidValue { return_code: self });
}
Ok(())
}
}

Expand Down
16 changes: 6 additions & 10 deletions ndarray-linalg/src/lapack/qr.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,14 +1,10 @@
//! QR decomposition

use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};
use num_traits::Zero;
use std::cmp::min;

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::into_result;

/// Wraps `*geqrf` and `*orgqr` (`*ungqr` for complex numbers)
pub trait QR_: Sized {
unsafe fn householder(l: MatrixLayout, a: &mut [Self]) -> Result<Vec<Self>>;
Expand All @@ -23,15 +19,15 @@ macro_rules! impl_qr {
let (row, col) = l.size();
let k = min(row, col);
let mut tau = vec![Self::zero(); k as usize];
let info = $qrf(l.lapacke_layout(), row, col, &mut a, l.lda(), &mut tau);
into_result(info, tau)
$qrf(l.lapacke_layout(), row, col, &mut a, l.lda(), &mut tau).as_lapack_result()?;
Ok(tau)
}

unsafe fn q(l: MatrixLayout, mut a: &mut [Self], tau: &[Self]) -> Result<()> {
let (row, col) = l.size();
let k = min(row, col);
let info = $gqr(l.lapacke_layout(), row, k, k, &mut a, l.lda(), &tau);
into_result(info, ())
$gqr(l.lapacke_layout(), row, k, k, &mut a, l.lda(), &tau).as_lapack_result()?;
Ok(())
}

unsafe fn qr(l: MatrixLayout, a: &mut [Self]) -> Result<Vec<Self>> {
Expand Down
31 changes: 14 additions & 17 deletions ndarray-linalg/src/lapack/solve.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,14 +1,9 @@
//! Solve linear problem using LU decomposition

use super::*;
use crate::{error::*, layout::MatrixLayout, types::*};
use num_traits::Zero;

use crate::error::*;
use crate::layout::MatrixLayout;
use crate::types::*;

use super::NormType;
use super::{into_result, Pivot, Transpose};

/// Wraps `*getrf`, `*getri`, and `*getrs`
pub trait Solve_: Scalar + Sized {
/// Computes the LU factorization of a general `m x n` matrix `a` using
Expand Down Expand Up @@ -41,29 +36,30 @@ macro_rules! impl_solve {
let (row, col) = l.size();
let k = ::std::cmp::min(row, col);
let mut ipiv = vec![0; k as usize];
let info = $getrf(l.lapacke_layout(), row, col, a, l.lda(), &mut ipiv);
into_result(info, ipiv)
$getrf(l.lapacke_layout(), row, col, a, l.lda(), &mut ipiv).as_lapack_result()?;
Ok(ipiv)
}

unsafe fn inv(l: MatrixLayout, a: &mut [Self], ipiv: &Pivot) -> Result<()> {
let (n, _) = l.size();
let info = $getri(l.lapacke_layout(), n, a, l.lda(), ipiv);
into_result(info, ())
$getri(l.lapacke_layout(), n, a, l.lda(), ipiv).as_lapack_result()?;
Ok(())
}

unsafe fn rcond(l: MatrixLayout, a: &[Self], anorm: Self::Real) -> Result<Self::Real> {
let (n, _) = l.size();
let mut rcond = Self::Real::zero();
let info = $gecon(
$gecon(
l.lapacke_layout(),
NormType::One as u8,
n,
a,
l.lda(),
anorm,
&mut rcond,
);
into_result(info, rcond)
)
.as_lapack_result()?;
Ok(rcond)
}

unsafe fn solve(
Expand All @@ -76,7 +72,7 @@ macro_rules! impl_solve {
let (n, _) = l.size();
let nrhs = 1;
let ldb = 1;
let info = $getrs(
$getrs(
l.lapacke_layout(),
t as u8,
n,
Expand All @@ -86,8 +82,9 @@ macro_rules! impl_solve {
ipiv,
b,
ldb,
);
into_result(info, ())
)
.as_lapack_result()?;
Ok(())
}
}
};
Expand Down
Loading