{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
module Numeric.LAPACK.Matrix.Square (
   Square,
   size,
   mapSize,
   toFull,
   toGeneral,
   fromGeneral,
   fromScalar,
   toScalar,
   fromList,
   autoFromList,

   transpose,
   adjoint,

   identity,
   identityFrom,
   identityFromWidth,
   identityFromHeight,
   diagonal,
   takeDiagonal,
   trace,

   stack, (|=|),

   multiply,
   square,
   power,
   congruence,
   congruenceAdjoint,

   solve,
   inverse,
   determinant,

   eigenvalues,
   schur,
   schurComplex,
   eigensystem,
   ComplexOf,
   ) where

import qualified Numeric.LAPACK.Matrix.Triangular as Triangular
import qualified Numeric.LAPACK.Matrix.Square.Eigen as Eigen
import qualified Numeric.LAPACK.Matrix.Square.Linear as Linear
import qualified Numeric.LAPACK.Matrix.Square.Basic as Basic
import qualified Numeric.LAPACK.Matrix.Basic as FullBasic

import qualified Numeric.LAPACK.Matrix.Array as ArrMatrix
import qualified Numeric.LAPACK.Matrix.Shape.Private as MatrixShape
import qualified Numeric.LAPACK.Matrix.Extent as Extent
import Numeric.LAPACK.Matrix.Array (Full, General, Square)
import Numeric.LAPACK.Matrix.Private (ShapeInt)
import Numeric.LAPACK.Vector (Vector)
import Numeric.LAPACK.Scalar (ComplexOf)

import qualified Numeric.Netlib.Class as Class

import qualified Data.Array.Comfort.Shape as Shape
import Data.Array.Comfort.Shape ((:+:))

import Foreign.Storable (Storable)

import Data.Tuple.HT (mapPair, mapSnd, mapTriple)
import Data.Complex (Complex)


size :: Square sh a -> sh
size :: Square sh a -> sh
size = Full Small Small sh sh -> sh
forall vert horiz height width.
(C vert, C horiz) =>
Full vert horiz height width -> height
MatrixShape.fullHeight (Full Small Small sh sh -> sh)
-> (Square sh a -> Full Small Small sh sh) -> Square sh a -> sh
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Full Small Small sh sh
forall sh a. ArrayMatrix sh a -> sh
ArrMatrix.shape

mapSize :: (sh0 -> sh1) -> Square sh0 a -> Square sh1 a
mapSize :: (sh0 -> sh1) -> Square sh0 a -> Square sh1 a
mapSize = (Array (Square sh0) a -> Array (Square sh1) a)
-> Square sh0 a -> Square sh1 a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Square sh0) a -> Array (Square sh1) a)
 -> Square sh0 a -> Square sh1 a)
-> ((sh0 -> sh1) -> Array (Square sh0) a -> Array (Square sh1) a)
-> (sh0 -> sh1)
-> Square sh0 a
-> Square sh1 a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (sh0 -> sh1) -> Array (Square sh0) a -> Array (Square sh1) a
forall sh0 sh1 a. (sh0 -> sh1) -> Square sh0 a -> Square sh1 a
Basic.mapSize

toGeneral :: Square sh a -> General sh sh a
toGeneral :: Square sh a -> General sh sh a
toGeneral = Square sh a -> General sh sh a
forall vert horiz sh a.
(C vert, C horiz) =>
Square sh a -> Full vert horiz sh sh a
toFull

toFull ::
   (Extent.C vert, Extent.C horiz) => Square sh a -> Full vert horiz sh sh a
toFull :: Square sh a -> Full vert horiz sh sh a
toFull = (Array (Square sh) a -> Array (Full vert horiz sh sh) a)
-> Square sh a -> Full vert horiz sh sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Full vert horiz sh sh) a
forall vert horiz sh a.
(C vert, C horiz) =>
Square sh a -> Full vert horiz sh sh a
Basic.toFull

fromGeneral :: (Eq sh) => General sh sh a -> Square sh a
fromGeneral :: General sh sh a -> Square sh a
fromGeneral = (Array (General sh sh) a -> Array (Square sh) a)
-> General sh sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (General sh sh) a -> Array (Square sh) a
forall sh a. Eq sh => General sh sh a -> Square sh a
Basic.fromGeneral


