{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE EmptyDataDecls #-}
module Numeric.LAPACK.Orthogonal.Private where

import qualified Numeric.LAPACK.Matrix.Divide as Divide
import qualified Numeric.LAPACK.Matrix.Multiply as Multiply
import qualified Numeric.LAPACK.Matrix.Type as Type
import qualified Numeric.LAPACK.Matrix.Array as ArrMatrix
import qualified Numeric.LAPACK.Matrix.Basic as Basic
import qualified Numeric.LAPACK.Matrix.Private as Matrix
import qualified Numeric.LAPACK.Matrix.Shape.Private as MatrixShape
import qualified Numeric.LAPACK.Matrix.Extent.Private as ExtentPriv
import qualified Numeric.LAPACK.Matrix.Extent as Extent
import qualified Numeric.LAPACK.Split as Split
import qualified Numeric.LAPACK.Shape as ExtShape
import Numeric.LAPACK.Output ((/+/))
import Numeric.LAPACK.Matrix.Plain.Format (formatArray)
import Numeric.LAPACK.Matrix.Type (FormatMatrix(formatMatrix))
import Numeric.LAPACK.Matrix.Triangular.Basic (Upper)
import Numeric.LAPACK.Matrix.Shape.Private
         (Order(RowMajor, ColumnMajor), sideSwapFromOrder)
import Numeric.LAPACK.Matrix.Extent.Private (Extent)
import Numeric.LAPACK.Matrix.Modifier
         (Transposition(NonTransposed, Transposed),
          Conjugation(NonConjugated, Conjugated))
import Numeric.LAPACK.Matrix.Private (Full)
import Numeric.LAPACK.Vector (Vector)
import Numeric.LAPACK.Scalar (RealOf, zero, isZero, absolute, conjugate)
import Numeric.LAPACK.Private
         (fill, copySubMatrix, copyBlock, conjugateToTemp, caseRealComplexFunc,
          withAutoWorkspaceInfo, errorCodeMsg)

import qualified Numeric.LAPACK.FFI.Generic as LapackGen
import qualified Numeric.Netlib.Utility as Call
import qualified Numeric.Netlib.Class as Class

import qualified Data.Array.Comfort.Storable.Unchecked.Monadic as ArrayIO
import qualified Data.Array.Comfort.Storable.Unchecked as Array
import qualified Data.Array.Comfort.Shape as Shape
import Data.Array.Comfort.Storable.Unchecked (Array(Array))

import Foreign.Marshal.Array (advancePtr)
import Foreign.ForeignPtr (ForeignPtr, withForeignPtr)

import Control.Monad.Trans.Cont (ContT(ContT), evalContT)
import Control.Monad.IO.Class (liftIO)
import Control.Monad (when)
import Control.Applicative (liftA2)

import qualified Data.List as List
import Data.Monoid ((<>))


data Hh vert horiz height width

data instance Type.Matrix (Hh vert horiz height width) a =
   Householder {
      Matrix (Hh vert horiz height width) a
-> Vector (Min height width) a
tau_ :: Vector (ExtShape.Min height width) a,
      Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_ ::
         Array
            (MatrixShape.Split MatrixShape.Reflector vert horiz height width) a
   } deriving (Int -> Matrix (Hh vert horiz height width) a -> ShowS
[Matrix (Hh vert horiz height width) a] -> ShowS
Matrix (Hh vert horiz height width) a -> String
(Int -> Matrix (Hh vert horiz height width) a -> ShowS)
-> (Matrix (Hh vert horiz height width) a -> String)
-> ([Matrix (Hh vert horiz height width) a] -> ShowS)
-> Show (Matrix (Hh vert horiz height width) a)
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
Int -> Matrix (Hh vert horiz height width) a -> ShowS
forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
[Matrix (Hh vert horiz height width) a] -> ShowS
forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
Matrix (Hh vert horiz height width) a -> String
showList :: [Matrix (Hh vert horiz height width) a] -> ShowS
$cshowList :: forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
[Matrix (Hh vert horiz height width) a] -> ShowS
show :: Matrix (Hh vert horiz height width) a -> String
$cshow :: forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
Matrix (Hh vert horiz height width) a -> String
showsPrec :: Int -> Matrix (Hh vert horiz height width) a -> ShowS
$cshowsPrec :: forall vert horiz height width a.
(C height, C width, Storable a, Show height, Show width, Show a,
 C vert, C horiz) =>
Int -> Matrix (Hh vert horiz height width) a -> ShowS
Show)

type Householder vert horiz height width =
         Type.Matrix (Hh vert horiz height width)
type General height width = Householder Extent.Big Extent.Big height width
type Tall height width = Householder Extent.Big Extent.Small height width
type Wide height width = Householder Extent.Small Extent.Big height width
type Square sh = Householder Extent.Small Extent.Small sh sh


extent_ ::
   Householder vert horiz height width a ->
   Extent vert horiz height width
extent_ :: Householder vert horiz height width a
-> Extent vert horiz height width
extent_ = Split Reflector vert horiz height width
-> Extent vert horiz height width
forall lower vert horiz height width.
Split lower vert horiz height width
-> Extent vert horiz height width
MatrixShape.splitExtent (Split Reflector vert horiz height width
 -> Extent vert horiz height width)
-> (Householder vert horiz height width a
    -> Split Reflector vert horiz height width)
-> Householder vert horiz height width a
-> Extent vert horiz height width
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Array (Split Reflector vert horiz height width) a
-> Split Reflector vert horiz height width
forall sh a. Array sh a -> sh
Array.shape (Array (Split Reflector vert horiz height width) a
 -> Split Reflector vert horiz height width)
-> (Householder vert horiz height width a
    -> Array (Split Reflector vert horiz height width) a)
-> Householder vert horiz height width a
-> Split Reflector vert horiz height width
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vert horiz height width a
-> Array (Split Reflector vert horiz height width) a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

mapExtent ::
   (Extent.C vertA, Extent.C horizA) =>
   (Extent.C vertB, Extent.C horizB) =>
   Extent.Map vertA horizA vertB horizB height width ->
   Householder vertA horizA height width a ->
   Householder vertB horizB height width a
mapExtent :: Map vertA horizA vertB horizB height width
-> Householder vertA horizA height width a
-> Householder vertB horizB height width a
mapExtent Map vertA horizA vertB horizB height width
f (Householder tau split) =
   Vector (Min height width) a
-> Array (Split Reflector vertB horizB height width) a
-> Householder vertB horizB height width a
forall vert horiz height width a.
Vector (Min height width) a
-> Array (Split Reflector vert horiz height width) a
-> Matrix (Hh vert horiz height width) a
Householder Vector (Min height width) a
tau (Array (Split Reflector vertB horizB height width) a
 -> Householder vertB horizB height width a)
-> Array (Split Reflector vertB horizB height width) a
-> Householder vertB horizB height width a
forall a b. (a -> b) -> a -> b
$ (Split Reflector vertA horizA height width
 -> Split Reflector vertB horizB height width)
-> Array (Split Reflector vertA horizA height width) a
-> Array (Split Reflector vertB horizB height width) a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape (Map vertA horizA vertB horizB height width
-> Split Reflector vertA horizA height width
-> Split Reflector vertB horizB height width
forall vertA horizA vertB horizB height width lower.
Map vertA horizA vertB horizB height width
-> Split lower vertA horizA height width
-> Split lower vertB horizB height width
MatrixShape.splitMapExtent Map vertA horizA vertB horizB height width
f) Array (Split Reflector vertA horizA height width) a
split

caseTallWide ::
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width) =>
   Householder vert horiz height width a ->
   Either (Tall height width a) (Wide height width a)
caseTallWide :: Householder vert horiz height width a
-> Either (Tall height width a) (Wide height width a)
caseTallWide (Householder tau (Array shape a)) =
   (Split Reflector Big Small height width
 -> Either (Tall height width a) (Wide height width a))
-> (Split Reflector Small Big height width
    -> Either (Tall height width a) (Wide height width a))
-> Either
     (Split Reflector Big Small height width)
     (Split Reflector Small Big height width)
