additional svd functions

1 files changed, 157 insertions, 68 deletions

diff --git a/lib/Numeric/LinearAlgebra/LAPACK.hs b/lib/Numeric/LinearAlgebra/LAPACK.hs
index 30a3d3b..2eae9b7 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK.hs
+++ b/lib/Numeric/LinearAlgebra/LAPACK.hs
@@ -1,4 +1,3 @@
-{-# OPTIONS_GHC #-}
 -----------------------------------------------------------------------------
 -- |
 -- Module      :  Numeric.LinearAlgebra.LAPACK
@@ -9,21 +8,34 @@
 -- Stability   :  provisional
 -- Portability :  portable (uses FFI)
 --
--- Wrappers for a few LAPACK functions (<http://www.netlib.org/lapack>).
+-- Functional interface to selected LAPACK functions (<http://www.netlib.org/lapack>).
 --
 -----------------------------------------------------------------------------
 
 module Numeric.LinearAlgebra.LAPACK (
+    -- * Matrix product
     multiplyR, multiplyC,
-    svdR, svdRdd, svdC,
-    eigC, eigR, eigS, eigH, eigS', eigH',
+    -- * Linear systems
     linearSolveR, linearSolveC,
+    lusR, lusC,
     linearSolveLSR, linearSolveLSC,
     linearSolveSVDR, linearSolveSVDC,
-    luR, luC, lusR, lusC,
+    -- * SVD
+    svR, svRd, svC, svCd,
+    svdR, svdRd, svdC, svdCd,
+    thinSVDR, thinSVDRd, thinSVDC, thinSVDCd,
+    rightSVR, rightSVC, leftSVR, leftSVC,
+    -- * Eigensystems
+    eigR, eigC, eigS, eigS', eigH, eigH',
+    -- * LU
+    luR, luC,
+    -- * Cholesky
     cholS, cholH,
+    -- * QR
     qrR, qrC,
+    -- * Hessenberg
     hessR, hessC,
+    -- * Schur
     schurR, schurC
 ) where
 
@@ -61,28 +73,27 @@ multiplyC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Compl
 multiplyC a b = multiplyAux zgemmc "zgemmc" a b
 
 -----------------------------------------------------------------------------
-foreign import ccall "LAPACK/lapack-aux.h svd_l_R" dgesvd :: TMMVM
-foreign import ccall "LAPACK/lapack-aux.h svd_l_C" zgesvd :: TCMCMVCM
-foreign import ccall "LAPACK/lapack-aux.h svd_l_Rdd" dgesdd :: TMMVM
+foreign import ccall "svd_l_R" dgesvd :: TMMVM
+foreign import ccall "svd_l_C" zgesvd :: TCMCMVCM
+foreign import ccall "svd_l_Rdd" dgesdd :: TMMVM
+foreign import ccall "svd_l_Cdd" zgesdd :: TCMCMVCM
 
--- | Wrapper for LAPACK's /dgesvd/, which computes the full svd decomposition of a real matrix.
---
--- @(u,s,v)=full svdR m@ so that @m=u \<\> s \<\> 'trans' v@.
+-- | Full SVD of a real matrix using LAPACK's /dgesvd/.
 svdR :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
 svdR = svdAux dgesvd "svdR" . fmat
 
--- | Wrapper for LAPACK's /dgesvd/, which computes the full svd decomposition of a real matrix.
---
--- @(u,s,v)=full svdRdd m@ so that @m=u \<\> s \<\> 'trans' v@.
-svdRdd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
-svdRdd = svdAux dgesdd "svdRdd" . fmat
+-- | Full SVD of a real matrix using LAPACK's /dgesdd/.
+svdRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
+svdRd = svdAux dgesdd "svdRdd" . fmat
 
--- | Wrapper for LAPACK's /zgesvd/, which computes the full svd decomposition of a complex matrix.
---
--- @(u,s,v)=full svdC m@ so that @m=u \<\> comp s \<\> 'trans' v@.
+-- | Full SVD of a complex matrix using LAPACK's /zgesvd/.
 svdC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
 svdC = svdAux zgesvd "svdC" . fmat
 
+-- | Full SVD of a complex matrix using LAPACK's /zgesdd/.
+svdCd :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
+svdCd = svdAux zgesdd "svdCdd" . fmat
+
 svdAux f st x = unsafePerformIO $ do
     u <- createMatrix ColumnMajor r r
     s <- createVector (min r c)
@@ -92,7 +103,104 @@ svdAux f st x = unsafePerformIO $ do
   where r = rows x
         c = cols x
 
+
+-- | Thin SVD of a real matrix, using LAPACK's /dgesvd/ with jobu == jobvt == \'S\'.
+thinSVDR :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
+thinSVDR = thinSVDAux dgesvd "thinSVDR" . fmat
+
+-- | Thin SVD of a complex matrix, using LAPACK's /zgesvd/ with jobu == jobvt == \'S\'.
+thinSVDC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
+thinSVDC = thinSVDAux zgesvd "thinSVDC" . fmat
+
+-- | Thin SVD of a real matrix, using LAPACK's /dgesdd/ with jobz == \'S\'.
+thinSVDRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
+thinSVDRd = thinSVDAux dgesdd "thinSVDRdd" . fmat
+
+-- | Thin SVD of a complex matrix, using LAPACK's /zgesdd/ with jobz == \'S\'.
+thinSVDCd :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
+thinSVDCd = thinSVDAux zgesdd "thinSVDCdd" . fmat
+
+thinSVDAux f st x = unsafePerformIO $ do
+    u <- createMatrix ColumnMajor r q
+    s <- createVector q
+    v <- createMatrix ColumnMajor q c
+    app4 f mat x mat u vec s mat v st
+    return (u,s,trans v)
+  where r = rows x
+        c = cols x
+        q = min r c
+
+
+-- | Singular values of a real matrix, using LAPACK's /dgesvd/ with jobu == jobvt == \'N\'.
+svR :: Matrix Double -> Vector Double
+svR = svAux dgesvd "svR" . fmat
+
+-- | Singular values of a complex matrix, using LAPACK's /zgesvd/ with jobu == jobvt == \'N\'.
+svC :: Matrix (Complex Double) -> Vector Double
+svC = svAux zgesvd "svC" . fmat
+
+-- | Singular values of a real matrix, using LAPACK's /dgesdd/ with jobz == \'N\'.
+svRd :: Matrix Double -> Vector Double
+svRd = svAux dgesdd "svRd" . fmat
+
+-- | Singular values of a complex matrix, using LAPACK's /zgesdd/ with jobz == \'N\'.
+svCd :: Matrix (Complex Double) -> Vector Double
+svCd = svAux zgesdd "svCd" . fmat
+
+svAux f st x = unsafePerformIO $ do
+    s <- createVector q
+    app2 g mat x vec s st
+    return s
+  where r = rows x
+        c = cols x
+        q = min r c
+        g ra ca pa nb pb = f ra ca pa 0 0 nullPtr nb pb 0 0 nullPtr
+
+
+-- | Singular values and all right singular vectors of a real matrix, using LAPACK's /dgesvd/ with jobu == \'N\' and jobvt == \'A\'.
+rightSVR :: Matrix Double -> (Vector Double, Matrix Double)
+rightSVR = rightSVAux dgesvd "rightSVR" . fmat
+
+-- | Singular values and all right singular vectors of a complex matrix, using LAPACK's /zgesvd/ with jobu == \'N\' and jobvt == \'A\'.
+rightSVC :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
+rightSVC = rightSVAux zgesvd "rightSVC" . fmat
+
+rightSVAux f st x = unsafePerformIO $ do
+    s <- createVector q
+    v <- createMatrix ColumnMajor c c
+    app3 g mat x vec s mat v st
+    return (s,trans v)
+  where r = rows x
+        c = cols x
+        q = min r c
+        g ra ca pa = f ra ca pa 0 0 nullPtr
+
+
+-- | Singular values and all left singular vectors of a real matrix, using LAPACK's /dgesvd/  with jobu == \'A\' and jobvt == \'N\'.
+leftSVR :: Matrix Double -> (Matrix Double, Vector Double)
+leftSVR = leftSVAux dgesvd "leftSVR" . fmat
+
+-- | Singular values and all left singular vectors of a complex matrix, using LAPACK's /zgesvd/ with jobu == \'A\' and jobvt == \'N\'.
+leftSVC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double)
+leftSVC = leftSVAux zgesvd "leftSVC" . fmat
+
+leftSVAux f st x = unsafePerformIO $ do
+    u <- createMatrix ColumnMajor r r
+    s <- createVector q
+    app3 g mat x mat u vec s st
+    return (u,s)
+  where r = rows x
+        c = cols x
+        q = min r c
+        g ra ca pa ru cu pu nb pb = f ra ca pa ru cu pu nb pb 0 0 nullPtr
+
 -----------------------------------------------------------------------------
+
+foreign import ccall "LAPACK/lapack-aux.h eig_l_R" dgeev :: TMMCVM
+foreign import ccall "LAPACK/lapack-aux.h eig_l_C" zgeev :: TCMCMCVCM
+foreign import ccall "LAPACK/lapack-aux.h eig_l_S" dsyev :: TMVM
+foreign import ccall "LAPACK/lapack-aux.h eig_l_H" zheev :: TCMVCM
+
 eigAux f st m
     | r == 1 = (fromList [flatten m `at` 0], singleton 1)
     | otherwise = unsafePerformIO $ do
@@ -104,28 +212,13 @@ eigAux f st m
   where r = rows m
 
 
-foreign import ccall "LAPACK/lapack-aux.h eig_l_C" zgeev :: TCMCMCVCM
-foreign import ccall "LAPACK/lapack-aux.h eig_l_R" dgeev :: TMMCVM
-foreign import ccall "LAPACK/lapack-aux.h eig_l_S" dsyev :: TMVM
-foreign import ccall "LAPACK/lapack-aux.h eig_l_H" zheev :: TCMVCM
-
--- | Wrapper for LAPACK's /zgeev/, which computes the eigenvalues and right eigenvectors of a general complex matrix:
---
--- if @(l,v)=eigC m@ then @m \<\> v = v \<\> diag l@.
---
--- The eigenvectors are the columns of v.
--- The eigenvalues are not sorted.
+-- | Eigenvalues and right eigenvectors of a general complex matrix, using LAPACK's /zgeev/.
+-- The eigenvectors are the columns of v. The eigenvalues are not sorted.
 eigC :: Matrix (Complex Double) -> (Vector (Complex Double), Matrix (Complex Double))
 eigC = eigAux zgeev "eigC" . fmat
 
------------------------------------------------------------------------------
-
--- | Wrapper for LAPACK's /dgeev/, which computes the eigenvalues and right eigenvectors of a general real matrix:
---
--- if @(l,v)=eigR m@ then @m \<\> v = v \<\> diag l@.
---
--- The eigenvectors are the columns of v.
--- The eigenvalues are not sorted.
+-- | Eigenvalues and right eigenvectors of a general real matrix, using LAPACK's /dgeev/.
+-- The eigenvectors are the columns of v. The eigenvalues are not sorted.
 eigR :: Matrix Double -> (Vector (Complex Double), Matrix (Complex Double))
 eigR m = (s', v'')
     where (s,v) = eigRaux (fmat m)
@@ -155,17 +248,16 @@ fixeig _ _ = error "fixeig with impossible inputs"
 
 -----------------------------------------------------------------------------
 
--- | Wrapper for LAPACK's /dsyev/, which computes the eigenvalues and right eigenvectors of a symmetric real matrix:
---
--- if @(l,v)=eigSl m@ then @m \<\> v = v \<\> diag l@.
---
+-- | Eigenvalues and right eigenvectors of a symmetric real matrix, using LAPACK's /dsyev/.
 -- The eigenvectors are the columns of v.
--- The eigenvalues are sorted in descending order (use eigS' for ascending order).
+-- The eigenvalues are sorted in descending order (use 'eigS'' for ascending order).
 eigS :: Matrix Double -> (Vector Double, Matrix Double)
 eigS m = (s', fliprl v)
     where (s,v) = eigS' (fmat m)
           s' = fromList . reverse . toList $  s
 
+-- | 'eigS' in ascending order
+eigS' :: Matrix Double -> (Vector Double, Matrix Double)
 eigS' m
     | r == 1 = (fromList [flatten m `at` 0], singleton 1)
     | otherwise = unsafePerformIO $ do
@@ -177,17 +269,16 @@ eigS' m
 
 -----------------------------------------------------------------------------
 
--- | Wrapper for LAPACK's /zheev/, which computes the eigenvalues and right eigenvectors of a hermitian complex matrix:
---
--- if @(l,v)=eigH m@ then @m \<\> s v = v \<\> diag l@.
---
+-- | Eigenvalues and right eigenvectors of a hermitian complex matrix, using LAPACK's /zheev/.
 -- The eigenvectors are the columns of v.
--- The eigenvalues are sorted in descending order (use eigH' for ascending order).
+-- The eigenvalues are sorted in descending order (use 'eigH'' for ascending order).
 eigH :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
 eigH m = (s', fliprl v)
     where (s,v) = eigH' (fmat m)
           s' = fromList . reverse . toList $  s
 
+-- | 'eigH' in ascending order
+eigH' :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
 eigH' m
     | r == 1 = (fromList [realPart (flatten m `at` 0)], singleton 1)
     | otherwise = unsafePerformIO $ do
@@ -212,11 +303,11 @@ linearSolveSQAux f st a b
         r  = rows b
         c  = cols b
 
--- | Wrapper for LAPACK's /dgesv/, which solves a general real linear system (for several right-hand sides) internally using the lu decomposition.
+-- | Solve a real linear system (for square coefficient matrix and several right-hand sides) using the LU decomposition, based on LAPACK's /dgesv/. For underconstrained or overconstrained systems use 'linearSolveLSR' or 'linearSolveSVDR'. See also 'lusR'.
 linearSolveR :: Matrix Double -> Matrix Double -> Matrix Double
 linearSolveR a b = linearSolveSQAux dgesv "linearSolveR" (fmat a) (fmat b)
 
--- | Wrapper for LAPACK's /zgesv/, which solves a general complex linear system (for several right-hand sides) internally using the lu decomposition.
+-- | Solve a complex linear system (for square coefficient matrix and several right-hand sides) using the LU decomposition, based on LAPACK's /zgesv/. For underconstrained or overconstrained systems use 'linearSolveLSC' or 'linearSolveSVDC'. See also 'lusC'.
 linearSolveC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
 linearSolveC a b = linearSolveSQAux zgesv "linearSolveC" (fmat a) (fmat b)
 
@@ -234,17 +325,17 @@ linearSolveAux f st a b = unsafePerformIO $ do
         n = cols a
         nrhs = cols b
 
--- | Wrapper for LAPACK's /dgels/, which obtains the least squared error solution of an overconstrained real linear system or the minimum norm solution of an underdetermined system, for several right-hand sides. For rank deficient systems use 'linearSolveSVDR'.
+-- | Least squared error solution of an overconstrained real linear system, or the minimum norm solution of an underconstrained system, using LAPACK's /dgels/. For rank-deficient systems use 'linearSolveSVDR'.
 linearSolveLSR :: Matrix Double -> Matrix Double -> Matrix Double
 linearSolveLSR a b = subMatrix (0,0) (cols a, cols b) $
                      linearSolveAux dgels "linearSolverLSR" (fmat a) (fmat b)
 
--- | Wrapper for LAPACK's /zgels/, which obtains the least squared error solution of an overconstrained complex linear system or the minimum norm solution of an underdetermined system, for several right-hand sides. For rank deficient systems use 'linearSolveSVDC'.
+-- | Least squared error solution of an overconstrained complex linear system, or the minimum norm solution of an underconstrained system, using LAPACK's /zgels/. For rank-deficient systems use 'linearSolveSVDC'.
 linearSolveLSC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
 linearSolveLSC a b = subMatrix (0,0) (cols a, cols b) $
                      linearSolveAux zgels "linearSolveLSC" (fmat a) (fmat b)
 
--- | Wrapper for LAPACK's /dgelss/, which obtains the minimum norm solution to a real linear least squares problem Ax=B using the svd, for several right-hand sides. Admits rank deficient systems but it is slower than 'linearSolveLSR'. The effective rank of A is determined by treating as zero those singular valures which are less than rcond times the largest singular value. If rcond == Nothing machine precision is used.
+-- | Minimum norm solution of a general real linear least squares problem Ax=B using the SVD, based on LAPACK's /dgelss/. Admits rank-deficient systems but it is slower than 'linearSolveLSR'. The effective rank of A is determined by treating as zero those singular valures which are less than rcond times the largest singular value. If rcond == Nothing machine precision is used.
 linearSolveSVDR :: Maybe Double   -- ^ rcond
                 -> Matrix Double  -- ^ coefficient matrix
                 -> Matrix Double  -- ^ right hand sides (as columns)
@@ -253,7 +344,7 @@ linearSolveSVDR (Just rcond) a b = subMatrix (0,0) (cols a, cols b) $
                                    linearSolveAux (dgelss rcond) "linearSolveSVDR" (fmat a) (fmat b)
 linearSolveSVDR Nothing a b = linearSolveSVDR (Just (-1)) (fmat a) (fmat b)
 
--- | Wrapper for LAPACK's /zgelss/, which obtains the minimum norm solution to a complex linear least squares problem Ax=B using the svd, for several right-hand sides. Admits rank deficient systems but it is slower than 'linearSolveLSC'. The effective rank of A is determined by treating as zero those singular valures which are less than rcond times the largest singular value. If rcond == Nothing machine precision is used.
+-- | Minimum norm solution of a general complex linear least squares problem Ax=B using the SVD, based on LAPACK's /zgelss/. Admits rank-deficient systems but it is slower than 'linearSolveLSC'. The effective rank of A is determined by treating as zero those singular valures which are less than rcond times the largest singular value. If rcond == Nothing machine precision is used.
 linearSolveSVDC :: Maybe Double            -- ^ rcond
                 -> Matrix (Complex Double) -- ^ coefficient matrix
                 -> Matrix (Complex Double) -- ^ right hand sides (as columns)
@@ -266,13 +357,11 @@ linearSolveSVDC Nothing a b = linearSolveSVDC (Just (-1)) (fmat a) (fmat b)
 foreign import ccall "LAPACK/lapack-aux.h chol_l_H" zpotrf :: TCMCM
 foreign import ccall "LAPACK/lapack-aux.h chol_l_S" dpotrf :: TMM
 
--- | Wrapper for LAPACK's /zpotrf/, which computes the Cholesky factorization of a
--- complex Hermitian positive definite matrix.
+-- | Cholesky factorization of a complex Hermitian positive definite matrix, using LAPACK's /zpotrf/.
 cholH :: Matrix (Complex Double) -> Matrix (Complex Double)
 cholH = cholAux zpotrf "cholH" . fmat
 
--- | Wrapper for LAPACK's /dpotrf/, which computes the Cholesky factorization of a
--- real symmetric positive definite matrix.
+-- | Cholesky factorization of a real symmetric positive definite matrix, using LAPACK's /dpotrf/.
 cholS :: Matrix Double -> Matrix Double
 cholS = cholAux dpotrf "cholS" . fmat
 
@@ -286,11 +375,11 @@ cholAux f st a = unsafePerformIO $ do
 foreign import ccall "LAPACK/lapack-aux.h qr_l_R" dgeqr2 :: TMVM
 foreign import ccall "LAPACK/lapack-aux.h qr_l_C" zgeqr2 :: TCMCVCM
 
--- | Wrapper for LAPACK's /dgeqr2/, which computes a QR factorization of a real matrix.
+-- | QR factorization of a real matrix, using LAPACK's /dgeqr2/.
 qrR :: Matrix Double -> (Matrix Double, Vector Double)
 qrR = qrAux dgeqr2 "qrR" . fmat
 
--- | Wrapper for LAPACK's /zgeqr2/, which computes a QR factorization of a complex matrix.
+-- | QR factorization of a complex matrix, using LAPACK's /zgeqr2/.
 qrC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector (Complex Double))
 qrC = qrAux zgeqr2 "qrC" . fmat
 
@@ -307,11 +396,11 @@ qrAux f st a = unsafePerformIO $ do
 foreign import ccall "LAPACK/lapack-aux.h hess_l_R" dgehrd :: TMVM
 foreign import ccall "LAPACK/lapack-aux.h hess_l_C" zgehrd :: TCMCVCM
 
--- | Wrapper for LAPACK's /dgehrd/, which computes a Hessenberg factorization of a square real matrix.
+-- | Hessenberg factorization of a square real matrix, using LAPACK's /dgehrd/.
 hessR :: Matrix Double -> (Matrix Double, Vector Double)
 hessR = hessAux dgehrd "hessR" . fmat
 
--- | Wrapper for LAPACK's /zgehrd/, which computes a Hessenberg factorization of a square complex matrix.
+-- | Hessenberg factorization of a square complex matrix, using LAPACK's /zgehrd/.
 hessC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector (Complex Double))
 hessC = hessAux zgehrd "hessC" . fmat
 
@@ -328,11 +417,11 @@ hessAux f st a = unsafePerformIO $ do
 foreign import ccall "LAPACK/lapack-aux.h schur_l_R" dgees :: TMMM
 foreign import ccall "LAPACK/lapack-aux.h schur_l_C" zgees :: TCMCMCM
 
--- | Wrapper for LAPACK's /dgees/, which computes a Schur factorization of a square real matrix.
+-- | Schur factorization of a square real matrix, using LAPACK's /dgees/.
 schurR :: Matrix Double -> (Matrix Double, Matrix Double)
 schurR = schurAux dgees "schurR" . fmat
 
--- | Wrapper for LAPACK's /zgees/, which computes a Schur factorization of a square complex matrix.
+-- | Schur factorization of a square complex matrix, using LAPACK's /zgees/.
 schurC :: Matrix (Complex Double) -> (Matrix (Complex Double), Matrix (Complex Double))
 schurC = schurAux zgees "schurC" . fmat
 
@@ -347,11 +436,11 @@ schurAux f st a = unsafePerformIO $ do
 foreign import ccall "LAPACK/lapack-aux.h lu_l_R" dgetrf :: TMVM
 foreign import ccall "LAPACK/lapack-aux.h lu_l_C" zgetrf :: TCMVCM
 
--- | Wrapper for LAPACK's /dgetrf/, which computes a LU factorization of a general real matrix.
+-- | LU factorization of a general real matrix, using LAPACK's /dgetrf/.
 luR :: Matrix Double -> (Matrix Double, [Int])
 luR = luAux dgetrf "luR" . fmat
 
--- | Wrapper for LAPACK's /zgees/, which computes a Schur factorization of a square complex matrix.
+-- | LU factorization of a general complex matrix, using LAPACK's /zgetrf/.
 luC :: Matrix (Complex Double) -> (Matrix (Complex Double), [Int])
 luC = luAux zgetrf "luC" . fmat
 
@@ -370,11 +459,11 @@ type TQ a = CInt -> CInt -> PC -> a
 foreign import ccall "LAPACK/lapack-aux.h luS_l_R" dgetrs :: TMVMM
 foreign import ccall "LAPACK/lapack-aux.h luS_l_C" zgetrs :: TQ (TW (TQ (TQ (IO CInt))))
 
--- | Wrapper for LAPACK's /dgetrs/, which solves a general real linear system (for several right-hand sides) from a precomputed LU decomposition.
+-- | Solve a real linear system from a precomputed LU decomposition ('luR'), using LAPACK's /dgetrs/.
 lusR :: Matrix Double -> [Int] -> Matrix Double -> Matrix Double
 lusR a piv b = lusAux dgetrs "lusR" (fmat a) piv (fmat b)
 
--- | Wrapper for LAPACK's /zgetrs/, which solves a general real linear system (for several right-hand sides) from a precomputed LU decomposition.
+-- | Solve a real linear system from a precomputed LU decomposition ('luC'), using LAPACK's /zgetrs/.
 lusC :: Matrix (Complex Double) -> [Int] -> Matrix (Complex Double) -> Matrix (Complex Double)
 lusC a piv b = lusAux zgetrs "lusC" (fmat a) piv (fmat b)