{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
{-# LANGUAGE ConstraintKinds #-}
module Numeric.LAPACK.Matrix (
   Matrix.Matrix,
   Full,
   General, Tall, Wide, Square.Square, ArrMatrix.LiberalSquare,
   Matrix.Quadratic,
   Triangular.Triangular, Triangular.Upper, Triangular.Lower,
   Triangular.Symmetric,
   Triangular.Hermitian,
   Triangular.HermitianPosDef,
   Triangular.HermitianPosSemidef,
   Triangular.FlexHermitian,
   Banded.Diagonal, Banded.FlexDiagonal,
   Banded.RectangularDiagonal,
   Banded.Banded,
   BandedHermitian.BandedHermitian,
   PermMatrix.Permutation,

   Matrix.Format(Matrix.format), Matrix.FormatExtra, Matrix.formatWithLayout,
   Matrix.Layout(Matrix.layout), Matrix.LayoutExtra,

   ShapeInt, shapeInt,
   Matrix.Transpose, Matrix.TransposeExtra,
   Matrix.transpose, MatrixClass.adjoint,
   Matrix.height, Matrix.width, Matrix.extent,
   Matrix.Box, Matrix.BoxExtra, Matrix.indices,
   ArrMatrix.reshape,
   ArrMatrix.mapShape,
   caseTallWide,
   fromScalar, toScalar,
   fromList,
   MatrixClass.Unpack, MatrixClass.UnpackExtra,
   MatrixClass.toFull, MatrixClass.unpack,
   Matrix.mapExtent, fromFull,
   asGeneral, asTall, asWide,
   tallFromGeneral, wideFromGeneral,
   generalizeTall, generalizeWide,
   MatrixClass.MapSquareSize, MatrixClass.MapSize,
   MatrixClass.mapHeight, MatrixClass.mapWidth,
   MatrixClass.mapSquareSize,
   Quadratic.identity,
   Quadratic.diagonal,
   fromRowsNonEmpty,    fromRowArray,    fromRows,
   fromRowsNonEmptyContainer,    fromRowContainer,
   fromColumnsNonEmpty, fromColumnArray, fromColumns,
   fromColumnsNonEmptyContainer, fromColumnContainer,
   singleRow,   singleColumn,
   flattenRow,  flattenColumn,
   liftRow,     liftColumn,
   unliftRow,   unliftColumn,
   toRows, toColumns,
   toRowArray, toColumnArray,
   toRowContainer, toColumnContainer,
   takeRow, takeColumn,
   takeRows, takeColumns, takeEqually,
   dropRows, dropColumns, dropEqually,
   takeTop, takeBottom,
   takeLeft, takeRight,
   takeRowArray, takeColumnArray,
   swapRows, swapColumns,
   reverseRows, reverseColumns,
   fromRowMajor, toRowMajor,
   ArrMatrix.forceOrder, ArrMatrix.adaptOrder,
   Basic.OrderBias, leftBias, rightBias, contiguousBias,
   (|||), beside,
   (===), above,
   stack,

   (|*-),
   tensorProduct,
   outer,
   kronecker,
   sumRank1,

   map,
   MatrixClass.Complex, MatrixClass.conjugate,
   MatrixClass.fromReal, MatrixClass.toComplex,
   MatrixClass.SquareShape, MatrixClass.SquareShapeExtra, MatrixClass.toSquare,
   OmniMatrix.identityFromShape,
   MatrixClass.identityFrom,
   MatrixClass.takeDiagonal, MatrixClass.trace,

   RealOf,
   rowSums, columnSums,
   rowArgAbsMaximums, columnArgAbsMaximums,
   scaleRows, scaleColumns,
   scaleRowsReal, scaleColumnsReal,
   (\*#), (#*\),
   (\\#), (#/\),
   Full.multiply,
   Full.multiplyVector,

   Matrix.ToQuadratic(..),

   MatrixClass.Homogeneous, MatrixClass.HomogeneousExtra,
   ArrMatrix.zero, MatrixClass.zeroFrom, MatrixClass.negate,
   MatrixClass.Scale, MatrixClass.ScaleExtra, MatrixClass.scale,
   MatrixClass.scaleReal, MatrixClass.scaleRealReal,
   (MatrixClass..*#),
   MatrixClass.Additive, MatrixClass.Subtractive,
   MatrixClass.AdditiveExtra, MatrixClass.SubtractiveExtra,
   MatrixClass.add, MatrixClass.sub,
   (MatrixClass.#+#), (MatrixClass.#-#),
   Multiply.Multiply, (Multiply.#*#),
   Multiply.MultiplyVector, Multiply.MultiplyVectorExtra,
   (Multiply.#*|), (Multiply.-*#),
   Multiply.MultiplySquare, Multiply.MultiplySquareExtra, multiplySquare,
   Matrix.MultiplySame(multiplySame), Matrix.MultiplySameExtra,
   Multiply.Power, Multiply.PowerExtra, Multiply.square,
   Multiply.power, Multiply.powers, Multiply.powers1,
   (Multiply.##*#), (Multiply.#*##),
   Indexed.Indexed, (Indexed.#!),

   Divide.Determinant, Divide.DeterminantExtra, Divide.determinant,
   Divide.Solve, Divide.SolveExtra,
   Divide.solve, Divide.solveLeft, Divide.solveRight,
   (Divide.##/#), (Divide.#\##),
   Divide.solveVector, (Divide.-/#), (Divide.#\|),
   Divide.Inverse, Divide.InverseExtra, Divide.inverse,
   Mod.Transposition(..),
   ) where

import qualified Numeric.LAPACK.Matrix.Permutation as PermMatrix
import qualified Numeric.LAPACK.Matrix.Array.Mosaic as Triangular
import qualified Numeric.LAPACK.Matrix.BandedHermitian as BandedHermitian
import qualified Numeric.LAPACK.Matrix.Banded as Banded
import qualified Numeric.LAPACK.Matrix.Quadratic as Quadratic
import qualified Numeric.LAPACK.Matrix.Square as Square
import qualified Numeric.LAPACK.Matrix.Full as Full

import qualified Numeric.LAPACK.Matrix.Extent as Extent
import qualified Numeric.LAPACK.Matrix.Basic as Basic
import qualified Numeric.LAPACK.Matrix.Array.Basic as OmniMatrix
import qualified Numeric.LAPACK.Matrix.Array.Private as ArrMatrix
import qualified Numeric.LAPACK.Matrix.Type.Private as Matrix
import qualified Numeric.LAPACK.Matrix.Plain as Plain
import qualified Numeric.LAPACK.Matrix.Modifier as Mod
import qualified Numeric.LAPACK.Matrix.Divide as Divide
import qualified Numeric.LAPACK.Matrix.Multiply as Multiply
import qualified Numeric.LAPACK.Matrix.Shape.Omni as Omni
import qualified Numeric.LAPACK.Matrix.Indexed as Indexed
import qualified Numeric.LAPACK.Matrix.Class as MatrixClass
import qualified Numeric.LAPACK.Matrix.Private as MatrixPriv
import qualified Numeric.LAPACK.Vector as Vector
import qualified Numeric.LAPACK.Shape as ExtShape
import Numeric.LAPACK.Matrix.Layout.Private (Order)
import Numeric.LAPACK.Matrix.Array.Private (Full, General, Tall, Wide)
import Numeric.LAPACK.Matrix.Private (ShapeInt, shapeInt)
import Numeric.LAPACK.Vector (Vector)
import Numeric.LAPACK.Scalar (RealOf)

import qualified Numeric.Netlib.Class as Class

import qualified Data.Array.Comfort.Storable as Array
import qualified Data.Array.Comfort.Boxed as BoxedArray
import qualified Data.Array.Comfort.Container as Container
import qualified Data.Array.Comfort.Shape as Shape
import Data.Array.Comfort.Storable.Unchecked (Array, (!))
import Data.Array.Comfort.Shape ((::+))

import Foreign.Storable (Storable)

import qualified Data.NonEmpty as NonEmpty
import qualified Data.Either.HT as EitherHT
import Data.Function.HT (Id)

import Prelude hiding (map)


fromFull ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz) =>
   Full meas vert horiz height width a -> General height width a
fromFull = ArrMatrix.lift1 Plain.fromFull

tallFromGeneral ::
   (Shape.C height, Shape.C width, Storable a) =>
   General height width a -> Tall height width a
tallFromGeneral = ArrMatrix.lift1 Plain.tallFromGeneral

wideFromGeneral ::
   (Shape.C height, Shape.C width, Storable a) =>
   General height width a -> Wide height width a
wideFromGeneral = ArrMatrix.lift1 Plain.wideFromGeneral

generalizeTall ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz) =>
   Full meas vert Extent.Small height width a ->
   Full Extent.Size vert horiz height width a
generalizeTall = Full.mapExtent Extent.generalizeTall

generalizeWide ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz) =>
   Full meas Extent.Small horiz height width a ->
   Full Extent.Size vert horiz height width a
generalizeWide = Full.mapExtent Extent.generalizeWide


asGeneral :: Id (General height width a)
asGeneral = id

asTall :: Id (Tall height width a)
asTall = id

asWide :: Id (Wide height width a)
asWide = id



fromScalar :: (Storable a) => a -> General () () a
fromScalar = fromFull . Square.fromScalar

toScalar :: (Storable a) => General () () a -> a
toScalar a =
   either id id (MatrixPriv.revealOrder (ArrMatrix.toVector a)) ! ((),())

fromList ::
   (Shape.C height, Shape.C width, Storable a) =>
   height -> width -> [a] -> General height width a
fromList height width = ArrMatrix.lift0 . Plain.fromList height width


{- |
Square matrices will be classified as 'Tall'.
-}
caseTallWide ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width) =>
   Full meas vert horiz height width a ->
   Either (Tall height width a) (Wide height width a)
caseTallWide =
   EitherHT.mapBoth ArrMatrix.lift0 ArrMatrix.lift0 .
   Basic.caseTallWide . ArrMatrix.toVector


singleRow :: Order -> Vector width a -> General () width a
singleRow order = ArrMatrix.lift0 . Basic.singleRow order

singleColumn :: Order -> Vector height a -> General height () a
singleColumn order = ArrMatrix.lift0 . Basic.singleColumn order

flattenRow :: General () width a -> Vector width a
flattenRow = Basic.flattenRow . ArrMatrix.toVector

flattenColumn :: General height () a -> Vector height a
flattenColumn = Basic.flattenColumn . ArrMatrix.toVector

liftRow ::
   Order ->
   (Vector height0 a -> Vector height1 b) ->
   General () height0 a -> General () height1 b
liftRow order = ArrMatrix.lift1 . Basic.liftRow order

liftColumn ::
   Order ->
   (Vector height0 a -> Vector height1 b) ->
   General height0 () a -> General height1 () b
liftColumn order = ArrMatrix.lift1 . Basic.liftColumn order

unliftRow ::
   Order ->
   (General () height0 a -> General () height1 b) ->
   Vector height0 a -> Vector height1 b
unliftRow order = Basic.unliftRow order .  ArrMatrix.unlift1

unliftColumn ::
   Order ->
   (General height0 () a -> General height1 () b) ->
   Vector height0 a -> Vector height1 b
unliftColumn order = Basic.unliftColumn order . ArrMatrix.unlift1


fromRowsNonEmpty ::
   (Shape.C width, Eq width, Storable a) =>
   NonEmpty.T [] (Vector width a) -> General ShapeInt width a
fromRowsNonEmpty = ArrMatrix.lift0 . Plain.fromRowsNonEmpty

fromRowArray ::
   (Shape.C height, Shape.C width, Eq width, Storable a) =>
   width -> BoxedArray.Array height (Vector width a) -> General height width a
fromRowArray width = ArrMatrix.lift0 . Plain.fromRowArray width

-- ToDo: generalize to a new NonEmpty.Head class
fromRowsNonEmptyContainer ::
   (f ~ NonEmpty.T g, Container.C g,
    Shape.C width, Eq width, Storable a) =>
   f (Vector width a) -> General (Container.Shape f) width a
fromRowsNonEmptyContainer = ArrMatrix.lift0 . Plain.fromRowsNonEmptyContainer

fromRowContainer ::
   (Container.C f, Shape.C width, Eq width, Storable a) =>
   width -> f (Vector width a) -> General (Container.Shape f) width a
fromRowContainer width = ArrMatrix.lift0 . Plain.fromRowContainer width

fromRows ::
   (Shape.C width, Eq width, Storable a) =>
   width -> [Vector width a] -> General ShapeInt width a
fromRows width = ArrMatrix.lift0 . Plain.fromRows width


fromColumnsNonEmpty ::
   (Shape.C height, Eq height, Storable a) =>
   NonEmpty.T [] (Vector height a) -> General height ShapeInt a
fromColumnsNonEmpty = ArrMatrix.lift0 . Plain.fromColumnsNonEmpty

fromColumnArray ::
   (Shape.C height, Eq height, Shape.C width, Storable a) =>
   height -> BoxedArray.Array width (Vector height a) -> General height width a
fromColumnArray height = ArrMatrix.lift0 . Plain.fromColumnArray height

fromColumnsNonEmptyContainer ::
   (f ~ NonEmpty.T g, Container.C g,
    Shape.C height, Eq height, Storable a) =>
   f (Vector height a) -> General height (Container.Shape f) a
fromColumnsNonEmptyContainer =
   ArrMatrix.lift0 . Plain.fromColumnsNonEmptyContainer

fromColumnContainer ::
   (Container.C f, Shape.C height, Eq height, Storable a) =>
   height -> f (Vector height a) -> General height (Container.Shape f) a
fromColumnContainer height = ArrMatrix.lift0 . Plain.fromColumnContainer height

fromColumns ::
   (Shape.C height, Eq height, Storable a) =>
   height -> [Vector height a] -> General height ShapeInt a
fromColumns height = ArrMatrix.lift0 . Plain.fromColumns height


toRows ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a -> [Vector width a]
toRows = Plain.toRows . ArrMatrix.toVector

toColumns ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a -> [Vector height a]
toColumns = Plain.toColumns . ArrMatrix.toVector

toRowArray ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   BoxedArray.Array height (Vector width a)
toRowArray = Plain.toRowArray . ArrMatrix.toVector

toColumnArray ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   BoxedArray.Array width (Vector height a)
toColumnArray = Plain.toColumnArray . ArrMatrix.toVector

toRowContainer ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Container.C f, Shape.C width, Class.Floating a) =>
   Full meas vert horiz (Container.Shape f) width a -> f (Vector width a)
toRowContainer = Plain.toRowContainer . ArrMatrix.toVector

toColumnContainer ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Container.C f, Class.Floating a) =>
   Full meas vert horiz height (Container.Shape f) a -> f (Vector height a)
toColumnContainer = Plain.toColumnContainer . ArrMatrix.toVector



{-
The parameter order is swapped with respect to 'takeRowArray'
but it is the order that is used most oftenly.
-}
takeRow ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.Indexed height, Shape.C width, Shape.Index height ~ ix,
    Class.Floating a) =>
   Full meas vert horiz height width a -> ix -> Vector width a
takeRow = Plain.takeRow . ArrMatrix.toVector

takeColumn ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.Indexed width, Shape.Index width ~ ix,
    Class.Floating a) =>
   Full meas vert horiz height width a -> ix -> Vector height a
takeColumn = Plain.takeColumn . ArrMatrix.toVector


takeTop ::
   (Extent.C vert, Shape.C height0, Shape.C height1, Shape.C width,
    Class.Floating a) =>
   Full Extent.Size vert Extent.Big (height0::+height1) width a ->
   Full Extent.Size vert Extent.Big height0 width a
takeTop = ArrMatrix.lift1 Basic.takeTop

takeBottom ::
   (Extent.C vert, Shape.C height0, Shape.C height1, Shape.C width,
    Class.Floating a) =>
   Full Extent.Size vert Extent.Big (height0::+height1) width a ->
   Full Extent.Size vert Extent.Big height1 width a
takeBottom = ArrMatrix.lift1 Basic.takeBottom

takeLeft ::
   (Extent.C vert, Shape.C height, Shape.C width0, Shape.C width1,
    Class.Floating a) =>
   Full Extent.Size Extent.Big vert height (width0::+width1) a ->
   Full Extent.Size Extent.Big vert height width0 a
takeLeft = ArrMatrix.lift1 Basic.takeLeft

takeRight ::
   (Extent.C vert, Shape.C height, Shape.C width0, Shape.C width1,
    Class.Floating a) =>
   Full Extent.Size Extent.Big vert height (width0::+width1) a ->
   Full Extent.Size Extent.Big vert height width1 a
takeRight = ArrMatrix.lift1 Basic.takeRight

takeRows, dropRows ::
   (Extent.C vert, Shape.C width, Class.Floating a) =>
   Int ->
   Full Extent.Size vert Extent.Big ShapeInt width a ->
   Full Extent.Size vert Extent.Big ShapeInt width a
takeRows = ArrMatrix.lift1 . Basic.takeRows
dropRows = ArrMatrix.lift1 . Basic.dropRows

takeColumns, dropColumns ::
   (Extent.C horiz, Shape.C height, Class.Floating a) =>
   Int ->
   Full Extent.Size Extent.Big horiz height ShapeInt a ->
   Full Extent.Size Extent.Big horiz height ShapeInt a
takeColumns = ArrMatrix.lift1 . Basic.takeColumns
dropColumns = ArrMatrix.lift1 . Basic.dropColumns


{- |
Take a left-top aligned square or as much as possible of it.
The advantange of this function is that it maintains the matrix size relation,
e.g. Square remains Square, Tall remains Tall.
-}
takeEqually ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz, Class.Floating a) =>
   Int ->
   Full meas vert horiz ShapeInt ShapeInt a ->
   Full meas vert horiz ShapeInt ShapeInt a
takeEqually = ArrMatrix.lift1 . Plain.takeEqually

{- |
Drop the same number of top-most rows and left-most columns.
The advantange of this function is that it maintains the matrix size relation,
e.g. Square remains Square, Tall remains Tall.
-}
dropEqually ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz, Class.Floating a) =>
   Int ->
   Full meas vert horiz ShapeInt ShapeInt a ->
   Full meas vert horiz ShapeInt ShapeInt a
dropEqually = ArrMatrix.lift1 . Plain.dropEqually


swapRows ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.Indexed height, Shape.C width, Class.Floating a) =>
   Shape.Index height -> Shape.Index height ->
   Full meas vert horiz height width a -> Full meas vert horiz height width a
swapRows i j = ArrMatrix.lift1 $ Plain.swapRows i j

swapColumns ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.Indexed width, Class.Floating a) =>
   Shape.Index width -> Shape.Index width ->
   Full meas vert horiz height width a -> Full meas vert horiz height width a