-> Either (Tall height width a) (Wide height width a)
forall a c b. (a -> c) -> (b -> c) -> Either a b -> c
either
      (Tall height width a
-> Either (Tall height width a) (Wide height width a)
forall a b. a -> Either a b
Left (Tall height width a
 -> Either (Tall height width a) (Wide height width a))
-> (Split Reflector Big Small height width -> Tall height width a)
-> Split Reflector Big Small height width
-> Either (Tall height width a) (Wide height width a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector (Min height width) a
-> Array (Split Reflector Big Small height width) a
-> Tall height width a
forall vert horiz height width a.
Vector (Min height width) a
-> Array (Split Reflector vert horiz height width) a
-> Matrix (Hh vert horiz height width) a
Householder Vector (Min height width) a
tau (Array (Split Reflector Big Small height width) a
 -> Tall height width a)
-> (Split Reflector Big Small height width
    -> Array (Split Reflector Big Small height width) a)
-> Split Reflector Big Small height width
-> Tall height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Split Reflector Big Small height width
 -> ForeignPtr a
 -> Array (Split Reflector Big Small height width) a)
-> ForeignPtr a
-> Split Reflector Big Small height width
-> Array (Split Reflector Big Small height width) a
forall a b c. (a -> b -> c) -> b -> a -> c
flip Split Reflector Big Small height width
-> ForeignPtr a -> Array (Split Reflector Big Small height width) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array ForeignPtr a
a)
      (Wide height width a
-> Either (Tall height width a) (Wide height width a)
forall a b. b -> Either a b
Right (Wide height width a
 -> Either (Tall height width a) (Wide height width a))
-> (Split Reflector Small Big height width -> Wide height width a)
-> Split Reflector Small Big height width
-> Either (Tall height width a) (Wide height width a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector (Min height width) a
-> Array (Split Reflector Small Big height width) a
-> Wide height width a
forall vert horiz height width a.
Vector (Min height width) a
-> Array (Split Reflector vert horiz height width) a
-> Matrix (Hh vert horiz height width) a
Householder Vector (Min height width) a
tau (Array (Split Reflector Small Big height width) a
 -> Wide height width a)
-> (Split Reflector Small Big height width
    -> Array (Split Reflector Small Big height width) a)
-> Split Reflector Small Big height width
-> Wide height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Split Reflector Small Big height width
 -> ForeignPtr a
 -> Array (Split Reflector Small Big height width) a)
-> ForeignPtr a
-> Split Reflector Small Big height width
-> Array (Split Reflector Small Big height width) a
forall a b c. (a -> b -> c) -> b -> a -> c
flip Split Reflector Small Big height width
-> ForeignPtr a -> Array (Split Reflector Small Big height width) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array ForeignPtr a
a) (Either
   (Split Reflector Big Small height width)
   (Split Reflector Small Big height width)
 -> Either (Tall height width a) (Wide height width a))
-> Either
     (Split Reflector Big Small height width)
     (Split Reflector Small Big height width)
-> Either (Tall height width a) (Wide height width a)
forall a b. (a -> b) -> a -> b
$
   Split Reflector vert horiz height width
-> Either
     (Split Reflector Big Small height width)
     (Split Reflector Small Big height width)
forall vert horiz height width lower.
(C vert, C horiz, C height, C width) =>
Split lower vert horiz height width
-> Either
     (Split lower Big Small height width)
     (Split lower Small Big height width)
MatrixShape.caseTallWideSplit Split Reflector vert horiz height width
shape


instance
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width) =>
      FormatMatrix (Hh vert horiz height width) where
   formatMatrix :: String -> Matrix (Hh vert horiz height width) a -> out
formatMatrix String
fmt (Householder tau m) =
      String -> Array (ZeroBased Int) a -> out
forall sh a out.
(FormatArray sh, Floating a, Output out) =>
String -> Array sh a -> out
formatArray String
fmt ((Min height width -> ZeroBased Int)
-> Array (Min height width) a -> Array (ZeroBased Int) a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape (Int -> ZeroBased Int
forall n. n -> ZeroBased n
Shape.ZeroBased (Int -> ZeroBased Int)
-> (Min height width -> Int) -> Min height width -> ZeroBased Int
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Min height width -> Int
forall sh. C sh => sh -> Int
Shape.size) Array (Min height width) a
tau)
      out -> out -> out
forall out. Output out => out -> out -> out
/+/
      String -> Array (Split Reflector vert horiz height width) a -> out
forall sh a out.
(FormatArray sh, Floating a, Output out) =>
String -> Array sh a -> out
formatArray String
fmt Array (Split Reflector vert horiz height width) a
m

fromMatrix ::
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width,
    Class.Floating a) =>
   Full vert horiz height width a ->
   Householder vert horiz height width a
fromMatrix :: Full vert horiz height width a
-> Householder vert horiz height width a
fromMatrix (Array shape :: Full vert horiz height width
shape@(MatrixShape.Full Order
order Extent vert horiz height width
extent) ForeignPtr a
a) =
   (Vector (Min height width) a
 -> Array (Split Reflector vert horiz height width) a
 -> Householder vert horiz height width a)
-> (Vector (Min height width) a,
    Array (Split Reflector vert horiz height width) a)
-> Householder vert horiz height width a
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry Vector (Min height width) a
-> Array (Split Reflector vert horiz height width) a
-> Householder vert horiz height width a
forall vert horiz height width a.
Vector (Min height width) a
-> Array (Split Reflector vert horiz height width) a
-> Matrix (Hh vert horiz height width) a
Householder ((Vector (Min height width) a,
  Array (Split Reflector vert horiz height width) a)
 -> Householder vert horiz height width a)
-> (Vector (Min height width) a,
    Array (Split Reflector vert horiz height width) a)
-> Householder vert horiz height width a
forall a b. (a -> b) -> a -> b
$
   Min height width
-> (Int
    -> Ptr a -> IO (Array (Split Reflector vert horiz height width) a))
-> (Vector (Min height width) a,
    Array (Split Reflector vert horiz height width) a)
forall sh a b.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO b) -> (Array sh a, b)
Array.unsafeCreateWithSizeAndResult
      ((height -> width -> Min height width)
-> (height, width) -> Min height width
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry height -> width -> Min height width
forall sh0 sh1. sh0 -> sh1 -> Min sh0 sh1
ExtShape.Min ((height, width) -> Min height width)
-> (height, width) -> Min height width
forall a b. (a -> b) -> a -> b
$ Extent vert horiz height width -> (height, width)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vert horiz height width
extent) ((Int
  -> Ptr a -> IO (Array (Split Reflector vert horiz height width) a))
 -> (Vector (Min height width) a,
     Array (Split Reflector vert horiz height width) a))
-> (Int
    -> Ptr a -> IO (Array (Split Reflector vert horiz height width) a))
-> (Vector (Min height width) a,
    Array (Split Reflector vert horiz height width) a)
forall a b. (a -> b) -> a -> b
$ \Int
_ Ptr a
tauPtr ->
   Split Reflector vert horiz height width
-> (Ptr a -> IO ())
-> IO (Array (Split Reflector vert horiz height width) a)
forall (m :: * -> *) sh a.
(PrimMonad m, C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> m (Array sh a)
ArrayIO.unsafeCreate
      (Reflector
-> Order
-> Extent vert horiz height width
-> Split Reflector vert horiz height width
forall lower vert horiz height width.
lower
-> Order
-> Extent vert horiz height width
-> Split lower vert horiz height width
MatrixShape.Split Reflector
MatrixShape.Reflector Order
order Extent vert horiz height width
extent) ((Ptr a -> IO ())
 -> IO (Array (Split Reflector vert horiz height width) a))
-> (Ptr a -> IO ())
-> IO (Array (Split Reflector vert horiz height width) a)
forall a b. (a -> b) -> a -> b
$ \Ptr a
qrPtr ->

   ContT () IO () -> IO ()
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT () IO () -> IO ()) -> ContT () IO () -> IO ()
forall a b. (a -> b) -> a -> b
$ do
      let (Int
m,Int
n) = Full vert horiz height width -> (Int, Int)
forall vert horiz height width.
(C vert, C horiz, C height, C width) =>
Full vert horiz height width -> (Int, Int)
MatrixShape.dimensions Full vert horiz height width
shape
      Ptr CInt
mPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
m
      Ptr CInt
nPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
a
      Ptr CInt
ldaPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
      IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
         Int -> Ptr a -> Ptr a -> IO ()
forall a. Floating a => Int -> Ptr a -> Ptr a -> IO ()
copyBlock (Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
n) Ptr a
aPtr Ptr a
qrPtr
         case Order
order of
            Order
RowMajor ->
               String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"gelqf" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
                  Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.gelqf Ptr CInt
mPtr Ptr CInt
nPtr Ptr a
qrPtr Ptr CInt
ldaPtr Ptr a
tauPtr
            Order