fromScalar :: (Storable a) => a -> Square () a
fromScalar :: a -> Square () a
fromScalar = Array (Square ()) a -> Square () a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0 (Array (Square ()) a -> Square () a)
-> (a -> Array (Square ()) a) -> a -> Square () a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> Array (Square ()) a
forall a. Storable a => a -> Square () a
Basic.fromScalar

toScalar :: (Storable a) => Square () a -> a
toScalar :: Square () a -> a
toScalar = Square () a -> a
forall a. Storable a => Square () a -> a
Basic.toScalar (Square () a -> a)
-> (Square () a -> Square () a) -> Square () a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square () a -> Square () a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector

fromList :: (Shape.C sh, Storable a) => sh -> [a] -> Square sh a
fromList :: sh -> [a] -> Square sh a
fromList sh
sh = Array (Square sh) a -> Square sh a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0 (Array (Square sh) a -> Square sh a)
-> ([a] -> Array (Square sh) a) -> [a] -> Square sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. sh -> [a] -> Array (Square sh) a
forall sh a. (C sh, Storable a) => sh -> [a] -> Square sh a
Basic.fromList sh
sh

autoFromList :: (Storable a) => [a] -> Square ShapeInt a
autoFromList :: [a] -> Square ShapeInt a
autoFromList = Array (Square ShapeInt) a -> Square ShapeInt a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0 (Array (Square ShapeInt) a -> Square ShapeInt a)
-> ([a] -> Array (Square ShapeInt) a) -> [a] -> Square ShapeInt a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [a] -> Array (Square ShapeInt) a
forall a. Storable a => [a] -> Square ShapeInt a
Basic.autoFromList

transpose :: Square sh a -> Square sh a
transpose :: Square sh a -> Square sh a
transpose = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Square sh) a
forall sh a. Square sh a -> Square sh a
Basic.transpose

{- |
conjugate transpose
-}
adjoint :: (Shape.C sh, Class.Floating a) => Square sh a -> Square sh a
adjoint :: Square sh a -> Square sh a
adjoint = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Square sh) a
forall sh a. (C sh, Floating a) => Square sh a -> Square sh a
Basic.adjoint

identity :: (Shape.C sh, Class.Floating a) => sh -> Square sh a
identity :: sh -> Square sh a
identity = Array (Square sh) a -> Square sh a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0 (Array (Square sh) a -> Square sh a)
-> (sh -> Array (Square sh) a) -> sh -> Square sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. sh -> Array (Square sh) a
forall sh a. (C sh, Floating a) => sh -> Square sh a
Basic.identity

identityFrom :: (Shape.C sh, Class.Floating a) => Square sh a -> Square sh a
identityFrom :: Square sh a -> Square sh a
identityFrom = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Square sh) a
forall sh a. (C sh, Floating a) => Square sh a -> Square sh a
Basic.identityFrom

identityFromWidth ::
   (Shape.C height, Shape.C width, Class.Floating a) =>
   General height width a -> Square width a
identityFromWidth :: General height width a -> Square width a
identityFromWidth = (Array (General height width) a -> Array (Square width) a)
-> General height width a -> Square width a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (General height width) a -> Array (Square width) a
forall height width a.
(C height, C width, Floating a) =>
General height width a -> Square width a
Basic.identityFromWidth

identityFromHeight ::
   (Shape.C height, Shape.C width, Class.Floating a) =>
   General height width a -> Square height a
identityFromHeight :: General height width a -> Square height a
identityFromHeight = (Array (General height width) a -> Array (Square height) a)
-> General height width a -> Square height a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (General height width) a -> Array (Square height) a
forall height width a.
(C height, C width, Floating a) =>
General height width a -> Square height a
Basic.identityFromHeight

diagonal :: (Shape.C sh, Class.Floating a) => Vector sh a -> Square sh a
diagonal :: Vector sh a -> Square sh a
diagonal = Array (Square sh) a -> Square sh a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0 (Array (Square sh) a -> Square sh a)
-> (Vector sh a -> Array (Square sh) a)
-> Vector sh a
-> Square sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector sh a -> Array (Square sh) a
forall sh a. (C sh, Floating a) => Vector sh a -> Square sh a
Basic.diagonal