swapColumns i j =  ArrMatrix.lift1 $ Plain.swapColumns i j


reverseRows ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    ExtShape.Permutable height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
reverseRows = ArrMatrix.lift1 Plain.reverseRows

reverseColumns ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, ExtShape.Permutable width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
reverseColumns = ArrMatrix.lift1 Plain.reverseColumns


{- |
The function is optimized for blocks of consecutive rows.
For scattered rows in column major order
the function has quite ugly memory access patterns.
-}
takeRowArray ::
   (Shape.Indexed height, Shape.C width, Shape.C sh, Class.Floating a) =>
   BoxedArray.Array sh (Shape.Index height) ->
   General height width a -> General sh width a
takeRowArray = ArrMatrix.lift1 . Plain.takeRowArray

takeColumnArray ::
   (Shape.C height, Shape.Indexed width, Shape.C sh, Class.Floating a) =>
   BoxedArray.Array sh (Shape.Index width) ->
   General height width a -> General height sh a
takeColumnArray = ArrMatrix.lift1 . Plain.takeColumnArray



fromRowMajor ::
   (Shape.C height, Shape.C width, Class.Floating a) =>
   Array (height,width) a -> General height width a
fromRowMajor = ArrMatrix.lift0 . Plain.fromRowMajor

toRowMajor ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a -> Array (height,width) a
toRowMajor = Plain.toRowMajor . ArrMatrix.toVector


