    {-# LANGUAGE TypeOperators #-}
    {-# LANGUAGE ViewPatterns #-}

    {-# OPTIONS_GHC -fno-warn-missing-signatures #-}

    -----------------------------------------------------------------------------
    -- |
    -- Module      :  Numeric.LinearAlgebra.LAPACK
    -- Copyright   :  (c) Alberto Ruiz 2006-14
    -- License     :  BSD3
    -- Maintainer  :  Alberto Ruiz
    -- Stability   :  provisional
    --
    -- Functional interface to selected LAPACK functions (<http://www.netlib.org/lapack>).
    --
    -----------------------------------------------------------------------------

    module Internal.LAPACK where

    import Data.Bifunctor (first)

    import Internal.Devel
    import Internal.Vector
    import Internal.Matrix hiding ((#), (#!))
    import Internal.Conversion
    import Internal.Element
    import Foreign.Ptr(nullPtr)
    import Foreign.C.Types
    import Control.Monad(when)
    import System.IO.Unsafe(unsafePerformIO)

    -----------------------------------------------------------------------------------

    infixr 1 #
    c
    a # :: c -> (b -> IO r) -> Trans c b -> IO r
    # b -> IO r
    b = c -> (b -> IO r) -> Trans c b -> IO r
    forall c b r. TransArray c => c -> (b -> IO r) -> Trans c b -> IO r
    forall b r. c -> (b -> IO r) -> Trans c b -> IO r
    apply c
    a b -> IO r
    b
    {-# INLINE (#) #-}

    c
    a #! :: c -> c -> Trans c (Trans c (IO r)) -> IO r
    #! c
    b = c
    a c -> (Trans c (IO r) -> IO r) -> Trans c (Trans c (IO r)) -> IO r
    forall {c} {b} {r}.
    TransArray c =>
    c -> (b -> IO r) -> Trans c b -> IO r
    # c
    b c -> (IO r -> IO r) -> Trans c (IO r) -> IO r
    forall {c} {b} {r}.
    TransArray c =>
    c -> (b -> IO r) -> Trans c b -> IO r
    # IO r -> IO r
    forall a. a -> a
    id
    {-# INLINE (#!) #-}

    -----------------------------------------------------------------------------------

    type TMMM t = t ::> t ::> t ::> Ok

    type F = Float
    type Q = Complex Float

    foreign import ccall unsafe "multiplyR" dgemmc :: CInt -> CInt -> TMMM R
    foreign import ccall unsafe "multiplyC" zgemmc :: CInt -> CInt -> TMMM C
    foreign import ccall unsafe "multiplyF" sgemmc :: CInt -> CInt -> TMMM F
    foreign import ccall unsafe "multiplyQ" cgemmc :: CInt -> CInt -> TMMM Q
    foreign import ccall unsafe "multiplyI" c_multiplyI :: I -> TMMM I
    foreign import ccall unsafe "multiplyL" c_multiplyL :: Z -> TMMM Z

    isT :: Matrix t -> a
    isT (Matrix t -> Bool
    forall t. Matrix t -> Bool
    rowOrder -> Bool
    False) = a
    0
    isT Matrix t
    _                   = a
    1

    tt :: Matrix t -> Matrix t
    tt x :: Matrix t
    x@(Matrix t -> Bool
    forall t. Matrix t -> Bool
    rowOrder -> Bool
    False) = Matrix t
    x
    tt Matrix t
    x                     = Matrix t -> Matrix t
    forall t. Matrix t -> Matrix t
    trans Matrix t
    x

    multiplyAux :: (t
     -> t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix t -> Matrix a
    multiplyAux t
    -> t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    a Matrix t
    b = IO (Matrix a) -> Matrix a
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a) -> Matrix a) -> IO (Matrix a) -> Matrix a
    forall a b. (a -> b) -> a -> b
    $ do
        Bool -> IO () -> IO ()
    forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
    when (Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    /= Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b) (IO () -> IO ()) -> IO () -> IO ()
    forall a b. (a -> b) -> a -> b
    $ String -> IO ()
    forall a. HasCallStack => String -> a
    error (String -> IO ()) -> String -> IO ()
    forall a b. (a -> b) -> a -> b
    $ String
    "inconsistent dimensions in matrix product "String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++
                                           (Int, Int) -> String
    forall a. Show a => a -> String
    show (Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a,Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a) String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " x " String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ (Int, Int) -> String
    forall a. Show a => a -> String
    show (Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b, Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    b)
        s <- MatrixOrder -> Int -> Int -> IO (Matrix a)
    forall a. Storable a => MatrixOrder -> Int -> Int -> IO (Matrix a)
    createMatrix MatrixOrder
    ColumnMajor (Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a) (Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    b)
        ((tt a) # (tt b) #! s) (f (isT a) (isT b)) #| st
        return s

    -- | Matrix product based on BLAS's /dgemm/.
    multiplyR :: Matrix Double -> Matrix Double -> Matrix Double
    multiplyR :: Matrix Double -> Matrix Double -> Matrix Double
    multiplyR Matrix Double
    a Matrix Double
    b = {-# SCC "multiplyR" #-} (CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> Matrix Double -> Matrix Double
    forall {a} {t} {t} {t} {t}.
    (Storable a, Storable t, Storable t, Num t, Num t) =>
    (t
     -> t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix t -> Matrix a
    multiplyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgemmc String
    "dgemmc" Matrix Double
    a Matrix Double
    b

    -- | Matrix product based on BLAS's /zgemm/.
    multiplyC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
    multiplyC :: Matrix (Complex Double)
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    multiplyC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {a} {t} {t} {t} {t}.
    (Storable a, Storable t, Storable t, Num t, Num t) =>
    (t
     -> t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix t -> Matrix a
    multiplyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgemmc String
    "zgemmc" Matrix (Complex Double)
    a Matrix (Complex Double)
    b

    -- | Matrix product based on BLAS's /sgemm/.
    multiplyF :: Matrix Float -> Matrix Float -> Matrix Float
    multiplyF :: Matrix Float -> Matrix Float -> Matrix Float
    multiplyF Matrix Float
    a Matrix Float
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Float
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Float
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Float
     -> IO CInt)
    -> String -> Matrix Float -> Matrix Float -> Matrix Float
    forall {a} {t} {t} {t} {t}.
    (Storable a, Storable t, Storable t, Num t, Num t) =>
    (t
     -> t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix t -> Matrix a
    multiplyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Float
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Float
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Float
    -> IO CInt
    sgemmc String
    "sgemmc" Matrix Float
    a Matrix Float
    b

    -- | Matrix product based on BLAS's /cgemm/.
    multiplyQ :: Matrix (Complex Float) -> Matrix (Complex Float) -> Matrix (Complex Float)
    multiplyQ :: Matrix (Complex Float)
    -> Matrix (Complex Float) -> Matrix (Complex Float)
    multiplyQ Matrix (Complex Float)
    a Matrix (Complex Float)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Float)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Float)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Float)
     -> IO CInt)
    -> String
    -> Matrix (Complex Float)
    -> Matrix (Complex Float)
    -> Matrix (Complex Float)
    forall {a} {t} {t} {t} {t}.
    (Storable a, Storable t, Storable t, Num t, Num t) =>
    (t
     -> t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix t -> Matrix a
    multiplyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Float)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Float)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Float)
    -> IO CInt
    cgemmc String
    "cgemmc" Matrix (Complex Float)
    a Matrix (Complex Float)
    b

    multiplyI :: I -> Matrix CInt -> Matrix CInt -> Matrix CInt
    multiplyI :: CInt -> Matrix CInt -> Matrix CInt -> Matrix CInt
    multiplyI CInt
    m Matrix CInt
    a Matrix CInt
    b = IO (Matrix CInt) -> Matrix CInt
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix CInt) -> Matrix CInt)
    -> IO (Matrix CInt) -> Matrix CInt
    forall a b. (a -> b) -> a -> b
    $ do
        Bool -> IO () -> IO ()
    forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
    when (Matrix CInt -> Int
    forall t. Matrix t -> Int
    cols Matrix CInt
    a Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    /= Matrix CInt -> Int
    forall t. Matrix t -> Int
    rows Matrix CInt
    b) (IO () -> IO ()) -> IO () -> IO ()
    forall a b. (a -> b) -> a -> b
    $ String -> IO ()
    forall a. HasCallStack => String -> a
    error (String -> IO ()) -> String -> IO ()
    forall a b. (a -> b) -> a -> b
    $
            String
    "inconsistent dimensions in matrix product "String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ Matrix CInt -> String
    forall t. Matrix t -> String
    shSize Matrix CInt
    a String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " x " String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ Matrix CInt -> String
    forall t. Matrix t -> String
    shSize Matrix CInt
    b
        s <- MatrixOrder -> Int -> Int -> IO (Matrix CInt)
    forall a. Storable a => MatrixOrder -> Int -> Int -> IO (Matrix a)
    createMatrix MatrixOrder
    ColumnMajor (Matrix CInt -> Int
    forall t. Matrix t -> Int
    rows Matrix CInt
    a) (Matrix CInt -> Int
    forall t. Matrix t -> Int
    cols Matrix CInt
    b)
        (a # b #! s) (c_multiplyI m) #|"c_multiplyI"
        return s

    multiplyL :: Z -> Matrix Z -> Matrix Z -> Matrix Z
    multiplyL :: Z -> Matrix Z -> Matrix Z -> Matrix Z
    multiplyL Z
    m Matrix Z
    a Matrix Z
    b = IO (Matrix Z) -> Matrix Z
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix Z) -> Matrix Z) -> IO (Matrix Z) -> Matrix Z
    forall a b. (a -> b) -> a -> b
    $ do
        Bool -> IO () -> IO ()
    forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
    when (Matrix Z -> Int
    forall t. Matrix t -> Int
    cols Matrix Z
    a Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    /= Matrix Z -> Int
    forall t. Matrix t -> Int
    rows Matrix Z
    b) (IO () -> IO ()) -> IO () -> IO ()
    forall a b. (a -> b) -> a -> b
    $ String -> IO ()
    forall a. HasCallStack => String -> a
    error (String -> IO ()) -> String -> IO ()
    forall a b. (a -> b) -> a -> b
    $
            String
    "inconsistent dimensions in matrix product "String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ Matrix Z -> String
    forall t. Matrix t -> String
    shSize Matrix Z
    a String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " x " String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ Matrix Z -> String
    forall t. Matrix t -> String
    shSize Matrix Z
    b
        s <- MatrixOrder -> Int -> Int -> IO (Matrix Z)
    forall a. Storable a => MatrixOrder -> Int -> Int -> IO (Matrix a)
    createMatrix MatrixOrder
    ColumnMajor (Matrix Z -> Int
    forall t. Matrix t -> Int
    rows Matrix Z
    a) (Matrix Z -> Int
    forall t. Matrix t -> Int
    cols Matrix Z
    b)
        (a # b #! s) (c_multiplyL m) #|"c_multiplyL"
        return s

    -----------------------------------------------------------------------------

    type TSVD t = t ::> t ::> R :> t ::> Ok

    foreign import ccall unsafe "svd_l_R" dgesvd :: TSVD R
    foreign import ccall unsafe "svd_l_C" zgesvd :: TSVD C
    foreign import ccall unsafe "svd_l_Rdd" dgesdd :: TSVD R
    foreign import ccall unsafe "svd_l_Cdd" zgesdd :: TSVD C

    -- | Full SVD of a real matrix using LAPACK's /dgesvd/.
    svdR :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    svdR :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    svdR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> (Matrix Double, Vector Double, Matrix Double)
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    svdAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesvd String
    "svdR"

    -- | Full SVD of a real matrix using LAPACK's /dgesdd/.
    svdRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    svdRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    svdRd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> (Matrix Double, Vector Double, Matrix Double)
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    svdAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesdd String
    "svdRdd"

    -- | Full SVD of a complex matrix using LAPACK's /zgesvd/.
    svdC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
    svdC :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    svdC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    svdAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesvd String
    "svdC"

    -- | Full SVD of a complex matrix using LAPACK's /zgesdd/.
    svdCd :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double, Matrix (Complex Double))
    svdCd :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    svdCd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    svdAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesdd String
    "svdCdd"

    svdAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    svdAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    x = IO (Matrix a, Vector a, Matrix a) -> (Matrix a, Vector a, Matrix a)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a, Vector a, Matrix a)
     -> (Matrix a, Vector a, Matrix a))
    -> IO (Matrix a, Vector a, Matrix a)
    -> (Matrix a, Vector a, Matrix a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    x
        u <- createMatrix ColumnMajor r r
        s <- createVector (min r c)
        v <- createMatrix ColumnMajor c c
        (a # u # s #! v) f #| st
        return (u,s,v)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    x
        c :: Int
    c = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    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 :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    thinSVDR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> (Matrix Double, Vector Double, Matrix Double)
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    thinSVDAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesvd String
    "thinSVDR"

    -- | 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 :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    thinSVDC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    thinSVDAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesvd String
    "thinSVDC"

    -- | Thin SVD of a real matrix, using LAPACK's /dgesdd/ with jobz == \'S\'.
    thinSVDRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    thinSVDRd :: Matrix Double -> (Matrix Double, Vector Double, Matrix Double)
    thinSVDRd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> (Matrix Double, Vector Double, Matrix Double)
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    thinSVDAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesdd String
    "thinSVDRdd"

    -- | 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 :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    thinSVDCd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double,
        Matrix (Complex Double))
    forall {t} {a} {a} {a}.
    (Element t, Storable a, Storable a, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    thinSVDAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesdd String
    "thinSVDCdd"

    thinSVDAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a, Matrix a)
    thinSVDAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    x = IO (Matrix a, Vector a, Matrix a) -> (Matrix a, Vector a, Matrix a)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a, Vector a, Matrix a)
     -> (Matrix a, Vector a, Matrix a))
    -> IO (Matrix a, Vector a, Matrix a)
    -> (Matrix a, Vector a, Matrix a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    x
        u <- createMatrix ColumnMajor r q
        s <- createVector q
        v <- createMatrix ColumnMajor q c
        (a # u # s #! v) f #| st
        return (u,s,v)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    x
        c :: Int
    c = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    x
        q :: Int
    q = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    r Int
    c

    -- | Singular values of a real matrix, using LAPACK's /dgesvd/ with jobu == jobvt == \'N\'.
    svR :: Matrix Double -> Vector Double
    svR :: Matrix Double -> Vector Double
    svR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> Vector Double
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    svAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesvd String
    "svR"

    -- | Singular values of a complex matrix, using LAPACK's /zgesvd/ with jobu == jobvt == \'N\'.
    svC :: Matrix (Complex Double) -> Vector Double
    svC :: Matrix (Complex Double) -> Vector Double
    svC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String -> Matrix (Complex Double) -> Vector Double
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    svAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesvd String
    "svC"

    -- | Singular values of a real matrix, using LAPACK's /dgesdd/ with jobz == \'N\'.
    svRd :: Matrix Double -> Vector Double
    svRd :: Matrix Double -> Vector Double
    svRd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> Vector Double
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    svAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesdd String
    "svRd"

    -- | Singular values of a complex matrix, using LAPACK's /zgesdd/ with jobz == \'N\'.
    svCd :: Matrix (Complex Double) -> Vector Double
    svCd :: Matrix (Complex Double) -> Vector Double
    svCd = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String -> Matrix (Complex Double) -> Vector Double
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    svAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesdd String
    "svCd"

    svAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    svAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    x = IO (Vector a) -> Vector a
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a) -> Vector a) -> IO (Vector a) -> Vector a
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    x
        s <- createVector q
        (a #! s) g #| st
        return s
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    x
        c :: Int
    c = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    x
        q :: Int
    q = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    r Int
    c
        g :: CInt -> CInt -> CInt -> CInt -> Ptr t -> CInt -> Ptr a -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa CInt
    nb Ptr a
    pb = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr CInt
    nb Ptr a
    pb t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    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 :: Matrix Double -> (Vector Double, Matrix Double)
    rightSVR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Vector Double, Matrix Double)
    forall {t} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Storable a, Storable a, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix a)
    rightSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesvd String
    "rightSVR"

    -- | 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 :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    rightSVC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Vector Double, Matrix (Complex Double))
    forall {t} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Storable a, Storable a, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix a)
    rightSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesvd String
    "rightSVC"

    rightSVAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix a)
    rightSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    x = IO (Vector a, Matrix a) -> (Vector a, Matrix a)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a, Matrix a) -> (Vector a, Matrix a))
    -> IO (Vector a, Matrix a) -> (Vector a, Matrix a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    x
        s <- createVector q
        v <- createMatrix ColumnMajor c c
        (a # s #! v) g #| st
        return (s,v)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    x
        c :: Int
    c = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    x
        q :: Int
    q = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    r Int
    c
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    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 :: Matrix Double -> (Matrix Double, Vector Double)
    leftSVR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, Vector Double)
    forall {t} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Storable a, Storable a, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a)
    leftSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesvd String
    "leftSVR"

    -- | 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 :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector Double)
    leftSVC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector Double)
    forall {t} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Storable a, Storable a, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a)
    leftSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesvd String
    "leftSVC"

    leftSVAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Vector a)
    leftSVAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    x = IO (Matrix a, Vector a) -> (Matrix a, Vector a)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a, Vector a) -> (Matrix a, Vector a))
    -> IO (Matrix a, Vector a) -> (Matrix a, Vector a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    x
        u <- createMatrix ColumnMajor r r
        s <- createVector q
        (a # u #! s) g #| st
        return (u,s)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    x
        c :: Int
    c = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    x
        q :: Int
    q = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    r Int
    c
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa CInt
    ru CInt
    cu CInt
    xru CInt
    xcu Ptr a
    pu CInt
    nb Ptr a
    pb = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa CInt
    ru CInt
    cu CInt
    xru CInt
    xcu Ptr a
    pu CInt
    nb Ptr a
    pb t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr

    -----------------------------------------------------------------------------

    foreign import ccall unsafe "eig_l_R" dgeev :: R ::> R ::> C :> R ::> Ok
    foreign import ccall unsafe "eig_l_G" dggev :: R ::> R ::> C :> R :> R ::> R ::> Ok
    foreign import ccall unsafe "eig_l_C" zgeev :: C ::> C ::> C :> C ::> Ok
    foreign import ccall unsafe "eig_l_GC" zggev :: C ::> C ::> C :> C :> C ::> C ::> Ok
    foreign import ccall unsafe "eig_l_S" dsyev :: CInt -> R :> R ::> Ok
    foreign import ccall unsafe "eig_l_H" zheev :: CInt -> R :> C ::> Ok

    eigAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix a)
    eigAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    m = IO (Vector a, Matrix a) -> (Vector a, Matrix a)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a, Matrix a) -> (Vector a, Matrix a))
    -> IO (Vector a, Matrix a) -> (Vector a, Matrix a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    m
        l <- createVector r
        v <- createMatrix ColumnMajor r r
        (a # l #! v) g #| st
        return (l,v)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    m
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr

    -- | 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 :: Matrix (Complex Double)
    -> (Vector (Complex Double), Matrix (Complex Double))
    eigC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Vector (Complex Double), Matrix (Complex Double))
    forall {t} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Storable a, Storable a, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix a)
    eigAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgeev String
    "eigC"

    eigOnlyAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    eigOnlyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    m = IO (Vector a) -> Vector a
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a) -> Vector a) -> IO (Vector a) -> Vector a
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    m
        l <- createVector r
        (a #! l) g #| st
        return l
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    m
        g :: CInt -> CInt -> CInt -> CInt -> Ptr t -> CInt -> Ptr a -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa CInt
    nl Ptr a
    pl = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr t
    pa t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr CInt
    nl Ptr a
    pl t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr

    -- | Eigenvalues of a general complex matrix, using LAPACK's /zgeev/ with jobz == \'N\'.
    -- The eigenvalues are not sorted.
    eigOnlyC :: Matrix (Complex Double) -> Vector (Complex Double)
    eigOnlyC :: Matrix (Complex Double) -> Vector (Complex Double)
    eigOnlyC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String -> Matrix (Complex Double) -> Vector (Complex Double)
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    eigOnlyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgeev String
    "eigOnlyC"

    -- | 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 :: Matrix Double -> (Vector (Complex Double), Matrix (Complex Double))
    eigR Matrix Double
    m = (Vector (Complex Double)
    s', Matrix (Complex Double)
    v'')
        where (Vector (Complex Double)
    s,Matrix Double
    v) = Matrix Double -> (Vector (Complex Double), Matrix Double)
    eigRaux Matrix Double
    m
              s' :: Vector (Complex Double)
    s' = Vector (Complex Double) -> Vector (Complex Double)
    forall {e}.
    RealElement e =>
    Vector (Complex e) -> Vector (Complex e)
    fixeig1 Vector (Complex Double)
    s
              v' :: [Vector Double]
    v' = Matrix Double -> [Vector Double]
    forall t. Element t => Matrix t -> [Vector t]
    toRows (Matrix Double -> [Vector Double])
    -> Matrix Double -> [Vector Double]
    forall a b. (a -> b) -> a -> b
    $ Matrix Double -> Matrix Double
    forall t. Matrix t -> Matrix t
    trans Matrix Double
    v
              v'' :: Matrix (Complex Double)
    v'' = [Vector (Complex Double)] -> Matrix (Complex Double)
    forall t. Element t => [Vector t] -> Matrix t
    fromColumns ([Vector (Complex Double)] -> Matrix (Complex Double))
    -> [Vector (Complex Double)] -> Matrix (Complex Double)
    forall a b. (a -> b) -> a -> b
    $ [Complex Double] -> [Vector Double] -> [Vector (Complex Double)]
    forall {e} {a}.
    (RealElement e, Num a, Eq a) =>
    [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeig (Vector (Complex Double) -> [Complex Double]
    forall a. Storable a => Vector a -> [a]
    toList Vector (Complex Double)
    s') [Vector Double]
    v'

    eigRaux :: Matrix Double -> (Vector (Complex Double), Matrix Double)
    eigRaux :: Matrix Double -> (Vector (Complex Double), Matrix Double)
    eigRaux Matrix Double
    m = IO (Vector (Complex Double), Matrix Double)
    -> (Vector (Complex Double), Matrix Double)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector (Complex Double), Matrix Double)
     -> (Vector (Complex Double), Matrix Double))
    -> IO (Vector (Complex Double), Matrix Double)
    -> (Vector (Complex Double), Matrix Double)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix Double -> IO (Matrix Double)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix Double
    m
        l <- createVector r
        v <- createMatrix ColumnMajor r r
        (a # l #! v) g #| "eigR"
        return (l,v)
      where
        r :: Int
    r = Matrix Double -> Int
    forall t. Matrix t -> Int
    rows Matrix Double
    m
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    g CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr Double
    pa = Double
    ::> (CInt
         -> CInt
         -> CInt
         -> CInt
         -> Ptr Double
         -> CInt
         -> Ptr (Complex Double)
         -> CInt
         -> CInt
         -> CInt
         -> CInt
         -> Ptr Double
         -> IO CInt)
    dgeev CInt
    ra CInt
    ca CInt
    xra CInt
    xca Ptr Double
    pa CInt
    0 CInt
    0 CInt
    0 CInt
    0 Ptr Double
    forall a. Ptr a
    nullPtr

    fixeig1 :: Vector (Complex e) -> Vector (Complex e)
    fixeig1 Vector (Complex e)
    s = (Vector e, Vector e) -> Vector (Complex e)
    forall e.
    RealElement e =>
    (Vector e, Vector e) -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    (c e, c e) -> c (Complex e)
    toComplex' (Int -> Int -> Vector e -> Vector e
    forall t. Storable t => Int -> Int -> Vector t -> Vector t
    subVector Int
    0 Int
    r (Vector (Complex e) -> Vector e
    forall a.
    (RealFloat a, Storable a) =>
    Vector (Complex a) -> Vector a
    asReal Vector (Complex e)
    s), Int -> Int -> Vector e -> Vector e
    forall t. Storable t => Int -> Int -> Vector t -> Vector t
    subVector Int
    r Int
    r (Vector (Complex e) -> Vector e
    forall a.
    (RealFloat a, Storable a) =>
    Vector (Complex a) -> Vector a
    asReal Vector (Complex e)
    s))
        where r :: Int
    r = Vector (Complex e) -> Int
    forall t. Storable t => Vector t -> Int
    dim Vector (Complex e)
    s

    fixeig :: [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeig  []  [Vector e]
    _ =  []
    fixeig [Complex a
    _] [Vector e
    v] = [Vector e -> Vector (Complex e)
    forall e. RealElement e => Vector e -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    c e -> c (Complex e)
    comp' Vector e
    v]
    fixeig ((a
    r1:+a
    i1):(a
    r2:+a
    i2):[Complex a]
    r) (Vector e
    v1:Vector e
    v2:[Vector e]
    vs)
        | a
    r1 a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == a
    r2 Bool -> Bool -> Bool
    && a
    i1 a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == (-a
    i2) = (Vector e, Vector e) -> Vector (Complex e)
    forall e.
    RealElement e =>
    (Vector e, Vector e) -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    (c e, c e) -> c (Complex e)
    toComplex' (Vector e
    v1,Vector e
    v2) Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : (Vector e, Vector e) -> Vector (Complex e)
    forall e.
    RealElement e =>
    (Vector e, Vector e) -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    (c e, c e) -> c (Complex e)
    toComplex' (Vector e
    v1, (e -> e) -> Vector e -> Vector e
    forall a b.
    (Storable a, Storable b) =>
    (a -> b) -> Vector a -> Vector b
    mapVector e -> e
    forall a. Num a => a -> a
    negate Vector e
    v2) Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeig [Complex a]
    r [Vector e]
    vs
        | Bool
    otherwise = Vector e -> Vector (Complex e)
    forall e. RealElement e => Vector e -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    c e -> c (Complex e)
    comp' Vector e
    v1 Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeig ((a
    r2a -> a -> Complex a
    forall a. a -> a -> Complex a
    :+a
    i2)Complex a -> [Complex a] -> [Complex a]
    forall a. a -> [a] -> [a]
    :[Complex a]
    r) (Vector e
    v2Vector e -> [Vector e] -> [Vector e]
    forall a. a -> [a] -> [a]
    :[Vector e]
    vs)
    fixeig [Complex a]
    _ [Vector e]
    _ = String -> [Vector (Complex e)]
    forall a. HasCallStack => String -> a
    error String
    "fixeig with impossible inputs"

    -- For dggev alpha(i) / beta(i), alpha(i+1) / beta(i+1) form a complex conjugate pair when Im alpha(i) != 0.
    -- However, this does not lead to Re alpha(i) == Re alpha(i+1), since beta(i) and beta(i+1)
    -- can be different. Therefore old 'fixeig' would fail for 'eigG'.
    fixeigG :: [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeigG  []  [Vector e]
    _  = []
    fixeigG [Complex a
    _] [Vector e
    v] = [Vector e -> Vector (Complex e)
    forall e. RealElement e => Vector e -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    c e -> c (Complex e)
    comp' Vector e
    v]
    fixeigG ((a
    _:+a
    ai1) : Complex a
    an : [Complex a]
    as) (Vector e
    v1:Vector e
    v2:[Vector e]
    vs)
        | a -> a
    forall a. Num a => a -> a
    abs a
    ai1 a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    > a
    1e-13 = (Vector e, Vector e) -> Vector (Complex e)
    forall e.
    RealElement e =>
    (Vector e, Vector e) -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    (c e, c e) -> c (Complex e)
    toComplex' (Vector e
    v1, Vector e
    v2) Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : (Vector e, Vector e) -> Vector (Complex e)
    forall e.
    RealElement e =>
    (Vector e, Vector e) -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    (c e, c e) -> c (Complex e)
    toComplex' (Vector e
    v1, (e -> e) -> Vector e -> Vector e
    forall a b.
    (Storable a, Storable b) =>
    (a -> b) -> Vector a -> Vector b
    mapVector e -> e
    forall a. Num a => a -> a
    negate Vector e
    v2) Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeigG [Complex a]
    as [Vector e]
    vs
        | Bool
    otherwise = Vector e -> Vector (Complex e)
    forall e. RealElement e => Vector e -> Vector (Complex e)
    forall (c :: * -> *) e.
    (Complexable c, RealElement e) =>
    c e -> c (Complex e)
    comp' Vector e
    v1 Vector (Complex e) -> [Vector (Complex e)] -> [Vector (Complex e)]
    forall a. a -> [a] -> [a]
    : [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeigG (Complex a
    anComplex a -> [Complex a] -> [Complex a]
    forall a. a -> [a] -> [a]
    :[Complex a]
    as) (Vector e
    v2Vector e -> [Vector e] -> [Vector e]
    forall a. a -> [a] -> [a]
    :[Vector e]
    vs)
    fixeigG [Complex a]
    _ [Vector e]
    _ = String -> [Vector (Complex e)]
    forall a. HasCallStack => String -> a
    error String
    "fixeigG with impossible inputs"

    -- | Eigenvalues of a general real matrix, using LAPACK's /dgeev/ with jobz == \'N\'.
    -- The eigenvalues are not sorted.
    eigOnlyR :: Matrix Double -> Vector (Complex Double)
    eigOnlyR :: Matrix Double -> Vector (Complex Double)
    eigOnlyR = Vector (Complex Double) -> Vector (Complex Double)
    forall {e}.
    RealElement e =>
    Vector (Complex e) -> Vector (Complex e)
    fixeig1 (Vector (Complex Double) -> Vector (Complex Double))
    -> (Matrix Double -> Vector (Complex Double))
    -> Matrix Double
    -> Vector (Complex Double)
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (Double
     ::> (CInt
          -> CInt
          -> CInt
          -> CInt
          -> Ptr Double
          -> CInt
          -> Ptr (Complex Double)
          -> CInt
          -> CInt
          -> CInt
          -> CInt
          -> Ptr Double
          -> IO CInt))
    -> String -> Matrix Double -> Vector (Complex Double)
    forall {t} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Storable a, Num t, Num t, Num t, Num t, Num t, Num t,
     Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Vector a
    eigOnlyAux Double
    ::> (CInt
         -> CInt
         -> CInt
         -> CInt
         -> Ptr Double
         -> CInt
         -> Ptr (Complex Double)
         -> CInt
         -> CInt
         -> CInt
         -> CInt
         -> Ptr Double
         -> IO CInt)
    dgeev String
    "eigOnlyR"

    -- | Generalized eigenvalues and right eigenvectors of a pair of real matrices, using LAPACK's /dggev/.
    -- The eigenvectors are the columns of v. The eigenvalues are represented as alphas / betas and not sorted.
    eigG :: Matrix Double -> Matrix Double -> (Vector (Complex Double), Vector Double, Matrix (Complex Double))
    eigG :: Matrix Double
    -> Matrix Double
    -> (Vector (Complex Double), Vector Double,
        Matrix (Complex Double))
    eigG Matrix Double
    a Matrix Double
    b = (Vector (Complex Double)
    alpha', Vector Double
    beta, Matrix (Complex Double)
    v'')
      where
        (Vector (Complex Double)
    alpha, Vector Double
    beta, Matrix Double
    v) = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> Matrix Double
    -> Matrix Double
    -> String
    -> (Vector (Complex Double), Vector Double, Matrix Double)
    forall {t} {t} {a} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Element t, Storable a, Storable a, Storable a, Num t,
     Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a, Matrix a)
    eigGaux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dggev Matrix Double
    a Matrix Double
    b String
    "eigG"
        alpha' :: Vector (Complex Double)
    alpha' = Vector (Complex Double) -> Vector (Complex Double)
    forall {e}.
    RealElement e =>
    Vector (Complex e) -> Vector (Complex e)
    fixeig1 Vector (Complex Double)
    alpha
        v' :: [Vector Double]
    v' = Matrix Double -> [Vector Double]
    forall t. Element t => Matrix t -> [Vector t]
    toRows (Matrix Double -> [Vector Double])
    -> Matrix Double -> [Vector Double]
    forall a b. (a -> b) -> a -> b
    $ Matrix Double -> Matrix Double
    forall t. Matrix t -> Matrix t
    trans Matrix Double
    v
        v'' :: Matrix (Complex Double)
    v'' = [Vector (Complex Double)] -> Matrix (Complex Double)
    forall t. Element t => [Vector t] -> Matrix t
    fromColumns ([Vector (Complex Double)] -> Matrix (Complex Double))
    -> [Vector (Complex Double)] -> Matrix (Complex Double)
    forall a b. (a -> b) -> a -> b
    $ [Complex Double] -> [Vector Double] -> [Vector (Complex Double)]
    forall {e} {a}.
    (RealElement e, Ord a, Fractional a) =>
    [Complex a] -> [Vector e] -> [Vector (Complex e)]
    fixeigG (Vector (Complex Double) -> [Complex Double]
    forall a. Storable a => Vector a -> [a]
    toList Vector (Complex Double)
    alpha') [Vector Double]
    v'

    eigGaux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a, Matrix a)
    eigGaux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f Matrix t
    ma Matrix t
    mb String
    st = IO (Vector a, Vector a, Matrix a) -> (Vector a, Vector a, Matrix a)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a, Vector a, Matrix a)
     -> (Vector a, Vector a, Matrix a))
    -> IO (Vector a, Vector a, Matrix a)
    -> (Vector a, Vector a, Matrix a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    ma
        b <- copy ColumnMajor mb
        alpha <- createVector r
        beta <- createVector r
        vr <- createMatrix ColumnMajor r r

        (a # b # alpha # beta #! vr) g #| st

        return (alpha, beta, vr)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    ma
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    g CInt
    ar CInt
    ac CInt
    xra CInt
    xca Ptr t
    pa CInt
    br CInt
    bc CInt
    xrb CInt
    xcb Ptr t
    pb CInt
    alphan Ptr a
    palpha CInt
    betan Ptr a
    pbeta = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f CInt
    ar CInt
    ac CInt
    xra CInt
    xca Ptr t
    pa CInt
    br CInt
    bc CInt
    xrb CInt
    xcb Ptr t
    pb CInt
    alphan Ptr a
    palpha CInt
    betan Ptr a
    pbeta t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr

    eigGOnlyAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a)
    eigGOnlyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f Matrix t
    ma Matrix t
    mb String
    st = IO (Vector a, Vector a) -> (Vector a, Vector a)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a, Vector a) -> (Vector a, Vector a))
    -> IO (Vector a, Vector a) -> (Vector a, Vector a)
    forall a b. (a -> b) -> a -> b
    $ do
        a <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    ma
        b <- copy ColumnMajor mb
        alpha <- createVector r
        beta <- createVector r

        (a # b # alpha #! beta) g #| st

        return (alpha, beta)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    ma
        g :: CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> IO CInt
    g CInt
    ar CInt
    ac CInt
    xra CInt
    xca Ptr t
    pa CInt
    br CInt
    bc CInt
    xrb CInt
    xcb Ptr t
    pb CInt
    alphan Ptr a
    palpha CInt
    betan Ptr a
    pbeta = CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr a
    -> CInt
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> t
    -> t
    -> t
    -> t
    -> Ptr a
    -> IO CInt
    f CInt
    ar CInt
    ac CInt
    xra CInt
    xca Ptr t
    pa CInt
    br CInt
    bc CInt
    xrb CInt
    xcb Ptr t
    pb CInt
    alphan Ptr a
    palpha CInt
    betan Ptr a
    pbeta t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr t
    0 t
    0 t
    0 t
    0 Ptr a
    forall a. Ptr a
    nullPtr

    -- | Generalized eigenvalues and right eigenvectors of a pair of complex matrices, using LAPACK's /zggev/.
    -- The eigenvectors are the columns of v. The eigenvalues are represented as alphas / betas and not sorted.
    eigGC :: Matrix (Complex Double) -> Matrix (Complex Double) -> (Vector (Complex Double), Vector (Complex Double), Matrix (Complex Double))
    eigGC :: Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> (Vector (Complex Double), Vector (Complex Double),
        Matrix (Complex Double))
    eigGC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> String
    -> (Vector (Complex Double), Vector (Complex Double),
        Matrix (Complex Double))
    forall {t} {t} {a} {a} {a} {t} {t} {t} {t} {a}.
    (Element t, Element t, Storable a, Storable a, Storable a, Num t,
     Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a, Matrix a)
    eigGaux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zggev Matrix (Complex Double)
    a Matrix (Complex Double)
    b String
    "eigGC"

    eigOnlyG :: Matrix Double -> Matrix Double -> (Vector (Complex Double), Vector Double)
    eigOnlyG :: Matrix Double
    -> Matrix Double -> (Vector (Complex Double), Vector Double)
    eigOnlyG Matrix Double
    a Matrix Double
    b = (Vector (Complex Double) -> Vector (Complex Double))
    -> (Vector (Complex Double), Vector Double)
    -> (Vector (Complex Double), Vector Double)
    forall a b c. (a -> b) -> (a, c) -> (b, c)
    forall (p :: * -> * -> *) a b c.
    Bifunctor p =>
    (a -> b) -> p a c -> p b c
    first Vector (Complex Double) -> Vector (Complex Double)
    forall {e}.
    RealElement e =>
    Vector (Complex e) -> Vector (Complex e)
    fixeig1 ((Vector (Complex Double), Vector Double)
     -> (Vector (Complex Double), Vector Double))
    -> (Vector (Complex Double), Vector Double)
    -> (Vector (Complex Double), Vector Double)
    forall a b. (a -> b) -> a -> b
    $ (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> Matrix Double
    -> Matrix Double
    -> String
    -> (Vector (Complex Double), Vector Double)
    forall {t} {t} {a} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Element t, Storable a, Storable a, Num t, Num t, Num t,
     Num t, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a)
    eigGOnlyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dggev Matrix Double
    a Matrix Double
    b String
    "eigOnlyG"

    eigOnlyGC :: Matrix (Complex Double) -> Matrix (Complex Double) -> (Vector (Complex Double), Vector (Complex Double))
    eigOnlyGC :: Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> (Vector (Complex Double), Vector (Complex Double))
    eigOnlyGC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> String
    -> (Vector (Complex Double), Vector (Complex Double))
    forall {t} {t} {a} {a} {t} {t} {t} {t} {t} {t} {t} {t} {a} {a}.
    (Element t, Element t, Storable a, Storable a, Num t, Num t, Num t,
     Num t, Num t, Num t, Num t, Num t) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr a
     -> CInt
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> t
     -> t
     -> t
     -> t
     -> Ptr a
     -> IO CInt)
    -> Matrix t -> Matrix t -> String -> (Vector a, Vector a)
    eigGOnlyAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zggev Matrix (Complex Double)
    a Matrix (Complex Double)
    b String
    "eigOnlyGC"

    -----------------------------------------------------------------------------

    eigSHAux :: (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix t)
    eigSHAux CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    m = IO (Vector a, Matrix t) -> (Vector a, Matrix t)
    forall a. IO a -> a
    unsafePerformIO (IO (Vector a, Matrix t) -> (Vector a, Matrix t))
    -> IO (Vector a, Matrix t) -> (Vector a, Matrix t)
    forall a b. (a -> b) -> a -> b
    $ do
        l <- Int -> IO (Vector a)
    forall a. Storable a => Int -> IO (Vector a)
    createVector Int
    r
        v <- copy ColumnMajor m
        (l #! v) f #| st
        return (l,v)
      where
        r :: Int
    r = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    m

    -- | 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).
    eigS :: Matrix Double -> (Vector Double, Matrix Double)
    eigS :: Matrix Double -> (Vector Double, Matrix Double)
    eigS Matrix Double
    m = (Vector Double
    s', Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fliprl Matrix Double
    v)
        where (Vector Double
    s,Matrix Double
    v) = Matrix Double -> (Vector Double, Matrix Double)
    eigS' Matrix Double
    m
              s' :: Vector Double
    s' = [Double] -> Vector Double
    forall a. Storable a => [a] -> Vector a
    fromList ([Double] -> Vector Double)
    -> (Vector Double -> [Double]) -> Vector Double -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . [Double] -> [Double]
    forall a. [a] -> [a]
    reverse ([Double] -> [Double])
    -> (Vector Double -> [Double]) -> Vector Double -> [Double]
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Vector Double -> [Double]
    forall a. Storable a => Vector a -> [a]
    toList (Vector Double -> Vector Double) -> Vector Double -> Vector Double
    forall a b. (a -> b) -> a -> b
    $  Vector Double
    s

    -- | 'eigS' in ascending order
    eigS' :: Matrix Double -> (Vector Double, Matrix Double)
    eigS' :: Matrix Double -> (Vector Double, Matrix Double)
    eigS' = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Vector Double, Matrix Double)
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix t)
    eigSHAux (CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dsyev CInt
    1) String
    "eigS'"

    -- | 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).
    eigH :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    eigH :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    eigH Matrix (Complex Double)
    m = (Vector Double
    s', Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fliprl Matrix (Complex Double)
    v)
      where
        (Vector Double
    s,Matrix (Complex Double)
    v) = Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    eigH' Matrix (Complex Double)
    m
        s' :: Vector Double
    s' = [Double] -> Vector Double
    forall a. Storable a => [a] -> Vector a
    fromList ([Double] -> Vector Double)
    -> (Vector Double -> [Double]) -> Vector Double -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . [Double] -> [Double]
    forall a. [a] -> [a]
    reverse ([Double] -> [Double])
    -> (Vector Double -> [Double]) -> Vector Double -> [Double]
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Vector Double -> [Double]
    forall a. Storable a => Vector a -> [a]
    toList (Vector Double -> Vector Double) -> Vector Double -> Vector Double
    forall a b. (a -> b) -> a -> b
    $  Vector Double
    s

    -- | 'eigH' in ascending order
    eigH' :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    eigH' :: Matrix (Complex Double) -> (Vector Double, Matrix (Complex Double))
    eigH' = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Vector Double, Matrix (Complex Double))
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix t)
    eigSHAux (CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zheev CInt
    1) String
    "eigH'"

    -- | Eigenvalues of a symmetric real matrix, using LAPACK's /dsyev/ with jobz == \'N\'.
    -- The eigenvalues are sorted in descending order.
    eigOnlyS :: Matrix Double -> Vector Double
    eigOnlyS :: Matrix Double -> Vector Double
    eigOnlyS = Vector Double -> Vector Double
    vrev (Vector Double -> Vector Double)
    -> (Matrix Double -> Vector Double)
    -> Matrix Double
    -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (Vector Double, Matrix Double) -> Vector Double
    forall a b. (a, b) -> a
    fst((Vector Double, Matrix Double) -> Vector Double)
    -> (Matrix Double -> (Vector Double, Matrix Double))
    -> Matrix Double
    -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Vector Double, Matrix Double)
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix t)
    eigSHAux (CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dsyev CInt
    0) String
    "eigS'"

    -- | Eigenvalues of a hermitian complex matrix, using LAPACK's /zheev/ with jobz == \'N\'.
    -- The eigenvalues are sorted in descending order.
    eigOnlyH :: Matrix (Complex Double) -> Vector Double
    eigOnlyH :: Matrix (Complex Double) -> Vector Double
    eigOnlyH = Vector Double -> Vector Double
    vrev (Vector Double -> Vector Double)
    -> (Matrix (Complex Double) -> Vector Double)
    -> Matrix (Complex Double)
    -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (Vector Double, Matrix (Complex Double)) -> Vector Double
    forall a b. (a, b) -> a
    fst((Vector Double, Matrix (Complex Double)) -> Vector Double)
    -> (Matrix (Complex Double)
        -> (Vector Double, Matrix (Complex Double)))
    -> Matrix (Complex Double)
    -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Vector Double, Matrix (Complex Double))
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Vector a, Matrix t)
    eigSHAux (CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zheev CInt
    0) String
    "eigH'"

    vrev :: Vector Double -> Vector Double
    vrev = Matrix Double -> Vector Double
    forall t. Element t => Matrix t -> Vector t
    flatten (Matrix Double -> Vector Double)
    -> (Vector Double -> Matrix Double)
    -> Vector Double
    -> Vector Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    flipud (Matrix Double -> Matrix Double)
    -> (Vector Double -> Matrix Double)
    -> Vector Double
    -> Matrix Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Int -> Vector Double -> Matrix Double
    forall t. Storable t => Int -> Vector t -> Matrix t
    reshape Int
    1

    -----------------------------------------------------------------------------
    foreign import ccall unsafe "linearSolveR_l" dgesv :: R ::> R ::> Ok
    foreign import ccall unsafe "linearSolveC_l" zgesv :: C ::> C ::> Ok

    linearSolveSQAux :: (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux IO (Matrix t) -> IO c
    g CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Matrix t
    a Matrix t
    b
        | Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    n2 Bool -> Bool -> Bool
    && Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    r = IO c -> c
    forall a. IO a -> a
    unsafePerformIO (IO c -> c) -> (IO (Matrix t) -> IO c) -> IO (Matrix t) -> c
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix t) -> IO c
    g (IO (Matrix t) -> c) -> IO (Matrix t) -> c
    forall a b. (a -> b) -> a -> b
    $ do
            a' <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
            s  <- copy ColumnMajor b
            (a' #! s) f #| st
            return s
        | Bool
    otherwise = String -> c
    forall a. HasCallStack => String -> a
    error (String -> c) -> String -> c
    forall a b. (a -> b) -> a -> b
    $ String
    st String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " of nonsquare matrix"
      where
        n1 :: Int
    n1 = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n2 :: Int
    n2 = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        r :: Int
    r  = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b

    -- | 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 :: Matrix Double -> Matrix Double -> Matrix Double
    linearSolveR Matrix Double
    a Matrix Double
    b = (IO (Matrix Double) -> IO (Matrix Double))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Matrix Double
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {c}.
    (Element t, Element t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux IO (Matrix Double) -> IO (Matrix Double)
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesv String
    "linearSolveR" Matrix Double
    a Matrix Double
    b

    mbLinearSolveR :: Matrix Double -> Matrix Double -> Maybe (Matrix Double)
    mbLinearSolveR :: Matrix Double -> Matrix Double -> Maybe (Matrix Double)
    mbLinearSolveR Matrix Double
    a Matrix Double
    b = (IO (Matrix Double) -> IO (Maybe (Matrix Double)))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Matrix Double
    -> Matrix Double
    -> Maybe (Matrix Double)
    forall {t} {t} {c}.
    (Element t, Element t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux IO (Matrix Double) -> IO (Maybe (Matrix Double))
    forall x. IO x -> IO (Maybe x)
    mbCatch CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgesv String
    "linearSolveR" Matrix Double
    a Matrix Double
    b

    -- | 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 :: Matrix (Complex Double)
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    linearSolveC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double)))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {c}.
    (Element t, Element t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double))
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesv String
    "linearSolveC" Matrix (Complex Double)
    a Matrix (Complex Double)
    b

    mbLinearSolveC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Maybe (Matrix (Complex Double))
    mbLinearSolveC :: Matrix (Complex Double)
    -> Matrix (Complex Double) -> Maybe (Matrix (Complex Double))
    mbLinearSolveC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (IO (Matrix (Complex Double))
     -> IO (Maybe (Matrix (Complex Double))))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Maybe (Matrix (Complex Double))
    forall {t} {t} {c}.
    (Element t, Element t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux IO (Matrix (Complex Double))
    -> IO (Maybe (Matrix (Complex Double)))
    forall x. IO x -> IO (Maybe x)
    mbCatch CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgesv String
    "linearSolveC" Matrix (Complex Double)
    a Matrix (Complex Double)
    b

    --------------------------------------------------------------------------------
    foreign import ccall unsafe "cholSolveR_l" dpotrs  :: R ::> R ::> Ok
    foreign import ccall unsafe "cholSolveC_l" zpotrs  :: C ::> C ::> Ok

    linearSolveSQAux2 :: (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux2 IO (Matrix t) -> IO c
    g CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Matrix t
    a Matrix t
    b
        | Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    n2 Bool -> Bool -> Bool
    && Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    r = IO c -> c
    forall a. IO a -> a
    unsafePerformIO (IO c -> c) -> (IO (Matrix t) -> IO c) -> IO (Matrix t) -> c
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix t) -> IO c
    g (IO (Matrix t) -> c) -> IO (Matrix t) -> c
    forall a b. (a -> b) -> a -> b
    $ do
            s <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    b
            (a #! s) f #| st
            return s
        | Bool
    otherwise = String -> c
    forall a. HasCallStack => String -> a
    error (String -> c) -> String -> c
    forall a b. (a -> b) -> a -> b
    $ String
    st String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " of nonsquare matrix"
      where
        n1 :: Int
    n1 = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n2 :: Int
    n2 = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        r :: Int
    r  = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b

    -- | Solves a symmetric positive definite system of linear equations using a precomputed Cholesky factorization obtained by 'cholS'.
    cholSolveR :: Matrix Double -> Matrix Double -> Matrix Double
    cholSolveR :: Matrix Double -> Matrix Double -> Matrix Double
    cholSolveR Matrix Double
    a Matrix Double
    b = (IO (Matrix Double) -> IO (Matrix Double))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Matrix Double
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux2 IO (Matrix Double) -> IO (Matrix Double)
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dpotrs String
    "cholSolveR" (Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix Double
    a) Matrix Double
    b

    -- | Solves a Hermitian positive definite system of linear equations using a precomputed Cholesky factorization obtained by 'cholH'.
    cholSolveC :: Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
    cholSolveC :: Matrix (Complex Double)
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    cholSolveC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double)))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveSQAux2 IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double))
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zpotrs String
    "cholSolveC" (Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix (Complex Double)
    a) Matrix (Complex Double)
    b

    --------------------------------------------------------------------------------
    foreign import ccall unsafe "triSolveR_l_u" dtrtrs_u  :: R ::> R ::> Ok
    foreign import ccall unsafe "triSolveC_l_u" ztrtrs_u  :: C ::> C ::> Ok
    foreign import ccall unsafe "triSolveR_l_l" dtrtrs_l  :: R ::> R ::> Ok
    foreign import ccall unsafe "triSolveC_l_l" ztrtrs_l  :: C ::> C ::> Ok

    linearSolveTRAux2 :: (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveTRAux2 IO (Matrix t) -> IO c
    g CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Matrix t
    a Matrix t
    b
        | Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    n2 Bool -> Bool -> Bool
    && Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    r = IO c -> c
    forall a. IO a -> a
    unsafePerformIO (IO c -> c) -> (IO (Matrix t) -> IO c) -> IO (Matrix t) -> c
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix t) -> IO c
    g (IO (Matrix t) -> c) -> IO (Matrix t) -> c
    forall a b. (a -> b) -> a -> b
    $ do
            s <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    b
            (a #! s) f #| st
            return s
        | Bool
    otherwise = String -> c
    forall a. HasCallStack => String -> a
    error (String -> c) -> String -> c
    forall a b. (a -> b) -> a -> b
    $ String
    st String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " of nonsquare matrix"
      where
        n1 :: Int
    n1 = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n2 :: Int
    n2 = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        r :: Int
    r  = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b

    data UpLo = Lower | Upper

    -- | Solves a triangular system of linear equations.
    triSolveR :: UpLo -> Matrix Double -> Matrix Double -> Matrix Double
    triSolveR :: UpLo -> Matrix Double -> Matrix Double -> Matrix Double
    triSolveR UpLo
    Lower Matrix Double
    a Matrix Double
    b = (IO (Matrix Double) -> IO (Matrix Double))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Matrix Double
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveTRAux2 IO (Matrix Double) -> IO (Matrix Double)
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dtrtrs_l String
    "triSolveR" (Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix Double
    a) Matrix Double
    b
    triSolveR UpLo
    Upper Matrix Double
    a Matrix Double
    b = (IO (Matrix Double) -> IO (Matrix Double))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Matrix Double
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveTRAux2 IO (Matrix Double) -> IO (Matrix Double)
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dtrtrs_u String
    "triSolveR" (Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix Double
    a) Matrix Double
    b

    -- | Solves a triangular system of linear equations.
    triSolveC :: UpLo -> Matrix (Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
    triSolveC :: UpLo
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    triSolveC UpLo
    Lower Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double)))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveTRAux2 IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double))
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    ztrtrs_l String
    "triSolveC" (Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix (Complex Double)
    a) Matrix (Complex Double)
    b
    triSolveC UpLo
    Upper Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double)))
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {c}.
    (Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Matrix t
    -> Matrix t
    -> c
    linearSolveTRAux2 IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double))
    forall a. a -> a
    id CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    ztrtrs_u String
    "triSolveC" (Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix (Complex Double)
    a) Matrix (Complex Double)
    b

    --------------------------------------------------------------------------------
    foreign import ccall unsafe "triDiagSolveR_l" dgttrs  :: R :> R :> R :> R ::> Ok
    foreign import ccall unsafe "triDiagSolveC_l" zgttrs  :: C :> C :> C :> C ::> Ok

    linearSolveGTAux2 :: (IO (Matrix t) -> IO c)
    -> (CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Vector t
    -> Vector t
    -> Vector t
    -> Matrix t
    -> c
    linearSolveGTAux2 IO (Matrix t) -> IO c
    g CInt
    -> Ptr t
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Vector t
    dl Vector t
    d Vector t
    du Matrix t
    b
        | Int
    ndl  Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    nd Int -> Int -> Int
    forall a. Num a => a -> a -> a
    - Int
    1 Bool -> Bool -> Bool
    &&
          Int
    ndu  Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    nd Int -> Int -> Int
    forall a. Num a => a -> a -> a
    - Int
    1 Bool -> Bool -> Bool
    &&
          Int
    nd   Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    r = IO c -> c
    forall a. IO a -> a
    unsafePerformIO (IO c -> c) -> (IO (Matrix t) -> IO c) -> IO (Matrix t) -> c
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix t) -> IO c
    g (IO (Matrix t) -> c) -> IO (Matrix t) -> c
    forall a b. (a -> b) -> a -> b
    $ do
            dl' <- [Vector t] -> Vector t
    forall a. HasCallStack => [a] -> a
    head ([Vector t] -> Vector t)
    -> (Matrix t -> [Vector t]) -> Matrix t -> Vector t
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Matrix t -> [Vector t]
    forall t. Element t => Matrix t -> [Vector t]
    toRows (Matrix t -> Vector t) -> IO (Matrix t) -> IO (Vector t)
    forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
    <$> MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor ([Vector t] -> Matrix t
    forall t. Element t => [Vector t] -> Matrix t
    fromRows [Vector t
    dl])
            du' <- head . toRows <$> copy ColumnMajor (fromRows [du])
            s <- copy ColumnMajor b
            (dl' # d # du' #! s) f #| st
            return s
        | Bool
    otherwise = String -> c
    forall a. HasCallStack => String -> a
    error (String -> c) -> String -> c
    forall a b. (a -> b) -> a -> b
    $ String
    st String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++ String
    " of nonsquare matrix"
      where
        ndl :: Int
    ndl  = Vector t -> Int
    forall t. Storable t => Vector t -> Int
    dim Vector t
    dl
        nd :: Int
    nd   = Vector t -> Int
    forall t. Storable t => Vector t -> Int
    dim Vector t
    d
        ndu :: Int
    ndu  = Vector t -> Int
    forall t. Storable t => Vector t -> Int
    dim Vector t
    du
        r :: Int
    r    = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b

    -- | Solves a tridiagonal system of linear equations.
    triDiagSolveR :: Vector Double
    -> Vector Double -> Vector Double -> Matrix Double -> Matrix Double
    triDiagSolveR Vector Double
    dl Vector Double
    d Vector Double
    du Matrix Double
    b = (IO (Matrix Double) -> IO (Matrix Double))
    -> (CInt
        -> Ptr Double
        -> CInt
        -> Ptr Double
        -> CInt
        -> Ptr Double
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr Double
        -> IO CInt)
    -> String
    -> Vector Double
    -> Vector Double
    -> Vector Double
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {t} {t} {c}.
    (Element t, Element t, Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Vector t
    -> Vector t
    -> Vector t
    -> Matrix t
    -> c
    linearSolveGTAux2 IO (Matrix Double) -> IO (Matrix Double)
    forall a. a -> a
    id CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgttrs String
    "triDiagSolveR" Vector Double
    dl Vector Double
    d Vector Double
    du Matrix Double
    b
    triDiagSolveC :: Vector (Complex Double)
    -> Vector (Complex Double)
    -> Vector (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    triDiagSolveC Vector (Complex Double)
    dl Vector (Complex Double)
    d Vector (Complex Double)
    du Matrix (Complex Double)
    b = (IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double)))
    -> (CInt
        -> Ptr (Complex Double)
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> Ptr (Complex Double)
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr (Complex Double)
        -> IO CInt)
    -> String
    -> Vector (Complex Double)
    -> Vector (Complex Double)
    -> Vector (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {t} {t} {c}.
    (Element t, Element t, Element t, Storable t) =>
    (IO (Matrix t) -> IO c)
    -> (CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> Ptr t
        -> CInt
        -> CInt
        -> CInt
        -> CInt
        -> Ptr t
        -> IO CInt)
    -> String
    -> Vector t
    -> Vector t
    -> Vector t
    -> Matrix t
    -> c
    linearSolveGTAux2 IO (Matrix (Complex Double)) -> IO (Matrix (Complex Double))
    forall a. a -> a
    id CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgttrs String
    "triDiagSolveC" Vector (Complex Double)
    dl Vector (Complex Double)
    d Vector (Complex Double)
    du Matrix (Complex Double)
    b

    -----------------------------------------------------------------------------------

    foreign import ccall unsafe "linearSolveLSR_l"   dgels ::           R ::> R ::> Ok
    foreign import ccall unsafe "linearSolveLSC_l"   zgels ::           C ::> C ::> Ok
    foreign import ccall unsafe "linearSolveSVDR_l" dgelss :: Double -> R ::> R ::> Ok
    foreign import ccall unsafe "linearSolveSVDC_l" zgelss :: Double -> C ::> C ::> Ok

    linearSolveAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix a -> Matrix a
    linearSolveAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> IO CInt
    f String
    st Matrix t
    a Matrix a
    b
        | Int
    m Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Matrix a -> Int
    forall t. Matrix t -> Int
    rows Matrix a
    b = IO (Matrix a) -> Matrix a
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a) -> Matrix a) -> IO (Matrix a) -> Matrix a
    forall a b. (a -> b) -> a -> b
    $ do
            a' <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
            r  <- createMatrix ColumnMajor (max m n) nrhs
            setRect 0 0 b r
            (a' #! r) f #| st
            return r
        | Bool
    otherwise = String -> Matrix a
    forall a. HasCallStack => String -> a
    error (String -> Matrix a) -> String -> Matrix a
    forall a b. (a -> b) -> a -> b
    $ String
    "different number of rows in linearSolve ("String -> String -> String
    forall a. [a] -> [a] -> [a]
    ++String
    stString -> String -> String
    forall a. [a] -> [a] -> [a]
    ++String
    ")"
      where
        m :: Int
    m = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        nrhs :: Int
    nrhs = Matrix a -> Int
    forall t. Matrix t -> Int
    cols Matrix a
    b

    -- | 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 :: Matrix Double -> Matrix Double -> Matrix Double
    linearSolveLSR Matrix Double
    a Matrix Double
    b = (Int, Int) -> (Int, Int) -> Matrix Double -> Matrix Double
    forall a.
    Element a =>
    (Int, Int) -> (Int, Int) -> Matrix a -> Matrix a
    subMatrix (Int
    0,Int
    0) (Matrix Double -> Int
    forall t. Matrix t -> Int
    cols Matrix Double
    a, Matrix Double -> Int
    forall t. Matrix t -> Int
    cols Matrix Double
    b) (Matrix Double -> Matrix Double) -> Matrix Double -> Matrix Double
    forall a b. (a -> b) -> a -> b
    $
                         (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> Matrix Double -> Matrix Double
    forall {t} {a}.
    (Element t, Element a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix a -> Matrix a
    linearSolveAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgels String
    "linearSolverLSR" Matrix Double
    a Matrix Double
    b

    -- | 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 :: Matrix (Complex Double)
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    linearSolveLSC Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (Int, Int)
    -> (Int, Int) -> Matrix (Complex Double) -> Matrix (Complex Double)
    forall a.
    Element a =>
    (Int, Int) -> (Int, Int) -> Matrix a -> Matrix a
    subMatrix (Int
    0,Int
    0) (Matrix (Complex Double) -> Int
    forall t. Matrix t -> Int
    cols Matrix (Complex Double)
    a, Matrix (Complex Double) -> Int
    forall t. Matrix t -> Int
    cols Matrix (Complex Double)
    b) (Matrix (Complex Double) -> Matrix (Complex Double))
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    forall a b. (a -> b) -> a -> b
    $
                         (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {a}.
    (Element t, Element a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix a -> Matrix a
    linearSolveAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgels String
    "linearSolveLSC" Matrix (Complex Double)
    a Matrix (Complex Double)
    b

    -- | 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)
                    -> Matrix Double  -- ^ solution vectors (as columns)
    linearSolveSVDR :: Maybe Double -> Matrix Double -> Matrix Double -> Matrix Double
    linearSolveSVDR (Just Double
    rcond) Matrix Double
    a Matrix Double
    b = (Int, Int) -> (Int, Int) -> Matrix Double -> Matrix Double
    forall a.
    Element a =>
    (Int, Int) -> (Int, Int) -> Matrix a -> Matrix a
    subMatrix (Int
    0,Int
    0) (Matrix Double -> Int
    forall t. Matrix t -> Int
    cols Matrix Double
    a, Matrix Double -> Int
    forall t. Matrix t -> Int
    cols Matrix Double
    b) (Matrix Double -> Matrix Double) -> Matrix Double -> Matrix Double
    forall a b. (a -> b) -> a -> b
    $
                                       (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> Matrix Double -> Matrix Double
    forall {t} {a}.
    (Element t, Element a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix a -> Matrix a
    linearSolveAux (Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgelss Double
    rcond) String
    "linearSolveSVDR" Matrix Double
    a Matrix Double
    b
    linearSolveSVDR Maybe Double
    Nothing Matrix Double
    a Matrix Double
    b = Maybe Double -> Matrix Double -> Matrix Double -> Matrix Double
    linearSolveSVDR (Double -> Maybe Double
    forall a. a -> Maybe a
    Just (-Double
    1)) Matrix Double
    a Matrix Double
    b

    -- | 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)
                    -> Matrix (Complex Double) -- ^ solution vectors (as columns)
    linearSolveSVDC :: Maybe Double
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    linearSolveSVDC (Just Double
    rcond) Matrix (Complex Double)
    a Matrix (Complex Double)
    b = (Int, Int)
    -> (Int, Int) -> Matrix (Complex Double) -> Matrix (Complex Double)
    forall a.
    Element a =>
    (Int, Int) -> (Int, Int) -> Matrix a -> Matrix a
    subMatrix (Int
    0,Int
    0) (Matrix (Complex Double) -> Int
    forall t. Matrix t -> Int
    cols Matrix (Complex Double)
    a, Matrix (Complex Double) -> Int
    forall t. Matrix t -> Int
    cols Matrix (Complex Double)
    b) (Matrix (Complex Double) -> Matrix (Complex Double))
    -> Matrix (Complex Double) -> Matrix (Complex Double)
    forall a b. (a -> b) -> a -> b
    $
                                       (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {a}.
    (Element t, Element a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> IO CInt)
    -> String -> Matrix t -> Matrix a -> Matrix a
    linearSolveAux (Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgelss Double
    rcond) String
    "linearSolveSVDC" Matrix (Complex Double)
    a Matrix (Complex Double)
    b
    linearSolveSVDC Maybe Double
    Nothing Matrix (Complex Double)
    a Matrix (Complex Double)
    b = Maybe Double
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    linearSolveSVDC (Double -> Maybe Double
    forall a. a -> Maybe a
    Just (-Double
    1)) Matrix (Complex Double)
    a Matrix (Complex Double)
    b

    -----------------------------------------------------------------------------------

    foreign import ccall unsafe "chol_l_H" zpotrf :: C ::> Ok
    foreign import ccall unsafe "chol_l_S" dpotrf :: R ::> Ok

    cholAux :: (CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> IO (Matrix t)
    cholAux CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    a = do
        r <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
        (r # id) f #| st
        return r

    -- | Cholesky factorization of a complex Hermitian positive definite matrix, using LAPACK's /zpotrf/.
    cholH :: Matrix (Complex Double) -> Matrix (Complex Double)
    cholH :: Matrix (Complex Double) -> Matrix (Complex Double)
    cholH = IO (Matrix (Complex Double)) -> Matrix (Complex Double)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix (Complex Double)) -> Matrix (Complex Double))
    -> (Matrix (Complex Double) -> IO (Matrix (Complex Double)))
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt -> CInt -> CInt -> CInt -> Ptr (Complex Double) -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> IO (Matrix (Complex Double))
    forall {t}.
    Element t =>
    (CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> IO (Matrix t)
    cholAux CInt -> CInt -> CInt -> CInt -> Ptr (Complex Double) -> IO CInt
    zpotrf String
    "cholH"

    -- | Cholesky factorization of a real symmetric positive definite matrix, using LAPACK's /dpotrf/.
    cholS :: Matrix Double -> Matrix Double
    cholS :: Matrix Double -> Matrix Double
    cholS =  IO (Matrix Double) -> Matrix Double
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix Double) -> Matrix Double)
    -> (Matrix Double -> IO (Matrix Double))
    -> Matrix Double
    -> Matrix Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt -> CInt -> CInt -> CInt -> Ptr Double -> IO CInt)
    -> String -> Matrix Double -> IO (Matrix Double)
    forall {t}.
    Element t =>
    (CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> IO (Matrix t)
    cholAux CInt -> CInt -> CInt -> CInt -> Ptr Double -> IO CInt
    dpotrf String
    "cholS"

    -- | Cholesky factorization of a complex Hermitian positive definite matrix, using LAPACK's /zpotrf/ ('Maybe' version).
    mbCholH :: Matrix (Complex Double) -> Maybe (Matrix (Complex Double))
    mbCholH :: Matrix (Complex Double) -> Maybe (Matrix (Complex Double))
    mbCholH = IO (Maybe (Matrix (Complex Double)))
    -> Maybe (Matrix (Complex Double))
    forall a. IO a -> a
    unsafePerformIO (IO (Maybe (Matrix (Complex Double)))
     -> Maybe (Matrix (Complex Double)))
    -> (Matrix (Complex Double)
        -> IO (Maybe (Matrix (Complex Double))))
    -> Matrix (Complex Double)
    -> Maybe (Matrix (Complex Double))
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix (Complex Double))
    -> IO (Maybe (Matrix (Complex Double)))
    forall x. IO x -> IO (Maybe x)
    mbCatch (IO (Matrix (Complex Double))
     -> IO (Maybe (Matrix (Complex Double))))
    -> (Matrix (Complex Double) -> IO (Matrix (Complex Double)))
    -> Matrix (Complex Double)
    -> IO (Maybe (Matrix (Complex Double)))
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt -> CInt -> CInt -> CInt -> Ptr (Complex Double) -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> IO (Matrix (Complex Double))
    forall {t}.
    Element t =>
    (CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> IO (Matrix t)
    cholAux CInt -> CInt -> CInt -> CInt -> Ptr (Complex Double) -> IO CInt
    zpotrf String
    "cholH"

    -- | Cholesky factorization of a real symmetric positive definite matrix, using LAPACK's /dpotrf/  ('Maybe' version).
    mbCholS :: Matrix Double -> Maybe (Matrix Double)
    mbCholS :: Matrix Double -> Maybe (Matrix Double)
    mbCholS =  IO (Maybe (Matrix Double)) -> Maybe (Matrix Double)
    forall a. IO a -> a
    unsafePerformIO (IO (Maybe (Matrix Double)) -> Maybe (Matrix Double))
    -> (Matrix Double -> IO (Maybe (Matrix Double)))
    -> Matrix Double
    -> Maybe (Matrix Double)
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . IO (Matrix Double) -> IO (Maybe (Matrix Double))
    forall x. IO x -> IO (Maybe x)
    mbCatch (IO (Matrix Double) -> IO (Maybe (Matrix Double)))
    -> (Matrix Double -> IO (Matrix Double))
    -> Matrix Double
    -> IO (Maybe (Matrix Double))
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (CInt -> CInt -> CInt -> CInt -> Ptr Double -> IO CInt)
    -> String -> Matrix Double -> IO (Matrix Double)
    forall {t}.
    Element t =>
    (CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> IO (Matrix t)
    cholAux CInt -> CInt -> CInt -> CInt -> Ptr Double -> IO CInt
    dpotrf String
    "cholS"

    -----------------------------------------------------------------------------------

    type TMVM t = t ::> t :> t ::> Ok

    foreign import ccall unsafe "qr_l_R" dgeqr2 :: R :> R ::> Ok
    foreign import ccall unsafe "qr_l_C" zgeqr2 :: C :> C ::> Ok

    -- | QR factorization of a real matrix, using LAPACK's /dgeqr2/.
    qrR :: Matrix Double -> (Matrix Double, Vector Double)
    qrR :: Matrix Double -> (Matrix Double, Vector Double)
    qrR = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, Vector Double)
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    qrAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgeqr2 String
    "qrR"

    -- | QR factorization of a complex matrix, using LAPACK's /zgeqr2/.
    qrC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector (Complex Double))
    qrC :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector (Complex Double))
    qrC = (CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector (Complex Double))
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    qrAux CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgeqr2 String
    "qrC"

    qrAux :: (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    qrAux CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    a = IO (Matrix t, Vector a) -> (Matrix t, Vector a)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t, Vector a) -> (Matrix t, Vector a))
    -> IO (Matrix t, Vector a) -> (Matrix t, Vector a)
    forall a b. (a -> b) -> a -> b
    $ do
        r <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
        tau <- createVector mn
        (tau #! r) f #| st
        return (r,tau)
      where
        m :: Int
    m = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        mn :: Int
    mn = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    m Int
    n

    foreign import ccall unsafe "c_dorgqr" dorgqr :: R :> R ::> Ok
    foreign import ccall unsafe "c_zungqr" zungqr :: C :> C ::> Ok

    -- | build rotation from reflectors
    qrgrR :: Int -> (Matrix Double, Vector Double) -> Matrix Double
    qrgrR :: Int -> (Matrix Double, Vector Double) -> Matrix Double
    qrgrR = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Int -> (Matrix Double, Vector Double) -> Matrix Double
    forall {t} {t}.
    (Element t, Element t, Num t) =>
    (CInt -> Ptr t -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Int -> (Matrix t, Vector t) -> Matrix t
    qrgrAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dorgqr String
    "qrgrR"
    -- | build rotation from reflectors
    qrgrC :: Int -> (Matrix (Complex Double), Vector (Complex Double)) -> Matrix (Complex Double)
    qrgrC :: Int
    -> (Matrix (Complex Double), Vector (Complex Double))
    -> Matrix (Complex Double)
    qrgrC = (CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Int
    -> (Matrix (Complex Double), Vector (Complex Double))
    -> Matrix (Complex Double)
    forall {t} {t}.
    (Element t, Element t, Num t) =>
    (CInt -> Ptr t -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Int -> (Matrix t, Vector t) -> Matrix t
    qrgrAux CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zungqr String
    "qrgrC"

    qrgrAux :: (CInt -> Ptr t -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Int -> (Matrix t, Vector t) -> Matrix t
    qrgrAux CInt -> Ptr t -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Int
    n (Matrix t
    a, Vector t
    tau) = IO (Matrix t) -> Matrix t
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t) -> Matrix t) -> IO (Matrix t) -> Matrix t
    forall a b. (a -> b) -> a -> b
    $ do
        res <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor ((Int, Int) -> (Int, Int) -> Matrix t -> Matrix t
    forall a.
    Element a =>
    (Int, Int) -> (Int, Int) -> Matrix a -> Matrix a
    subMatrix (Int
    0,Int
    0) (Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a,Int
    n) Matrix t
    a)
        ((subVector 0 n tau') #! res) f #| st
        return res
      where
        tau' :: Vector t
    tau' = [Vector t] -> Vector t
    forall t. Storable t => [Vector t] -> Vector t
    vjoin [Vector t
    tau, t -> Int -> Vector t
    forall a. Element a => a -> Int -> Vector a
    constantD t
    0 Int
    n]

    -----------------------------------------------------------------------------------
    foreign import ccall unsafe "hess_l_R" dgehrd :: R :> R ::> Ok
    foreign import ccall unsafe "hess_l_C" zgehrd :: C :> C ::> Ok

    -- | Hessenberg factorization of a square real matrix, using LAPACK's /dgehrd/.
    hessR :: Matrix Double -> (Matrix Double, Vector Double)
    hessR :: Matrix Double -> (Matrix Double, Vector Double)
    hessR = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, Vector Double)
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    hessAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgehrd String
    "hessR"

    -- | Hessenberg factorization of a square complex matrix, using LAPACK's /zgehrd/.
    hessC :: Matrix (Complex Double) -> (Matrix (Complex Double), Vector (Complex Double))
    hessC :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector (Complex Double))
    hessC = (CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Vector (Complex Double))
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    hessAux CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgehrd String
    "hessC"

    hessAux :: (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, Vector a)
    hessAux CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    a = IO (Matrix t, Vector a) -> (Matrix t, Vector a)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t, Vector a) -> (Matrix t, Vector a))
    -> IO (Matrix t, Vector a) -> (Matrix t, Vector a)
    forall a b. (a -> b) -> a -> b
    $ do
        r <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
        tau <- createVector (mn-1)
        (tau #! r) f #| st
        return (r,tau)
      where
        m :: Int
    m = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        mn :: Int
    mn = Int -> Int -> Int
    forall a. Ord a => a -> a -> a
    min Int
    m Int
    n

    -----------------------------------------------------------------------------------
    foreign import ccall unsafe "schur_l_R" dgees :: R ::> R ::> Ok
    foreign import ccall unsafe "schur_l_C" zgees :: C ::> C ::> Ok

    -- | Schur factorization of a square real matrix, using LAPACK's /dgees/.
    schurR :: Matrix Double -> (Matrix Double, Matrix Double)
    schurR :: Matrix Double -> (Matrix Double, Matrix Double)
    schurR = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, Matrix Double)
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Matrix t)
    schurAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgees String
    "schurR"

    -- | Schur factorization of a square complex matrix, using LAPACK's /zgees/.
    schurC :: Matrix (Complex Double) -> (Matrix (Complex Double), Matrix (Complex Double))
    schurC :: Matrix (Complex Double)
    -> (Matrix (Complex Double), Matrix (Complex Double))
    schurC = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), Matrix (Complex Double))
    forall {t} {a}.
    (Element t, Storable a) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Matrix t)
    schurAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgees String
    "schurC"

    schurAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr a
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> (Matrix a, Matrix t)
    schurAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr a
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Matrix t
    a = IO (Matrix a, Matrix t) -> (Matrix a, Matrix t)
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix a, Matrix t) -> (Matrix a, Matrix t))
    -> IO (Matrix a, Matrix t) -> (Matrix a, Matrix t)
    forall a b. (a -> b) -> a -> b
    $ do
        u <- MatrixOrder -> Int -> Int -> IO (Matrix a)
    forall a. Storable a => MatrixOrder -> Int -> Int -> IO (Matrix a)
    createMatrix MatrixOrder
    ColumnMajor Int
    n Int
    n
        s <- copy ColumnMajor a
        (u #! s) f #| st
        return (u,s)
      where
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a

    -----------------------------------------------------------------------------------
    foreign import ccall unsafe "lu_l_R" dgetrf :: R :> R ::> Ok
    foreign import ccall unsafe "lu_l_C" zgetrf :: R :> C ::> Ok

    -- | LU factorization of a general real matrix, using LAPACK's /dgetrf/.
    luR :: Matrix Double -> (Matrix Double, [Int])
    luR :: Matrix Double -> (Matrix Double, [Int])
    luR = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, [Int])
    forall {t} {a} {b}.
    (Element t, Storable a, RealFrac a, Integral b) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    luAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgetrf String
    "luR"

    -- | LU factorization of a general complex matrix, using LAPACK's /zgetrf/.
    luC :: Matrix (Complex Double) -> (Matrix (Complex Double), [Int])
    luC :: Matrix (Complex Double) -> (Matrix (Complex Double), [Int])
    luC = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), [Int])
    forall {t} {a} {b}.
    (Element t, Storable a, RealFrac a, Integral b) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    luAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgetrf String
    "luC"

    luAux :: (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    luAux CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    a = IO (Matrix t, [b]) -> (Matrix t, [b])
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t, [b]) -> (Matrix t, [b]))
    -> IO (Matrix t, [b]) -> (Matrix t, [b])
    forall a b. (a -> b) -> a -> b
    $ do
        lu <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
        piv <- createVector (min n m)
        (piv #! lu) f #| st
        return (lu, map (pred.round) (toList piv))
      where
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        m :: Int
    m = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a

    -----------------------------------------------------------------------------------

    foreign import ccall unsafe "luS_l_R" dgetrs :: R ::> R :> R ::> Ok
    foreign import ccall unsafe "luS_l_C" zgetrs :: C ::> R :> C ::> Ok

    -- | Solve a real linear system from a precomputed LU decomposition ('luR'), using LAPACK's /dgetrs/.
    lusR :: Matrix Double -> [Int] -> Matrix Double -> Matrix Double
    lusR :: Matrix Double -> [Int] -> Matrix Double -> Matrix Double
    lusR Matrix Double
    a [Int]
    piv Matrix Double
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> [Int]
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {b}.
    (Element t, Storable t, Integral b) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> [b] -> Matrix t -> Matrix t
    lusAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dgetrs String
    "lusR" (Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix Double
    a) [Int]
    piv Matrix Double
    b

    -- | Solve a complex linear system from a precomputed LU decomposition ('luC'), using LAPACK's /zgetrs/.
    lusC :: Matrix (Complex Double) -> [Int] -> Matrix (Complex Double) -> Matrix (Complex Double)
    lusC :: Matrix (Complex Double)
    -> [Int] -> Matrix (Complex Double) -> Matrix (Complex Double)
    lusC Matrix (Complex Double)
    a [Int]
    piv Matrix (Complex Double)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> [Int]
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {b}.
    (Element t, Storable t, Integral b) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> [b] -> Matrix t -> Matrix t
    lusAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zgetrs String
    "lusC" (Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix (Complex Double)
    a) [Int]
    piv Matrix (Complex Double)
    b

    lusAux :: (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> [b] -> Matrix t -> Matrix t
    lusAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr t
    -> IO CInt
    f String
    st Matrix t
    a [b]
    piv Matrix t
    b
        | Int
    n1Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    n2 Bool -> Bool -> Bool
    && Int
    n2Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    ==Int
    n =IO (Matrix t) -> Matrix t
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t) -> Matrix t) -> IO (Matrix t) -> Matrix t
    forall a b. (a -> b) -> a -> b
    $ do
             x <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    b
             (a # piv' #! x) f #| st
             return x
        | Bool
    otherwise = String -> Matrix t
    forall a. HasCallStack => String -> a
    error String
    st
      where
        n1 :: Int
    n1 = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    a
        n2 :: Int
    n2 = Matrix t -> Int
    forall t. Matrix t -> Int
    cols Matrix t
    a
        n :: Int
    n = Matrix t -> Int
    forall t. Matrix t -> Int
    rows Matrix t
    b
        piv' :: Vector Double
    piv' = [Double] -> Vector Double
    forall a. Storable a => [a] -> Vector a
    fromList ((b -> Double) -> [b] -> [Double]
    forall a b. (a -> b) -> [a] -> [b]
    map (b -> Double
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral(b -> Double) -> (b -> b) -> b -> Double
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    .b -> b
    forall a. Enum a => a -> a
    succ) [b]
    piv) :: Vector Double

    -----------------------------------------------------------------------------------
    foreign import ccall unsafe "ldl_R" dsytrf :: R :> R ::> Ok
    foreign import ccall unsafe "ldl_C" zhetrf :: R :> C ::> Ok

    -- | LDL factorization of a symmetric real matrix, using LAPACK's /dsytrf/.
    ldlR :: Matrix Double -> (Matrix Double, [Int])
    ldlR :: Matrix Double -> (Matrix Double, [Int])
    ldlR = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String -> Matrix Double -> (Matrix Double, [Int])
    forall {t} {a} {b}.
    (Element t, Storable a, RealFrac a, Integral b) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    ldlAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dsytrf String
    "ldlR"

    -- | LDL factorization of a hermitian complex matrix, using LAPACK's /zhetrf/.
    ldlC :: Matrix (Complex Double) -> (Matrix (Complex Double), [Int])
    ldlC :: Matrix (Complex Double) -> (Matrix (Complex Double), [Int])
    ldlC = (CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> (Matrix (Complex Double), [Int])
    forall {t} {a} {b}.
    (Element t, Storable a, RealFrac a, Integral b) =>
    (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    ldlAux CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zhetrf String
    "ldlC"

    ldlAux :: (CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt)
    -> String -> Matrix t -> (Matrix t, [b])
    ldlAux CInt -> Ptr a -> CInt -> CInt -> CInt -> CInt -> Ptr t -> IO CInt
    f String
    st Matrix t
    a = IO (Matrix t, [b]) -> (Matrix t, [b])
    forall a. IO a -> a
    unsafePerformIO (IO (Matrix t, [b]) -> (Matrix t, [b]))
    -> IO (Matrix t, [b]) -> (Matrix t, [b])
    forall a b. (a -> b) -> a -> b
    $ do
        ldl <- MatrixOrder -> Matrix t -> IO (Matrix t)
    forall t. Element t => MatrixOrder -> Matrix t -> IO (Matrix t)
    copy MatrixOrder
    ColumnMajor Matrix t
    a
        piv <- createVector (rows a)
        (piv #! ldl) f #| st
        return (ldl, map (pred.round) (toList piv))

    -----------------------------------------------------------------------------------

    foreign import ccall unsafe "ldl_S_R" dsytrs :: R ::> R :> R ::> Ok
    foreign import ccall unsafe "ldl_S_C" zsytrs :: C ::> R :> C ::> Ok

    -- | Solve a real linear system from a precomputed LDL decomposition ('ldlR'), using LAPACK's /dsytrs/.
    ldlsR :: Matrix Double -> [Int] -> Matrix Double -> Matrix Double
    ldlsR :: Matrix Double -> [Int] -> Matrix Double -> Matrix Double
    ldlsR Matrix Double
    a [Int]
    piv Matrix Double
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr Double
     -> IO CInt)
    -> String
    -> Matrix Double
    -> [Int]
    -> Matrix Double
    -> Matrix Double
    forall {t} {t} {b}.
    (Element t, Storable t, Integral b) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> [b] -> Matrix t -> Matrix t
    lusAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr Double
    -> IO CInt
    dsytrs String
    "ldlsR" (Matrix Double -> Matrix Double
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix Double
    a) [Int]
    piv Matrix Double
    b

    -- | Solve a complex linear system from a precomputed LDL decomposition ('ldlC'), using LAPACK's /zsytrs/.
    ldlsC :: Matrix (Complex Double) -> [Int] -> Matrix (Complex Double) -> Matrix (Complex Double)
    ldlsC :: Matrix (Complex Double)
    -> [Int] -> Matrix (Complex Double) -> Matrix (Complex Double)
    ldlsC Matrix (Complex Double)
    a [Int]
    piv Matrix (Complex Double)
    b = (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr (Complex Double)
     -> IO CInt)
    -> String
    -> Matrix (Complex Double)
    -> [Int]
    -> Matrix (Complex Double)
    -> Matrix (Complex Double)
    forall {t} {t} {b}.
    (Element t, Storable t, Integral b) =>
    (CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> CInt
     -> Ptr Double
     -> CInt
     -> CInt
     -> CInt
     -> CInt
     -> Ptr t
     -> IO CInt)
    -> String -> Matrix t -> [b] -> Matrix t -> Matrix t
    lusAux CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> CInt
    -> Ptr Double
    -> CInt
    -> CInt
    -> CInt
    -> CInt
    -> Ptr (Complex Double)
    -> IO CInt
    zsytrs String
    "ldlsC" (Matrix (Complex Double) -> Matrix (Complex Double)
    forall t. Element t => Matrix t -> Matrix t
    fmat Matrix (Complex Double)
    a) [Int]
    piv Matrix (Complex Double)
    b