{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
module Numeric.LAPACK.Matrix.Banded.Basic (
   Banded,
   General,
   Square,
   Upper,
   Lower,
   Diagonal,
   fromList,
   squareFromList,
   lowerFromList,
   upperFromList,
   mapExtent,
   mapHeight,
   mapWidth,
   identityFatOrder,
   diagonal,
   fromDiagonal,
   takeDiagonal,
   toFull,
   toLowerTriangular,
   toUpperTriangular,
   transpose,
   adjoint,
   multiplyVector,
   multiply,
   multiplyFull,
   ) where

import qualified Numeric.LAPACK.Matrix.Shape.Private as MatrixShape
import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent
import qualified Numeric.LAPACK.Matrix.Triangular.Private as TriangularPriv
import qualified Numeric.LAPACK.Matrix.Triangular.Basic as Triangular
import qualified Numeric.LAPACK.Matrix.RowMajor as RowMajor
import qualified Numeric.LAPACK.Matrix.Private as Matrix
import qualified Numeric.LAPACK.Vector as Vector
import qualified Numeric.LAPACK.Private as Private
import qualified Numeric.LAPACK.ShapeStatic as ShapeStatic
import Numeric.LAPACK.Matrix.Shape.Private
         (Order(RowMajor,ColumnMajor), transposeFromOrder, swapOnRowMajor,
          UnaryProxy, addOffDiagonals)
import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated))
import Numeric.LAPACK.Vector (Vector)
import Numeric.LAPACK.Scalar (zero, one)
import Numeric.LAPACK.Private
         (fill, pointerSeq, pokeCInt, copySubMatrix, copySubTrapezoid)

import qualified Numeric.BLAS.FFI.Generic as BlasGen
import qualified Numeric.Netlib.Utility as Call
import qualified Numeric.Netlib.Class as Class

import qualified Type.Data.Num.Unary.Literal as TypeNum
import qualified Type.Data.Num.Unary.Proof as Proof
import qualified Type.Data.Num.Unary as Unary
import Type.Data.Num.Unary ((:+:))
import Type.Data.Num (integralFromProxy)
import Type.Base.Proxy (Proxy(Proxy))

import qualified Data.Array.Comfort.Storable.Unchecked as Array
import qualified Data.Array.Comfort.Storable as CheckedArray
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 Foreign.Ptr (Ptr)
import Foreign.Storable (Storable)

import qualified Control.Monad.Trans.Maybe as MM
import qualified Control.Monad.Trans.Reader as MR
import Control.Monad.Trans.Cont (ContT(ContT), evalContT)
import Control.Monad.IO.Class (liftIO)
import Control.Monad (mzero, void)

import Data.Foldable (forM_)
import Data.Tuple.HT (swap)
import Data.Ord.HT (limit)


type Banded sub super vert horiz height width =
      Array (MatrixShape.Banded sub super vert horiz height width)

type General sub super height width =
      Array (MatrixShape.BandedGeneral sub super height width)

type Square sub super size =
      Array (MatrixShape.BandedSquare sub super size)

type Lower sub size = Square sub TypeNum.U0 size
type Upper super size = Square TypeNum.U0 super size

type Diagonal size = Square TypeNum.U0 TypeNum.U0 size


fromList ::
   (Unary.Natural sub, Unary.Natural super,
    Shape.C height, Shape.C width, Storable a) =>
   (UnaryProxy sub, UnaryProxy super) -> Order -> height -> width -> [a] ->
   General sub super height width a