infixr 3 |||
infixr 2 ===

(|||) ::
   (Extent.C vertA, Extent.C vertB, Extent.C vertC,
    Extent.Append vertA vertB ~ vertC,
    Shape.C height, Eq height, Shape.C widthA, Shape.C widthB,
    Class.Floating a) =>
   Full Extent.Size vertA Extent.Big height widthA a ->
   Full Extent.Size vertB Extent.Big height widthB a ->
   Full Extent.Size vertC Extent.Big height (widthA::+widthB) a
(|||) = beside rightBias Extent.appendAny

(===) ::
   (Extent.C horizA, Extent.C horizB, Extent.C horizC,
    Extent.Append horizA horizB ~ horizC,
    Shape.C width, Eq width, Shape.C heightA, Shape.C heightB,
    Class.Floating a) =>
   Full Extent.Size Extent.Big horizA heightA width a ->
   Full Extent.Size Extent.Big horizB heightB width a ->
   Full Extent.Size Extent.Big horizC (heightA::+heightB) width a
(===) = above rightBias Extent.appendAny

beside ::
   (Extent.C vertA, Extent.C vertB, Extent.C vertC,
    Shape.C height, Eq height, Shape.C widthA, Shape.C widthB,
    Class.Floating a) =>
   Basic.OrderBias ->
   Extent.AppendMode vertA vertB vertC height widthA widthB ->
   Full Extent.Size vertA Extent.Big height widthA a ->
   Full Extent.Size vertB Extent.Big height widthB a ->
   Full Extent.Size vertC Extent.Big height (widthA::+widthB) a