takeDiagonal :: (Shape.C sh, Class.Floating a) => Square sh a -> Vector sh a
takeDiagonal :: Square sh a -> Vector sh a
takeDiagonal = Square sh a -> Vector sh a
forall sh a. (C sh, Floating a) => Square sh a -> Vector sh a
Basic.takeDiagonal (Square sh a -> Vector sh a)
-> (Square sh a -> Square sh a) -> Square sh a -> Vector sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Square sh a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector

trace :: (Shape.C sh, Class.Floating a) => Square sh a -> a
trace :: Square sh a -> a
trace = Square sh a -> a
forall sh a. (C sh, Floating a) => Square sh a -> a
Basic.trace (Square sh a -> a)
-> (Square sh a -> Square sh a) -> Square sh a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Square sh a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector

infix 3 |=|

(|=|) ::
   (Extent.C vert, Extent.C horiz,
    Shape.C sizeA, Eq sizeA, Shape.C sizeB, Eq sizeB, Class.Floating a) =>
   (Square sizeA a, Full vert horiz sizeA sizeB a) ->
   (Full horiz vert sizeB sizeA a, Square sizeB a) ->
   Square (sizeA:+:sizeB) a
(Square sizeA a
a,Full vert horiz sizeA sizeB a
b) |=| :: (Square sizeA a, Full vert horiz sizeA sizeB a)
-> (Full horiz vert sizeB sizeA a, Square sizeB a)
-> Square (sizeA :+: sizeB) a
|=| (Full horiz vert sizeB sizeA a
c,Square sizeB a
d)  =  Square sizeA a
-> Full vert horiz sizeA sizeB a
-> Full horiz vert sizeB sizeA a
-> Square sizeB a
-> Square (sizeA :+: sizeB) a
forall vert horiz sizeA sizeB a.
(C vert, C horiz, C sizeA, Eq sizeA, C sizeB, Eq sizeB,
 Floating a) =>
Square sizeA a
-> Full vert horiz sizeA sizeB a
-> Full horiz vert sizeB sizeA a
-> Square sizeB a
-> Square (sizeA :+: sizeB) a
stack Square sizeA a
a Full vert horiz sizeA sizeB a
b Full horiz vert sizeB sizeA a
c Square sizeB a
d

stack ::
   (Extent.C vert, Extent.C horiz,
    Shape.C sizeA, Eq sizeA, Shape.C sizeB, Eq sizeB, Class.Floating a) =>
   Square sizeA a -> Full vert horiz sizeA sizeB a ->
   Full horiz vert sizeB sizeA a -> Square sizeB a ->
   Square (sizeA:+:sizeB) a
stack :: Square sizeA a
-> Full vert horiz sizeA sizeB a
-> Full horiz vert sizeB sizeA a
-> Square sizeB a
-> Square (sizeA :+: sizeB) a
stack = (Array (Square sizeA) a
 -> Array (Full vert horiz sizeA sizeB) a
 -> Array (Full horiz vert sizeB sizeA) a
 -> Array (Square sizeB) a
 -> Array (Square (sizeA :+: sizeB)) a)
-> Square sizeA a
-> Full vert horiz sizeA sizeB a
-> Full horiz vert sizeB sizeA a
-> Square sizeB a
-> Square (sizeA :+: sizeB) a
forall shA a shB b shC c shD d shE e.
(Array shA a
 -> Array shB b -> Array shC c -> Array shD d -> Array shE e)
-> ArrayMatrix shA a
-> ArrayMatrix shB b
-> ArrayMatrix shC c
-> ArrayMatrix shD d
-> ArrayMatrix shE e
ArrMatrix.lift4 Array (Square sizeA) a
-> Array (Full vert horiz sizeA sizeB) a
-> Array (Full horiz vert sizeB sizeA) a
-> Array (Square sizeB) a
-> Array (Square (sizeA :+: sizeB)) a
forall vert horiz sizeA sizeB a.
(C vert, C horiz, C sizeA, Eq sizeA, C sizeB, Eq sizeB,
 Floating a) =>
Square sizeA a
-> Full vert horiz sizeA sizeB a
-> Full horiz vert sizeB sizeA a
-> Square sizeB a
-> Square (sizeA :+: sizeB) a
Basic.stack

multiply ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   Square sh a -> Square sh a -> Square sh a
multiply :: Square sh a -> Square sh a -> Square sh a
multiply = (Array (Full Small Small sh sh) a
 -> Array (Full Small Small sh sh) a
 -> Array (Full Small Small sh sh) a)