fromList :: (UnaryProxy sub, UnaryProxy super)
-> Order
-> height
-> width
-> [a]
-> General sub super height width a
fromList (UnaryProxy sub, UnaryProxy super)
offDiag Order
order height
height width
width =
   (UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent Big Big height width
-> [a]
-> General sub super height width a
forall sub super vert horiz height width a.
(Natural sub, Natural super, C vert, C horiz, C height, C width,
 Storable a) =>
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> [a]
-> Banded sub super vert horiz height width a
fromListGen (UnaryProxy sub, UnaryProxy super)
offDiag Order
order (height -> width -> Extent Big Big height width
forall height width. height -> width -> General height width
Extent.general height
height width
width)

squareFromList ::
   (Unary.Natural sub, Unary.Natural super, Shape.C size, Storable a) =>
   (UnaryProxy sub, UnaryProxy super) -> Order -> size -> [a] ->
   Square sub super size a
squareFromList :: (UnaryProxy sub, UnaryProxy super)
-> Order -> size -> [a] -> Square sub super size a
squareFromList (UnaryProxy sub, UnaryProxy super)
offDiag Order
order size
size =
   (UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent Small Small size size
-> [a]
-> Square sub super size a
forall sub super vert horiz height width a.
(Natural sub, Natural super, C vert, C horiz, C height, C width,
 Storable a) =>
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> [a]
-> Banded sub super vert horiz height width a
fromListGen (UnaryProxy sub, UnaryProxy super)
offDiag Order
order (size -> Extent Small Small size size
forall sh. sh -> Square sh
Extent.square size
size)

lowerFromList ::
   (Unary.Natural sub, Shape.C size, Storable a) =>
   UnaryProxy sub -> Order -> size -> [a] -> Lower sub size a
lowerFromList :: UnaryProxy sub -> Order -> size -> [a] -> Lower sub size a
lowerFromList UnaryProxy sub
numOff Order
order size
size =
   (UnaryProxy sub, UnaryProxy U0)
-> Order -> Extent Small Small size size -> [a] -> Lower sub size a
forall sub super vert horiz height width a.
(Natural sub, Natural super, C vert, C horiz, C height, C width,
 Storable a) =>
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> [a]
-> Banded sub super vert horiz height width a
fromListGen (UnaryProxy sub
numOff,UnaryProxy U0
forall a. Proxy a
Proxy) Order
order (size -> Extent Small Small size size
forall sh. sh -> Square sh
Extent.square size
size)

upperFromList ::
   (Unary.Natural super, Shape.C size, Storable a) =>
   UnaryProxy super -> Order -> size -> [a] -> Upper super size a
upperFromList :: UnaryProxy super -> Order -> size -> [a] -> Upper super size a
upperFromList UnaryProxy super
numOff Order
order size
size =
   (UnaryProxy U0, UnaryProxy super)
-> Order
-> Extent Small Small size size
-> [a]
-> Upper super size a
forall sub super vert horiz height width a.
(Natural sub, Natural super, C vert, C horiz, C height, C width,
 Storable a) =>
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> [a]
-> Banded sub super vert horiz height width a
fromListGen (UnaryProxy U0
forall a. Proxy a
Proxy,UnaryProxy super
numOff) Order
order (size -> Extent Small Small size size
forall sh. sh -> Square sh
Extent.square size
size)

fromListGen ::
   (Unary.Natural sub, Unary.Natural super,
    Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Storable a) =>
   (UnaryProxy sub, UnaryProxy super) -> Order ->
   Extent.Extent vert horiz height width -> [a] ->
   Banded sub super vert horiz height width a
fromListGen :: (UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> [a]
-> Banded sub super vert horiz height width a
fromListGen (UnaryProxy sub, UnaryProxy super)
offDiag Order
order Extent vert horiz height width
extent =
   Banded sub super vert horiz height width
-> [a] -> Banded sub super vert horiz height width a
forall sh a. (C sh, Storable a) => sh -> [a] -> Array sh a
CheckedArray.fromList ((UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy sub, UnaryProxy super)
offDiag Order
order Extent vert horiz height width
extent)


mapExtent ::
   (Extent.C vertA, Extent.C horizA) =>
   (Extent.C vertB, Extent.C horizB) =>
   Extent.Map vertA horizA vertB horizB height width ->
   Banded super sub vertA horizA height width a ->
   Banded super sub vertB horizB height width a
mapExtent :: Map vertA horizA vertB horizB height width
-> Banded super sub vertA horizA height width a
-> Banded super sub vertB horizB height width a
mapExtent Map vertA horizA vertB horizB height width
f = (Banded super sub vertA horizA height width
 -> Banded super sub vertB horizB height width)
-> Banded super sub vertA horizA height width a
-> Banded super sub vertB horizB height width a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape ((Banded super sub vertA horizA height width
  -> Banded super sub vertB horizB height width)
 -> Banded super sub vertA horizA height width a
 -> Banded super sub vertB horizB height width a)
-> (Banded super sub vertA horizA height width
    -> Banded super sub vertB horizB height width)
-> Banded super sub vertA horizA height width a
-> Banded super sub vertB horizB height width a
forall a b. (a -> b) -> a -> b
$ Map vertA horizA vertB horizB height width
-> Banded super sub vertA horizA height width
-> Banded super sub vertB horizB height width
forall vertA horizA vertB horizB height width sub super.
Map vertA horizA vertB horizB height width
-> Banded sub super vertA horizA height width
-> Banded sub super vertB horizB height width
MatrixShape.bandedMapExtent Map vertA horizA vertB horizB height width
f

mapHeight ::
   (Extent.GeneralTallWide vert horiz,
    Extent.GeneralTallWide horiz vert) =>
   (heightA -> heightB) ->
   Banded super sub vert horiz heightA width a ->
   Banded super sub vert horiz heightB width a
mapHeight :: (heightA -> heightB)
-> Banded super sub vert horiz heightA width a
-> Banded super sub vert horiz heightB width a
mapHeight heightA -> heightB
f =
   (Banded super sub vert horiz heightA width
 -> Banded super sub vert horiz heightB width)
-> Banded super sub vert horiz heightA width a
-> Banded super sub vert horiz heightB width a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape
      (\(MatrixShape.Banded (UnaryProxy super, UnaryProxy sub)
offDiag Order
order Extent vert horiz heightA width
extent) ->
         (UnaryProxy super, UnaryProxy sub)
-> Order
-> Extent vert horiz heightB width
-> Banded super sub vert horiz heightB width
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy super, UnaryProxy sub)
offDiag Order
order (Extent vert horiz heightB width
 -> Banded super sub vert horiz heightB width)
-> Extent vert horiz heightB width
-> Banded super sub vert horiz heightB width
forall a b. (a -> b) -> a -> b
$ (heightA -> heightB)
-> Extent vert horiz heightA width
-> Extent vert horiz heightB width
forall vert horiz heightA heightB width.
GeneralTallWide vert horiz =>
(heightA -> heightB)
-> Extent vert horiz heightA width
-> Extent vert horiz heightB width
Extent.mapHeight heightA -> heightB
f Extent vert horiz heightA width
extent)

mapWidth ::
   (Extent.GeneralTallWide vert horiz,
    Extent.GeneralTallWide horiz vert) =>
   (widthA -> widthB) ->
   Banded super sub vert horiz height widthA a ->
   Banded super sub vert horiz height widthB a
mapWidth :: (widthA -> widthB)
-> Banded super sub vert horiz height widthA a
-> Banded super sub vert horiz height widthB a
mapWidth widthA -> widthB
f =
   (Banded super sub vert horiz height widthA
 -> Banded super sub vert horiz height widthB)
-> Banded super sub vert horiz height widthA a
-> Banded super sub vert horiz height widthB a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape
      (\(MatrixShape.Banded (UnaryProxy super, UnaryProxy sub)
offDiag Order
order Extent vert horiz height widthA
extent) ->
         (UnaryProxy super, UnaryProxy sub)
-> Order
-> Extent vert horiz height widthB
-> Banded super sub vert horiz height widthB
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy super, UnaryProxy sub)
offDiag Order
order (Extent vert horiz height widthB
 -> Banded super sub vert horiz height widthB)
-> Extent vert horiz height widthB
-> Banded super sub vert horiz height widthB
forall a b. (a -> b) -> a -> b
$ (widthA -> widthB)
-> Extent vert horiz height widthA
-> Extent vert horiz height widthB
forall vert horiz widthA widthB height.
GeneralTallWide vert horiz =>
(widthA -> widthB)
-> Extent vert horiz height widthA
-> Extent vert horiz height widthB
Extent.mapWidth widthA -> widthB
f Extent vert horiz height widthA
extent)

transpose ::
   (Extent.C vert, Extent.C horiz) =>
   Banded sub super vert horiz height width a ->
   Banded super sub horiz vert width height a
transpose :: Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
transpose = (Banded sub super vert horiz height width
 -> Banded super sub horiz vert width height)
-> Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape Banded sub super vert horiz height width
-> Banded super sub horiz vert width height
forall vert horiz sub super height width.
(C vert, C horiz) =>
Banded sub super vert horiz height width
-> Banded super sub horiz vert width height
MatrixShape.bandedTranspose

adjoint ::
   (Unary.Natural sub, Unary.Natural super, Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Class.Floating a) =>
   Banded sub super vert horiz height width a ->
   Banded super sub horiz vert width height a
adjoint :: Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
adjoint = Banded super sub horiz vert width height a
-> Banded super sub horiz vert width height a
forall sh a. (C sh, Floating a) => Vector sh a -> Vector sh a
Vector.conjugate (Banded super sub horiz vert width height a
 -> Banded super sub horiz vert width height a)
-> (Banded sub super vert horiz height width a
    -> Banded super sub horiz vert width height a)
-> Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
forall vert horiz sub super height width a.
(C vert, C horiz) =>
Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
transpose


identityFatOrder ::
   (Unary.Natural sub, Unary.Natural super, Shape.C size, Class.Floating a) =>
   Order -> size -> Square sub super size a
identityFatOrder :: Order -> size -> Square sub super size a
identityFatOrder Order
order =
   case Order
order of
      Order
RowMajor -> (size -> Extent Small Small size size)
-> Matrix size (Section sub super) a -> Square sub super size a
forall sub super height vert horiz width a.
(Natural sub, Natural super, C height) =>
(height -> Extent vert horiz height width)
-> Matrix height (Section sub super) a
-> Banded sub super vert horiz height width a
fromRowMajor size -> Extent Small Small size size
forall sh. sh -> Square sh
Extent.square (Matrix size (Section sub super) a -> Square sub super size a)
-> (size -> Matrix size (Section sub super) a)
-> size
-> Square sub super size a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. size -> Matrix size (Section sub super) a
forall sub super height a.
(Natural sub, Natural super, C height, Floating a) =>
height -> Matrix height (Section sub super) a
identityStripe
      Order
ColumnMajor -> (size -> Extent Small Small size size)
-> Matrix size (Section super sub) a -> Square sub super size a
forall sub super width vert horiz height a.
(Natural sub, Natural super, C width) =>
(width -> Extent vert horiz height width)
-> Matrix width (Section super sub) a
-> Banded sub super vert horiz height width a
fromColumnMajor size -> Extent Small Small size size
forall sh. sh -> Square sh
Extent.square (Matrix size (Section super sub) a -> Square sub super size a)
-> (size -> Matrix size (Section super sub) a)
-> size
-> Square sub super size a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. size -> Matrix size (Section super sub) a
forall sub super height a.
(Natural sub, Natural super, C height, Floating a) =>
height -> Matrix height (Section sub super) a
identityStripe

type Section sub super =
      ShapeStatic.ZeroBased sub Shape.:+: ()
         Shape.:+: ShapeStatic.ZeroBased super

identityStripe ::
   (Unary.Natural sub, Unary.Natural super, Shape.C height, Class.Floating a) =>
   height -> RowMajor.Matrix height (Section sub super) a
identityStripe :: height -> Matrix height (Section sub super) a
identityStripe height
sh =
   Either Conjugation Conjugation
-> Vector height a
-> Vector (Section sub super) a
-> Matrix height (Section sub super) a
forall height width a.
(C height, C width, Floating a) =>
Either Conjugation Conjugation
-> Vector height a -> Vector width a -> Matrix height width a
RowMajor.tensorProduct (Conjugation -> Either Conjugation Conjugation
forall a b. a -> Either a b
Left Conjugation
NonConjugated)
      (height -> Vector height a
forall sh a. (C sh, Floating a) => sh -> Vector sh a
Vector.one height
sh) (Section sub super
-> Index (Section sub super) -> Vector (Section sub super) a
forall sh a.
(Indexed sh, Floating a) =>
sh -> Index sh -> Vector sh a
Vector.unit Section sub super
forall sh. Static sh => sh
Shape.static (Index (Section sub super) -> Vector (Section sub super) a)
-> Index (Section sub super) -> Vector (Section sub super) a
forall a b. (a -> b) -> a -> b
$ Either () (Index super)
-> Either (Index sub) (Either () (Index super))
forall a b. b -> Either a b
Right (() -> Either () (Index super)
forall a b. a -> Either a b
Left ()))

fromRowMajor ::
   (Unary.Natural sub, Unary.Natural super, Shape.C height) =>
   (height -> Extent.Extent vert horiz height width) ->
   RowMajor.Matrix height (Section sub super) a ->
   Banded sub super vert horiz height width a
fromRowMajor :: (height -> Extent vert horiz height width)
-> Matrix height (Section sub super) a
-> Banded sub super vert horiz height width a
fromRowMajor height -> Extent vert horiz height width
extent =
   ((height, Section sub super)
 -> Banded sub super vert horiz height width)
-> Matrix height (Section sub super) a
-> Banded sub super vert horiz height width a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape
      (\(height
size,
         ShapeStatic.ZeroBased UnaryProxy sub
kl Shape.:+: ()
            Shape.:+: ShapeStatic.ZeroBased UnaryProxy super
ku) ->
         (UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy sub
kl,UnaryProxy super
ku) Order
RowMajor (Extent vert horiz height width
 -> Banded sub super vert horiz height width)
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
forall a b. (a -> b) -> a -> b
$ height -> Extent vert horiz height width
extent height
size)

fromColumnMajor ::
   (Unary.Natural sub, Unary.Natural super, Shape.C width) =>
   (width -> Extent.Extent vert horiz height width) ->
   RowMajor.Matrix width (Section super sub) a ->
   Banded sub super vert horiz height width a
fromColumnMajor :: (width -> Extent vert horiz height width)
-> Matrix width (Section super sub) a
-> Banded sub super vert horiz height width a
fromColumnMajor width -> Extent vert horiz height width
extent =
   ((width, Section super sub)
 -> Banded sub super vert horiz height width)
-> Matrix width (Section super sub) a
-> Banded sub super vert horiz height width a
forall sh0 sh1 a. (sh0 -> sh1) -> Array sh0 a -> Array sh1 a
Array.mapShape
      (\(width
size,
         ShapeStatic.ZeroBased UnaryProxy super
ku Shape.:+: ()
            Shape.:+: ShapeStatic.ZeroBased UnaryProxy sub
kl) ->
         (UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy sub
kl,UnaryProxy super
ku) Order
ColumnMajor (Extent vert horiz height width
 -> Banded sub super vert horiz height width)
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
forall a b. (a -> b) -> a -> b
$ width -> Extent vert horiz height width
extent width
size)


diagonal ::
   (Shape.C sh, Class.Floating a) => Order -> Vector sh a -> Diagonal sh a
diagonal :: Order -> Vector sh a -> Diagonal sh a
diagonal Order
order (Array sh
sh ForeignPtr a
x) =
   Banded U0 U0 Small Small sh sh -> ForeignPtr a -> Diagonal sh a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array ((UnaryProxy U0, UnaryProxy U0)
-> Order -> sh -> Banded U0 U0 Small Small sh sh
forall sub super size.
(UnaryProxy sub, UnaryProxy super)
-> Order -> size -> Banded sub super Small Small size size
MatrixShape.bandedSquare (UnaryProxy U0
forall a. Proxy a
Proxy,UnaryProxy U0
forall a. Proxy a
Proxy) Order
order sh
sh) ForeignPtr a
x

fromDiagonal ::
   (Shape.C sh, Class.Floating a) =>
   TriangularPriv.FlexDiagonal diag sh a -> Diagonal sh a
fromDiagonal :: FlexDiagonal diag sh a -> Diagonal sh a
fromDiagonal (Array (MatrixShape.Triangular diag
_diag (Empty, Empty)
_uplo Order
order sh
sh) ForeignPtr a
x) =
   Banded U0 U0 Small Small sh sh -> ForeignPtr a -> Diagonal sh a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array ((UnaryProxy U0, UnaryProxy U0)
-> Order -> sh -> Banded U0 U0 Small Small sh sh
forall sub super size.
(UnaryProxy sub, UnaryProxy super)
-> Order -> size -> Banded sub super Small Small size size
MatrixShape.bandedSquare (UnaryProxy U0
forall a. Proxy a
Proxy,UnaryProxy U0
forall a. Proxy a
Proxy) Order
order sh
sh) ForeignPtr a
x

takeDiagonal ::
   (Unary.Natural sub, Unary.Natural super, Shape.C sh, Class.Floating a) =>
   Square sub super sh a -> Vector sh a
takeDiagonal :: Square sub super sh a -> Vector sh a
takeDiagonal (Array (MatrixShape.Banded (UnaryProxy sub
sub,UnaryProxy super
super) Order
order Extent Small Small sh sh
extent) ForeignPtr a
x) =
   let size :: sh
size = Extent Small Small sh sh -> sh
forall height width. Extent Small Small height width -> height
Extent.squareSize Extent Small Small sh sh
extent
       kl :: Int
kl = UnaryProxy sub -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy sub
sub
       ku :: Int
ku = UnaryProxy super -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy super
super
   in if (Int
kl,Int
ku) (Int, Int) -> (Int, Int) -> Bool
forall a. Eq a => a -> a -> Bool
== (Int
0,Int
0)
        then sh -> ForeignPtr a -> Vector sh a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array sh
size ForeignPtr a
x
        else
            sh -> (Int -> Ptr a -> IO ()) -> Vector sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize sh
size ((Int -> Ptr a -> IO ()) -> Vector sh a)
-> (Int -> Ptr a -> IO ()) -> Vector sh a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
yPtr -> 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
nPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
               Ptr a
xPtr <- ((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
x
               let k :: Int
k =
                     case Order
order of
                        Order
RowMajor -> Int
kl
                        Order
ColumnMajor -> Int
ku
               Ptr CInt
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kuInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
               Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
               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 CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
BlasGen.copy Ptr CInt
nPtr (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
xPtr Int
k) Ptr CInt
incxPtr Ptr a
yPtr Ptr CInt
incyPtr


multiplyVector ::
   (Unary.Natural sub, Unary.Natural super,
    Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width, Eq width,
    Class.Floating a) =>
   Banded sub super vert horiz height width a ->
   Vector width a -> Vector height a
multiplyVector :: Banded sub super vert horiz height width a
-> Vector width a -> Vector height a
multiplyVector
   (Array (MatrixShape.Banded (UnaryProxy sub, UnaryProxy super)
numOff Order
order Extent vert horiz height width
extent) ForeignPtr a
a) (Array width
width ForeignPtr a
x) =
      let height :: height
height = 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
      in height -> (Ptr a -> IO ()) -> Vector height a
forall sh a.
(C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> Array sh a
Array.unsafeCreate height
height ((Ptr a -> IO ()) -> Vector height a)
-> (Ptr a -> IO ()) -> Vector height a
forall a b. (a -> b) -> a -> b
$ \Ptr a
yPtr -> do

   String -> Bool -> IO ()
Call.assert String
"Banded.multiplyVector: shapes mismatch"
      (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 width -> width -> Bool
forall a. Eq a => a -> a -> Bool
== width
width)
   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 -> (Int, Int))
-> Full vert horiz height width -> (Int, Int)
forall a b. (a -> b) -> 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
extent
   let (Int
kl,Int
ku) = Order -> (UnaryProxy sub, UnaryProxy super) -> (Int, Int)
forall sub super.
(Natural sub, Natural super) =>
Order -> (UnaryProxy sub, UnaryProxy super) -> (Int, Int)
MatrixShape.numOffDiagonals Order
order (UnaryProxy sub, UnaryProxy super)
numOff
   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
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
$ Order -> Char
transposeFromOrder Order
order
      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
klPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
kl
      Ptr CInt
kuPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
ku
      Ptr a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
one
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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 -> FortranIO () (Ptr CInt)) -> Int -> FortranIO () (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
      Ptr a
xPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
x
      Ptr CInt
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr a
betaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
zero
      Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      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 CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
Private.gbmv Ptr CChar
transPtr Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
klPtr Ptr CInt
kuPtr
            Ptr a
alphaPtr Ptr a
aPtr Ptr CInt
ldaPtr Ptr a
xPtr Ptr CInt
incxPtr Ptr a
betaPtr Ptr a
yPtr Ptr CInt
incyPtr


multiply ::
   (Unary.Natural subA, Unary.Natural superA,
    Unary.Natural subB, Unary.Natural superB,
    (subA :+: subB) ~ subC,
    (superA :+: superB) ~ superC,
    Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Shape.C fuse, Eq fuse,
    Class.Floating a) =>
   Banded subA superA vert horiz height fuse a ->
   Banded subB superB vert horiz fuse width a ->
   Banded subC superC vert horiz height width a
multiply :: Banded subA superA vert horiz height fuse a
-> Banded subB superB vert horiz fuse width a
-> Banded subC superC vert horiz height width a
multiply
      (Array (MatrixShape.Banded (UnaryProxy subA, UnaryProxy superA)
numOffA Order
orderA Extent vert horiz height fuse
extentA) ForeignPtr a
a)
      (Array (MatrixShape.Banded (UnaryProxy subB, UnaryProxy superB)
numOffB Order
orderB Extent vert horiz fuse width
extentB) ForeignPtr a
b) =
   case ((UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> ((Nat subC, Nat superC), (UnaryProxy subC, UnaryProxy superC))
forall subA superA subB superB subC superC.
(Natural subA, Natural superA, Natural subB, Natural superB,
 (subA :+: subB) ~ subC, (superA :+: superB) ~ superC) =>
(UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> ((Nat subC, Nat superC), (UnaryProxy subC, UnaryProxy superC))
addOffDiagonals (UnaryProxy subA, UnaryProxy superA)
numOffA (UnaryProxy subB, UnaryProxy superB)
numOffB, 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 horiz height fuse
extentA Extent vert horiz fuse width
extentB) of
      (((Nat subC, Nat superC), (UnaryProxy subC, UnaryProxy superC))
_, Maybe (Extent vert horiz height width)
Nothing) -> String -> Banded subC superC vert horiz height width a
forall a. HasCallStack => String -> a
error String
"Banded.multiply: shapes mismatch"
      (((Nat subC
Proof.Nat, Nat superC
Proof.Nat), (UnaryProxy subC, UnaryProxy superC)
numOffC), Just Extent vert horiz height width
extent) ->
         Banded subC superC vert horiz height width
-> (Ptr a -> IO ()) -> Banded subC superC vert horiz height width a
forall sh a.
(C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> Array sh a
Array.unsafeCreate
               ((UnaryProxy subC, UnaryProxy superC)
-> Order
-> Extent vert horiz height width
-> Banded subC superC vert horiz height width
forall sub super vert horiz height width.
(UnaryProxy sub, UnaryProxy super)
-> Order
-> Extent vert horiz height width
-> Banded sub super vert horiz height width
MatrixShape.Banded (UnaryProxy subC, UnaryProxy superC)
numOffC Order
orderB Extent vert horiz height width
extent) ((Ptr a -> IO ()) -> Banded subC superC vert horiz height width a)
-> (Ptr a -> IO ()) -> Banded subC superC vert horiz height width a
forall a b. (a -> b) -> a -> b
$ \Ptr a
cPtr ->
            let (height
height,fuse
fuse) = Extent vert horiz height fuse -> (height, fuse)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vert horiz height fuse
extentA
                width :: width
width = 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
            in case (Order
orderA,Order
orderB) of
                  (Order
ColumnMajor,Order
ColumnMajor) ->
                     Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall subA superA subB superB height width fuse a.
(Natural subA, Natural superA, Natural subB, Natural superB,
 C height, C width, C fuse, Floating a) =>
Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnMajor Order
ColumnMajor
                        (UnaryProxy subA, UnaryProxy superA)
numOffA (UnaryProxy subB, UnaryProxy superB)
numOffB (height
height,fuse
fuse,width
width) ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr
                  (Order
RowMajor,Order
ColumnMajor) ->
                     Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall subA superA subB superB height width fuse a.
(Natural subA, Natural superA, Natural subB, Natural superB,
 C height, C width, C fuse, Floating a) =>
Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnMajor Order
RowMajor
                        (UnaryProxy subA, UnaryProxy superA)
numOffA (UnaryProxy subB, UnaryProxy superB)
numOffB (height
height,fuse
fuse,width
width) ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr
                  (Order
ColumnMajor,Order
RowMajor) ->
                     (UnaryProxy superB, UnaryProxy subB)
-> (UnaryProxy superA, UnaryProxy subA)
-> (width, fuse, height)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall subA superA subB superB height width fuse a.
(Natural subA, Natural superA, Natural subB, Natural superB,
 C height, C width, C fuse, Floating a) =>
(UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnRowMajor
                        ((UnaryProxy subB, UnaryProxy superB)
-> (UnaryProxy superB, UnaryProxy subB)
forall a b. (a, b) -> (b, a)
swap (UnaryProxy subB, UnaryProxy superB)
numOffB) ((UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy superA, UnaryProxy subA)
forall a b. (a, b) -> (b, a)
swap (UnaryProxy subA, UnaryProxy superA)
numOffA)
                        (width
width,fuse
fuse,height
height) ForeignPtr a
b ForeignPtr a
a Ptr a
cPtr
                  (Order
RowMajor,Order
RowMajor) ->
                     Order
-> (UnaryProxy superB, UnaryProxy subB)
-> (UnaryProxy superA, UnaryProxy subA)
-> (width, fuse, height)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall subA superA subB superB height width fuse a.
(Natural subA, Natural superA, Natural subB, Natural superB,
 C height, C width, C fuse, Floating a) =>
Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnMajor Order
ColumnMajor
                        ((UnaryProxy subB, UnaryProxy superB)
-> (UnaryProxy superB, UnaryProxy subB)
forall a b. (a, b) -> (b, a)
swap (UnaryProxy subB, UnaryProxy superB)
numOffB) ((UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy superA, UnaryProxy subA)
forall a b. (a, b) -> (b, a)
swap (UnaryProxy subA, UnaryProxy superA)
numOffA)
                        (width
width,fuse
fuse,height
height) ForeignPtr a
b ForeignPtr a
a Ptr a
cPtr

multiplyColumnMajor ::
   (Unary.Natural subA, Unary.Natural superA,
    Unary.Natural subB, Unary.Natural superB,
    Shape.C height, Shape.C width, Shape.C fuse,
    Class.Floating a) =>
   Order ->
   (UnaryProxy subA, UnaryProxy superA) ->
   (UnaryProxy subB, UnaryProxy superB) ->
   (height, fuse, width) ->
   ForeignPtr a -> ForeignPtr a -> Ptr a -> IO ()
multiplyColumnMajor :: Order
-> (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnMajor Order
orderA (UnaryProxy subA
subA,UnaryProxy superA
superA) (UnaryProxy subB
subB,UnaryProxy superB
superB)
      (height
height,fuse
fuse,width
width) ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr = do
   let m :: Int
m = height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height
   let k :: Int
k = fuse -> Int
forall sh. C sh => sh -> Int
Shape.size fuse
fuse
   let n :: Int
n = width -> Int
forall sh. C sh => sh -> Int
Shape.size width
width
   let (Int
kla,Int
kua) = (UnaryProxy subA -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy subA
subA, UnaryProxy superA -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy superA
superA)
   let (Int
klb,Int
kub) = (UnaryProxy subB -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy subB
subB, UnaryProxy superB -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy superB
superB)
   let ku :: Int
ku = Int
kuaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kub
   let kl :: Int
kl = Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klb
   let lda0 :: Int
lda0 = Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kua
   let ldb0 :: Int
ldb0 = Int
klbInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kub
   let ldc0 :: Int
ldc0 = Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ldb0
   let lda :: Int
lda = Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
   let ldc :: Int
ldc = Int
ldc0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
   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
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
$ Order -> Char
transposeFromOrder Order
orderA
      Ptr CInt
mPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
nPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
klPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
kuPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      let ((Ptr CInt
miPtr,Ptr CInt
kliPtr),(Ptr CInt
niPtr,Ptr CInt
kuiPtr)) =
            Order
-> ((Ptr CInt, Ptr CInt), (Ptr CInt, Ptr CInt))
-> ((Ptr CInt, Ptr CInt), (Ptr CInt, Ptr CInt))
forall a. Order -> (a, a) -> (a, a)
swapOnRowMajor Order
orderA ((Ptr CInt
mPtr,Ptr CInt
klPtr),(Ptr CInt
nPtr,Ptr CInt
kuPtr))
      Ptr a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
one
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
lda
      Ptr a
bPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr a
betaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
zero
      Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      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] -> (Int -> IO ()) -> IO ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ (Int -> [Int] -> [Int]
forall a. Int -> [a] -> [a]
take Int
n [Int
0..]) ((Int -> IO ()) -> IO ()) -> (Int -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
            let top :: Int
top = Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
ku)
            let bottom :: Int
bottom = Int -> Int -> Int
forall a. Ord a => a -> a -> a
min Int
m (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
            let left :: Int
left = Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
kub)
            let right :: Int
right = Int -> Int -> Int
forall a. Ord a => a -> a -> a
min Int
k (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klbInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
miPtr (Int -> IO ()) -> Int -> IO ()
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
bottomInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
top
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
niPtr (Int -> IO ()) -> Int -> IO ()
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
rightInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
left
            let d :: Int
d = Int
topInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
left; kli :: Int
kli = Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
d; kui :: Int
kui = Int
kuaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
d
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
kuiPtr Int
kui
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
kliPtr Int
kli
            let j0 :: Int
j0 = Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldc
            let j1 :: Int
j1 = Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldc0 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
topInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
            let j2 :: Int
j2 = Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldc0 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
bottomInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
            a -> Int -> Ptr a -> IO ()
forall a. Floating a => a -> Int -> Ptr a -> IO ()
fill a
forall a. Floating a => a
zero (Int
j1Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
j0) (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr Int
j0)
            let aOffset :: Int
aOffset =
                  case Order
orderA of
                     Order
ColumnMajor -> Int
left
                     Order
RowMajor -> Int
top
            Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
Private.gbmv Ptr CChar
transPtr Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
klPtr Ptr CInt
kuPtr
               Ptr a
alphaPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
aOffsetInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda)) Ptr CInt
ldaPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldb0 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
leftInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kub)) Ptr CInt
incxPtr
               Ptr a
betaPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr Int
j1) Ptr CInt
incyPtr
            a -> Int -> Ptr a -> IO ()
forall a. Floating a => a -> Int -> Ptr a -> IO ()
fill a
forall a. Floating a => a
zero (Int
j0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ldcInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
j2) (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr Int
j2)

multiplyColumnRowMajor ::
   (Unary.Natural subA, Unary.Natural superA,
    Unary.Natural subB, Unary.Natural superB,
    Shape.C height, Shape.C width, Shape.C fuse,
    Class.Floating a) =>
   (UnaryProxy subA, UnaryProxy superA) ->
   (UnaryProxy subB, UnaryProxy superB) ->
   (height, fuse, width) ->
   ForeignPtr a -> ForeignPtr a -> Ptr a -> IO ()
multiplyColumnRowMajor :: (UnaryProxy subA, UnaryProxy superA)
-> (UnaryProxy subB, UnaryProxy superB)
-> (height, fuse, width)
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyColumnRowMajor (UnaryProxy subA
subA,UnaryProxy superA
superA) (UnaryProxy subB
subB,UnaryProxy superB
superB)
      (height
height,fuse
fuse,width
width) ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr = do
   let m :: Int
m = height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height
   let k :: Int
k = fuse -> Int
forall sh. C sh => sh -> Int
Shape.size fuse
fuse
   let n :: Int
n = width -> Int
forall sh. C sh => sh -> Int
Shape.size width
width
   let (Int
kla,Int
kua) = (UnaryProxy subA -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy subA
subA, UnaryProxy superA -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy superA
superA)
   let (Int
klb,Int
kub) = (UnaryProxy subB -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy subB
subB, UnaryProxy superB -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy superB
superB)
   let ku :: Int
ku = Int
kuaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kub
   let kl :: Int
kl = Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klb
   let lda0 :: Int
lda0 = Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kua
   let ldb0 :: Int
ldb0 = Int
klbInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kub
   let ldc0 :: Int
ldc0 = Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
   let ldc :: Int
ldc = Int
ldc0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
   a -> Int -> Ptr a -> IO ()
forall a. Floating a => a -> Int -> Ptr a -> IO ()
fill a
forall a. Floating a => a
zero (Int
ldcInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
n) 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 CInt
mPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
nPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
one
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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 a
bPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr CInt
ldc0Ptr <- 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
$ Int
ldc0 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ if Int
ldb0Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
==Int
0 then Int
1 else Int
0
      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] -> (Int -> IO ()) -> IO ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ (Int -> [Int] -> [Int]
forall a. Int -> [a] -> [a]
take Int
k [Int
0..]) ((Int -> IO ()) -> IO ()) -> (Int -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
            let top :: Int
top = Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
kua)
            let bottom :: Int
bottom = Int -> Int -> Int
forall a. Ord a => a -> a -> a
min Int
m (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klaInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
            let left :: Int
left = Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
klb)
            let right :: Int
right = Int -> Int -> Int
forall a. Ord a => a -> a -> a
min Int
n (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kubInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
mPtr (Int -> IO ()) -> Int -> IO ()
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
bottomInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
top
            Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
nPtr (Int -> IO ()) -> Int -> IO ()
forall a b. (a -> b) -> a -> b
$ Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
rightInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
left
            Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> IO ()
BlasGen.geru Ptr CInt
mPtr Ptr CInt
nPtr Ptr a
alphaPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
topInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kua)) Ptr CInt
incxPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldb0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
leftInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
klb)) Ptr CInt
incyPtr
               (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
cPtr (Int
leftInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldc0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
topInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku)) Ptr CInt
ldc0Ptr