beside orderBias = ArrMatrix.lift2 . Basic.beside orderBias

above ::
   (Extent.C horizA, Extent.C horizB, Extent.C horizC,
    Shape.C width, Eq width, Shape.C heightA, Shape.C heightB,
    Class.Floating a) =>
   Basic.OrderBias ->
   Extent.AppendMode horizA horizB horizC width heightA heightB ->
   Full Extent.Size Extent.Big horizA heightA width a ->
   Full Extent.Size Extent.Big horizB heightB width a ->
   Full Extent.Size Extent.Big horizC (heightA::+heightB) width a
above orderBias = ArrMatrix.lift2 . Basic.above orderBias

{- |
Use the element order of the first operand.
-}
leftBias :: Basic.OrderBias
leftBias = Basic.LeftBias

{- |
Use the element order of the second operand.
-}
rightBias :: Basic.OrderBias
rightBias = Basic.RightBias

{- |
Choose element order such that, if possible,
one part can be copied as one block.
For 'above' this means that 'RowMajor' is chosen
whenever at least one operand is 'RowMajor'
and 'ColumnMajor' is chosen when both operands are 'ColumnMajor'.
-}
contiguousBias :: Basic.OrderBias
contiguousBias = Basic.ContiguousBias


stack ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C heightA, Eq heightA, Shape.C heightB, Eq heightB,
    Shape.C widthA, Eq widthA, Shape.C widthB, Eq widthB, Class.Floating a) =>
   Full meas vert horiz heightA widthA a -> General heightA widthB a ->
   General heightB widthA a -> Full meas vert horiz heightB widthB a ->
   Full meas vert horiz (heightA::+heightB) (widthA::+widthB) a