ColumnMajor ->
               String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"geqrf" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
                  Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.geqrf Ptr CInt
mPtr Ptr CInt
nPtr Ptr a
qrPtr Ptr CInt
ldaPtr Ptr a
tauPtr

determinantR ::
   (Extent.C vert, Shape.C height, Shape.C width, Class.Floating a) =>
   Householder vert Extent.Small height width a -> a
determinantR :: Householder vert Small height width a -> a
determinantR = Split Reflector vert Small height width a -> a
forall vert height width a lower.
(C vert, C height, C width, Floating a) =>
Split lower vert Small height width a -> a
Split.determinantR (Split Reflector vert Small height width a -> a)
-> (Householder vert Small height width a
    -> Split Reflector vert Small height width a)
-> Householder vert Small height width a
-> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vert Small height width a
-> Split Reflector vert Small height width a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

{-
For complex numbers LAPACK uses not exactly reflections,
i.e. the determinants of the primitive transformations are not necessarily -1.

It holds: det(I-tau*v*v^H) = 1-tau*v^H*v
   because of https://en.wikipedia.org/wiki/Sylvester's_determinant_theorem
   simple proof from: https://en.wikipedia.org/wiki/Matrix_determinant_lemma
   I  0 . I+u*vt u .  I  0  =  I+u*vt     u      .  I  0 = I u
   vt 1     0    1   -vt 1     vt+vt*u*vt vt*u+1   -vt 1   0 vt*u+1

We already know:
   v^H*v is real and greater or equal to 1, because v[i] = 1,
   and determinant has absolute value 1.

Let k = v^H*v.
For which real k lies 1-tau*k on the unit circle?

   (1-taur*k)^2 + (taui*k)^2 = 1
   1-2*taur*k+(taur^2+taui^2)*k^2 = 1
   (taur^2 + taui^2)*k^2 - 2*taur*k = 0   (k/=0)
   (taur^2 + taui^2)*k - 2*taur = 0
   k = 2*taur / (taur^2 + taui^2)

   1-tau*k
      = (taur^2 + taui^2 - tau*2*taur) / (taur^2 + taui^2)
      = (taur^2 + taui^2 - 2*(taur+i*taui)*taur) / (taur^2 + taui^2)
      = (-taur^2 + taui^2 - 2*(i*taui)*taur) / (taur^2 + taui^2)
      = -(taur + i*taui)^2 / (taur^2 + taui^2)
-}
determinant ::
   (Shape.C sh, Class.Floating a) =>
   Square sh a -> a
determinant :: Square sh a -> a
determinant (Householder tau split) =
   (a -> a -> a) -> a -> [a] -> a
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
List.foldl' a -> a -> a
forall a. Num a => a -> a -> a
(*) (Split Reflector Small Small sh sh a -> a
forall vert height width a lower.
(C vert, C height, C width, Floating a) =>
Split lower vert Small height width a -> a
Split.determinantR Split Reflector Small Small sh sh a
split) ([a] -> a) -> [a] -> a
forall a b. (a -> b) -> a -> b
$
   (case Split Reflector Small Small sh sh -> Order
forall lower vert horiz height width.
Split lower vert horiz height width -> Order
MatrixShape.splitOrder (Split Reflector Small Small sh sh -> Order)
-> Split Reflector Small Small sh sh -> Order
forall a b. (a -> b) -> a -> b
$ Split Reflector Small Small sh sh a
-> Split Reflector Small Small sh sh
forall sh a. Array sh a -> sh
Array.shape Split Reflector Small Small sh sh a
split of
      Order
RowMajor -> (a -> a) -> [a] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map a -> a
forall a. Floating a => a -> a
conjugate
      Order
ColumnMajor -> [a] -> [a]
forall a. a -> a
id) ([a] -> [a]) -> [a] -> [a]
forall a b. (a -> b) -> a -> b
$
   (a -> a) -> [a] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map (a -> a
forall a. Num a => a -> a
negate(a -> a) -> (a -> a) -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.(a -> Int -> a
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
2::Int))(a -> a) -> (a -> a) -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.a -> a
forall a. Num a => a -> a
signum) ([a] -> [a]) -> [a] -> [a]
forall a b. (a -> b) -> a -> b
$
   (a -> Bool) -> [a] -> [a]
forall a. (a -> Bool) -> [a] -> [a]
filter (Bool -> Bool
not (Bool -> Bool) -> (a -> Bool) -> a -> Bool
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> Bool
forall a. Floating a => a -> Bool
isZero) ([a] -> [a]) -> [a] -> [a]
forall a b. (a -> b) -> a -> b
$ Array (Min sh sh) a -> [a]
forall sh a. (C sh, Storable a) => Array sh a -> [a]
Array.toList Array (Min sh sh) a
tau

determinantAbsolute ::
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width,
    Class.Floating a) =>
   Householder vert horiz height width a -> RealOf a
determinantAbsolute :: Householder vert horiz height width a -> RealOf a
determinantAbsolute =
   a -> RealOf a
forall a. Floating a => a -> RealOf a
absolute (a -> RealOf a)
-> (Householder vert horiz height width a -> a)
-> Householder vert horiz height width a
-> RealOf a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Householder Big Small height width a -> a)
-> (Wide height width a -> a)
-> Either
     (Householder Big Small height width a) (Wide height width a)
-> a
forall a c b. (a -> c) -> (b -> c) -> Either a b -> c
either Householder Big Small height width a -> a
forall vert height width a.
(C vert, C height, C width, Floating a) =>
Householder vert Small height width a -> a
determinantR (a -> Wide height width a -> a
forall a b. a -> b -> a
const a
forall a. Floating a => a
zero) (Either
   (Householder Big Small height width a) (Wide height width a)
 -> a)
-> (Householder vert horiz height width a
    -> Either
         (Householder Big Small height width a) (Wide height width a))
-> Householder vert horiz height width a
-> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vert horiz height width a
-> Either
     (Householder Big Small height width a) (Wide height width a)
forall vert horiz height width a.
(C vert, C horiz, C height, C width) =>
Householder vert horiz height width a
-> Either (Tall height width a) (Wide height width a)
caseTallWide


leastSquares ::
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Eq width, Shape.C nrhs,
    Class.Floating a) =>
   Householder horiz Extent.Small height width a ->
   Full vert horiz height nrhs a ->
   Full vert horiz width nrhs a
leastSquares :: Householder horiz Small height width a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
leastSquares Householder horiz Small height width a
qr =
   Transposition
-> Conjugation
-> Householder horiz Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
forall vertA vert horiz height width nrhs a.
(C vertA, C vert, C horiz, C height, C width, Eq width, C nrhs,
 Floating a) =>
Transposition
-> Conjugation
-> Householder vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
tallSolveR Transposition
NonTransposed Conjugation
NonConjugated Householder horiz Small height width a
qr (Full vert horiz width nrhs a -> Full vert horiz width nrhs a)
-> (Full vert horiz height nrhs a -> Full vert horiz width nrhs a)
-> Full vert horiz height nrhs a
-> Full vert horiz width nrhs a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder horiz Small height width a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
forall vert horiz height width fuse a.
(C vert, C horiz, C height, C width, C fuse, Eq fuse,
 Floating a) =>
Householder horiz Small fuse height a
-> Full vert horiz fuse width a -> Full vert horiz height width a
tallMultiplyQAdjoint Householder horiz Small height width a
qr

minimumNorm ::
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Eq width, Shape.C nrhs,
    Class.Floating a) =>
   Householder vert Extent.Small width height a ->
   Full vert horiz height nrhs a ->
   Full vert horiz width nrhs a
minimumNorm :: Householder vert Small width height a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
minimumNorm Householder vert Small width height a
qr = Householder vert Small width height a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
forall vert horiz height width fuse a.
(C vert, C horiz, C height, Eq height, C width, C fuse, Eq fuse,
 Floating a) =>
Householder vert Small height fuse a
-> Full vert horiz fuse width a -> Full vert horiz height width a
tallMultiplyQ Householder vert Small width height a
qr (Full vert horiz height nrhs a -> Full vert horiz width nrhs a)
-> (Full vert horiz height nrhs a -> Full vert horiz height nrhs a)
-> Full vert horiz height nrhs a
-> Full vert horiz width nrhs a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Transposition
-> Conjugation
-> Householder vert Small width height a
-> Full vert horiz height nrhs a
-> Full vert horiz height nrhs a
forall vertA vert horiz height width nrhs a.
(C vertA, C vert, C horiz, C height, C width, Eq width, C nrhs,
 Floating a) =>
