module Numeric.LAPACK.Matrix.HermitianPositiveDefinite (
   solve,
   solveDecomposed,
   inverse,
   decompose,
   determinant,
   ) where

import qualified Numeric.LAPACK.Matrix.HermitianPositiveDefinite.Linear
                                                                  as Linear
import qualified Numeric.LAPACK.Matrix.Array as ArrMatrix
import qualified Numeric.LAPACK.Matrix.Extent as Extent
import Numeric.LAPACK.Matrix.Array.Triangular (Hermitian, Upper)
import Numeric.LAPACK.Matrix.Array (Full)
import Numeric.LAPACK.Scalar (RealOf)

import qualified Numeric.Netlib.Class as Class

import qualified Data.Array.Comfort.Shape as Shape


solve ::
   (Extent.C vert, Extent.C horiz,
    Shape.C sh, Eq sh, Shape.C nrhs, Class.Floating a) =>
   Hermitian sh a -> Full vert horiz sh nrhs a -> Full vert horiz sh nrhs a
solve = ArrMatrix.lift2 Linear.solve

{- |
> solve a b == solveDecomposed (decompose a) b
> solve (gramian u) b == solveDecomposed u b
-}
solveDecomposed ::
   (Extent.C vert, Extent.C horiz,
    Shape.C sh, Eq sh, Shape.C nrhs, Class.Floating a) =>
   Upper sh a -> Full vert horiz sh nrhs a -> Full vert horiz sh nrhs a
solveDecomposed = ArrMatrix.lift2 Linear.solveDecomposed

inverse :: (Shape.C sh, Class.Floating a) => Hermitian sh a -> Hermitian sh a
inverse = ArrMatrix.lift1 Linear.inverse

{- |
Cholesky decomposition
-}
decompose :: (Shape.C sh, Class.Floating a) => Hermitian sh a -> Upper sh a
decompose = ArrMatrix.lift1 Linear.decompose

determinant :: (Shape.C sh, Class.Floating a) => Hermitian sh a -> RealOf a
determinant = Linear.determinant . ArrMatrix.toVector