stack = ArrMatrix.lift4 Basic.stack



rowSums ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a -> Vector height a
rowSums = Plain.rowSums . ArrMatrix.toVector

columnSums ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Full meas vert horiz height width a -> Vector width a
columnSums = Plain.columnSums . ArrMatrix.toVector


rowArgAbsMaximums ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.InvIndexed width, Shape.Index width ~ ix, Storable ix,
    Class.Floating a) =>
   Full meas vert horiz height width a -> (Vector height ix, Vector height a)
rowArgAbsMaximums = Plain.rowArgAbsMaximums . ArrMatrix.toVector

columnArgAbsMaximums ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.InvIndexed height, Shape.C width,
    Shape.Index height ~ ix, Storable ix,
    Class.Floating a) =>
   Full meas vert horiz height width a -> (Vector width ix, Vector width a)
columnArgAbsMaximums = Plain.columnArgAbsMaximums . ArrMatrix.toVector


map ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Storable a, Storable b) =>
   (a -> b) ->
   Full meas vert horiz height width a -> Full meas vert horiz height width b
map = ArrMatrix.lift1 . Array.map


infixl 7 |*-

(|*-) ::
   (Shape.C height, Eq height, Shape.C width, Eq width, Class.Floating a) =>
   Vector height a -> Vector width a -> General height width a