-> Square sh a -> Square sh a -> Square sh a
forall shA a shB b shC c.
(Array shA a -> Array shB b -> Array shC c)
-> ArrayMatrix shA a -> ArrayMatrix shB b -> ArrayMatrix shC c
ArrMatrix.lift2 Array (Full Small Small sh sh) a
-> Array (Full Small Small sh sh) a
-> Array (Full Small Small sh sh) a
forall vert horiz height fuse width a.
(C vert, C horiz, C height, C fuse, Eq fuse, C width,
 Floating a) =>
Full vert horiz height fuse a
-> Full vert horiz fuse width a -> Full vert horiz height width a
FullBasic.multiply

square :: (Shape.C sh, Class.Floating a) => Square sh a -> Square sh a
square :: Square sh a -> Square sh a
square = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Square sh) a
forall sh a. (C sh, Floating a) => Square sh a -> Square sh a
Basic.square

power ::
   (Shape.C sh, Class.Floating a) =>
   Integer -> Square sh a -> Square sh a
power :: Integer -> Square sh a -> Square sh a
power = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Square sh) a -> Array (Square sh) a)
 -> Square sh a -> Square sh a)
-> (Integer -> Array (Square sh) a -> Array (Square sh) a)
-> Integer
-> Square sh a
-> Square sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Integer -> Array (Square sh) a -> Array (Square sh) a
forall sh a.
(C sh, Floating a) =>
Integer -> Square sh a -> Square sh a
Basic.power

{- |
congruence B A = A^H * B * A


The meaning and order of matrix factors of these functions is consistent:

* 'Numeric.LAPACK.Matrix.Square.congruence'
* 'Numeric.LAPACK.Matrix.Hermitian.gramian'
* 'Numeric.LAPACK.Matrix.Hermitian.anticommutator'
* 'Numeric.LAPACK.Matrix.Hermitian.congruence'
* 'Numeric.LAPACK.Matrix.Hermitian.congruenceDiagonal'
-}
congruence ::
   (Shape.C height, Eq height, Shape.C width, Class.Floating a) =>
   Square height a -> General height width a -> Square width a
congruence :: Square height a -> General height width a -> Square width a
congruence = (Array (Square height) a
 -> Array (General height width) a -> Array (Square width) a)
-> Square height a -> General height width a -> Square width a
forall shA a shB b shC c.
(Array shA a -> Array shB b -> Array shC c)
-> ArrayMatrix shA a -> ArrayMatrix shB b -> ArrayMatrix shC c
ArrMatrix.lift2 Array (Square height) a
-> Array (General height width) a -> Array (Square width) a
forall height width a.
(C height, Eq height, C width, Floating a) =>
Square height a -> General height width a -> Square width a
Basic.congruence

{- |
congruenceAdjoint A B = A * B * A^H
-}
congruenceAdjoint ::
   (Shape.C height, Shape.C width, Eq width, Class.Floating a) =>
   General height width a -> Square width a -> Square height a
congruenceAdjoint :: General height width a -> Square width a -> Square height a
congruenceAdjoint = (Array (General height width) a
 -> Array (Square width) a -> Array (Square height) a)
-> General height width a -> Square width a -> Square height a
forall shA a shB b shC c.
(Array shA a -> Array shB b -> Array shC c)
-> ArrayMatrix shA a -> ArrayMatrix shB b -> ArrayMatrix shC c
ArrMatrix.lift2 Array (General height width) a
-> Array (Square width) a -> Array (Square height) a
forall height width a.
(C height, C width, Eq width, Floating a) =>
General height width a -> Square width a -> Square height a
Basic.congruenceAdjoint



solve ::
   (Extent.C vert, Extent.C horiz,
    Shape.C sh, Eq sh, Shape.C nrhs, Class.Floating a) =>
   Square sh a -> Full vert horiz sh nrhs a -> Full vert horiz sh nrhs a
solve :: Square sh a
-> Full vert horiz sh nrhs a -> Full vert horiz sh nrhs a
solve = (Array (Square sh) a
 -> Array (Full vert horiz sh nrhs) a
 -> Array (Full vert horiz sh nrhs) a)