multiplyFull ::
   (Unary.Natural sub, Unary.Natural super,
    Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Shape.C fuse, Eq fuse,
    Class.Floating a) =>
   Banded sub super vert horiz height fuse a ->
   Matrix.Full vert horiz fuse width a -> Matrix.Full vert horiz height width a
multiplyFull :: Banded sub super vert horiz height fuse a
-> Full vert horiz fuse width a -> Full vert horiz height width a
multiplyFull
      (Array (MatrixShape.Banded (UnaryProxy sub, UnaryProxy super)
numOff Order
orderA Extent vert horiz height fuse
extentA) ForeignPtr a
a)
      (Array (MatrixShape.Full Order
orderB 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 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
"Banded.multiplyFull: shapes mismatch"
      Just Extent vert horiz height width
extent ->
         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
orderB 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
cPtr ->
            let (height
height,fuse
fuse) = Extent vert horiz height fuse -> (height, fuse)
forall vert horiz height width.
(C vert, C horiz) =>
Extent vert horiz height width -> (height, width)
Extent.dimensions Extent vert horiz height fuse
extentA
                width :: width
width = 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
            in case Order
orderB of
                  Order
ColumnMajor ->
                     (UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> Extent vert horiz height fuse
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall sub super vert horiz height width fuse a.
(Natural sub, Natural super, C vert, C horiz, C height, C width,
 C fuse, Floating a) =>
(UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> Extent vert horiz height fuse
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyFullColumnMajor
                        (UnaryProxy sub, UnaryProxy super)
numOff (height
height,fuse
fuse,width
width) Order
orderA Extent vert horiz height fuse
extentA ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr
                  Order
RowMajor ->
                     (UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
forall sub super height width fuse a.
(Natural sub, Natural super, C height, C width, C fuse,
 Floating a) =>
(UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyFullRowMajor
                        (UnaryProxy sub, UnaryProxy super)
numOff (height
height,fuse
fuse,width
width) Order
orderA ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr

multiplyFullColumnMajor ::
   (Unary.Natural sub, Unary.Natural super,
    Extent.C vert, Extent.C horiz,
    Shape.C height, Shape.C width, Shape.C fuse,
    Class.Floating a) =>
   (UnaryProxy sub, UnaryProxy super) ->
   (height, fuse, width) ->
   Order -> Extent.Extent vert horiz height fuse ->
   ForeignPtr a -> ForeignPtr a -> Ptr a -> IO ()
multiplyFullColumnMajor :: (UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> Extent vert horiz height fuse
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyFullColumnMajor (UnaryProxy sub, UnaryProxy super)
numOff (height
height,fuse
fuse,width
width) Order
orderA Extent vert horiz height fuse
extentA ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr = do
   let (Int
m,Int
n) = Full vert horiz height fuse -> (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 fuse -> (Int, Int))
-> Full vert horiz height fuse -> (Int, Int)
forall a b. (a -> b) -> a -> b
$ Order
-> Extent vert horiz height fuse -> Full vert horiz height fuse
forall vert horiz height width.
Order
-> Extent vert horiz height width -> Full vert horiz height width
MatrixShape.Full Order
orderA Extent vert horiz height fuse
extentA
   let k :: Int
k = width -> Int
forall sh. C sh => sh -> Int
Shape.size width
width
   let (Int
kl,Int
ku) = Order -> (UnaryProxy sub, UnaryProxy super) -> (Int, Int)
forall sub super.
(Natural sub, Natural super) =>
Order -> (UnaryProxy sub, UnaryProxy super) -> (Int, Int)
MatrixShape.numOffDiagonals Order
orderA (UnaryProxy sub, UnaryProxy super)
numOff
   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
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
$ Order -> Char
transposeFromOrder Order
orderA
      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
klPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
kl
      Ptr CInt
kuPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
ku
      Ptr a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
one
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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 -> FortranIO () (Ptr CInt)) -> Int -> FortranIO () (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
      Ptr a
bPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr a
betaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
zero
      Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      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 a, Ptr a)] -> ((Ptr a, Ptr a) -> IO ()) -> IO ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ (Int -> [(Ptr a, Ptr a)] -> [(Ptr a, Ptr a)]
forall a. Int -> [a] -> [a]
take Int
k ([(Ptr a, Ptr a)] -> [(Ptr a, Ptr a)])
-> [(Ptr a, Ptr a)] -> [(Ptr a, Ptr a)]
forall a b. (a -> b) -> a -> b
$
                [Ptr a] -> [Ptr a] -> [(Ptr a, Ptr a)]
forall a b. [a] -> [b] -> [(a, b)]
zip (Int -> Ptr a -> [Ptr a]
forall a. Storable a => Int -> Ptr a -> [Ptr a]
pointerSeq (fuse -> Int
forall sh. C sh => sh -> Int
Shape.size fuse
fuse) Ptr a
bPtr)
                    (Int -> Ptr a -> [Ptr a]
forall a. Storable a => Int -> Ptr a -> [Ptr a]
pointerSeq (height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height) Ptr a
cPtr)) (((Ptr a, Ptr a) -> IO ()) -> IO ())
-> ((Ptr a, Ptr a) -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
            \(Ptr a
xPtr,Ptr a
yPtr) ->
               Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
Private.gbmv Ptr CChar
transPtr Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
klPtr Ptr CInt
kuPtr
                  Ptr a
alphaPtr Ptr a
aPtr Ptr CInt
ldaPtr Ptr a
xPtr Ptr CInt
incxPtr
                  Ptr a
betaPtr Ptr a
yPtr Ptr CInt
incyPtr

multiplyFullRowMajor ::
   (Unary.Natural sub, Unary.Natural super,
    Shape.C height, Shape.C width, Shape.C fuse,
    Class.Floating a) =>
   (UnaryProxy sub, UnaryProxy super) ->
   (height, fuse, width) ->
   Order -> ForeignPtr a -> ForeignPtr a -> Ptr a -> IO ()
multiplyFullRowMajor :: (UnaryProxy sub, UnaryProxy super)
-> (height, fuse, width)
-> Order
-> ForeignPtr a
-> ForeignPtr a
-> Ptr a
-> IO ()
multiplyFullRowMajor (UnaryProxy sub
sub,UnaryProxy super
super) (height
height,fuse
fuse,width
width) Order
orderA ForeignPtr a
a ForeignPtr a
b Ptr a
cPtr = do
   let m :: Int
m = height -> Int
forall sh. C sh => sh -> Int
Shape.size height
height
   let n :: Int
n = fuse -> Int
forall sh. C sh => sh -> Int
Shape.size fuse
fuse
   let k :: Int
k = width -> Int
forall sh. C sh => sh -> Int
Shape.size width
width
   let kl :: Int
kl = UnaryProxy sub -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy sub
sub
   let ku :: Int
ku = UnaryProxy super -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy super
super
   let lda0 :: Int
lda0 = Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
   let lda :: Int
lda = Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
   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
transPtr <- Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'N'
      Ptr CInt
kPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
k
      Ptr CInt
dPtr <- FortranIO () (Ptr CInt)
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
one
      Ptr a
aPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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 a
bPtr <- ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a)
forall k (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO ()) -> IO ()) -> FortranIO () (Ptr a))
-> ((Ptr a -> IO ()) -> IO ()) -> FortranIO () (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
ldbPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
k
      Ptr CInt
incxPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> FortranIO () (Ptr CInt)) -> Int -> FortranIO () (Ptr CInt)
forall a b. (a -> b) -> a -> b
$
         case Order
orderA of
            Order
RowMajor -> Int
1
            Order
ColumnMajor -> Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
1 Int
lda0
      Ptr a
betaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
zero
      Ptr CInt
incyPtr <- Int -> FortranIO () (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      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))]
-> ((Int, (Ptr a, Ptr a)) -> IO ()) -> IO ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ (Int -> [(Int, (Ptr a, Ptr a))] -> [(Int, (Ptr a, Ptr a))]
forall a. Int -> [a] -> [a]
take Int
m ([(Int, (Ptr a, Ptr a))] -> [(Int, (Ptr a, Ptr a))])
-> [(Int, (Ptr a, Ptr a))] -> [(Int, (Ptr a, Ptr a))]
forall a b. (a -> b) -> a -> b
$ [Int] -> [(Ptr a, Ptr a)] -> [(Int, (Ptr a, Ptr a))]
forall a b. [a] -> [b] -> [(a, b)]
zip [Int
0..] ([(Ptr a, Ptr a)] -> [(Int, (Ptr a, Ptr a))])
-> [(Ptr a, Ptr a)] -> [(Int, (Ptr a, Ptr a))]
forall a b. (a -> b) -> a -> b
$
                [Ptr a] -> [Ptr a] -> [(Ptr a, Ptr a)]
forall a b. [a] -> [b] -> [(a, b)]
zip (Int -> Ptr a -> [Ptr a]
forall a. Storable a => Int -> Ptr a -> [Ptr a]
pointerSeq Int
lda Ptr a
aPtr) (Int -> Ptr a -> [Ptr a]
forall a. Storable a => Int -> Ptr a -> [Ptr a]
pointerSeq Int
k Ptr a
cPtr)) (((Int, (Ptr a, Ptr a)) -> IO ()) -> IO ())
-> ((Int, (Ptr a, Ptr a)) -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$
            \(Int
i,(Ptr a
xPtr,Ptr a
yPtr)) -> do
               let firstRow :: Int
firstRow = (Int, Int) -> Int -> Int
forall a. Ord a => (a, a) -> a -> a
limit (Int
0,Int
n) (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
kl)
               let last1Row :: Int
last1Row = (Int, Int) -> Int -> Int
forall a. Ord a => (a, a) -> a -> a
limit (Int
0,Int
n) (Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kuInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
               let biPtr :: Ptr a
biPtr = Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
firstRowInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k)
               let xOffset :: Int
xOffset =
                     case Order
orderA of
                        Order
RowMajor -> Int
firstRowInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kl
                        Order
ColumnMajor -> (Int
firstRowInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
i)Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
               let xiPtr :: Ptr a
xiPtr = Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
xPtr Int
xOffset
               Ptr CInt -> Int -> IO ()
pokeCInt Ptr CInt
dPtr (Int -> IO ()) -> Int -> IO ()
forall a b. (a -> b) -> a -> b
$ Int
last1Row Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
firstRow
               Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
Private.gemv Ptr CChar
transPtr Ptr CInt
kPtr Ptr CInt
dPtr
                  Ptr a
alphaPtr Ptr a
biPtr Ptr CInt
ldbPtr Ptr a
xiPtr Ptr CInt
incxPtr
                  Ptr a
betaPtr Ptr a
yPtr Ptr CInt
incyPtr


toLowerTriangular ::
   (Unary.Natural sub, Shape.C sh, Class.Floating a) =>
   Lower sub sh a -> Triangular.Lower sh a
toLowerTriangular :: Lower sub sh a -> Lower sh a
toLowerTriangular =
   Triangular Empty NonUnit Filled sh a -> Lower sh a
forall lo up diag sh a.
(Content lo, Content up, TriDiag diag) =>
Triangular lo diag up sh a -> Triangular up diag lo sh a
Triangular.transpose (Triangular Empty NonUnit Filled sh a -> Lower sh a)
-> (Lower sub sh a -> Triangular Empty NonUnit Filled sh a)
-> Lower sub sh a
-> Lower sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Upper sub sh a -> Triangular Empty NonUnit Filled sh a
forall super sh a.
(Natural super, C sh, Floating a) =>
Upper super sh a -> Upper sh a
toUpperTriangular (Upper sub sh a -> Triangular Empty NonUnit Filled sh a)
-> (Lower sub sh a -> Upper sub sh a)
-> Lower sub sh a
-> Triangular Empty NonUnit Filled sh a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Lower sub sh a -> Upper sub sh a
forall vert horiz sub super height width a.
(C vert, C horiz) =>
Banded sub super vert horiz height width a
-> Banded super sub horiz vert width height a
transpose

toUpperTriangular ::
   (Unary.Natural super, Shape.C sh, Class.Floating a) =>
   Upper super sh a -> Triangular.Upper sh a
toUpperTriangular :: Upper super sh a -> Upper sh a
toUpperTriangular (Array (MatrixShape.Banded (UnaryProxy U0
_sub,UnaryProxy super
super) Order
order Extent Small Small sh sh
extent) ForeignPtr a
a) =
   let size :: sh
size = Extent Small Small sh sh -> sh
forall height width. Extent Small Small height width -> height
Extent.squareSize Extent Small Small sh sh
extent
   in Triangular Empty NonUnit Filled sh
-> (Int -> Ptr a -> IO ()) -> Upper sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize
         (NonUnit
-> (Empty, Filled)
-> Order
-> sh
-> Triangular Empty NonUnit Filled sh
forall lo diag up size.
diag -> (lo, up) -> Order -> size -> Triangular lo diag up size
MatrixShape.Triangular NonUnit
MatrixShape.NonUnit (Empty, Filled)
MatrixShape.upper
            Order
order sh
size) ((Int -> Ptr a -> IO ()) -> Upper sh a)
-> (Int -> Ptr a -> IO ()) -> Upper sh a
forall a b. (a -> b) -> a -> b
$
      Int -> Order -> Int -> ForeignPtr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Order -> Int -> ForeignPtr a -> Int -> Ptr a -> IO ()
TriangularPriv.fromBanded
         (UnaryProxy super -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy super
super) Order
order (sh -> Int
forall sh. C sh => sh -> Int
Shape.size sh
size) ForeignPtr a
a

toFull ::
   (Unary.Natural sub, Unary.Natural super,
    Extent.C vert, Extent.C horiz, Shape.C height, Shape.C width,
    Class.Floating a) =>
   Banded sub super vert horiz height width a ->
   Matrix.Full vert horiz height width a
toFull :: Banded sub super vert horiz height width a
-> Full vert horiz height width a
toFull (Array (MatrixShape.Banded (UnaryProxy sub
sub,UnaryProxy super
super) Order
order Extent vert horiz height width
extent) ForeignPtr a
a) =
   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
extent) ((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
bSize Ptr a
bPtr ->
   ForeignPtr a -> (Ptr a -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
a ((Ptr a -> IO ()) -> IO ()) -> (Ptr a -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \Ptr a
aPtr -> 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
      a -> Int -> Ptr a -> IO ()
forall a. Floating a => a -> Int -> Ptr a -> IO ()
fill a
forall a. Floating a => a
zero Int
bSize Ptr a
bPtr
      case Order
order of
         Order
ColumnMajor -> (UnaryProxy sub, UnaryProxy super)
-> (height, width) -> Ptr a -> Ptr a -> IO ()
forall sub super height width a.
(Natural sub, Natural super, C height, C width, Floating a) =>
(UnaryProxy sub, UnaryProxy super)
-> (height, width) -> Ptr a -> Ptr a -> IO ()
toFullColumnMajor (UnaryProxy sub
sub,UnaryProxy super
super) (height
height,width
width) Ptr a
aPtr Ptr a
bPtr
         Order
RowMajor -> (UnaryProxy super, UnaryProxy sub)
-> (width, height) -> Ptr a -> Ptr a -> IO ()
forall sub super height width a.
(Natural sub, Natural super, C height, C width, Floating a) =>
(UnaryProxy sub, UnaryProxy super)
-> (height, width) -> Ptr a -> Ptr a -> IO ()
toFullColumnMajor (UnaryProxy super
super,UnaryProxy sub
sub) (width
width,height
height) Ptr a
aPtr Ptr a
bPtr

toFullColumnMajor ::
   (Unary.Natural sub, Unary.Natural super, Shape.C height, Shape.C width,
    Class.Floating a) =>
   (UnaryProxy sub, UnaryProxy super) -> (height,width) ->
   Ptr a -> Ptr a -> IO ()
toFullColumnMajor :: (UnaryProxy sub, UnaryProxy super)
-> (height, width) -> Ptr a -> Ptr a -> IO ()
toFullColumnMajor (UnaryProxy sub
sub,UnaryProxy super
super) (height
height,width
width) Ptr a
aPtr Ptr a
bPtr = do
   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 kl :: Int
kl = UnaryProxy sub -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy sub
sub
   let ku :: Int
ku = UnaryProxy super -> Int
forall x y. (Integer x, Num y) => Proxy x -> y
integralFromProxy UnaryProxy super
super
   let lda0 :: Int
lda0 = Int
klInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
   let lda :: Int
lda = Int
lda0Int -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1

   IO (Maybe ()) -> IO ()
forall (f :: * -> *) a. Functor f => f a -> f ()
void (IO (Maybe ()) -> IO ()) -> IO (Maybe ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ MaybeT IO () -> IO (Maybe ())
forall (m :: * -> *) a. MaybeT m a -> m (Maybe a)
MM.runMaybeT (MaybeT IO () -> IO (Maybe ())) -> MaybeT IO () -> IO (Maybe ())
forall a b. (a -> b) -> a -> b
$ (ReaderT Int (MaybeT IO) () -> Int -> MaybeT IO ())
-> Int -> ReaderT Int (MaybeT IO) () -> MaybeT IO ()
forall a b c. (a -> b -> c) -> b -> a -> c
flip ReaderT Int (MaybeT IO) () -> Int -> MaybeT IO ()
forall r (m :: * -> *) a. ReaderT r m a -> r -> m a
MR.runReaderT Int
n (ReaderT Int (MaybeT IO) () -> MaybeT IO ())
-> ReaderT Int (MaybeT IO) () -> MaybeT IO ()
forall a b. (a -> b) -> a -> b
$
      if Int
m Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
> Int
lda0
         then do -- diagonal stripe
            let col0 :: Int
col0 = Int
ku
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col0 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copyUpperTrapezoid (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kl) Int
col Int
lda0 (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr Int
ku) Int
m Ptr a
bPtr
            let col1 :: Int
col1 = Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
kl
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col1 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubMatrix Int
lda (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
col0)
                  Int
lda (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
col0Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda))
                  (Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
col0Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
m))
            let col2 :: Int
col2 = Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col2 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubTrapezoid Char
'L' Int
lda0 (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
col1)
                  Int
lda0 (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
col1Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda))
                  Int
m (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
col1Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
lda0))
         else do -- full block in the middle
            let col0 :: Int
col0 = Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
kl
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col0 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copyUpperTrapezoid (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
kl) Int
col Int
lda0 (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr Int
ku) Int
m Ptr a
bPtr
            let col1 :: Int
col1 = Int
ku
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col1 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubMatrix Int
m (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
col0)
                  Int
lda0 (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
col0Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
ldaInt -> Int -> Int
forall a. Num a => a -> a -> a
+(Int
col1Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
col0)))
                  Int
m (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
col0Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
m))
            let col2 :: Int
