{-# LANGUAGE TypeFamilies #-}
module Test.Multiply (
   multiplySquare,
   squareSquare,
   power,

   multiplyIdentityVector,
   multiplyIdentityFull,
   multiplyVector,
   multiplyFull,
   multiplyVectorLeft,
   multiplyVectorRight,
   multiplyLeft,
   multiplyRight,
   ) where

import qualified Test.Utility as Util
import Test.Utility (approxArray, approxMatrix, approxVector)

import qualified Numeric.LAPACK.Matrix.Square as Square
import qualified Numeric.LAPACK.Matrix.Array as ArrMatrix
import qualified Numeric.LAPACK.Matrix.Shape as MatrixShape
import qualified Numeric.LAPACK.Matrix as Matrix
import qualified Numeric.LAPACK.Vector as Vector
import Numeric.LAPACK.Matrix (ShapeInt, (-*#), (##*#), (#*##), (#*|))
import Numeric.LAPACK.Vector (Vector)
import Numeric.LAPACK.Scalar (RealOf)

import qualified Numeric.Netlib.Class as Class



multiplySquare ::
   (Matrix.Power typ xl xu, Matrix.MultiplySquare typ xl xu,
    Matrix.SquareShape typ, Matrix.ToQuadratic typ,
    MatrixShape.PowerStrip lower, MatrixShape.PowerStrip upper,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   Matrix.Quadratic typ xl xu lower upper ShapeInt a -> Bool
multiplySquare a =
   Util.approxArray -- (Scalar.selectReal 1e-1 1e-5)
      (Matrix.toSquare $ Matrix.square a)
      (a #*## Matrix.toSquare a)

squareSquare ::
   (Matrix.Power typ xl xu, Matrix.MultiplySquare typ xl xu,
    Matrix.SquareShape typ, Matrix.ToQuadratic typ,
    MatrixShape.PowerStrip lower, MatrixShape.PowerStrip upper,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   Matrix.Quadratic typ xl xu lower upper ShapeInt a -> Bool
squareSquare a =
   Util.approxArray -- (Scalar.selectReal 1e-1 1e-5)
      (Matrix.toSquare $ Matrix.square a)
      (Square.square $ Matrix.toSquare a)

power ::
   (Matrix.Power typ xl xu, Matrix.MultiplySquare typ xl xu,
    Matrix.SquareShape typ, Matrix.ToQuadratic typ,
    MatrixShape.PowerStrip lower, MatrixShape.PowerStrip upper,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   Int -> Matrix.Quadratic typ xl xu lower upper ShapeInt a -> Bool
power n0 a =
   let n = fromIntegral n0
       b = Matrix.toSquare (Matrix.power (n+1) a)
       c = a #*## Matrix.toSquare (Matrix.power n a)
       normInf1 = Vector.normInf1 . ArrMatrix.toVector
   in Util.approxArrayTol (1e-6 * (normInf1 b + normInf1 c)) b c



multiplyIdentityVector ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a, Vector ShapeInt a) -> Bool
multiplyIdentityVector (eye,a) = approxVector a (eye #*| a)

multiplyIdentityFull ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a,
    Matrix.General ShapeInt ShapeInt a) ->
   Bool
multiplyIdentityFull (eye,a) = approxArray a (eye #*## a)

multiplyVector ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a, Vector ShapeInt a) -> Bool
multiplyVector (a,x) = approxVector (Matrix.toSquare a #*| x) (a #*| x)

multiplyFull ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a,
    Matrix.General ShapeInt ShapeInt a) ->
   Bool
multiplyFull (a,b) = approxArray (Matrix.toSquare a #*## b) (a #*## b)


multiplyVectorLeft ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (Vector ShapeInt a, ArrMatrix.Quadratic pack prop lo up ShapeInt a) -> Bool
multiplyVectorLeft (x,a) =
   approxVector (x -*# Matrix.toSquare a) (x -*# a)

multiplyVectorRight ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a, Vector ShapeInt a) -> Bool
multiplyVectorRight (a,x) =
   approxVector (Matrix.toSquare a #*| x) (a #*| x)


multiplyLeft ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (Matrix.General ShapeInt ShapeInt a,
    ArrMatrix.Quadratic pack prop lo up ShapeInt a) -> Bool
multiplyLeft (a,b) =
   approxMatrix 1e-5 (a ##*# Matrix.toSquare b) (a ##*# b)

multiplyRight ::
   (MatrixShape.Packing pack) =>
   (MatrixShape.Property prop, MatrixShape.PowerStrip lo, MatrixShape.PowerStrip up,
    Class.Floating a, RealOf a ~ ar, Class.Real ar) =>
   (ArrMatrix.Quadratic pack prop lo up ShapeInt a,
    Matrix.General ShapeInt ShapeInt a) -> Bool
multiplyRight (a,b) =
   approxArray (Matrix.toSquare a #*## b) (a #*## b)