-> Square sh a
-> Full vert horiz sh nrhs a
-> Full vert horiz sh nrhs a
forall shA a shB b shC c.
(Array shA a -> Array shB b -> Array shC c)
-> ArrayMatrix shA a -> ArrayMatrix shB b -> ArrayMatrix shC c
ArrMatrix.lift2 Array (Square sh) a
-> Array (Full vert horiz sh nrhs) a
-> Array (Full vert horiz sh nrhs) a
forall vert horiz sh nrhs a.
(C vert, C horiz, C sh, Eq sh, C nrhs, Floating a) =>
Square sh a
-> Full vert horiz sh nrhs a -> Full vert horiz sh nrhs a
Linear.solve

inverse :: (Shape.C sh, Class.Floating a) => Square sh a -> Square sh a
inverse :: Square sh a -> Square sh a
inverse = (Array (Square sh) a -> Array (Square sh) a)
-> Square sh a -> Square sh a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 Array (Square sh) a -> Array (Square sh) a
forall sh a. (C sh, Floating a) => Square sh a -> Square sh a
Linear.inverse

determinant :: (Shape.C sh, Class.Floating a) => Square sh a -> a
determinant :: Square sh a -> a
determinant = Square sh a -> a
forall sh a. (C sh, Floating a) => Square sh a -> a
Linear.determinant (Square sh a -> a)
-> (Square sh a -> Square sh a) -> Square sh a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Square sh a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector



eigenvalues ::
   (Shape.C sh, Class.Floating a) =>
   Square sh a -> Vector sh (ComplexOf a)
