Safe Haskell | None |
---|---|
Language | Haskell2010 |
- qr :: (Epsilon a, Floating a, Real a) => SpMatrix a -> (SpMatrix a, SpMatrix a)
- lu :: (Epsilon a, Fractional a, Real a) => SpMatrix a -> (SpMatrix a, SpMatrix a)
- chol :: (Epsilon a, Real a, Floating a) => SpMatrix a -> SpMatrix a
- conditionNumberSM :: (Epsilon a, RealFloat a) => SpMatrix a -> a
- hhMat :: Num a => a -> SpVector a -> SpMatrix a
- hhRefl :: Num a => SpVector a -> SpMatrix a
- givens :: (Floating a, Epsilon a, Ord a) => SpMatrix a -> IxRow -> IxCol -> SpMatrix a
- arnoldi :: (Floating a, Eq a) => SpMatrix a -> Int -> (SpMatrix a, SpMatrix a)
- eigsQR :: (Epsilon a, Real a, Floating a) => Int -> SpMatrix a -> SpVector a
- eigRayleigh :: Int -> SpMatrix Double -> (SpVector Double, Double) -> (SpVector Double, Double)
- linSolve :: LinSolveMethod -> SpMatrix Double -> SpVector Double -> SpVector Double
- data LinSolveMethod
- (<\>) :: SpMatrix Double -> SpVector Double -> SpVector Double
- luSolve :: (Fractional a, Eq a, Epsilon a) => SpMatrix a -> SpMatrix a -> SpVector a -> SpVector a
- cgne :: SpMatrix Double -> SpVector Double -> SpVector Double -> CGNE
- tfqmr :: SpMatrix Double -> SpVector Double -> SpVector Double -> TFQMR
- bicgstab :: SpMatrix Double -> SpVector Double -> SpVector Double -> SpVector Double -> BICGSTAB
- cgs :: SpMatrix Double -> SpVector Double -> SpVector Double -> SpVector Double -> CGS
- bcg :: SpMatrix Double -> SpVector Double -> SpVector Double -> BCG
- _xCgne :: CGNE -> SpVector Double
- _xTfq :: TFQMR -> SpVector Double
- _xBicgstab :: BICGSTAB -> SpVector Double
- _x :: CGS -> SpVector Double
- _xBcg :: BCG -> SpVector Double
- cgsStep :: SpMatrix Double -> SpVector Double -> CGS -> CGS
- bicgstabStep :: SpMatrix Double -> SpVector Double -> BICGSTAB -> BICGSTAB
- data CGNE
- data TFQMR
- data BICGSTAB
- data CGS
- data BCG
- ilu0 :: (Epsilon a, Real a, Fractional a) => SpMatrix a -> (SpMatrix a, SpMatrix a)
- mSsor :: Fractional a => SpMatrix a -> a -> (SpMatrix a, SpMatrix a)
- diagPartitions :: SpMatrix a -> (SpMatrix a, SpMatrix a, SpMatrix a)
- randArray :: PrimMonad m => Int -> Double -> Double -> m [Double]
- randMat :: PrimMonad m => Int -> m (SpMatrix Double)
- randVec :: PrimMonad m => Int -> m (SpVector Double)
- randSpMat :: Int -> Int -> IO (SpMatrix Double)
- randSpVec :: Int -> Int -> IO (SpVector Double)
- sparsifySV :: Epsilon a => SpVector a -> SpVector a
- modifyInspectN :: MonadState s m => Int -> ([s] -> Bool) -> (s -> s) -> m s
- runAppendN' :: (t -> t) -> Int -> t -> [t]
- untilConverged :: MonadState a m => (a -> SpVector Double) -> (a -> a) -> m a
- diffSqL :: Floating a => [a] -> a
Matrix factorizations
qr :: (Epsilon a, Floating a, Real a) => SpMatrix a -> (SpMatrix a, SpMatrix a) Source #
Given a matrix A, returns a pair of matrices (Q, R) such that Q R = A, Q is orthogonal and R is upper triangular. Applies Givens rotation iteratively to zero out sub-diagonal elements
lu :: (Epsilon a, Fractional a, Real a) => SpMatrix a -> (SpMatrix a, SpMatrix a) Source #
Given a matrix A, returns a pair of matrices (L, U) such that L U = A
chol :: (Epsilon a, Real a, Floating a) => SpMatrix a -> SpMatrix a Source #
Given a positive semidefinite matrix A, returns a lower-triangular matrix L such that L L^T = A
Condition number
conditionNumberSM :: (Epsilon a, RealFloat a) => SpMatrix a -> a Source #
uses the R matrix from the QR factorization
Householder reflection
hhRefl :: Num a => SpVector a -> SpMatrix a Source #
a vector x
uniquely defines an orthogonal plane; the Householder operator reflects any point v
with respect to this plane:
v' = (I - 2 x >< x) v
Givens' rotation
givens :: (Floating a, Epsilon a, Ord a) => SpMatrix a -> IxRow -> IxCol -> SpMatrix a Source #
Givens method, row version: choose other row index i' s.t. i' is : * below the diagonal * corresponding element is nonzero
QR.C1 ) To zero out entry A(i, j) we must find row k such that A(k, j) is non-zero but A has zeros in row k for all columns less than j.
NB: the current version is quite inefficient in that:
1. the Givens' matrix G_i
is different from Identity only in 4 entries
2. at each iteration i
we multiply G_i
by the previous partial result M
. Since this corresponds to a rotation, and the givensCoef
function already computes the value of the resulting non-zero component (output r
), `G_i ## M` can be simplified by just changing two entries of M
(i.e. zeroing one out and changing the other into r
).
Arnoldi iteration
arnoldi :: (Floating a, Eq a) => SpMatrix a -> Int -> (SpMatrix a, SpMatrix a) Source #
Given a matrix A and a positive integer n
, this procedure finds the basis of an order n
Krylov subspace (as the columns of matrix Q), along with an upper Hessenberg matrix H, such that A = Q^T H Q.
At the i`th iteration, it finds (i + 1) coefficients (the i`th column of the Hessenberg matrix H) and the (i + 1)`th Krylov vector.
Eigensolvers
eigsQR :: (Epsilon a, Real a, Floating a) => Int -> SpMatrix a -> SpVector a Source #
`eigsQR n mm` performs n
iterations of the QR algorithm on matrix mm
, and returns a SpVector containing all eigenvalues
eigRayleigh :: Int -> SpMatrix Double -> (SpVector Double, Double) -> (SpVector Double, Double) Source #
`eigsRayleigh n mm` performs n
iterations of the Rayleigh algorithm on matrix mm
and returns the eigenpair closest to the initialization. It displays cubic-order convergence, but it also requires an educated guess on the initial eigenpair
eigRayleigh :: Int -- max # iterations
-> SpMatrix Double -- matrix
-> (SpVector Double, Double) -- initial guess of (eigenvector, eigenvalue)
-> (SpVector Double, Double) -- final estimate of (eigenvector, eigenvalue)
Linear solvers
linSolve :: LinSolveMethod -> SpMatrix Double -> SpVector Double -> SpVector Double Source #
Linear solve with _deterministic_ starting vector (every component at 0.1)
(<\>) :: SpMatrix Double -> SpVector Double -> SpVector Double Source #
\ : linSolve using the BiCGSTAB method as default
Direct methods
luSolve :: (Fractional a, Eq a, Epsilon a) => SpMatrix a -> SpMatrix a -> SpVector a -> SpVector a Source #
Direct solver based on a triangular factorization of the system matrix.
Iterative methods
bicgstab :: SpMatrix Double -> SpVector Double -> SpVector Double -> SpVector Double -> BICGSTAB Source #
iterate solver until convergence or until max # of iterations is reached
cgs :: SpMatrix Double -> SpVector Double -> SpVector Double -> SpVector Double -> CGS Source #
iterate solver until convergence or until max # of iterations is reached
bicgstabStep :: SpMatrix Double -> SpVector Double -> BICGSTAB -> BICGSTAB Source #
one step of BiCGSTAB
Preconditioners
ilu0 :: (Epsilon a, Real a, Fractional a) => SpMatrix a -> (SpMatrix a, SpMatrix a) Source #
used for Incomplete LU : remove entries in m
corresponding to zero entries in m2
mSsor :: Fractional a => SpMatrix a -> a -> (SpMatrix a, SpMatrix a) Source #
`mSsor aa omega` : if `omega = 1` it returns the symmetric Gauss-Seidel preconditioner
Matrix partitioning
diagPartitions :: SpMatrix a -> (SpMatrix a, SpMatrix a, SpMatrix a) Source #
Partition a matrix into strictly subdiagonal, diagonal and strictly superdiagonal parts
Random arrays
Random matrices and vectors
Sparse "
Sparsify data
Iteration combinators
modifyInspectN :: MonadState s m => Int -> ([s] -> Bool) -> (s -> s) -> m s Source #
Keep a moving window buffer (length 2) of state x
to assess convergence, stop when either a condition on that list is satisfied or when max # of iterations is reached (i.e. same thing as loopUntilAcc
but this one runs in the State monad)
runAppendN' :: (t -> t) -> Int -> t -> [t] Source #
", NO convergence check
untilConverged :: MonadState a m => (a -> SpVector Double) -> (a -> a) -> m a Source #
iterate until convergence is verified or we run out of a fixed iteration budget