Transposition
-> Conjugation
-> Householder vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
tallSolveR Transposition
Transposed Conjugation
Conjugated Householder vert Small width height a
qr

takeRows ::
   (Extent.C vert, Extent.C horiz,
    Eq fuse, Shape.C fuse, Shape.C height, Shape.C width, Class.Floating a) =>
   Extent Extent.Small horiz height fuse ->
   Full vert horiz fuse width a ->
   Full vert horiz height width a
takeRows :: Extent Small horiz height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
takeRows Extent Small horiz height fuse
extentA (Array (MatrixShape.Full Order
order Extent vert horiz fuse width
extentB) ForeignPtr a
b) =
   case Extent vert horiz height fuse
-> Extent vert horiz fuse width
-> Maybe (Extent vert horiz height width)
forall vert horiz fuse height width.
(C vert, C horiz, Eq fuse) =>
Extent vert horiz height fuse
-> Extent vert horiz fuse width
-> Maybe (Extent vert horiz height width)
Extent.fuse (Extent Small horiz height fuse -> Extent vert horiz height fuse
forall vert horiz height width.
(C vert, C horiz) =>
Extent Small horiz height width -> Extent vert horiz height width
ExtentPriv.generalizeWide Extent Small horiz height fuse
extentA) Extent vert horiz fuse width
extentB of
      Maybe (Extent vert horiz height width)
Nothing -> String -> Full vert horiz height width a
forall a. HasCallStack => String -> a
error String
"Householder.takeRows: heights mismatch"
      Just Extent vert horiz height width
extentC ->
         fuse
-> Int
-> ForeignPtr a
-> Full vert horiz height width
-> Full vert horiz height width a
forall vert horiz heightA height width a.
(C vert, C horiz, C heightA, C height, C width, Floating a) =>
heightA
-> Int
-> ForeignPtr a
-> Full vert horiz height width
-> Full vert horiz height width a
Basic.takeSub
            (Extent vert horiz fuse width -> fuse
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> height
Extent.height Extent vert horiz fuse width
extentB) Int
0 ForeignPtr a
b (Order
-> Extent vert horiz height width -> Full vert horiz height width
forall vert horiz height width.
Order
-> Extent vert horiz height width -> Full vert horiz height width
MatrixShape.Full Order
order Extent vert horiz height width
extentC)

addRows ::
   (Extent.C vert, Extent.C horiz,
    Eq fuse, Shape.C fuse, Shape.C height, Shape.C width, Class.Floating a) =>
   Extent vert Extent.Small height fuse ->
   Full vert horiz fuse width a ->
   Full vert horiz height width a
addRows :: Extent vert Small height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
addRows Extent vert Small height fuse
extentA (Array shapeB :: Full vert horiz fuse width
shapeB@(MatrixShape.Full Order
order Extent vert horiz fuse width
extentB) ForeignPtr a
b) =
   case Extent vert horiz height fuse
-> Extent vert horiz fuse width
-> Maybe (Extent vert horiz height width)
forall vert horiz fuse height width.
(C vert, C horiz, Eq fuse) =>
Extent vert horiz height fuse
-> Extent vert horiz fuse width
-> Maybe (Extent vert horiz height width)
Extent.fuse (Extent vert Small height fuse -> Extent vert horiz height fuse
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert Small height width -> Extent vert horiz height width
ExtentPriv.generalizeTall Extent vert Small height fuse
extentA) Extent vert horiz fuse width
extentB of
      Maybe (Extent vert horiz height width)
Nothing -> String -> Full vert horiz height width a
forall a. HasCallStack => String -> a
error String
"Householder.addRows: heights mismatch"
      Just Extent vert horiz height width
extentC ->
         Full vert horiz height width
-> (Int -> Ptr a -> IO ()) -> Full vert horiz height width a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (Order
-> Extent vert horiz height width -> Full vert horiz height width
forall vert horiz height width.
Order
-> Extent vert horiz height width -> Full vert horiz height width
MatrixShape.Full Order
order Extent vert horiz height width
extentC) ((Int -> Ptr a -> IO ()) -> Full vert horiz height width a)
-> (Int -> Ptr a -> IO ()) -> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$
            \Int
cSize Ptr a
cPtr ->
         ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
b ((Ptr a -> IO ()) -> IO ()) -> (Ptr a -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \Ptr a
bPtr ->
         case Order
order of
            Order
RowMajor -> do
               let bSize :: Int
bSize = Full vert horiz fuse width -> Int
forall sh. C sh => sh -> Int
Shape.size Full vert horiz fuse width
shapeB
               Int -> Ptr a -> Ptr a -> IO ()
forall a. Floating a => Int -> Ptr a -> Ptr a -> IO ()
copyBlock Int
bSize Ptr a
bPtr Ptr a
cPtr
               a -> Int -> Ptr a -> IO ()
forall a. Floating a => a -> Int -> Ptr a -> IO ()
fill a
forall a. Floating a => a
zero (Int
cSize Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
bSize) (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr Int
bSize)
            Order
ColumnMajor -> do
               let n :: Int
n  = width -> Int
forall sh. C sh => sh -> Int
Shape.size (width -> Int) -> width -> Int
forall a b. (a -> b) -> a -> b
$ Extent vert horiz fuse width -> width
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> width
Extent.width  Extent vert horiz fuse width
extentB
                   mb :: Int
mb = fuse -> Int
forall sh. C sh => sh -> Int
Shape.size (fuse -> Int) -> fuse -> Int
forall a b. (a -> b) -> a -> b
$ Extent vert horiz fuse width -> fuse
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> height
Extent.height Extent vert horiz fuse width
extentB
                   mc :: Int
mc = height -> Int
forall sh. C sh => sh -> Int
Shape.size (height -> Int) -> height -> Int
forall a b. (a -> b) -> a -> b
$ Extent vert horiz height width -> height
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> height
Extent.height Extent vert horiz height width
extentC
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubMatrix Int
mb Int
n Int
mb Ptr a
bPtr Int
mc Ptr a
cPtr
               ContT () IO () -> IO ()
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT () IO () -> IO ()) -> ContT () IO () -> IO ()
forall a b. (a -> b) -> a -> b
$ do
                  Ptr CChar
uploPtr <- Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'A'
                  Ptr CInt
mPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int
mcInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
mb)
                  Ptr CInt
nPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
                  Ptr CInt
ldcPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
mc
                  Ptr a
zPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
zero
                  IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$
                     Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
LapackGen.laset Ptr CChar
uploPtr Ptr CInt
mPtr Ptr CInt
nPtr Ptr a
zPtr Ptr a
zPtr
                        (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr Int
mb) Ptr CInt
ldcPtr


extractQ ::
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width,
    Class.Floating a) =>
   Householder vert horiz height width a -> Matrix.Square height a
extractQ :: Householder vert horiz height width a -> Square height a
extractQ
   (Householder tau (Array (MatrixShape.Split _ order extent) qr)) =
      Vector (Min height width) a
-> width
-> Order
-> Extent Small Small height height
-> ForeignPtr a
-> Square height a
forall vert horiz height width widthQR a.
(C vert, C horiz, C height, C width, C widthQR, Floating a) =>
Vector (Min height widthQR) a
-> widthQR
-> Order
-> Extent vert horiz height width
-> ForeignPtr a
-> Full vert horiz height width a
extractQAux Vector (Min height width) a
tau (Extent vert horiz height width -> width
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> width
Extent.width Extent vert horiz height width
extent) Order
order
         (height -> Extent Small Small height height
forall sh. sh -> Square sh
Extent.square (height -> Extent Small Small height height)
-> height -> Extent Small Small height height
forall a b. (a -> b) -> a -> b
$ Extent vert horiz height width -> height
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> height
Extent.height Extent vert horiz height width
extent) ForeignPtr a
qr

tallExtractQ ::
   (Extent.C vert, Shape.C height, Shape.C width, Class.Floating a) =>
   Householder vert Extent.Small height width a ->
   Full vert Extent.Small height width a
tallExtractQ :: Householder vert Small height width a
-> Full vert Small height width a
tallExtractQ
   (Householder tau (Array (MatrixShape.Split _ order extent) qr)) =
      Vector (Min height width) a
