{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE GADTs #-}
module Numeric.LAPACK.Matrix.Array.Hermitian (
   Semidefinite,
   AnyHermitian,
   AnyHermitianP,
   assureFullRank,
   assureAnyRank,
   relaxSemidefinite,
   relaxIndefinite,
   assurePositiveDefiniteness,
   relaxDefiniteness,
   asUnknownDefiniteness,
   ) where

import qualified Numeric.LAPACK.Matrix.Shape.Omni as Omni
import Numeric.LAPACK.Matrix.Array.Private (Matrix(Array), Quadratic)
import Numeric.LAPACK.Matrix.Layout.Private (Packed)

import qualified Type.Data.Bool as TBool
import Type.Data.Bool (False, True)

import qualified Data.Array.Comfort.Storable.Unchecked as Array
import Data.Function.HT (Id)


class (TBool.C neg, TBool.C pos) => Semidefinite neg pos where
instance Semidefinite False True where
instance Semidefinite True False where

type AnyHermitian neg zero pos bands sh =
      AnyHermitianP Packed neg zero pos bands sh
type AnyHermitianP pack neg zero pos bands sh =
      Quadratic pack (Omni.Hermitian neg zero pos) bands bands sh

{-
If a semidefinite matrix has full rank, then it is definite.
You can verify this by watching the eigenvalue decomposition of the matrix.
-}
assureFullRank ::
   (Semidefinite neg pos, TBool.C zero) =>
   AnyHermitianP pack neg zero pos bands sh a ->
   AnyHermitianP pack neg False pos bands sh a
assureFullRank = alterDefiniteness

assureAnyRank ::
   (Semidefinite neg pos, TBool.C zero) =>
   AnyHermitianP pack neg True pos bands sh a ->
   AnyHermitianP pack neg zero pos bands sh a
assureAnyRank = alterDefiniteness

relaxSemidefinite ::
   (TBool.C neg, TBool.C zero, TBool.C pos) =>
   AnyHermitianP pack neg False pos bands sh a ->
   AnyHermitianP pack neg zero pos bands sh a
relaxSemidefinite = alterDefiniteness

relaxIndefinite ::
   (TBool.C neg, TBool.C zero, TBool.C pos) =>
   AnyHermitianP pack neg zero pos bands sh a ->
   Quadratic pack Omni.HermitianUnknownDefiniteness bands bands sh a
relaxIndefinite = alterDefiniteness

assurePositiveDefiniteness ::
   (TBool.C neg, TBool.C zero, TBool.C pos) =>
   AnyHermitianP pack neg zero pos bands sh a ->
   Quadratic pack Omni.HermitianPositiveDefinite bands bands sh a
assurePositiveDefiniteness = alterDefiniteness

relaxDefiniteness ::
   (TBool.C neg, TBool.C zero, TBool.C pos) =>
   Quadratic pack Omni.HermitianPositiveDefinite bands bands sh a ->
   AnyHermitianP pack neg zero pos bands sh a
relaxDefiniteness = alterDefiniteness

alterDefiniteness ::
   (TBool.C negA, TBool.C zeroA, TBool.C posA,
    TBool.C negB, TBool.C zeroB, TBool.C posB) =>
   AnyHermitianP pack negA zeroA posA bands sh a ->
   AnyHermitianP pack negB zeroB posB bands sh a
alterDefiniteness (Array a) =
   Array $ flip Array.mapShape a $ \omni ->
      case omni of
         Omni.Full sh -> Omni.Full sh
         Omni.Hermitian sh -> Omni.Hermitian sh
         Omni.BandedHermitian sh -> Omni.BandedHermitian sh


asUnknownDefiniteness ::
   Id (Quadratic pack Omni.HermitianUnknownDefiniteness bands bands sh a)
asUnknownDefiniteness = id