col2 = Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
ku
            Int -> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a. Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col2 ((Int -> IO ()) -> ReaderT Int (MaybeT IO) ())
-> (Int -> IO ()) -> ReaderT Int (MaybeT IO) ()
forall a b. (a -> b) -> a -> b
$ \Int
col ->
               Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubTrapezoid Char
'L' Int
m (Int
colInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
col1)
                  Int
lda0 (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr (Int
kuInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
lda))
                  Int
m (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr (Int
kuInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
m))

withRightBound ::
   Int -> (Int -> IO a) -> MR.ReaderT Int (MM.MaybeT IO) a
withRightBound :: Int -> (Int -> IO a) -> ReaderT Int (MaybeT IO) a
withRightBound Int
col Int -> IO a
act = do
   Int
n <- ReaderT Int (MaybeT IO) Int
forall (m :: * -> *) r. Monad m => ReaderT r m r
MR.ask
   if Int
nInt -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<=Int
col
     then IO a -> ReaderT Int (MaybeT IO) a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (Int -> IO a
act Int
n) ReaderT Int (MaybeT IO) a
-> ReaderT Int (MaybeT IO) a -> ReaderT Int (MaybeT IO) a
forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> ReaderT Int (MaybeT IO) a
forall (m :: * -> *) a. MonadPlus m => m a
mzero
     else IO a -> ReaderT Int (MaybeT IO) a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (Int -> IO a
act Int
col)

copyUpperTrapezoid ::
   (Class.Floating a) =>
   Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copyUpperTrapezoid :: Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copyUpperTrapezoid Int
m Int
n Int
lda Ptr a
aPtr Int
ldb Ptr a
bPtr = do
   let d :: Int
d = Int
mInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
n
   Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubMatrix Int
d Int
n Int
lda Ptr a
aPtr Int
ldb Ptr a
bPtr
   Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
forall a.
Floating a =>
Char -> Int -> Int -> Int -> Ptr a -> Int -> Ptr a -> IO ()
copySubTrapezoid Char
'U' Int
n Int
n
      Int
lda (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
aPtr Int
d)
      Int
ldb (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
bPtr Int
d)