-> width
-> Order
-> Extent vert Small height width
-> ForeignPtr a
-> Full vert Small height width a
forall vert horiz height width widthQR a.
(C vert, C horiz, C height, C width, C widthQR, Floating a) =>
Vector (Min height widthQR) a
-> widthQR
-> Order
-> Extent vert horiz height width
-> ForeignPtr a
-> Full vert horiz height width a
extractQAux Vector (Min height width) a
tau (Extent vert Small height width -> width
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> width
Extent.width Extent vert Small height width
extent) Order
order Extent vert Small height width
extent ForeignPtr a
qr


extractQAux ::
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Shape.C widthQR,
    Class.Floating a) =>
   Vector (ExtShape.Min height widthQR) a -> widthQR ->
   Order -> Extent vert horiz height width -> ForeignPtr a ->
   Full vert horiz height width a
extractQAux :: Vector (Min height widthQR) a
-> widthQR
-> Order
-> Extent vert horiz height width
-> ForeignPtr a
-> Full vert horiz height width a
extractQAux (Array Min height widthQR
widthTau ForeignPtr a
tau) widthQR
widthQR Order
order Extent vert horiz height width
extent ForeignPtr a
qr =
   Full vert horiz height width
-> (Ptr a -> IO ()) -> Full vert horiz height width a
forall sh a.
(C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> Array sh a
Array.unsafeCreate (Order
-> Extent vert horiz height width -> Full vert horiz height width
forall vert horiz height width.
Order
-> Extent vert horiz height width -> Full vert horiz height width
MatrixShape.Full Order
order Extent vert horiz height width
extent) ((Ptr a -> IO ()) -> Full vert horiz height width a)
-> (Ptr a -> IO ()) -> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$ \Ptr a
qPtr -> do

   let (height
height,width
width) = Extent vert horiz height width -> (height, width)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vert horiz height width
extent
   let m :: Int
m = height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height
   let n :: Int
n = width -> Int
forall sh. C sh => sh -> Int
Shape.size width
width
   let k :: Int
k = Min height widthQR -> Int
forall sh. C sh => sh -> Int
Shape.size Min height widthQR
widthTau
   ContT () IO () -> IO ()
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT () IO () -> IO ()) -> ContT () IO () -> IO ()
forall a b. (a -> b) -> a -> b
$ do
      Ptr CInt
mPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
m
      Ptr CInt
nPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      Ptr CInt
kPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
k
      Ptr a
qrPtr <- ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
qr
      Ptr a
tauPtr <- ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
tau
      case Order
order of
         Order
RowMajor -> do
            Ptr CInt
ldaPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
n
            IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubMatrix Int
k Int
m (widthQR -> Int
forall sh. C sh => sh -> Int
Shape.size widthQR
widthQR) Ptr a
qrPtr Int
n Ptr a
qPtr
               String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"unglq" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
                  Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.unglq Ptr CInt
nPtr Ptr CInt
mPtr Ptr CInt
kPtr Ptr a
qPtr Ptr CInt
ldaPtr Ptr a
tauPtr
         Order
ColumnMajor -> do
            Ptr CInt
ldaPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
            IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
               Int -> Ptr a -> Ptr a -> IO ()
forall a. Floating a => Int -> Ptr a -> Ptr a -> IO ()
copyBlock (Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k) Ptr a
qrPtr Ptr a
qPtr
               String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"ungqr" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
                  Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.ungqr Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
kPtr Ptr a
qPtr Ptr CInt
ldaPtr Ptr a
tauPtr


tallMultiplyQ ::
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Shape.C fuse, Eq fuse,
    Class.Floating a) =>
   Householder vert Extent.Small height fuse a ->
   Full vert horiz fuse width a ->
   Full vert horiz height width a
tallMultiplyQ :: Householder vert Small height fuse a
-> Full vert horiz fuse width a -> Full vert horiz height width a
tallMultiplyQ Householder vert Small height fuse a
qr =
   Transposition
-> Conjugation
-> Householder vert Small height fuse a
-> Full vert horiz height width a
-> Full vert horiz height width a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
NonTransposed Conjugation
NonConjugated Householder vert Small height fuse a
qr (Full vert horiz height width a -> Full vert horiz height width a)
-> (Full vert horiz fuse width a -> Full vert horiz height width a)
-> Full vert horiz fuse width a
-> Full vert horiz height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Extent vert Small height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
forall vert horiz fuse height width a.
(C vert, C horiz, Eq fuse, C fuse, C height, C width,
 Floating a) =>
Extent vert Small height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
addRows (Householder vert Small height fuse a
-> Extent vert Small height fuse
forall vert horiz height width a.
Householder vert horiz height width a
-> Extent vert horiz height width
extent_ Householder vert Small height fuse a
qr)

tallMultiplyQAdjoint ::
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Shape.C fuse, Eq fuse, Class.Floating a) =>
   Householder horiz Extent.Small fuse height a ->
   Full vert horiz fuse width a ->
   Full vert horiz height width a
tallMultiplyQAdjoint :: Householder horiz Small fuse height a
-> Full vert horiz fuse width a -> Full vert horiz height width a
tallMultiplyQAdjoint Householder horiz Small fuse height a
qr =
   Extent Small horiz height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
forall vert horiz fuse height width a.
(C vert, C horiz, Eq fuse, C fuse, C height, C width,
 Floating a) =>
Extent Small horiz height fuse
-> Full vert horiz fuse width a -> Full vert horiz height width a
takeRows (Extent horiz Small fuse height -> Extent Small horiz height fuse
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> Extent horiz vert width height
Extent.transpose (Extent horiz Small fuse height -> Extent Small horiz height fuse)
-> Extent horiz Small fuse height -> Extent Small horiz height fuse
forall a b. (a -> b) -> a -> b
$ Householder horiz Small fuse height a
-> Extent horiz Small fuse height
forall vert horiz height width a.
Householder vert horiz height width a
-> Extent vert horiz height width
extent_ Householder horiz Small fuse height a
qr) (Full vert horiz fuse width a -> Full vert horiz height width a)
-> (Full vert horiz fuse width a -> Full vert horiz fuse width a)
-> Full vert horiz fuse width a
-> Full vert horiz height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   Transposition
-> Conjugation
-> Householder horiz Small fuse height a
-> Full vert horiz fuse width a
-> Full vert horiz fuse width a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
Transposed Conjugation
Conjugated Householder horiz Small fuse height a
qr


multiplyQ ::
   (Extent.C vertA, Extent.C horizA, Shape.C widthA,
    Extent.C vertB, Extent.C horizB, Shape.C widthB,
    Shape.C height, Eq height, Class.Floating a) =>
   Transposition -> Conjugation ->
   Householder vertA horizA height widthA a ->
   Full vertB horizB height widthB a ->
   Full vertB horizB height widthB a
multiplyQ :: Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
transposed Conjugation
conjugated
   (Householder
      (Array widthTau tau)
      (Array shapeA@(MatrixShape.Split _ orderA extentA) qr))
   (Array shapeB :: Full vertB horizB height widthB
shapeB@(MatrixShape.Full Order
orderB Extent vertB horizB height widthB
extentB) ForeignPtr a
b) =

   Full vertB horizB height widthB
-> (Int -> Ptr a -> IO ()) -> Full vertB horizB height widthB a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize Full vertB horizB height widthB
shapeB ((Int -> Ptr a -> IO ()) -> Full vertB horizB height widthB a)
-> (Int -> Ptr a -> IO ()) -> Full vertB horizB height widthB a
forall a b. (a -> b) -> a -> b
$ \Int
cSize Ptr a
cPtr -> do

   let (height
heightA,widthA
widthA) = Extent vertA horizA height widthA -> (height, widthA)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vertA horizA height widthA
extentA
   let (height
height,widthB
width) = Extent vertB horizB height widthB -> (height, widthB)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vertB horizB height widthB
extentB
   String -> Bool -> IO ()
Call.assert String
"Householder.multiplyQ: height shapes mismatch"
      (height
heightA height -> height -> Bool
forall a. Eq a => a -> a -> Bool
== height
height)

   let (Char
side,(Int
m,Int
n)) =
         Order -> (Int, Int) -> (Char, (Int, Int))
forall a. Order -> (a, a) -> (Char, (a, a))
sideSwapFromOrder Order
orderB (height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height, widthB -> Int
forall sh. C sh => sh -> Int
Shape.size widthB
width)

   ContT () IO () -> IO ()
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT () IO () -> IO ()) -> ContT () IO () -> IO ()
forall a b. (a -> b) -> a -> b
$ do
      Ptr CChar
sidePtr <- Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
side
      Ptr CInt
mPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
m
      Ptr CInt
nPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      let k :: Int
k = Min height widthA -> Int
forall sh. C sh => sh -> Int
Shape.size Min height widthA
widthTau
      Ptr CInt
kPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
k
      Ptr CChar
transPtr <-
         Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char (Char -> FortranIO () (Ptr CChar))
-> Char -> FortranIO () (Ptr CChar)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> Transposition -> Char
forall a (f :: * -> *). Floating a => f a -> Transposition -> Char
adjointFromTranspose ForeignPtr a
qr (Transposition -> Char) -> Transposition -> Char
forall a b. (a -> b) -> a -> b
$
         Transposition
transposed Transposition -> Transposition -> Transposition
forall a. Semigroup a => a -> a -> a
<> if Order
orderAOrder -> Order -> Bool
forall a. Eq a => a -> a -> Bool
==Order
orderB then Transposition
NonTransposed else Transposition
Transposed
      (Ptr a
qrPtr,Ptr a
tauPtr) <-
         if (Order
orderAOrder -> Order -> Bool
forall a. Eq a => a -> a -> Bool
==Order
orderB)
            Bool -> Bool -> Bool
forall a. Eq a => a -> a -> Bool
==
            (Transposition
transposedTransposition -> Transposition -> Bool
forall a. Eq a => a -> a -> Bool
==Transposition
NonTransposed Bool -> Bool -> Bool
&& Conjugation
conjugatedConjugation -> Conjugation -> Bool
forall a. Eq a => a -> a -> Bool
==Conjugation
NonConjugated
             Bool -> Bool -> Bool
||
             Transposition
transposedTransposition -> Transposition -> Bool
forall a. Eq a => a -> a -> Bool
==Transposition
Transposed Bool -> Bool -> Bool
&& Conjugation
conjugatedConjugation -> Conjugation -> Bool
forall a. Eq a => a -> a -> Bool
==Conjugation
Conjugated)
            then
               (Ptr a -> Ptr a -> (Ptr a, Ptr a))
-> ContT () IO (Ptr a)
-> ContT () IO (Ptr a)
-> ContT () IO (Ptr a, Ptr a)
forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 (,)
                  (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
qr)
                  (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
tau)
            else
               (Ptr a -> Ptr a -> (Ptr a, Ptr a))
-> ContT () IO (Ptr a)
-> ContT () IO (Ptr a)
-> ContT () IO (Ptr a, Ptr a)
forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 (,)
                  (Int -> ForeignPtr a -> ContT () IO (Ptr a)
forall a r. Floating a => Int -> ForeignPtr a -> ContT r IO (Ptr a)
conjugateToTemp (Split Reflector vertA horizA height widthA -> Int
forall sh. C sh => sh -> Int
Shape.size Split Reflector vertA horizA height widthA
shapeA) ForeignPtr a
qr)
                  (Int -> ForeignPtr a -> ContT () IO (Ptr a)
forall a r. Floating a => Int -> ForeignPtr a -> ContT r IO (Ptr a)
conjugateToTemp Int
k ForeignPtr a
tau)
      Ptr a