x|*-y = ArrMatrix.lift0 $ x Plain.|*- y

{- |
> tensorProduct order x y = singleColumn order x #*# singleRow order y
-}
tensorProduct ::
   (Shape.C height, Eq height, Shape.C width, Eq width, Class.Floating a) =>
   Order -> Vector height a -> Vector width a -> General height width a
tensorProduct order x y = ArrMatrix.lift0 $ Plain.tensorProduct order x y

{- |
> outer order x y = tensorProduct order x (Vector.conjugate y)
-}
outer ::
   (Shape.C height, Eq height, Shape.C width, Eq width, Class.Floating a) =>
   Order -> Vector height a -> Vector width a -> General height width a
outer order x y = ArrMatrix.lift0 $ Plain.outer order x y

kronecker ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C heightA, Shape.C widthA, Shape.C heightB, Shape.C widthB,
    Class.Floating a) =>
   Full meas vert horiz heightA widthA a ->
   Full meas vert horiz heightB widthB a ->
   Full meas vert horiz (heightA,heightB) (widthA,widthB) a
kronecker = ArrMatrix.lift2 Plain.kronecker

sumRank1 ::
   (Shape.C height, Eq height, Shape.C width, Eq width, Class.Floating a) =>
   (height,width) ->
   [(a, (Vector height a, Vector width a))] -> General height width a
