Split LAPACK error into Computational failure and invalid value

Author: termoshtt

ndarray-linalg/src/error.rs

@@ -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)]

ndarray-linalg/src/lapack/cholesky.rs

@@ -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`
@@ -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(
@@ -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(())
             }
         }
     };

ndarray-linalg/src/lapack/eig.rs

@@ -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(
@@ -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,
@@ -42,8 +38,9 @@ macro_rules! impl_eig_complex {
                     n,
                     &mut vr,
                     n,
-                );
-                into_result(info, (w, vr))
+                )
+                .as_lapack_result()?;
+                Ok((w, vr))
             }
         }
     };
@@ -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], ... ].
                 //
@@ -121,7 +119,8 @@ macro_rules! impl_eig_real {
                     })
                     .collect();
 
-                into_result(info, (w, v))
+                info.as_lapack_result()?;
+                Ok((w, v))
             }
         }
     };

ndarray-linalg/src/lapack/eigh.rs

@@ -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(
@@ -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(
@@ -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,
@@ -62,8 +59,9 @@ macro_rules! impl_eigh {
                     &mut b,
                     n,
                     &mut w,
-                );
-                into_result(info, w)
+                )
+                .as_lapack_result()?;
+                Ok(w)
             }
         }
     };

ndarray-linalg/src/lapack/least_squares.rs

@@ -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
@@ -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,
@@ -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(
@@ -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,
@@ -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,
+                })
             }
         }
     };

ndarray-linalg/src/lapack/mod.rs

@@ -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(())
     }
 }
 

ndarray-linalg/src/lapack/qr.rs

@@ -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>>;
@@ -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>> {

ndarray-linalg/src/lapack/solve.rs

@@ -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
@@ -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(
@@ -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,
@@ -86,8 +82,9 @@ macro_rules! impl_solve {
                     ipiv,
                     b,
                     ldb,
-                );
-                into_result(info, ())
+                )
+                .as_lapack_result()?;
+                Ok(())
             }
         }
     };