bPtr <- ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> ContT () IO (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
b
      Ptr CInt
ldcPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
      IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ Int -> Ptr a -> Ptr a -> IO ()
forall a. Floating a => Int -> Ptr a -> Ptr a -> IO ()
copyBlock Int
cSize Ptr a
bPtr Ptr a
cPtr
      case Order
orderA of
         Order
ColumnMajor -> do
            Ptr CInt
ldaPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim (Int -> FortranIO () (Ptr CInt)) -> Int -> FortranIO () (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ height -> Int
forall sh. C sh => sh -> Int
Shape.size height
heightA
            IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"unmqr" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
               Ptr CChar
-> Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.unmqr Ptr CChar
sidePtr Ptr CChar
transPtr
                  Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
kPtr Ptr a
qrPtr Ptr CInt
ldaPtr Ptr a
tauPtr Ptr a
cPtr Ptr CInt
ldcPtr
         Order
RowMajor -> do
            Ptr CInt
ldaPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim (Int -> FortranIO () (Ptr CInt)) -> Int -> FortranIO () (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ widthA -> Int
forall sh. C sh => sh -> Int
Shape.size widthA
widthA
            -- work-around for https://github.com/Reference-LAPACK/lapack/issues/260
            IO () -> ContT () IO ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ Bool -> IO () -> IO ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
when (Int
kInt -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>Int
0) (IO () -> IO ()) -> IO () -> IO ()
forall a b. (a -> b) -> a -> b
$
               String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a.
Floating a =>
String
-> String -> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
withAutoWorkspaceInfo String
errorCodeMsg String
"unmlq" ((Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ())
-> (Ptr a -> Ptr CInt -> Ptr CInt -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
               Ptr CChar
-> Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr CInt
-> IO ()
LapackGen.unmlq Ptr CChar
sidePtr Ptr CChar
transPtr
                  Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
kPtr Ptr a
qrPtr Ptr CInt
ldaPtr Ptr a
tauPtr Ptr a
cPtr Ptr CInt
ldcPtr

adjointFromTranspose :: (Class.Floating a) => f a -> Transposition -> Char
adjointFromTranspose :: f a -> Transposition -> Char
adjointFromTranspose f a
qr Transposition
Transposed = f a -> Char
forall a (f :: * -> *). Floating a => f a -> Char
invChar f a
qr
adjointFromTranspose f a
_ Transposition
NonTransposed = Char
'N'

invChar :: (Class.Floating a) => f a -> Char
invChar :: f a -> Char
invChar f a
f = f a -> Char -> Char -> Char
forall a (f :: * -> *) b. Floating a => f a -> b -> b -> b
caseRealComplexFunc f a
f Char
'T' Char
'C'


extractR ::
   (Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width,
    Class.Floating a) =>
   Householder vert horiz height width a ->
   Full vert horiz height width a
extractR :: Householder vert horiz height width a
-> Full vert horiz height width a
extractR = Either Reflector Triangle
-> Split Reflector vert horiz height width a
-> Full vert horiz height width a
forall vert horiz height width a lower.
(C vert, C horiz, C height, C width, Floating a) =>
Either lower Triangle
-> Split lower vert horiz height width a
-> Full vert horiz height width a
Split.extractTriangle (Triangle -> Either Reflector Triangle
forall a b. b -> Either a b
Right Triangle
MatrixShape.Triangle) (Split Reflector vert horiz height width a
 -> Full vert horiz height width a)
-> (Householder vert horiz height width a
    -> Split Reflector vert horiz height width a)
-> Householder vert horiz height width a
-> Full vert horiz height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vert horiz height width a
-> Split Reflector vert horiz height width a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

tallExtractR ::
   (Extent.C vert, Shape.C height, Shape.C width, Class.Floating a) =>
   Householder vert Extent.Small height width a -> Upper width a
tallExtractR :: Householder vert Small height width a -> Upper width a
tallExtractR = Split Reflector vert Small height width a -> Upper width a
forall vert height width a lower.
(C vert, C height, C width, Floating a) =>
Split lower vert Small height width a -> Upper width a
Split.tallExtractR (Split Reflector vert Small height width a -> Upper width a)
-> (Householder vert Small height width a
    -> Split Reflector vert Small height width a)
-> Householder vert Small height width a
-> Upper width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vert Small height width a
-> Split Reflector vert Small height width a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

tallMultiplyR ::
   (Extent.C vertA, Extent.C vert, Extent.C horiz, Shape.C height, Eq height,
    Shape.C heightA, Shape.C widthB, Class.Floating a) =>
   Transposition ->
   Householder vertA Extent.Small heightA height a ->
   Full vert horiz height widthB a ->
   Full vert horiz height widthB a
tallMultiplyR :: Transposition
-> Householder vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
tallMultiplyR Transposition
transposed = Transposition
-> Split Reflector vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
forall vertA vert horiz height heightA widthB a lower.
(C vertA, C vert, C horiz, C height, Eq height, C heightA,
 C widthB, Floating a) =>
Transposition
-> Split lower vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
Split.tallMultiplyR Transposition
transposed (Split Reflector vertA Small heightA height a
 -> Full vert horiz height widthB a
 -> Full vert horiz height widthB a)
-> (Householder vertA Small heightA height a
    -> Split Reflector vertA Small heightA height a)
-> Householder vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vertA Small heightA height a
-> Split Reflector vertA Small heightA height a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

tallSolveR ::
   (Extent.C vertA, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Eq width, Shape.C nrhs, Class.Floating a) =>
   Transposition -> Conjugation ->
   Householder vertA Extent.Small height width a ->
   Full vert horiz width nrhs a -> Full vert horiz width nrhs a
tallSolveR :: Transposition
-> Conjugation
-> Householder vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
tallSolveR Transposition
transposed Conjugation
conjugated =
   Transposition
-> Conjugation
-> Split Reflector vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
forall vertA vert horiz height width nrhs a lower.
(C vertA, C vert, C horiz, C height, C width, Eq width, C nrhs,
 Floating a) =>
Transposition
-> Conjugation
-> Split lower vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
Split.tallSolveR Transposition
transposed Conjugation
conjugated (Split Reflector vertA Small height width a
 -> Full vert horiz width nrhs a -> Full vert horiz width nrhs a)
-> (Householder vertA Small height width a
    -> Split Reflector vertA Small height width a)
-> Householder vertA Small height width a
-> Full vert horiz width nrhs a
-> Full vert horiz width nrhs a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder vertA Small height width a
-> Split Reflector vertA Small height width a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_


instance
   (Extent.C vert, Extent.C horiz) =>
      Type.Box (Hh vert horiz height width) where
   type HeightOf (Hh vert horiz height width) = height
   type WidthOf (Hh vert horiz height width) = width
   height :: Matrix (Hh vert horiz height width) a
-> HeightOf (Hh vert horiz height width)
height = Split Reflector vert horiz height width -> height
forall vert horiz lower height width.
(C vert, C horiz) =>
Split lower vert horiz height width -> height
MatrixShape.splitHeight (Split Reflector vert horiz height width -> height)
-> (Matrix (Hh vert horiz height width) a
    -> Split Reflector vert horiz height width)
-> Matrix (Hh vert horiz height width) a
-> height
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Array (Split Reflector vert horiz height width) a
-> Split Reflector vert horiz height width
forall sh a. Array sh a -> sh
Array.shape (Array (Split Reflector vert horiz height width) a
 -> Split Reflector vert horiz height width)
-> (Matrix (Hh vert horiz height width) a
    -> Array (Split Reflector vert horiz height width) a)
-> Matrix (Hh vert horiz height width) a
-> Split Reflector vert horiz height width
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_
   width :: Matrix (Hh vert horiz height width) a
-> WidthOf (Hh vert horiz height width)
width = Split Reflector vert horiz height width -> width
forall vert horiz lower height width.
(C vert, C horiz) =>
Split lower vert horiz height width -> width
MatrixShape.splitWidth (Split Reflector vert horiz height width -> width)
-> (Matrix (Hh vert horiz height width) a
    -> Split Reflector vert horiz height width)
-> Matrix (Hh vert horiz height width) a
-> width
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Array (Split Reflector vert horiz height width) a
-> Split Reflector vert horiz height width
forall sh a. Array sh a -> sh
Array.shape (Array (Split Reflector vert horiz height width) a
 -> Split Reflector vert horiz height width)
-> (Matrix (Hh vert horiz height width) a
    -> Array (Split Reflector vert horiz height width) a)
-> Matrix (Hh vert horiz height width) a
-> Split Reflector vert horiz height width
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
forall vert horiz height width a.
Matrix (Hh vert horiz height width) a
-> Array (Split Reflector vert horiz height width) a
split_

instance
   (Extent.C vert, Extent.C horiz,
    Shape.C height, Eq height, Shape.C width, Eq width) =>
      Multiply.MultiplyVector (Hh vert horiz height width) where
   matrixVector :: Matrix (Hh vert horiz height width) a
-> Vector width a
-> Vector (HeightOf (Hh vert horiz height width)) a
matrixVector Matrix (Hh vert horiz height width) a
qr Vector width a
x =
      Order
-> (General height () a -> General height () a)
-> Vector height a
-> Vector height a
forall height0 a height1 b.
Order
-> (General height0 () a -> General height1 () b)
-> Vector height0 a
-> Vector height1 b
Basic.unliftColumn Order
MatrixShape.ColumnMajor
         (Transposition
-> Conjugation
-> Matrix (Hh vert horiz height width) a
-> General height () a
-> General height () a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
NonTransposed Conjugation
NonConjugated Matrix (Hh vert horiz height width) a
qr) (Vector height a -> Vector height a)
-> Vector height a -> Vector height a
forall a b. (a -> b) -> a -> b
$
      Full vert horiz height width a -> Vector width a -> Vector height a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Eq width, Floating a) =>
Full vert horiz height width a -> Vector width a -> Vector height a
Basic.multiplyVector (Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Floating a) =>
Householder vert horiz height width a
-> Full vert horiz height width a
extractR Matrix (Hh vert horiz height width) a
qr) Vector width a
Vector width a
x
   vectorMatrix :: Vector height a
-> Matrix (Hh vert horiz height width) a
-> Vector (WidthOf (Hh vert horiz height width)) a
vectorMatrix Vector height a
x Matrix (Hh vert horiz height width) a
qr =
      Full horiz vert width height a -> Vector height a -> Vector width a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Eq width, Floating a) =>
Full vert horiz height width a -> Vector width a -> Vector height a
Basic.multiplyVector (Full vert horiz height width a -> Full horiz vert width height a
forall vert horiz height width a.
(C vert, C horiz) =>
Full vert horiz height width a -> Full horiz vert width height a
Basic.transpose (Full vert horiz height width a -> Full horiz vert width height a)
-> Full vert horiz height width a -> Full horiz vert width height a
forall a b. (a -> b) -> a -> b
$ Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Floating a) =>
Householder vert horiz height width a
-> Full vert horiz height width a
extractR Matrix (Hh vert horiz height width) a
qr) (Vector height a -> Vector width a)
-> Vector height a -> Vector width a
forall a b. (a -> b) -> a -> b
$
      Order
-> (General height () a -> General height () a)
-> Vector height a
-> Vector height a
forall height0 a height1 b.
Order
-> (General height0 () a -> General height1 () b)
-> Vector height0 a
-> Vector height1 b
Basic.unliftColumn Order
MatrixShape.ColumnMajor
         (Transposition
-> Conjugation
-> Matrix (Hh vert horiz height width) a
-> General height () a
-> General height () a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
Transposed Conjugation
NonConjugated Matrix (Hh vert horiz height width) a
qr) Vector height a
Vector height a
x

instance
   (vert ~ Extent.Small, horiz ~ Extent.Small,
    Shape.C height, height ~ width) =>
      Multiply.MultiplySquare (Hh vert horiz height width) where

   squareFull :: Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a -> Full vert horiz height width a
squareFull Matrix (Hh vert horiz height width) a
qr =
      (Array (Full vert horiz width width) a
 -> Array (Full vert horiz height width) a)
-> ArrayMatrix (Full vert horiz width width) a
-> ArrayMatrix (Full vert horiz height width) a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Full vert horiz width width) a
  -> Array (Full vert horiz height width) a)
 -> ArrayMatrix (Full vert horiz width width) a
 -> ArrayMatrix (Full vert horiz height width) a)
-> (Array (Full vert horiz width width) a
    -> Array (Full vert horiz height width) a)
-> ArrayMatrix (Full vert horiz width width) a
-> ArrayMatrix (Full vert horiz height width) a
forall a b. (a -> b) -> a -> b
$
         Transposition
-> Conjugation
-> Matrix (Hh vert horiz height width) a
-> Array (Full vert horiz height width) a
-> Array (Full vert horiz height width) a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
NonTransposed Conjugation
NonConjugated Matrix (Hh vert horiz height width) a
qr (Array (Full vert horiz height width) a
 -> Array (Full vert horiz height width) a)
-> (Array (Full vert horiz width width) a
    -> Array (Full vert horiz height width) a)
-> Array (Full vert horiz width width) a
-> Array (Full vert horiz height width) a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Transposition
-> Householder vert Small height width a
-> Array (Full vert horiz width width) a
-> Array (Full vert horiz width width) a
forall vertA vert horiz height heightA widthB a.
(C vertA, C vert, C horiz, C height, Eq height, C heightA,
 C widthB, Floating a) =>
Transposition
-> Householder vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
tallMultiplyR Transposition
NonTransposed Matrix (Hh vert horiz height width) a
Householder vert Small height width a
qr

   fullSquare :: Full vert horiz height width a
-> Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a
fullSquare = (Householder vert Small width width a
 -> Full vert horiz height width a
 -> Full vert horiz height width a)
-> Full vert horiz height width a
-> Householder vert Small width width a
-> Full vert horiz height width a
forall a b c. (a -> b -> c) -> b -> a -> c
flip ((Householder vert Small width width a
  -> Full vert horiz height width a
  -> Full vert horiz height width a)
 -> Full vert horiz height width a
 -> Householder vert Small width width a
 -> Full vert horiz height width a)