sumRank1 dims = ArrMatrix.lift0 . Plain.sumRank1 dims


infixl 7 #*\, #/\
infixr 7 \*#, \\#

scaleRows, (\*#) ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Class.Floating a) =>
   Vector height a ->
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
scaleRows = ArrMatrix.lift1 . Basic.scaleRows
(\*#) = scaleRows

scaleColumns ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Eq width, Class.Floating a) =>
   Vector width a ->
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
scaleColumns = ArrMatrix.lift1 . Basic.scaleColumns

(#*\) ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Eq width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   Vector width a ->
   Full meas vert horiz height width a
(#*\) = flip scaleColumns

scaleRowsReal ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Class.Floating a) =>
   Vector height (RealOf a) ->
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
scaleRowsReal = ArrMatrix.lift1 . Basic.scaleRowsReal

scaleColumnsReal ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Eq width, Class.Floating a) =>
   Vector width (RealOf a) ->
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
scaleColumnsReal = ArrMatrix.lift1 . Basic.scaleColumnsReal

(\\#) ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Class.Floating a) =>
   Vector height a ->
   Full meas vert horiz height width a ->
   Full meas vert horiz height width a
(\\#) = scaleRows . Vector.recip

(#/\) ::
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Eq width, Class.Floating a) =>
   Full meas vert horiz height width a ->
   Vector width a ->
   Full meas vert horiz height width a
(#/\) a x = scaleColumns (Vector.recip x) a


multiplySquare ::
   (Multiply.MultiplySquare typ) =>
   (Multiply.MultiplySquareExtra typ xl, Multiply.MultiplySquareExtra typ xu) =>
   (Omni.Strip lower, Omni.Strip upper) =>
   (Extent.Measure meas, Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Class.Floating a) =>
   Mod.Transposition ->
   Matrix.Quadratic typ xl xu lower upper height a ->
   Full meas vert horiz height width a -> Full meas vert horiz height width a
multiplySquare = Multiply.transposableSquare