eigenvalues :: Square sh a -> Vector sh (ComplexOf a)
eigenvalues = Square sh a -> Vector sh (ComplexOf a)
forall sh a.
(C sh, Floating a) =>
Square sh a -> Vector sh (ComplexOf a)
Eigen.values (Square sh a -> Vector sh (ComplexOf a))
-> (Square sh a -> Square sh a)
-> Square sh a
-> Vector sh (ComplexOf a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Square sh a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector

{- |
If @(q,r) = schur a@, then @a = q \<\> r \<\> adjoint q@,
where @q@ is unitary (orthogonal)
and @r@ is a right-upper triangular matrix for complex @a@
and a 1x1-or-2x2-block upper triangular matrix for real @a@.
With @takeDiagonal r@ you get all eigenvalues of @a@ if @a@ is complex
and the real parts of the eigenvalues if @a@ is real.
Complex conjugated eigenvalues of a real matrix @a@
are encoded as 2x2 blocks along the diagonal.


The meaning and order of matrix factors of these functions is consistent:

* 'Numeric.LAPACK.Matrix.Square.schur'
* 'Numeric.LAPACK.Matrix.Square.schurComplex'
* 'Numeric.LAPACK.Matrix.Hermitian.eigensystem'
* 'Numeric.LAPACK.Matrix.BandedHermitian.eigensystem'
* 'Numeric.LAPACK.Matrix.Square.congruenceAdjoint'
* 'Numeric.LAPACK.Matrix.Hermitian.gramianAdjoint'
* 'Numeric.LAPACK.Matrix.Hermitian.anticommutatorAdjoint'
* 'Numeric.LAPACK.Matrix.Hermitian.congruenceAdjoint'
* 'Numeric.LAPACK.Matrix.Hermitian.congruenceDiagonalAdjoint'
-}
schur ::
   (Shape.C sh, Class.Floating a) =>
   Square sh a -> (Square sh a, Square sh a)
schur :: Square sh a -> (Square sh a, Square sh a)
schur =
   (Array (Square sh) a -> Square sh a,
 Array (Square sh) a -> Square sh a)
-> (Array (Square sh) a, Array (Square sh) a)
-> (Square sh a, Square sh a)
forall a c b d. (a -> c, b -> d) -> (a, b) -> (c, d)
mapPair (Array (Square sh) a -> Square sh a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0, Array (Square sh) a -> Square sh a
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0) ((Array (Square sh) a, Array (Square sh) a)
 -> (Square sh a, Square sh a))
-> (Square sh a -> (Array (Square sh) a, Array (Square sh) a))
-> Square sh a
-> (Square sh a, Square sh a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Array (Square sh) a -> (Array (Square sh) a, Array (Square sh) a)
forall sh a.
(C sh, Floating a) =>
Square sh a -> (Square sh a, Square sh a)
Eigen.schur (Array (Square sh) a -> (Array (Square sh) a, Array (Square sh) a))
-> (Square sh a -> Array (Square sh) a)
-> Square sh a
-> (Array (Square sh) a, Array (Square sh) a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Array (Square sh) a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector

schurComplex ::
   (Shape.C sh, Class.Real a, Complex a ~ ac) =>
   Square sh ac -> (Square sh ac, Triangular.Upper sh ac)
schurComplex :: Square sh ac -> (Square sh ac, Upper sh ac)
schurComplex = (Square sh ac -> Upper sh ac)
-> (Square sh ac, Square sh ac) -> (Square sh ac, Upper sh ac)
forall b c a. (b -> c) -> (a, b) -> (a, c)
mapSnd Square sh ac -> Upper sh ac
forall vert height width a.
(C vert, C height, C width, Floating a) =>
Full vert Small height width a -> Upper width a
Triangular.takeUpper ((Square sh ac, Square sh ac) -> (Square sh ac, Upper sh ac))
-> (Square sh ac -> (Square sh ac, Square sh ac))
-> Square sh ac
-> (Square sh ac, Upper sh ac)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh ac -> (Square sh ac, Square sh ac)
forall sh a.
(C sh, Floating a) =>
Square sh a -> (Square sh a, Square sh a)
schur


{- |
@(vr,d,vlAdj) = eigensystem a@

Counterintuitively, @vr@ contains the right eigenvectors as columns
and @vlAdj@ contains the left conjugated eigenvectors as rows.
The idea is to provide a decomposition of @a@.
If @a@ is diagonalizable, then @vr@ and @vlAdj@
are almost inverse to each other.
More precisely, @vlAdj \<\> vr@ is a diagonal matrix,
but not necessarily an identity matrix.
This is because all eigenvectors are normalized to Euclidean norm 1.
With the following scaling, the decomposition becomes perfect:

> let scal = takeDiagonal $ vlAdj <> vr
> a == vr #*\ Vector.divide d scal ##*# vlAdj

If @a@ is non-diagonalizable
then some columns of @vr@ and corresponding rows of @vlAdj@ are left zero
and the above property does not hold.


The meaning and order of result matrices of these functions is consistent:

* 'Numeric.LAPACK.Matrix.Square.eigensystem'
* 'Numeric.LAPACK.Matrix.Triangular.eigensystem'
* 'Numeric.LAPACK.Singular.decompose'
* 'Numeric.LAPACK.Singular.decomposeTall'
* 'Numeric.LAPACK.Singular.decomposeWide'
-}
eigensystem ::
   (Shape.C sh, Class.Floating a, ComplexOf a ~ ac) =>
   Square sh a -> (Square sh ac, Vector sh ac, Square sh ac)
eigensystem :: Square sh a -> (Square sh ac, Vector sh ac, Square sh ac)
eigensystem =
   (Array (Square sh) ac -> Square sh ac,
 Vector sh ac -> Vector sh ac, Array (Square sh) ac -> Square sh ac)
-> (Array (Square sh) ac, Vector sh ac, Array (Square sh) ac)
-> (Square sh ac, Vector sh ac, Square sh ac)
forall a d b e c f.
(a -> d, b -> e, c -> f) -> (a, b, c) -> (d, e, f)
mapTriple (Array (Square sh) ac -> Square sh ac
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0, Vector sh ac -> Vector sh ac
forall a. a -> a
id, Array (Square sh) ac -> Square sh ac
forall shA a. Array shA a -> ArrayMatrix shA a
ArrMatrix.lift0) ((Array (Square sh) ac, Vector sh ac, Array (Square sh) ac)
 -> (Square sh ac, Vector sh ac, Square sh ac))
-> (Square sh a
    -> (Array (Square sh) ac, Vector sh ac, Array (Square sh) ac))
-> Square sh a
-> (Square sh ac, Vector sh ac, Square sh ac)
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   Square sh a
-> (Array (Square sh) ac, Vector sh ac, Array (Square sh) ac)
forall sh a ac.
(C sh, Floating a, ComplexOf a ~ ac) =>
Square sh a -> (Square sh ac, Vector sh ac, Square sh ac)
Eigen.decompose (Square sh a
 -> (Array (Square sh) ac, Vector sh ac, Array (Square sh) ac))
-> (Square sh a -> Square sh a)
-> Square sh a
-> (Array (Square sh) ac, Vector sh ac, Array (Square sh) ac)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Square sh a -> Square sh a
forall sh a. ArrayMatrix sh a -> Array sh a
ArrMatrix.toVector