-> (Householder vert Small width width a
    -> Full vert horiz height width a
    -> Full vert horiz height width a)
-> Full vert horiz height width a
-> Householder vert Small width width a
-> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$ \Householder vert Small width width a
qr ->
      (Array (Full vert horiz height width) a
 -> Array (Full vert horiz height width) a)
-> Full vert horiz height width a -> Full vert horiz height width a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Full vert horiz height width) a
  -> Array (Full vert horiz height width) a)
 -> Full vert horiz height width a
 -> Full vert horiz height width a)
-> (Array (Full vert horiz height width) a
    -> Array (Full vert horiz height width) a)
-> Full vert horiz height width a
-> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$
         Full horiz vert width height a
-> Array (Full vert horiz height width) a
forall vert horiz height width a.
(C vert, C horiz) =>
Full vert horiz height width a -> Full horiz vert width height a
Basic.transpose (Full horiz vert width height a
 -> Array (Full vert horiz height width) a)
-> (Array (Full vert horiz height width) a
    -> Full horiz vert width height a)
-> Array (Full vert horiz height width) a
-> Array (Full vert horiz height width) a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Transposition
-> Householder vert Small width width a
-> Full horiz vert width height a
-> Full horiz vert width height a
forall vertA vert horiz height heightA widthB a.
(C vertA, C vert, C horiz, C height, Eq height, C heightA,
 C widthB, Floating a) =>
Transposition
-> Householder vertA Small heightA height a
-> Full vert horiz height widthB a
-> Full vert horiz height widthB a
tallMultiplyR Transposition
Transposed Householder vert Small width width a
qr (Full horiz vert width height a -> Full horiz vert width height a)
-> (Array (Full vert horiz height width) a
    -> Full horiz vert width height a)
-> Array (Full vert horiz height width) a
-> Full horiz vert width height a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Transposition
-> Conjugation
-> Householder vert Small width width a
-> Full horiz vert width height a
-> Full horiz vert width height a
forall vertA horizA widthA vertB horizB widthB height a.
(C vertA, C horizA, C widthA, C vertB, C horizB, C widthB,
 C height, Eq height, Floating a) =>
Transposition
-> Conjugation
-> Householder vertA horizA height widthA a
-> Full vertB horizB height widthB a
-> Full vertB horizB height widthB a
multiplyQ Transposition
Transposed Conjugation
NonConjugated Householder vert Small width width a
qr (Full horiz vert width height a -> Full horiz vert width height a)
-> (Array (Full vert horiz height width) a
    -> Full horiz vert width height a)
-> Array (Full vert horiz height width) a
-> Full horiz vert width height a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Array (Full vert horiz height width) a
-> Full horiz vert width height a
forall vert horiz height width a.
(C vert, C horiz) =>
Full vert horiz height width a -> Full horiz vert width height a
Basic.transpose

instance
   (vert ~ Extent.Small, horiz ~ Extent.Small,
    Shape.C height, height ~ width) =>
      Divide.Determinant (Hh vert horiz height width) where
   determinant :: Matrix (Hh vert horiz height width) a -> a
determinant = Matrix (Hh vert horiz height width) a -> a
forall sh a. (C sh, Floating a) => Square sh a -> a
determinant

instance
   (vert ~ Extent.Small, horiz ~ Extent.Small,
    Shape.C height, height ~ width) =>
      Divide.Solve (Hh vert horiz height width) where
   solveRight :: Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a -> Full vert horiz height width a
solveRight = (Array (Full vert horiz height width) a
 -> Array (Full vert horiz height width) a)
-> Full vert horiz height width a -> Full vert horiz height width a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Full vert horiz height width) a
  -> Array (Full vert horiz height width) a)
 -> Full vert horiz height width a
 -> Full vert horiz height width a)
-> (Householder Small Small height height a
    -> Array (Full vert horiz height width) a
    -> Array (Full vert horiz height width) a)
-> Householder Small Small height height a
-> Full vert horiz height width a
-> Full vert horiz height width a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Householder horiz Small height height a
-> Array (Full vert horiz height width) a
-> Array (Full vert horiz height width) a
forall vert horiz height width nrhs a.
(C vert, C horiz, C height, Eq height, C width, Eq width, C nrhs,
 Floating a) =>
Householder horiz Small height width a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
leastSquares (Householder horiz Small height height a
 -> Array (Full vert horiz height width) a
 -> Array (Full vert horiz height width) a)
-> (Householder Small Small height height a
    -> Householder horiz Small height height a)
-> Householder Small Small height height a
-> Array (Full vert horiz height width) a
-> Array (Full vert horiz height width) a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Map Small Small horiz Small height height
-> Householder Small Small height height a
-> Householder horiz Small height height a
forall vertA horizA vertB horizB height width a.
(C vertA, C horizA, C vertB, C horizB) =>
Map vertA horizA vertB horizB height width
-> Householder vertA horizA height width a
-> Householder vertB horizB height width a
mapExtent Map Small Small horiz Small height height
forall vert horiz size.
(C vert, C horiz) =>
Map Small Small vert horiz size size
Extent.fromSquare
   solveLeft :: Full vert horiz height width a
-> Matrix (Hh vert horiz height width) a
-> Full vert horiz height width a
solveLeft =
      (Householder Small Small width width a
 -> Full vert horiz height width a
 -> Full vert horiz height width a)
-> Full vert horiz height width a
-> Householder Small Small width width a
-> Full vert horiz height width a
forall a b c. (a -> b -> c) -> b -> a -> c
flip ((Householder Small Small width width a
  -> Full vert horiz height width a
  -> Full vert horiz height width a)
 -> Full vert horiz height width a
 -> Householder Small Small width width a
 -> Full vert horiz height width a)
-> (Householder Small Small width width a
    -> Full vert horiz height width a
    -> Full vert horiz height width a)
-> Full vert horiz height width a
-> Householder Small Small width width a
-> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$ \Householder Small Small width width a
a -> (Array (Full vert horiz height width) a
 -> Array (Full vert horiz height width) a)
-> Full vert horiz height width a -> Full vert horiz height width a
forall shA a shB b.
(Array shA a -> Array shB b)
-> ArrayMatrix shA a -> ArrayMatrix shB b
ArrMatrix.lift1 ((Array (Full vert horiz height width) a
  -> Array (Full vert horiz height width) a)
 -> Full vert horiz height width a
 -> Full vert horiz height width a)
-> (Array (Full vert horiz height width) a
    -> Array (Full vert horiz height width) a)
-> Full vert horiz height width a
-> Full vert horiz height width a
forall a b. (a -> b) -> a -> b
$
         Full horiz vert width height a
-> Array (Full vert horiz height width) a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Floating a) =>
Full vert horiz height width a -> Full horiz vert width height a
Basic.adjoint (Full horiz vert width height a
 -> Array (Full vert horiz height width) a)
-> (Array (Full vert horiz height width) a
    -> Full horiz vert width height a)
-> Array (Full vert horiz height width) a
-> Array (Full vert horiz height width) a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Householder horiz Small width width a
-> Full horiz vert width height a -> Full horiz vert width height a
forall vert horiz height width nrhs a.
(C vert, C horiz, C height, Eq height, C width, Eq width, C nrhs,
 Floating a) =>
Householder vert Small width height a
-> Full vert horiz height nrhs a -> Full vert horiz width nrhs a
minimumNorm (Map Small Small horiz Small width width
-> Householder Small Small width width a
-> Householder horiz Small width width a
forall vertA horizA vertB horizB height width a.
(C vertA, C horizA, C vertB, C horizB) =>
Map vertA horizA vertB horizB height width
-> Householder vertA horizA height width a
-> Householder vertB horizB height width a
mapExtent Map Small Small horiz Small width width
forall vert horiz size.
(C vert, C horiz) =>
Map Small Small vert horiz size size
Extent.fromSquare Householder Small Small width width a
a) (Full horiz vert width height a -> Full horiz vert width height a)
-> (Array (Full vert horiz height width) a
    -> Full horiz vert width height a)
-> Array (Full vert horiz height width) a
-> Full horiz vert width height a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
         Array (Full vert horiz height width) a
-> Full horiz vert width height a
forall vert horiz height width a.
(C vert, C horiz, C height, C width, Floating a) =>
Full vert horiz height width a -> Full horiz vert width height a
Basic.adjoint