{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
module Numeric.BLAS.Vector.Slice (
   T,
   shape,
   Vector,
   RealOf,
   ComplexOf,
   slice,
   fromStorableVector,
   fromVector,
   toVector,
   extract,
   access,
   replicate,
   append,
   Chunk, chunk,
   concat,
   dot, inner,
   sum,
   absSum,
   norm2,
   norm2Squared,
   normInf,
   normInf1,
   argAbsMaximum,
   argAbs1Maximum,
   product,
   scale, scaleReal,
   add, sub,
   negate, raise,
   mac,
   mul, mulConj,
   minimum, argMinimum,
   maximum, argMaximum,
   limits, argLimits,

   conjugate,
   fromReal,
   toComplex,
   realFromComplexVector,
   realPart,
   imaginaryPart,
   zipComplex,
   unzipComplex,
   ) where

import qualified Numeric.BLAS.Slice as Slice

import qualified Numeric.BLAS.Scalar as Scalar
import qualified Numeric.BLAS.Private as Private
import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated, Conjugated))
import Numeric.BLAS.Scalar (ComplexOf, RealOf)
import Numeric.BLAS.Private (ComplexShape, copyConjugate, realPtr)

import qualified Numeric.BLAS.FFI.Generic as Blas
import qualified Numeric.BLAS.FFI.Complex as BlasComplex
import qualified Numeric.BLAS.FFI.Real as BlasReal
import qualified Numeric.Netlib.Utility as Call
import qualified Numeric.Netlib.Class as Class

import qualified Foreign.Marshal.Array.Guarded as ForeignArray
import Foreign.Marshal.Array (advancePtr)
import Foreign.ForeignPtr (withForeignPtr, castForeignPtr)
import Foreign.Ptr (Ptr)
import Foreign.Storable (Storable, peek, peekElemOff)
import Foreign.C.Types (CInt)

import System.IO.Unsafe (unsafePerformIO)

import Control.Monad.Trans.Cont (ContT(ContT), evalContT)
import Control.Monad.IO.Class (liftIO)
import Control.Monad (fmap, return, (=<<))
import Control.Applicative (liftA2, pure, (<$>), (<*>))

import qualified Data.Array.Comfort.Storable.Unchecked as Array
import qualified Data.Array.Comfort.Shape as Shape
import qualified Data.StorableVector.Base as SVB
import qualified Data.Foldable as Fold
import qualified Data.NonEmpty as NonEmpty
import qualified Data.List as List
import Data.Array.Comfort.Storable.Unchecked (Array(Array))
import Data.Array.Comfort.Shape ((::+)((::+)))

import Text.Show (Show)
import Data.Function (($), (.))
import Data.Complex (Complex)
import Data.Maybe (Maybe(Nothing,Just), maybe)
import Data.Tuple.HT (mapFst, uncurry3)
import Data.Tuple (fst, snd, uncurry)
import Data.Ord ((>=))
import Data.Eq (Eq, (==))

import Prelude (Int, fromIntegral, (-), (+), (*), error, IO)


{- $setup
>>> :set -XTypeOperators
>>> :set -XGADTs
>>> import qualified Test.Slice as TestSlice
>>> import Test.NumberModule.Numeric.BLAS.Vector
>>>    (maxElem, maxDim, genVector, number_, real_, complex_)
>>> import Test.NumberModule.Type (Number_)
>>> import Test.Generator (genNumber)
>>> import Test.Utility (approx, approxReal)
>>>
>>> import qualified Numeric.BLAS.Matrix.RowMajor as Matrix
>>> import qualified Numeric.BLAS.Vector.Slice as VectorSlice
>>> import qualified Numeric.BLAS.Vector as Vector
>>> import qualified Numeric.BLAS.Slice as Slice
>>> import qualified Numeric.BLAS.Scalar as Scalar
>>> import qualified Numeric.Netlib.Class as Class
>>> import qualified Data.Array.Comfort.Storable as Array
>>> import qualified Data.Array.Comfort.Shape as Shape
>>> import qualified Data.List as List
>>> import Numeric.BLAS.Vector ((+++), (|+|))
>>> import Numeric.BLAS.Scalar (RealOf, absolute, realPart, minusOne)
>>> import Data.Array.Comfort.Shape ((::+)((::+)))
>>> import Data.Tuple.HT (mapPair)
>>> import Data.Complex (Complex)
>>> import Control.Applicative (liftA2)
>>>
>>> import qualified Test.QuickCheck as QC
>>> import Test.QuickCheck ((==>))
>>>
>>> type Real_ = RealOf Number_
>>> type Complex_ = Complex Real_
>>>
>>> maxDim1 :: Int
>>> maxDim1 = 10
>>>
>>> type ShapeInt = Shape.ZeroBased Int
>>> type Shape = ShapeInt ::+ (ShapeInt, ShapeInt) ::+ ShapeInt
>>> type Vector = Vector.Vector Shape
>>> type Sliced = VectorSlice.T Shape ShapeInt
>>>
>>> genDim :: QC.Gen Int
>>> genDim = QC.choose (0,maxDim)
>>>
>>> genShapeDim :: Int -> QC.Gen Shape
>>> genShapeDim numRows = do
>>>    left <- fmap Shape.ZeroBased $ QC.choose (0,maxDim)
>>>    right <- fmap Shape.ZeroBased $ QC.choose (0,maxDim)
>>>    columns <- fmap Shape.ZeroBased $ QC.choose (1,maxDim1)
>>>    return (left ::+ (Shape.ZeroBased numRows, columns) ::+ right)
>>>
>>> genShape :: QC.Gen Shape
>>> genShape = genShapeDim =<< QC.choose (0,maxDim1)
>>>
>>> forAll_ :: (Show a) => QC.Gen a -> (a -> QC.Property) -> QC.Property
>>> forAll_ = QC.forAll
>>>
>>> isNonEmpty :: Shape.C sh => VectorSlice.T shA sh a -> Bool
>>> isNonEmpty xs = Shape.size (VectorSlice.shape xs) > 0
>>>
>>> takeColumn ::
>>>    (Shape.Indexed sh1, Shape.C sh, Shape.C sh0, Shape.C sh2) =>
>>>    Shape.Index sh1 ->
>>>    Vector.Vector (sh0 ::+ (sh, sh1) ::+ sh2) a ->
>>>    VectorSlice.T (sh0 ::+ (sh, sh1) ::+ sh2) sh a
>>> takeColumn c =
>>>    VectorSlice.slice (Slice.column c . Slice.left . Slice.right) .
>>>    VectorSlice.fromVector
>>>
>>> listFromSlice ::
>>>    (Shape.C sh, Class.Floating a) => VectorSlice.T shA sh a -> [a]
>>> listFromSlice = Vector.toList . VectorSlice.toVector
>>>
>>> genSlicedDim ::
>>>    (Class.Floating a) =>
>>>    Int -> QC.Gen a -> QC.Gen (Int, Vector a)
>>> genSlicedDim numRows genElem = do
>>>    shape@(_::+(_rows,columns)::+_) <- genShapeDim numRows
>>>    c <- QC.elements (Shape.indices columns)
>>>    fmap ((,) c) $ genVector shape genElem
>>>
>>> genSliced ::
>>>    (Class.Floating a) =>
>>>    QC.Gen a -> QC.Gen (Int, Vector a)
>>> genSliced genElem = flip genSlicedDim genElem =<< genDim
>>>
>>> shrinkSliced ::
>>>    (Shape.C sh0, Shape.Indexed sh1, QC.Arbitrary a, Class.Floating a) =>
>>>    (Shape.Index sh1,
>>>     Vector.Vector (ShapeInt ::+ (sh0, sh1) ::+ ShapeInt) a) ->
>>>    [(Shape.Index sh1,
>>>      Vector.Vector (ShapeInt ::+ (sh0, sh1) ::+ ShapeInt) a)]
>>> shrinkSliced (c,xs) =
>>>    let xs0 = Vector.takeLeft xs in
>>>    let xs1 = Vector.takeRight xs in
>>>    let xs10 = Vector.takeLeft xs1 in
>>>    let xs11 = Vector.takeRight xs1 in
>>>    map (\(ysl,ysr) ->
>>>          (c,
>>>           Vector.autoFromList ysl +++ xs10 +++ Vector.autoFromList ysr)) $
>>>    QC.shrink (Vector.toList xs0, Vector.toList xs11)
>>>
>>> forSliced ::
>>>    (QC.Testable prop, QC.Arbitrary a, Class.Floating a, Show a) =>
>>>    QC.Gen a -> (Sliced a -> prop) -> QC.Property
>>> forSliced genElem prop =
>>>    QC.forAllShrink (genSliced genElem) shrinkSliced
>>>       (prop . uncurry takeColumn)
>>>
>>> genSliced2 ::
>>>    (Class.Floating a) =>
>>>    QC.Gen a -> QC.Gen ((Int, Vector a), (Int, Vector a))
>>> genSliced2 genElem = do
>>>    dim <- genDim
>>>    liftA2 (,) (genSlicedDim dim genElem) (genSlicedDim dim genElem)
>>>
>>> forSliced2 ::
>>>    (QC.Testable prop, Class.Floating a, Show a) =>
>>>    QC.Gen a -> (Sliced a -> Sliced a -> prop) -> QC.Property
>>> forSliced2 genElem prop =
>>>    QC.forAll (genSliced2 genElem)
>>>       (uncurry prop . mapPair (uncurry takeColumn, uncurry takeColumn))
-}


type Vector = Array


shape :: T sh slice a -> slice
shape :: forall sh slice a. T sh slice a -> slice
shape (Cons (Slice.Cons Int
_s Int
_k slice
slc) Vector sh a
_arr) = slice
slc

mapShape :: (slice0 -> slice1) -> T sh slice0 a -> T sh slice1 a
mapShape :: forall slice0 slice1 sh a.
(slice0 -> slice1) -> T sh slice0 a -> T sh slice1 a
mapShape slice0 -> slice1
f (Cons (Slice.Cons Int
s Int
k slice0
slc) Vector sh a
arr) =
   T slice1 -> Vector sh a -> T sh slice1 a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons (Int -> Int -> slice1 -> T slice1
forall sh. Int -> Int -> sh -> T sh
Slice.Cons Int
s Int
k (slice0 -> slice1
f slice0
slc)) Vector sh a
arr


increment :: T sh slice a -> Int
increment :: forall sh slice a. T sh slice a -> Int
increment (Cons (Slice.Cons Int
_s Int
k slice
_slc) Vector sh a
_arr) = Int
k


startArg ::
   (Storable a) =>
   T sh slice a -> Call.FortranIO r (Ptr a)
startArg :: forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg (Cons (Slice.Cons Int
s Int
_k slice
_slice) (Array sh
_sh ForeignPtr a
x)) = do
   Ptr a
sxPtr <- ((Ptr a -> IO r) -> IO r) -> FortranIO r (Ptr a)
forall {k} (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr a -> IO r) -> IO r) -> FortranIO r (Ptr a))
-> ((Ptr a -> IO r) -> IO r) -> FortranIO r (Ptr a)
forall a b. (a -> b) -> a -> b
$ ForeignPtr a -> (Ptr a -> IO r) -> IO r
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr a
x
   Ptr a -> FortranIO r (Ptr a)
forall a. a -> ContT r IO a
forall (m :: * -> *) a. Monad m => a -> m a
return (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
sxPtr Int
s)

sliceArg ::
   (Storable a) =>
   T sh slice a -> Call.FortranIO r (Ptr a, Ptr CInt)
sliceArg :: forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T sh slice a
x =
   (Ptr a -> Ptr CInt -> (Ptr a, Ptr CInt))
-> ContT r IO (Ptr a)
-> ContT r IO (Ptr CInt)
-> ContT r IO (Ptr a, Ptr CInt)
forall a b c.
(a -> b -> c) -> ContT r IO a -> ContT r IO b -> ContT r IO c
forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 (,) (T sh slice a -> ContT r IO (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T sh slice a
x) (Int -> ContT r IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> ContT r IO (Ptr CInt)) -> Int -> ContT r IO (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ T sh slice a -> Int
forall sh slice a. T sh slice a -> Int
increment T sh slice a
x)

sizeSliceArg ::
   (Shape.C sh, Storable a) =>
   T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg :: forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
x =
   (Ptr CInt -> (Ptr a, Ptr CInt) -> (Ptr CInt, Ptr a, Ptr CInt))
-> ContT r IO (Ptr CInt)
-> ContT r IO (Ptr a, Ptr CInt)
-> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
forall a b c.
(a -> b -> c) -> ContT r IO a -> ContT r IO b -> ContT r IO c
forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2
      (\Ptr CInt
nPtr (Ptr a
xPtr,Ptr CInt
incxPtr) -> (Ptr CInt
nPtr, Ptr a
xPtr,Ptr CInt
incxPtr))
      (Int -> ContT r IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> ContT r IO (Ptr CInt)) -> Int -> ContT r IO (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x)
      (T shA sh a -> ContT r IO (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x)

infixl 4 <*|>

(<*|>) ::
   (Storable a) =>
   Call.FortranIO r (Ptr a -> Ptr CInt -> b) ->
   T sh slice a ->
   Call.FortranIO r b
FortranIO r (Ptr a -> Ptr CInt -> b)
f <*|> :: forall a r b sh slice.
Storable a =>
FortranIO r (Ptr a -> Ptr CInt -> b)
-> T sh slice a -> FortranIO r b
<*|> T sh slice a
x  =  ((Ptr a -> Ptr CInt -> b) -> (Ptr a, Ptr CInt) -> b)
-> FortranIO r (Ptr a -> Ptr CInt -> b)
-> ContT r IO ((Ptr a, Ptr CInt) -> b)
forall a b. (a -> b) -> ContT r IO a -> ContT r IO b
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (Ptr a -> Ptr CInt -> b) -> (Ptr a, Ptr CInt) -> b
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry FortranIO r (Ptr a -> Ptr CInt -> b)
f ContT r IO ((Ptr a, Ptr CInt) -> b)
-> ContT r IO (Ptr a, Ptr CInt) -> ContT r IO b
forall a b. ContT r IO (a -> b) -> ContT r IO a -> ContT r IO b
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> T sh slice a -> ContT r IO (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T sh slice a
x


data T sh slice a = Cons (Slice.T slice) (Vector sh a)
   deriving (Int -> T sh slice a -> ShowS
[T sh slice a] -> ShowS
T sh slice a -> String
(Int -> T sh slice a -> ShowS)
-> (T sh slice a -> String)
-> ([T sh slice a] -> ShowS)
-> Show (T sh slice a)
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
Int -> T sh slice a -> ShowS
forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
[T sh slice a] -> ShowS
forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
T sh slice a -> String
$cshowsPrec :: forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
Int -> T sh slice a -> ShowS
showsPrec :: Int -> T sh slice a -> ShowS
$cshow :: forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
T sh slice a -> String
show :: T sh slice a -> String
$cshowList :: forall sh slice a.
(C sh, Storable a, Show slice, Show sh, Show a) =>
[T sh slice a] -> ShowS
showList :: [T sh slice a] -> ShowS
Show)

toVector ::
   (Shape.C slice, Class.Floating a) =>
   T sh slice a -> Vector slice a
toVector :: forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector T sh slice a
x =
   slice -> (Int -> Ptr a -> IO ()) -> Array slice a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T sh slice a -> slice
forall sh slice a. T sh slice a -> slice
shape T sh slice a
x) ((Int -> Ptr a -> IO ()) -> Array slice a)
-> (Int -> Ptr a -> IO ()) -> Array slice a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
syPtr ->
   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
$ FortranIO () (IO ()) -> ContT () IO ()
forall r a. FortranIO r (IO a) -> FortranIO r a
Call.run (FortranIO () (IO ()) -> ContT () IO ())
-> FortranIO () (IO ()) -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$
      (Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ())
-> ContT
     () IO (Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ())
forall a. a -> ContT () IO a
forall (f :: * -> *) a. Applicative f => a -> f a
pure 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 ()
Blas.copy
         ContT
  () IO (Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ())
-> ContT () IO (Ptr CInt)
-> ContT () IO (Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ())
forall a b. ContT () IO (a -> b) -> ContT () IO a -> ContT () IO b
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
         ContT () IO (Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ())
-> T sh slice a -> FortranIO () (Ptr a -> Ptr CInt -> IO ())
forall a r b sh slice.
Storable a =>
FortranIO r (Ptr a -> Ptr CInt -> b)
-> T sh slice a -> FortranIO r b
<*|> T sh slice a
x
         FortranIO () (Ptr a -> Ptr CInt -> IO ())
-> ContT () IO (Ptr a) -> ContT () IO (Ptr CInt -> IO ())
forall a b. ContT () IO (a -> b) -> ContT () IO a -> ContT () IO b
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Ptr a -> ContT () IO (Ptr a)
forall a. a -> ContT () IO a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Ptr a
syPtr
         ContT () IO (Ptr CInt -> IO ())
-> ContT () IO (Ptr CInt) -> FortranIO () (IO ())
forall a b. ContT () IO (a -> b) -> ContT () IO a -> ContT () IO b
forall (f :: * -> *) a b. Applicative f => f (a -> b) -> f a -> f b
<*> Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1

{- |
Non-copying conversion from @StorableVector@.
-}
fromStorableVector :: SVB.Vector a -> T ShapeInt ShapeInt a
fromStorableVector :: forall a. Vector a -> T ShapeInt ShapeInt a
fromStorableVector Vector a
xs =
   case Vector a -> (ForeignPtr a, Int, Int)
forall a. Vector a -> (ForeignPtr a, Int, Int)
SVB.toForeignPtr Vector a
xs of
      (ForeignPtr a
buffer, Int
s, Int
l) ->
         T ShapeInt -> Vector ShapeInt a -> T ShapeInt ShapeInt a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons
            (Int -> Int -> ShapeInt -> T ShapeInt
forall sh. Int -> Int -> sh -> T sh
Slice.Cons Int
s Int
1 (Int -> ShapeInt
forall n. n -> ZeroBased n
Shape.ZeroBased Int
l))
            (ShapeInt -> ForeignPtr a -> Vector ShapeInt a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array (Int -> ShapeInt
forall n. n -> ZeroBased n
Shape.ZeroBased (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
l)) ForeignPtr a
buffer)

fromVector :: (Shape.C sh) => Vector sh a -> T sh sh a
fromVector :: forall sh a. C sh => Vector sh a -> T sh sh a
fromVector Vector sh a
xs = T sh -> Vector sh a -> T sh sh a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons (sh -> T sh
forall sh. C sh => sh -> T sh
Slice.fromShape (sh -> T sh) -> sh -> T sh
forall a b. (a -> b) -> a -> b
$ Vector sh a -> sh
forall sh a. Array sh a -> sh
Array.shape Vector sh a
xs) Vector sh a
xs

slice :: (Slice.T shA -> Slice.T shB) -> T sh shA a -> T sh shB a
slice :: forall shA shB sh a. (T shA -> T shB) -> T sh shA a -> T sh shB a
slice T shA -> T shB
f (Cons T shA
slc Vector sh a
xs) = T shB -> Vector sh a -> T sh shB a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons (T shA -> T shB
f T shA
slc) Vector sh a
xs

{- |
prop> QC.forAll genShape $ \shape@(_::+(_rows,columns)::+_) -> QC.forAll (QC.elements (Shape.indices columns)) $ \c -> QC.forAll (genVector shape number_) $ \xs -> VectorSlice.extract (Slice.column c . Slice.left . Slice.right) xs == Matrix.takeColumn c (Vector.takeLeft (Vector.takeRight (xs :: Vector Number_)))

prop> forAll_ (TestSlice.genShapeSelect 4 100) $ \(TestSlice.ShapeSelect sh select) -> QC.forAll (genVector sh number_) $ \xs -> case TestSlice.instantiate sh select of TestSlice.Extraction slice cut -> VectorSlice.extract slice xs == cut xs
-}
extract ::
   (Shape.C slice, Shape.C sh, Class.Floating a) =>
   (Slice.T sh -> Slice.T slice) -> Vector sh a -> Vector slice a
extract :: forall slice sh a.
(C slice, C sh, Floating a) =>
(T sh -> T slice) -> Vector sh a -> Vector slice a
extract T sh -> T slice
slc Vector sh a
xs = T sh slice a -> Vector slice a
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector (T sh slice a -> Vector slice a) -> T sh slice a -> Vector slice a
forall a b. (a -> b) -> a -> b
$ (T sh -> T slice) -> T sh sh a -> T sh slice a
forall shA shB sh a. (T shA -> T shB) -> T sh shA a -> T sh shB a
slice T sh -> T slice
slc (T sh sh a -> T sh slice a) -> T sh sh a -> T sh slice a
forall a b. (a -> b) -> a -> b
$ Vector sh a -> T sh sh a
forall sh a. C sh => Vector sh a -> T sh sh a
fromVector Vector sh a
xs


access, (!) ::
   (Shape.C shA, Shape.Indexed sh, Storable a) =>
   T shA sh a -> Shape.Index sh -> a
access :: forall shA sh a.
(C shA, Indexed sh, Storable a) =>
T shA sh a -> Index sh -> a
access (Cons (Slice.Cons Int
s Int
k sh
ssh) (Array shA
sh ForeignPtr a
x)) Index sh
ix =
   Deferred shA -> ForeignPtr a -> Array (Deferred shA) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array (shA -> Deferred shA
forall sh. sh -> Deferred sh
Shape.Deferred shA
sh) ForeignPtr a
x
   Array (Deferred shA) a -> Index (Deferred shA) -> a
forall sh a.
(Indexed sh, Storable a) =>
Array sh a -> Index sh -> a
Array.!
   Int -> DeferredIndex shA
forall sh. Int -> DeferredIndex sh
Shape.DeferredIndex (Int
s Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
k Int -> Int -> Int
forall a. Num a => a -> a -> a
* sh -> Index sh -> Int
forall sh. Indexed sh => sh -> Index sh -> Int
Shape.offset sh
ssh Index sh
ix)
! :: forall shA sh a.
(C shA, Indexed sh, Storable a) =>
T shA sh a -> Index sh -> a
(!) = T shA sh a -> Index sh -> a
forall shA sh a.
(C shA, Indexed sh, Storable a) =>
T shA sh a -> Index sh -> a
access



{- |
prop> QC.forAll genShape $ \shape a -> VectorSlice.toVector (VectorSlice.replicate shape a) == Vector.constant shape (a :: Number_)
-}
replicate :: (Shape.C sh, Storable a) => sh -> a -> T () sh a
replicate :: forall sh a. (C sh, Storable a) => sh -> a -> T () sh a
replicate sh
sh a
a =
   T sh -> Vector () a -> T () sh a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons
      (Int -> Int -> sh -> T sh
forall sh. Int -> Int -> sh -> T sh
Slice.Cons Int
0 Int
0 sh
sh)
      (a -> () -> [(Index (), a)] -> Vector () a
forall sh a.
(Indexed sh, Storable a) =>
a -> sh -> [(Index sh, a)] -> Array sh a
Array.fromAssociations a
a () [])


type ShapeInt = Shape.ZeroBased Int

{- |
>>> List.map realPart $ Vector.toList $ VectorSlice.append (VectorSlice.fromVector $ Vector.autoFromList [3,1,4,1]) (VectorSlice.replicate (Shape.ZeroBased (3::Int)) (0::Number_))
[3.0,1.0,4.0,1.0,0.0,0.0,0.0]
-}
append ::
   (Shape.C sliceA, Shape.C sliceB, Class.Floating a) =>
   T shA sliceA a -> T shB sliceB a -> Vector (sliceA::+sliceB) a
append :: forall sliceA sliceB a shA shB.
(C sliceA, C sliceB, Floating a) =>
T shA sliceA a -> T shB sliceB a -> Vector (sliceA ::+ sliceB) a
append T shA sliceA a
xs T shB sliceB a
ys =
   (sliceA ::+ sliceB)
-> (Ptr a -> IO ()) -> Array (sliceA ::+ sliceB) a
forall sh a.
(C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> Array sh a
Array.unsafeCreate (T shA sliceA a -> sliceA
forall sh slice a. T sh slice a -> slice
shape T shA sliceA a
xs sliceA -> sliceB -> sliceA ::+ sliceB
forall sh0 sh1. sh0 -> sh1 -> sh0 ::+ sh1
::+ T shB sliceB a -> sliceB
forall sh slice a. T sh slice a -> slice
shape T shB sliceB a
ys) ((Ptr a -> IO ()) -> Array (sliceA ::+ sliceB) a)
-> (Ptr a -> IO ()) -> Array (sliceA ::+ sliceB) a
forall a b. (a -> b) -> a -> b
$ \Ptr a
zPtr ->
   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 chunkX :: Chunk a
chunkX = T shA ShapeInt a -> Chunk a
forall a sh. Storable a => T sh ShapeInt a -> Chunk a
chunk (T shA ShapeInt a -> Chunk a) -> T shA ShapeInt a -> Chunk a
forall a b. (a -> b) -> a -> b
$ (sliceA -> ShapeInt) -> T shA sliceA a -> T shA ShapeInt a
forall slice0 slice1 sh a.
(slice0 -> slice1) -> T sh slice0 a -> T sh slice1 a
mapShape (Int -> ShapeInt
forall n. n -> ZeroBased n
Shape.ZeroBased (Int -> ShapeInt) -> (sliceA -> Int) -> sliceA -> ShapeInt
forall b c a. (b -> c) -> (a -> b) -> a -> c
. sliceA -> Int
forall sh. C sh => sh -> Int
Shape.size) T shA sliceA a
xs
      let chunkY :: Chunk a
chunkY = T shB ShapeInt a -> Chunk a
forall a sh. Storable a => T sh ShapeInt a -> Chunk a
chunk (T shB ShapeInt a -> Chunk a) -> T shB ShapeInt a -> Chunk a
forall a b. (a -> b) -> a -> b
$ (sliceB -> ShapeInt) -> T shB sliceB a -> T shB ShapeInt a
forall slice0 slice1 sh a.
(slice0 -> slice1) -> T sh slice0 a -> T sh slice1 a
mapShape (Int -> ShapeInt
forall n. n -> ZeroBased n
Shape.ZeroBased (Int -> ShapeInt) -> (sliceB -> Int) -> sliceB -> ShapeInt
forall b c a. (b -> c) -> (a -> b) -> a -> c
. sliceB -> Int
forall sh. C sh => sh -> Int
Shape.size) T shB sliceB a
ys
      Ptr CInt
nxPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> ContT () IO (Ptr CInt)) -> Int -> ContT () IO (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ Chunk a -> Int
forall a. Chunk a -> Int
chunkSize Chunk a
chunkX
      Ptr CInt
nyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> ContT () IO (Ptr CInt)) -> Int -> ContT () IO (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ Chunk a -> Int
forall a. Chunk a -> Int
chunkSize Chunk a
chunkY
      (Ptr a
xPtr, Ptr CInt
incxPtr) <- Chunk a -> ContT () IO (Ptr a, Ptr CInt)
forall a. Chunk a -> ContT () IO (Ptr a, Ptr CInt)
chunkData Chunk a
chunkX
      (Ptr a
yPtr, Ptr CInt
incyPtr) <- Chunk a -> ContT () IO (Ptr a, Ptr CInt)
forall a. Chunk a -> ContT () IO (Ptr a, Ptr CInt)
chunkData Chunk a
chunkY
      Ptr CInt
inczPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 ()
Blas.copy Ptr CInt
nxPtr Ptr a
xPtr Ptr CInt
incxPtr Ptr a
zPtr Ptr CInt
inczPtr
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 ()
Blas.copy Ptr CInt
nyPtr Ptr a
yPtr Ptr CInt
incyPtr
                     (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
zPtr (Int -> Ptr a) -> Int -> Ptr a
forall a b. (a -> b) -> a -> b
$ Chunk a -> Int
forall a. Chunk a -> Int
chunkSize Chunk a
chunkX) Ptr CInt
inczPtr

data Chunk a =
   Chunk {
      forall a. Chunk a -> Int
chunkSize :: Int,
      forall a. Chunk a -> ContT () IO (Ptr a, Ptr CInt)
chunkData :: ContT () IO (Ptr a, Ptr CInt)
   }

chunk :: (Storable a) => T sh ShapeInt a -> Chunk a
chunk :: forall a sh. Storable a => T sh ShapeInt a -> Chunk a
chunk x :: T sh ShapeInt a
x@(Cons T ShapeInt
slc Vector sh a
_dat) = Int -> ContT () IO (Ptr a, Ptr CInt) -> Chunk a
forall a. Int -> ContT () IO (Ptr a, Ptr CInt) -> Chunk a
Chunk (ShapeInt -> Int
forall sh. C sh => sh -> Int
Shape.size (ShapeInt -> Int) -> ShapeInt -> Int
forall a b. (a -> b) -> a -> b
$ T ShapeInt -> ShapeInt
forall sh. T sh -> sh
Slice.shape T ShapeInt
slc) (T sh ShapeInt a -> ContT () IO (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T sh ShapeInt a
x)

{- |
>>> List.map realPart $ Vector.toList $ let matrix = Vector.fromList (Shape.ZeroBased (4::Int), Shape.ZeroBased (3::Int)) $ map (\x -> fromInteger x :: Number_) [0 ..] in VectorSlice.concat $ map (\k -> VectorSlice.chunk (VectorSlice.slice (Slice.column k) $ VectorSlice.fromVector matrix)) [0,1,2]
[0.0,3.0,6.0,9.0,1.0,4.0,7.0,10.0,2.0,5.0,8.0,11.0]
-}
concat :: (Class.Floating a) => [Chunk a] -> Vector ShapeInt a
concat :: forall a. Floating a => [Chunk a] -> Vector ShapeInt a
concat [Chunk a]
xs =
   let offsets :: T [] Int
offsets = (Int -> Int -> Int) -> Int -> [Int] -> T [] Int
forall (f :: * -> *) b a.
Traversable f =>
(b -> a -> b) -> b -> f a -> T f b
NonEmpty.scanl Int -> Int -> Int
forall a. Num a => a -> a -> a
(+) Int
0 ([Int] -> T [] Int) -> [Int] -> T [] Int
forall a b. (a -> b) -> a -> b
$ (Chunk a -> Int) -> [Chunk a] -> [Int]
forall a b. (a -> b) -> [a] -> [b]
List.map Chunk a -> Int
forall a. Chunk a -> Int
chunkSize [Chunk a]
xs
   in ShapeInt -> (Ptr a -> IO ()) -> Array ShapeInt a
forall sh a.
(C sh, Storable a) =>
sh -> (Ptr a -> IO ()) -> Array sh a
Array.unsafeCreate (Int -> ShapeInt
forall n. n -> ZeroBased n
Shape.ZeroBased (Int -> ShapeInt) -> Int -> ShapeInt
forall a b. (a -> b) -> a -> b
$ T [] Int -> Int
forall (f :: * -> *) a. Foldable f => T f a -> a
NonEmpty.last T [] Int
offsets) ((Ptr a -> IO ()) -> Array ShapeInt a)
-> (Ptr a -> IO ()) -> Array ShapeInt a
forall a b. (a -> b) -> a -> b
$ \Ptr a
yPtr ->
      [(Int, Chunk a)] -> ((Int, Chunk a) -> IO ()) -> IO ()
forall (t :: * -> *) (f :: * -> *) a b.
(Foldable t, Applicative f) =>
t a -> (a -> f b) -> f ()
Fold.for_ ([Int] -> [Chunk a] -> [(Int, Chunk a)]
forall a b. [a] -> [b] -> [(a, b)]
List.zip (T [] Int -> [Int]
forall (f :: * -> *) a. Cons f => T f a -> f a
NonEmpty.flatten T [] Int
offsets) [Chunk a]
xs) (((Int, Chunk a) -> IO ()) -> IO ())
-> ((Int, Chunk a) -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \(Int
offset, Chunk a
x) ->
      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 -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> ContT () IO (Ptr CInt)) -> Int -> ContT () IO (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ Chunk a -> Int
forall a. Chunk a -> Int
chunkSize Chunk a
x
         (Ptr a
xPtr, Ptr CInt
incxPtr) <- Chunk a -> ContT () IO (Ptr a, Ptr CInt)
forall a. Chunk a -> ContT () IO (Ptr a, Ptr CInt)
chunkData Chunk a
x
         Ptr CInt
incyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
         IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 ()
Blas.copy Ptr CInt
nPtr Ptr a
xPtr Ptr CInt
incxPtr (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
yPtr Int
offset) Ptr CInt
incyPtr


dot ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   T shA sh a -> T shB sh a -> a
dot :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> a
dot =
   Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a
forall (f :: * -> *) (g :: * -> *) a. Dot f g a -> f a -> g a -> a
runDot (Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a)
-> Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a
forall a b. (a -> b) -> a -> b
$
   Dot (T shA sh) (T shB sh) Float
-> Dot (T shA sh) (T shB sh) Double
-> Dot (T shA sh) (T shB sh) (Complex Float)
-> Dot (T shA sh) (T shB sh) (Complex Double)
-> Dot (T shA sh) (T shB sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> T shB sh Float -> Float)
-> Dot (T shA sh) (T shB sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh Float -> T shB sh Float -> Float
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> a
dotReal)
      ((T shA sh Double -> T shB sh Double -> Double)
-> Dot (T shA sh) (T shB sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh Double -> T shB sh Double -> Double
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> a
dotReal)
      ((T shA sh (Complex Float)
 -> T shB sh (Complex Float) -> Complex Float)
-> Dot (T shA sh) (T shB sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh (Complex Float)
-> T shB sh (Complex Float) -> Complex Float
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh (Complex a) -> T shB sh (Complex a) -> Complex a
dotComplex)
      ((T shA sh (Complex Double)
 -> T shB sh (Complex Double) -> Complex Double)
-> Dot (T shA sh) (T shB sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh (Complex Double)
-> T shB sh (Complex Double) -> Complex Double
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh (Complex a) -> T shB sh (Complex a) -> Complex a
dotComplex)


{- |
prop> forSliced2 number_ $ \xs ys -> VectorSlice.inner xs ys == Vector.dot (VectorSlice.conjugate xs) (VectorSlice.toVector ys)
prop> forSliced number_ $ \xs -> VectorSlice.inner xs xs == Scalar.fromReal (VectorSlice.norm2Squared xs)
-}
inner ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   T shA sh a -> T shB sh a -> a
inner :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> a
inner =
   Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a
forall (f :: * -> *) (g :: * -> *) a. Dot f g a -> f a -> g a -> a
runDot (Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a)
-> Dot (T shA sh) (T shB sh) a -> T shA sh a -> T shB sh a -> a
forall a b. (a -> b) -> a -> b
$
   Dot (T shA sh) (T shB sh) Float
-> Dot (T shA sh) (T shB sh) Double
-> Dot (T shA sh) (T shB sh) (Complex Float)
-> Dot (T shA sh) (T shB sh) (Complex Double)
-> Dot (T shA sh) (T shB sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> T shB sh Float -> Float)
-> Dot (T shA sh) (T shB sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh Float -> T shB sh Float -> Float
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> a
dotReal)
      ((T shA sh Double -> T shB sh Double -> Double)
-> Dot (T shA sh) (T shB sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot T shA sh Double -> T shB sh Double -> Double
forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> a
dotReal)
      ((T shA sh (Complex Float)
 -> T shB sh (Complex Float) -> Complex Float)
-> Dot (T shA sh) (T shB sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot ((T shA sh (Complex Float)
  -> T shB sh (Complex Float) -> Complex Float)
 -> Dot (T shA sh) (T shB sh) (Complex Float))
-> (T shA sh (Complex Float)
    -> T shB sh (Complex Float) -> Complex Float)
-> Dot (T shA sh) (T shB sh) (Complex Float)
forall a b. (a -> b) -> a -> b
$ Vector sh (Complex Float)
-> T shB sh (Complex Float) -> Complex Float
forall sh a shB.
(C sh, Eq sh, Real a) =>
Vector sh (Complex a) -> T shB sh (Complex a) -> Complex a
innerComplex (Vector sh (Complex Float)
 -> T shB sh (Complex Float) -> Complex Float)
-> (T shA sh (Complex Float) -> Vector sh (Complex Float))
-> T shA sh (Complex Float)
-> T shB sh (Complex Float)
-> Complex Float
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh (Complex Float) -> Vector sh (Complex Float)
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh (Complex Double)
 -> T shB sh (Complex Double) -> Complex Double)
-> Dot (T shA sh) (T shB sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a -> a) -> Dot f g a
Dot ((T shA sh (Complex Double)
  -> T shB sh (Complex Double) -> Complex Double)
 -> Dot (T shA sh) (T shB sh) (Complex Double))
-> (T shA sh (Complex Double)
    -> T shB sh (Complex Double) -> Complex Double)
-> Dot (T shA sh) (T shB sh) (Complex Double)
forall a b. (a -> b) -> a -> b
$ Vector sh (Complex Double)
-> T shB sh (Complex Double) -> Complex Double
forall sh a shB.
(C sh, Eq sh, Real a) =>
Vector sh (Complex a) -> T shB sh (Complex a) -> Complex a
innerComplex (Vector sh (Complex Double)
 -> T shB sh (Complex Double) -> Complex Double)
-> (T shA sh (Complex Double) -> Vector sh (Complex Double))
-> T shA sh (Complex Double)
-> T shB sh (Complex Double)
-> Complex Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh (Complex Double) -> Vector sh (Complex Double)
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)


newtype Dot f g a = Dot {forall (f :: * -> *) (g :: * -> *) a. Dot f g a -> f a -> g a -> a
runDot :: f a -> g a -> a}

dotReal ::
   (Shape.C sh, Eq sh, Class.Real a) =>
   T shA sh a -> T shB sh a -> a
dotReal :: forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> a
dotReal T shA sh a
x T shB sh a
y = IO a -> a
forall a. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$ do
   let shX :: sh
shX = T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x
   let shY :: sh
shY = T shB sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh a
y
   String -> Bool -> IO ()
Call.assert String
"dot: shapes mismatch" (sh
shX sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== sh
shY)
   ContT a IO a -> IO a
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT a IO a -> IO a) -> ContT a IO a -> IO a
forall a b. (a -> b) -> a -> b
$ do
      Ptr CInt
nPtr <- Int -> FortranIO a (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> FortranIO a (Ptr CInt)) -> Int -> FortranIO a (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ sh -> Int
forall sh. C sh => sh -> Int
Shape.size sh
shX
      (Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> FortranIO a (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x
      (Ptr a
syPtr, Ptr CInt
incyPtr) <- T shB sh a -> FortranIO a (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shB sh a
y
      IO a -> ContT a IO a
forall a. IO a -> ContT a IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO a -> ContT a IO a) -> IO a -> ContT a IO a
forall a b. (a -> b) -> a -> b
$ Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
forall a.
Real a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.dot Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr Ptr a
syPtr Ptr CInt
incyPtr

{-
We cannot use 'cdot' because Haskell's FFI
does not support Complex numbers as return values.
-}
dotComplex ::
   (Shape.C sh, Eq sh, Class.Real a) =>
   T shA sh (Complex a) -> T shB sh (Complex a) -> Complex a
dotComplex :: forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh (Complex a) -> T shB sh (Complex a) -> Complex a
dotComplex T shA sh (Complex a)
x T shB sh (Complex a)
y = IO (Complex a) -> Complex a
forall a. IO a -> a
unsafePerformIO (IO (Complex a) -> Complex a) -> IO (Complex a) -> Complex a
forall a b. (a -> b) -> a -> b
$ do
   String -> Bool -> IO ()
Call.assert String
"dot: shapes mismatch" (T shA sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh (Complex a)
x sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== T shB sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh (Complex a)
y)
   ContT (Complex a) IO (Complex a) -> IO (Complex a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (Complex a) IO (Complex a) -> IO (Complex a))
-> ContT (Complex a) IO (Complex a) -> IO (Complex a)
forall a b. (a -> b) -> a -> b
$ do
      Ptr CChar
transPtr <- Char -> FortranIO (Complex a) (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'N'
      Ptr CInt
mPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr CInt
nPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> FortranIO (Complex a) (Ptr CInt))
-> Int -> FortranIO (Complex a) (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh (Complex a)
x
      Ptr (Complex a)
alphaPtr <- Complex a -> FortranIO (Complex a) (Ptr (Complex a))
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number Complex a
forall a. Floating a => a
Scalar.one
      (Ptr (Complex a)
xPtr, Ptr CInt
ldxPtr) <- T shA sh (Complex a)
-> FortranIO (Complex a) (Ptr (Complex a), Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh (Complex a)
x
      (Ptr (Complex a)
yPtr, Ptr CInt
incyPtr) <- T shB sh (Complex a)
-> FortranIO (Complex a) (Ptr (Complex a), Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shB sh (Complex a)
y
      Ptr (Complex a)
betaPtr <- Complex a -> FortranIO (Complex a) (Ptr (Complex a))
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number Complex a
forall a. Floating a => a
Scalar.zero
      Ptr (Complex a)
zPtr <- FortranIO (Complex a) (Ptr (Complex a))
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
inczPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT (Complex a) IO ()
forall a. IO a -> ContT (Complex a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT (Complex a) IO ())
-> IO () -> ContT (Complex a) IO ()
forall a b. (a -> b) -> a -> b
$
         Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr (Complex a)
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr (Complex 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
mPtr Ptr CInt
nPtr Ptr (Complex a)
alphaPtr Ptr (Complex a)
xPtr Ptr CInt
ldxPtr
            Ptr (Complex a)
yPtr Ptr CInt
incyPtr Ptr (Complex a)
betaPtr Ptr (Complex a)
zPtr Ptr CInt
inczPtr
      IO (Complex a) -> ContT (Complex a) IO (Complex a)
forall a. IO a -> ContT (Complex a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Complex a) -> ContT (Complex a) IO (Complex a))
-> IO (Complex a) -> ContT (Complex a) IO (Complex a)
forall a b. (a -> b) -> a -> b
$ Ptr (Complex a) -> IO (Complex a)
forall a. Storable a => Ptr a -> IO a
peek Ptr (Complex a)
zPtr

{-
In contrast to Vector.inner
we cannot use gemv with 'C' because we need leading dimension > 1.
-}
innerComplex ::
   (Shape.C sh, Eq sh, Class.Real a) =>
   Vector sh (Complex a) -> T shB sh (Complex a) -> Complex a
innerComplex :: forall sh a shB.
(C sh, Eq sh, Real a) =>
Vector sh (Complex a) -> T shB sh (Complex a) -> Complex a
innerComplex (Array sh
shX ForeignPtr (Complex a)
x) T shB sh (Complex a)
y = IO (Complex a) -> Complex a
forall a. IO a -> a
unsafePerformIO (IO (Complex a) -> Complex a) -> IO (Complex a) -> Complex a
forall a b. (a -> b) -> a -> b
$ do
   String -> Bool -> IO ()
Call.assert String
"dot: shapes mismatch" (sh
shX sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== T shB sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh (Complex a)
y)
   ContT (Complex a) IO (Complex a) -> IO (Complex a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (Complex a) IO (Complex a) -> IO (Complex a))
-> ContT (Complex a) IO (Complex a) -> IO (Complex a)
forall a b. (a -> b) -> a -> b
$ do
      let m :: Int
m = sh -> Int
forall sh. C sh => sh -> Int
Shape.size sh
shX
      Ptr CChar
transPtr <- Char -> FortranIO (Complex a) (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'C'
      Ptr CInt
mPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
m
      Ptr CInt
nPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      Ptr (Complex a)
alphaPtr <- Complex a -> FortranIO (Complex a) (Ptr (Complex a))
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number Complex a
forall a. Floating a => a
Scalar.one
      Ptr (Complex a)
xPtr <- ((Ptr (Complex a) -> IO (Complex a)) -> IO (Complex a))
-> FortranIO (Complex a) (Ptr (Complex a))
forall {k} (r :: k) (m :: k -> *) a.
((a -> m r) -> m r) -> ContT r m a
ContT (((Ptr (Complex a) -> IO (Complex a)) -> IO (Complex a))
 -> FortranIO (Complex a) (Ptr (Complex a)))
-> ((Ptr (Complex a) -> IO (Complex a)) -> IO (Complex a))
-> FortranIO (Complex a) (Ptr (Complex a))
forall a b. (a -> b) -> a -> b
$ ForeignPtr (Complex a)
-> (Ptr (Complex a) -> IO (Complex a)) -> IO (Complex a)
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr (Complex a)
x
      Ptr CInt
ldxPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
      (Ptr (Complex a)
yPtr, Ptr CInt
incyPtr) <- T shB sh (Complex a)
-> FortranIO (Complex a) (Ptr (Complex a), Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shB sh (Complex a)
y
      Ptr (Complex a)
betaPtr <- Complex a -> FortranIO (Complex a) (Ptr (Complex a))
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number Complex a
forall a. Floating a => a
Scalar.zero
      Ptr (Complex a)
zPtr <- FortranIO (Complex a) (Ptr (Complex a))
forall a r. Storable a => FortranIO r (Ptr a)
Call.alloca
      Ptr CInt
inczPtr <- Int -> FortranIO (Complex a) (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT (Complex a) IO ()
forall a. IO a -> ContT (Complex a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT (Complex a) IO ())
-> IO () -> ContT (Complex a) IO ()
forall a b. (a -> b) -> a -> b
$
         Ptr CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr (Complex a)
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr (Complex 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
mPtr Ptr CInt
nPtr Ptr (Complex a)
alphaPtr Ptr (Complex a)
xPtr Ptr CInt
ldxPtr
            Ptr (Complex a)
yPtr Ptr CInt
incyPtr Ptr (Complex a)
betaPtr Ptr (Complex a)
zPtr Ptr CInt
inczPtr
      IO (Complex a) -> ContT (Complex a) IO (Complex a)
forall a. IO a -> ContT (Complex a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Complex a) -> ContT (Complex a) IO (Complex a))
-> IO (Complex a) -> ContT (Complex a) IO (Complex a)
forall a b. (a -> b) -> a -> b
$ Ptr (Complex a) -> IO (Complex a)
forall a. Storable a => Ptr a -> IO a
peek Ptr (Complex a)
zPtr


{- |
prop> forSliced number_ $ \xs -> VectorSlice.sum xs == List.sum (listFromSlice xs)
-}
sum :: (Shape.C sh, Class.Floating a) => T shA sh a -> a
sum :: forall sh a shA. (C sh, Floating a) => T shA sh a -> a
sum T shA sh a
x = IO a -> a
forall a. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$ ContT a IO a -> IO a
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT a IO a -> IO a) -> ContT a IO a -> IO a
forall a b. (a -> b) -> a -> b
$ do
   Ptr a
xPtr <- T shA sh a -> FortranIO a (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shA sh a
x
   IO a -> ContT a IO a
forall a. IO a -> ContT a IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO a -> ContT a IO a) -> IO a -> ContT a IO a
forall a b. (a -> b) -> a -> b
$ Int -> Ptr a -> Int -> IO a
forall a. Floating a => Int -> Ptr a -> Int -> IO a
Private.sum (sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) Ptr a
xPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x)


{- |
Sum of the absolute values of real numbers or components of complex numbers.
For real numbers it is equivalent to 'Numeric.LAPACK.Vector.norm1'.
-}
absSum :: (Shape.C sh, Class.Floating a) => T shA sh a -> RealOf a
absSum :: forall sh a shA. (C sh, Floating a) => T shA sh a -> RealOf a
absSum T shA sh a
arr = IO (RealOf a) -> RealOf a
forall a. IO a -> a
unsafePerformIO (IO (RealOf a) -> RealOf a) -> IO (RealOf a) -> RealOf a
forall a b. (a -> b) -> a -> b
$
   ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (RealOf a) IO (RealOf a) -> IO (RealOf a))
-> ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$ IO (RealOf a) -> ContT (RealOf a) IO (RealOf a)
forall a. IO a -> ContT (RealOf a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (RealOf a) -> ContT (RealOf a) IO (RealOf a))
-> ((Ptr CInt, Ptr a, Ptr CInt) -> IO (RealOf a))
-> (Ptr CInt, Ptr a, Ptr CInt)
-> ContT (RealOf a) IO (RealOf a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a))
-> (Ptr CInt, Ptr a, Ptr CInt) -> IO (RealOf a)
forall a b c d. (a -> b -> c -> d) -> (a, b, c) -> d
uncurry3 Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
asum ((Ptr CInt, Ptr a, Ptr CInt) -> ContT (RealOf a) IO (RealOf a))
-> ContT (RealOf a) IO (Ptr CInt, Ptr a, Ptr CInt)
-> ContT (RealOf a) IO (RealOf a)
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< T shA sh a -> ContT (RealOf a) IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
arr

asum :: Class.Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
asum :: forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
asum =
   Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a. Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
getNrm (Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a))
-> Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$
   Nrm Float
-> Nrm Double
-> Nrm (Complex Float)
-> Nrm (Complex Double)
-> Nrm a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((Ptr CInt -> Ptr Float -> Ptr CInt -> IO (RealOf Float))
-> Nrm Float
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr Float -> Ptr CInt -> IO Float
Ptr CInt -> Ptr Float -> Ptr CInt -> IO (RealOf Float)
forall a. Real a => Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.asum) ((Ptr CInt -> Ptr Double -> Ptr CInt -> IO (RealOf Double))
-> Nrm Double
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr Double -> Ptr CInt -> IO Double
Ptr CInt -> Ptr Double -> Ptr CInt -> IO (RealOf Double)
forall a. Real a => Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.asum)
      ((Ptr CInt
 -> Ptr (Complex Float) -> Ptr CInt -> IO (RealOf (Complex Float)))
-> Nrm (Complex Float)
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr (Complex Float) -> Ptr CInt -> IO Float
Ptr CInt
-> Ptr (Complex Float) -> Ptr CInt -> IO (RealOf (Complex Float))
forall a. Real a => Ptr CInt -> Ptr (Complex a) -> Ptr CInt -> IO a
BlasComplex.casum) ((Ptr CInt
 -> Ptr (Complex Double)
 -> Ptr CInt
 -> IO (RealOf (Complex Double)))
-> Nrm (Complex Double)
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr (Complex Double) -> Ptr CInt -> IO Double
Ptr CInt
-> Ptr (Complex Double) -> Ptr CInt -> IO (RealOf (Complex Double))
forall a. Real a => Ptr CInt -> Ptr (Complex a) -> Ptr CInt -> IO a
BlasComplex.casum)


{- |
Euclidean norm of a vector or Frobenius norm of a matrix.
-}
norm2 :: (Shape.C sh, Class.Floating a) => T shA sh a -> RealOf a
norm2 :: forall sh a shA. (C sh, Floating a) => T shA sh a -> RealOf a
norm2 T shA sh a
arr = IO (RealOf a) -> RealOf a
forall a. IO a -> a
unsafePerformIO (IO (RealOf a) -> RealOf a) -> IO (RealOf a) -> RealOf a
forall a b. (a -> b) -> a -> b
$
   ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (RealOf a) IO (RealOf a) -> IO (RealOf a))
-> ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$ IO (RealOf a) -> ContT (RealOf a) IO (RealOf a)
forall a. IO a -> ContT (RealOf a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (RealOf a) -> ContT (RealOf a) IO (RealOf a))
-> ((Ptr CInt, Ptr a, Ptr CInt) -> IO (RealOf a))
-> (Ptr CInt, Ptr a, Ptr CInt)
-> ContT (RealOf a) IO (RealOf a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a))
-> (Ptr CInt, Ptr a, Ptr CInt) -> IO (RealOf a)
forall a b c d. (a -> b -> c -> d) -> (a, b, c) -> d
uncurry3 Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
nrm2 ((Ptr CInt, Ptr a, Ptr CInt) -> ContT (RealOf a) IO (RealOf a))
-> ContT (RealOf a) IO (Ptr CInt, Ptr a, Ptr CInt)
-> ContT (RealOf a) IO (RealOf a)
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< T shA sh a -> ContT (RealOf a) IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
arr

nrm2 :: Class.Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
nrm2 :: forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
nrm2 =
   Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a. Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
getNrm (Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a))
-> Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$
   Nrm Float
-> Nrm Double
-> Nrm (Complex Float)
-> Nrm (Complex Double)
-> Nrm a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((Ptr CInt -> Ptr Float -> Ptr CInt -> IO (RealOf Float))
-> Nrm Float
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr Float -> Ptr CInt -> IO Float
Ptr CInt -> Ptr Float -> Ptr CInt -> IO (RealOf Float)
forall a. Real a => Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.nrm2) ((Ptr CInt -> Ptr Double -> Ptr CInt -> IO (RealOf Double))
-> Nrm Double
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr Double -> Ptr CInt -> IO Double
Ptr CInt -> Ptr Double -> Ptr CInt -> IO (RealOf Double)
forall a. Real a => Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.nrm2)
      ((Ptr CInt
 -> Ptr (Complex Float) -> Ptr CInt -> IO (RealOf (Complex Float)))
-> Nrm (Complex Float)
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr (Complex Float) -> Ptr CInt -> IO Float
Ptr CInt
-> Ptr (Complex Float) -> Ptr CInt -> IO (RealOf (Complex Float))
forall a. Real a => Ptr CInt -> Ptr (Complex a) -> Ptr CInt -> IO a
BlasComplex.cnrm2) ((Ptr CInt
 -> Ptr (Complex Double)
 -> Ptr CInt
 -> IO (RealOf (Complex Double)))
-> Nrm (Complex Double)
forall a. (Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)) -> Nrm a
Nrm Ptr CInt -> Ptr (Complex Double) -> Ptr CInt -> IO Double
Ptr CInt
-> Ptr (Complex Double) -> Ptr CInt -> IO (RealOf (Complex Double))
forall a. Real a => Ptr CInt -> Ptr (Complex a) -> Ptr CInt -> IO a
BlasComplex.cnrm2)

newtype Nrm a = Nrm {forall a. Nrm a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)
getNrm :: Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)}


newtype Norm f a = Norm {forall (f :: * -> *) a. Norm f a -> f a -> RealOf a
getNorm :: f a -> RealOf a}

norm2Squared :: (Shape.C sh, Class.Floating a) => T shA sh a -> RealOf a
norm2Squared :: forall sh a shA. (C sh, Floating a) => T shA sh a -> RealOf a
norm2Squared =
   Norm (T shA sh) a -> T shA sh a -> RealOf a
forall (f :: * -> *) a. Norm f a -> f a -> RealOf a
getNorm (Norm (T shA sh) a -> T shA sh a -> RealOf a)
-> Norm (T shA sh) a -> T shA sh a -> RealOf a
forall a b. (a -> b) -> a -> b
$
   Norm (T shA sh) Float
-> Norm (T shA sh) Double
-> Norm (T shA sh) (Complex Float)
-> Norm (T shA sh) (Complex Double)
-> Norm (T shA sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> RealOf Float) -> Norm (T shA sh) Float
forall (f :: * -> *) a. (f a -> RealOf a) -> Norm f a
Norm T shA sh Float -> Float
T shA sh Float -> RealOf Float
forall sh a shA. (C sh, Real a) => T shA sh a -> a
norm2SquaredReal)
      ((T shA sh Double -> RealOf Double) -> Norm (T shA sh) Double
forall (f :: * -> *) a. (f a -> RealOf a) -> Norm f a
Norm T shA sh Double -> Double
T shA sh Double -> RealOf Double
forall sh a shA. (C sh, Real a) => T shA sh a -> a
norm2SquaredReal)
      ((T shA sh (Complex Float) -> RealOf (Complex Float))
-> Norm (T shA sh) (Complex Float)
forall (f :: * -> *) a. (f a -> RealOf a) -> Norm f a
Norm T shA sh (Complex Float) -> Float
T shA sh (Complex Float) -> RealOf (Complex Float)
forall sh a shA. (C sh, Real a) => T shA sh (Complex a) -> a
norm2SquaredComplex)
      ((T shA sh (Complex Double) -> RealOf (Complex Double))
-> Norm (T shA sh) (Complex Double)
forall (f :: * -> *) a. (f a -> RealOf a) -> Norm f a
Norm T shA sh (Complex Double) -> Double
T shA sh (Complex Double) -> RealOf (Complex Double)
forall sh a shA. (C sh, Real a) => T shA sh (Complex a) -> a
norm2SquaredComplex)

norm2SquaredReal :: (Shape.C sh, Class.Real a) => T shA sh a -> a
norm2SquaredReal :: forall sh a shA. (C sh, Real a) => T shA sh a -> a
norm2SquaredReal T shA sh a
x =
   IO a -> a
forall a. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$ ContT a IO a -> IO a
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT a IO a -> IO a) -> ContT a IO a -> IO a
forall a b. (a -> b) -> a -> b
$ do
      (Ptr CInt
nPtr, Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> ContT a IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
x
      IO a -> ContT a IO a
forall a. IO a -> ContT a IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO a -> ContT a IO a) -> IO a -> ContT a IO a
forall a b. (a -> b) -> a -> b
$ Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
forall a.
Real a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.dot Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr Ptr a
sxPtr Ptr CInt
incxPtr

norm2SquaredComplex :: (Shape.C sh, Class.Real a) => T shA sh (Complex a) -> a
norm2SquaredComplex :: forall sh a shA. (C sh, Real a) => T shA sh (Complex a) -> a
norm2SquaredComplex T shA sh (Complex a)
x =
   IO a -> a
forall a. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$ ContT a IO a -> IO a
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT a IO a -> IO a) -> ContT a IO a -> IO a
forall a b. (a -> b) -> a -> b
$ do
      Ptr CInt
nPtr <- Int -> FortranIO a (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (Int -> FortranIO a (Ptr CInt)) -> Int -> FortranIO a (Ptr CInt)
forall a b. (a -> b) -> a -> b
$ sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh (Complex a)
x
      Ptr (Complex a)
xPtr <- T shA sh (Complex a) -> FortranIO a (Ptr (Complex a))
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shA sh (Complex a)
x
      let xrPtr :: Ptr (RealOf (Complex a))
xrPtr = Ptr (Complex a) -> Ptr (RealOf (Complex a))
forall a. Ptr a -> Ptr (RealOf a)
realPtr Ptr (Complex a)
xPtr
      let xiPtr :: Ptr a
xiPtr = Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
Ptr (RealOf (Complex a))
xrPtr Int
1
      Ptr CInt
incxPtr <- Int -> FortranIO a (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint (T shA sh (Complex a) -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh (Complex a)
x Int -> Int -> Int
forall a. Num a => a -> a -> a
* Int
2)
      IO a -> ContT a IO a
forall a. IO a -> ContT a IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO a -> ContT a IO a) -> IO a -> ContT a IO a
forall a b. (a -> b) -> a -> b
$
         (a -> a -> a) -> IO a -> IO a -> IO a
forall a b c. (a -> b -> c) -> IO a -> IO b -> IO c
forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 a -> a -> a
forall a. Num a => a -> a -> a
(+)
            (Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
forall a.
Real a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.dot Ptr CInt
nPtr Ptr a
Ptr (RealOf (Complex a))
xrPtr Ptr CInt
incxPtr Ptr a
Ptr (RealOf (Complex a))
xrPtr Ptr CInt
incxPtr)
            (Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
forall a.
Real a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO a
BlasReal.dot Ptr CInt
nPtr Ptr a
xiPtr Ptr CInt
incxPtr Ptr a
xiPtr Ptr CInt
incxPtr)


{- |
prop> forSliced number_ $ \xs -> VectorSlice.normInf xs == List.maximum (0 : List.map absolute (listFromSlice xs))
-}
normInf :: (Shape.C sh, Class.Floating a) => T shA sh a -> RealOf a
normInf :: forall sh a shA. (C sh, Floating a) => T shA sh a -> RealOf a
normInf T shA sh a
arr = IO (RealOf a) -> RealOf a
forall a. IO a -> a
unsafePerformIO (IO (RealOf a) -> RealOf a) -> IO (RealOf a) -> RealOf a
forall a b. (a -> b) -> a -> b
$
   (Maybe (Int, a) -> RealOf a)
-> IO (Maybe (Int, a)) -> IO (RealOf a)
forall a b. (a -> b) -> IO a -> IO b
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (a -> RealOf a
forall a. Floating a => a -> RealOf a
Scalar.absolute (a -> RealOf a)
-> (Maybe (Int, a) -> a) -> Maybe (Int, a) -> RealOf a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> ((Int, a) -> a) -> Maybe (Int, a) -> a
forall b a. b -> (a -> b) -> Maybe a -> b
maybe a
forall a. Floating a => a
Scalar.zero (Int, a) -> a
forall a b. (a, b) -> b
snd) (IO (Maybe (Int, a)) -> IO (RealOf a))
-> IO (Maybe (Int, a)) -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$ T shA sh a -> IO (Maybe (Int, a))
forall sh a shA.
(C sh, Floating a) =>
T shA sh a -> IO (Maybe (Int, a))
absMax T shA sh a
arr

{- |
Computes (almost) the infinity norm of the vector.
For complex numbers every element is replaced
by the sum of the absolute component values first.

prop> forSliced number_ $ \xs -> VectorSlice.normInf1 xs == List.maximum (0 : List.map Scalar.norm1 (listFromSlice xs))
-}
normInf1 :: (Shape.C sh, Class.Floating a) => T shA sh a -> RealOf a
normInf1 :: forall sh a shA. (C sh, Floating a) => T shA sh a -> RealOf a
normInf1 T shA sh a
x = IO (RealOf a) -> RealOf a
forall a. IO a -> a
unsafePerformIO (IO (RealOf a) -> RealOf a) -> IO (RealOf a) -> RealOf a
forall a b. (a -> b) -> a -> b
$
   ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (RealOf a) IO (RealOf a) -> IO (RealOf a))
-> ContT (RealOf a) IO (RealOf a) -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$ do
      (Ptr CInt
nPtr, Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> ContT (RealOf a) IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
x
      IO (RealOf a) -> ContT (RealOf a) IO (RealOf a)
forall a. IO a -> ContT (RealOf a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (RealOf a) -> ContT (RealOf a) IO (RealOf a))
-> IO (RealOf a) -> ContT (RealOf a) IO (RealOf a)
forall a b. (a -> b) -> a -> b
$
         (Maybe (Int, a) -> RealOf a)
-> IO (Maybe (Int, a)) -> IO (RealOf a)
forall a b. (a -> b) -> IO a -> IO b
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (a -> RealOf a
forall a. Floating a => a -> RealOf a
Scalar.norm1 (a -> RealOf a)
-> (Maybe (Int, a) -> a) -> Maybe (Int, a) -> RealOf a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> ((Int, a) -> a) -> Maybe (Int, a) -> a
forall b a. b -> (a -> b) -> Maybe a -> b
maybe a
forall a. Floating a => a
Scalar.zero (Int, a) -> a
forall a b. (a, b) -> b
snd) (IO (Maybe (Int, a)) -> IO (RealOf a))
-> IO (Maybe (Int, a)) -> IO (RealOf a)
forall a b. (a -> b) -> a -> b
$
         Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
forall a. Storable a => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 Ptr a
sxPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x) (CInt -> IO (Maybe (Int, a))) -> IO CInt -> IO (Maybe (Int, a))
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
forall a. Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
Blas.iamax Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr


{- |
Returns the index and value of the element with the maximal absolute value.
Caution: It actually returns the value of the element, not its absolute value!

prop> forSliced number_ $ \xs -> isNonEmpty xs ==> let (xi,xm) = VectorSlice.argAbsMaximum xs in VectorSlice.access xs xi == xm
prop> forSliced number_ $ \xs -> isNonEmpty xs ==> let (_xi,xm) = VectorSlice.argAbsMaximum xs in List.all (\x -> absolute x <= absolute xm) $ listFromSlice xs
prop> forSliced number_ $ \xs -> forSliced number_ $ \ys -> isNonEmpty xs && isNonEmpty ys ==> let (_xi,xm) = VectorSlice.argAbsMaximum xs; (_yi,ym) = VectorSlice.argAbsMaximum ys; (zi,zm) = Vector.argAbsMaximum (VectorSlice.toVector xs +++ VectorSlice.toVector ys) in case zi of Left _ -> xm==zm && absolute xm >= absolute ym; Right _ -> ym==zm && absolute xm < absolute ym
-}
argAbsMaximum ::
   (Shape.InvIndexed sh, Class.Floating a) =>
   T shA sh a -> (Shape.Index sh, a)
argAbsMaximum :: forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbsMaximum T shA sh a
arr = IO (Index sh, a) -> (Index sh, a)
forall a. IO a -> a
unsafePerformIO (IO (Index sh, a) -> (Index sh, a))
-> IO (Index sh, a) -> (Index sh, a)
forall a b. (a -> b) -> a -> b
$
   (Maybe (Int, a) -> (Index sh, a))
-> IO (Maybe (Int, a)) -> IO (Index sh, a)
forall a b. (a -> b) -> IO a -> IO b
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap
      ((Index sh, a)
-> ((Int, a) -> (Index sh, a)) -> Maybe (Int, a) -> (Index sh, a)
forall b a. b -> (a -> b) -> Maybe a -> b
maybe
         (String -> (Index sh, a)
forall a. HasCallStack => String -> a
error String
"Vector.argAbsMaximum: empty vector")
         ((Int -> Index sh) -> (Int, a) -> (Index sh, a)
forall a c b. (a -> c) -> (a, b) -> (c, b)
mapFst (sh -> Int -> Index sh
forall sh. InvIndexed sh => sh -> Int -> Index sh
Shape.uncheckedIndexFromOffset (sh -> Int -> Index sh) -> sh -> Int -> Index sh
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
arr))) (IO (Maybe (Int, a)) -> IO (Index sh, a))
-> IO (Maybe (Int, a)) -> IO (Index sh, a)
forall a b. (a -> b) -> a -> b
$
   T shA sh a -> IO (Maybe (Int, a))
forall sh a shA.
(C sh, Floating a) =>
T shA sh a -> IO (Maybe (Int, a))
absMax T shA sh a
arr

absMax ::
   (Shape.C sh, Class.Floating a) =>
   T shA sh a -> IO (Maybe (Int, a))
absMax :: forall sh a shA.
(C sh, Floating a) =>
T shA sh a -> IO (Maybe (Int, a))
absMax T shA sh a
x =
   case T shA sh a -> ComplexSingleton a
forall a (f :: * -> *). Floating a => f a -> ComplexSingleton a
Scalar.complexSingletonOfFunctor T shA sh a
x of
      ComplexSingleton a
Scalar.Real -> ContT (Maybe (Int, a)) IO (Maybe (Int, a)) -> IO (Maybe (Int, a))
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (Maybe (Int, a)) IO (Maybe (Int, a)) -> IO (Maybe (Int, a)))
-> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
-> IO (Maybe (Int, a))
forall a b. (a -> b) -> a -> b
$ do
         (Ptr CInt
nPtr, Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> ContT (Maybe (Int, a)) IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
x
         IO (Maybe (Int, a)) -> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
forall a. IO a -> ContT (Maybe (Int, a)) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Maybe (Int, a)) -> ContT (Maybe (Int, a)) IO (Maybe (Int, a)))
-> IO (Maybe (Int, a))
-> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
forall a b. (a -> b) -> a -> b
$
            Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
forall a. Storable a => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 Ptr a
sxPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x) (CInt -> IO (Maybe (Int, a))) -> IO CInt -> IO (Maybe (Int, a))
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
forall a. Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
Blas.iamax Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr
      ComplexSingleton a
Scalar.Complex -> ContT (Maybe (Int, a)) IO (Maybe (Int, a)) -> IO (Maybe (Int, a))
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (Maybe (Int, a)) IO (Maybe (Int, a)) -> IO (Maybe (Int, a)))
-> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
-> IO (Maybe (Int, a))
forall a b. (a -> b) -> a -> b
$ do
         let n :: Int
n = sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x
         Ptr a
sxPtr <- T shA sh a -> FortranIO (Maybe (Int, a)) (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shA sh a
x
         let incx :: Int
incx = T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x
         IO (Maybe (Int, a)) -> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
forall a. IO a -> ContT (Maybe (Int, a)) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Maybe (Int, a)) -> ContT (Maybe (Int, a)) IO (Maybe (Int, a)))
-> IO (Maybe (Int, a))
-> ContT (Maybe (Int, a)) IO (Maybe (Int, a))
forall a b. (a -> b) -> a -> b
$ Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
forall a. Storable a => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 Ptr a
sxPtr Int
incx (CInt -> IO (Maybe (Int, a))) -> IO CInt -> IO (Maybe (Int, a))
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Int -> Ptr (Complex a1) -> Int -> IO CInt
forall a. Real a => Int -> Ptr (Complex a) -> Int -> IO CInt
absMaxComplex Int
n Ptr a
Ptr (Complex a1)
sxPtr Int
incx

absMaxComplex :: (Class.Real a) => Int -> Ptr (Complex a) -> Int -> IO CInt
absMaxComplex :: forall a. Real a => Int -> Ptr (Complex a) -> Int -> IO CInt
absMaxComplex Int
n Ptr (Complex a)
sxPtr Int
incx =
   Int -> (Ptr a -> IO CInt) -> IO CInt
forall a b. Storable a => Int -> (Ptr a -> IO b) -> IO b
ForeignArray.alloca Int
n ((Ptr a -> IO CInt) -> IO CInt) -> (Ptr a -> IO CInt) -> IO CInt
forall a b. (a -> b) -> a -> b
$ \Ptr a
syPtr -> do
      let xrPtr :: Ptr (RealOf (Complex a))
xrPtr = Ptr (Complex a) -> Ptr (RealOf (Complex a))
forall a. Ptr a -> Ptr (RealOf a)
realPtr Ptr (Complex a)
sxPtr
      let incx2 :: Int
incx2 = Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
incx
      Conjugation
-> Int -> Ptr a -> Int -> Ptr a -> Int -> Ptr a -> Int -> IO ()
forall a.
Floating a =>
Conjugation
-> Int -> Ptr a -> Int -> Ptr a -> Int -> Ptr a -> Int -> IO ()
Private.mul    Conjugation
NonConjugated Int
n Ptr a
Ptr (RealOf (Complex a))
xrPtr Int
incx2 Ptr a
Ptr (RealOf (Complex a))
xrPtr Int
incx2 Ptr a
syPtr Int
1
      let xiPtr :: Ptr a
xiPtr = Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
Ptr (RealOf (Complex a))
xrPtr Int
1
      Conjugation
-> Int
-> Ptr a
-> Int
-> Ptr a
-> Int
-> a
-> Ptr a
-> Int
-> IO ()
forall a.
Floating a =>
Conjugation
-> Int
-> Ptr a
-> Int
-> Ptr a
-> Int
-> a
-> Ptr a
-> Int
-> IO ()
Private.mulAdd Conjugation
NonConjugated Int
n Ptr a
xiPtr Int
incx2 Ptr a
xiPtr Int
incx2 a
forall a. Floating a => a
Scalar.one Ptr a
syPtr Int
1
      ContT CInt IO CInt -> IO CInt
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT CInt IO CInt -> IO CInt) -> ContT CInt IO CInt -> IO CInt
forall a b. (a -> b) -> a -> b
$ do
         Ptr CInt
nPtr <- Int -> FortranIO CInt (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
         Ptr CInt
incyPtr <- Int -> FortranIO CInt (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
         IO CInt -> ContT CInt IO CInt
forall a. IO a -> ContT CInt IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO CInt -> ContT CInt IO CInt) -> IO CInt -> ContT CInt IO CInt
forall a b. (a -> b) -> a -> b
$ Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
forall a. Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
Blas.iamax Ptr CInt
nPtr Ptr a
syPtr Ptr CInt
incyPtr


{- |
Returns the index and value of the element with the maximal absolute value.
The function does not strictly compare the absolute value of a complex number
but the sum of the absolute complex components.
Caution: It actually returns the value of the element, not its absolute value!

prop> forSliced real_ $ \xs -> isNonEmpty xs ==> VectorSlice.argAbsMaximum xs == VectorSlice.argAbs1Maximum xs
-}
argAbs1Maximum ::
   (Shape.InvIndexed sh, Class.Floating a) =>
   T shA sh a -> (Shape.Index sh, a)
argAbs1Maximum :: forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbs1Maximum T shA sh a
x = IO (Index sh, a) -> (Index sh, a)
forall a. IO a -> a
unsafePerformIO (IO (Index sh, a) -> (Index sh, a))
-> IO (Index sh, a) -> (Index sh, a)
forall a b. (a -> b) -> a -> b
$
   ContT (Index sh, a) IO (Index sh, a) -> IO (Index sh, a)
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT (Index sh, a) IO (Index sh, a) -> IO (Index sh, a))
-> ContT (Index sh, a) IO (Index sh, a) -> IO (Index sh, a)
forall a b. (a -> b) -> a -> b
$ do
      (Ptr CInt
nPtr, Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> ContT (Index sh, a) IO (Ptr CInt, Ptr a, Ptr CInt)
forall sh a shA r.
(C sh, Storable a) =>
T shA sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt)
sizeSliceArg T shA sh a
x
      IO (Index sh, a) -> ContT (Index sh, a) IO (Index sh, a)
forall a. IO a -> ContT (Index sh, a) IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Index sh, a) -> ContT (Index sh, a) IO (Index sh, a))
-> IO (Index sh, a) -> ContT (Index sh, a) IO (Index sh, a)
forall a b. (a -> b) -> a -> b
$
         (Maybe (Int, a) -> (Index sh, a))
-> IO (Maybe (Int, a)) -> IO (Index sh, a)
forall a b. (a -> b) -> IO a -> IO b
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap
            ((Index sh, a)
-> ((Int, a) -> (Index sh, a)) -> Maybe (Int, a) -> (Index sh, a)
forall b a. b -> (a -> b) -> Maybe a -> b
maybe
               (String -> (Index sh, a)
forall a. HasCallStack => String -> a
error String
"Vector.argAbs1Maximum: empty vector")
               ((Int -> Index sh) -> (Int, a) -> (Index sh, a)
forall a c b. (a -> c) -> (a, b) -> (c, b)
mapFst (sh -> Int -> Index sh
forall sh. InvIndexed sh => sh -> Int -> Index sh
Shape.uncheckedIndexFromOffset (sh -> Int -> Index sh) -> sh -> Int -> Index sh
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x))) (IO (Maybe (Int, a)) -> IO (Index sh, a))
-> IO (Maybe (Int, a)) -> IO (Index sh, a)
forall a b. (a -> b) -> a -> b
$
         Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
forall a. Storable a => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 Ptr a
sxPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x) (CInt -> IO (Maybe (Int, a))) -> IO CInt -> IO (Maybe (Int, a))
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
forall a. Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO CInt
Blas.iamax Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr

peekElemOff1 :: (Storable a) => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 :: forall a. Storable a => Ptr a -> Int -> CInt -> IO (Maybe (Int, a))
peekElemOff1 Ptr a
ptr Int
inc CInt
k1 =
   let k1i :: Int
k1i = CInt -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral CInt
k1
       ki :: Int
ki = Int
k1iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1
   in if Int
k1i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0
         then Maybe (Int, a) -> IO (Maybe (Int, a))
forall a. a -> IO a
forall (m :: * -> *) a. Monad m => a -> m a
return Maybe (Int, a)
forall a. Maybe a
Nothing
         else (Int, a) -> Maybe (Int, a)
forall a. a -> Maybe a
Just ((Int, a) -> Maybe (Int, a))
-> (a -> (Int, a)) -> a -> Maybe (Int, a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (,) Int
ki (a -> Maybe (Int, a)) -> IO a -> IO (Maybe (Int, a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Ptr a -> Int -> IO a
forall a. Storable a => Ptr a -> Int -> IO a
peekElemOff Ptr a
ptr (Int
kiInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
inc)


{- |
prop> QC.forAll genShape $ \sh@(_::+(_rows,columns)::+_) -> QC.forAll (QC.elements (Shape.indices columns)) $ \c -> QC.forAll (genVector sh $ genNumber 3) $ \xt -> let xs = takeColumn c xt in approx 1e-2 (VectorSlice.product xs) (List.product (listFromSlice (xs :: Sliced Number_)))
-}
product :: (Shape.C sh, Class.Floating a) => T shA sh a -> a
product :: forall sh a shA. (C sh, Floating a) => T shA sh a -> a
product T shA sh a
x = IO a -> a
forall a. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$ ContT a IO a -> IO a
forall (m :: * -> *) r. Monad m => ContT r m r -> m r
evalContT (ContT a IO a -> IO a) -> ContT a IO a -> IO a
forall a b. (a -> b) -> a -> b
$ do
   Ptr a
xPtr <- T shA sh a -> FortranIO a (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shA sh a
x
   IO a -> ContT a IO a
forall a. IO a -> ContT a IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO a -> ContT a IO a) -> IO a -> ContT a IO a
forall a b. (a -> b) -> a -> b
$ Int -> Ptr a -> Int -> IO a
forall a. Floating a => Int -> Ptr a -> Int -> IO a
Private.product (sh -> Int
forall sh. C sh => sh -> Int
Shape.size (sh -> Int) -> sh -> Int
forall a b. (a -> b) -> a -> b
$ T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) Ptr a
xPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
x)


{- |
For restrictions see 'limits'.

prop> forSliced real_ $ \xs -> isNonEmpty xs ==> VectorSlice.minimum xs == List.minimum (listFromSlice xs)
prop> forSliced real_ $ \xs -> isNonEmpty xs ==> VectorSlice.maximum xs == List.maximum (listFromSlice xs)
-}
minimum, maximum :: (Shape.C shA, Shape.C sh, Class.Real a) => T shA sh a -> a
minimum :: forall shA sh a. (C shA, C sh, Real a) => T shA sh a -> a
minimum = (a, a) -> a
forall a b. (a, b) -> a
fst ((a, a) -> a) -> (T shA sh a -> (a, a)) -> T shA sh a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh a -> (a, a)
forall shA sh a. (C shA, C sh, Real a) => T shA sh a -> (a, a)
limits
maximum :: forall shA sh a. (C shA, C sh, Real a) => T shA sh a -> a
maximum = (a, a) -> a
forall a b. (a, b) -> b
snd ((a, a) -> a) -> (T shA sh a -> (a, a)) -> T shA sh a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh a -> (a, a)
forall shA sh a. (C shA, C sh, Real a) => T shA sh a -> (a, a)
limits

{- |
For restrictions see 'limits'.
-}
argMinimum, argMaximum ::
   (Shape.C shA, Shape.InvIndexed sh, Shape.Index sh ~ ix, Class.Real a) =>
   T shA sh a -> (ix,a)
argMinimum :: forall shA sh ix a.
(C shA, InvIndexed sh, Index sh ~ ix, Real a) =>
T shA sh a -> (ix, a)
argMinimum = ((ix, a), (ix, a)) -> (ix, a)
forall a b. (a, b) -> a
fst (((ix, a), (ix, a)) -> (ix, a))
-> (T shA sh a -> ((ix, a), (ix, a))) -> T shA sh a -> (ix, a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh a -> ((ix, a), (ix, a))
forall shA sh ix a.
(C shA, InvIndexed sh, Index sh ~ ix, Real a) =>
T shA sh a -> ((ix, a), (ix, a))
argLimits
argMaximum :: forall shA sh ix a.
(C shA, InvIndexed sh, Index sh ~ ix, Real a) =>
T shA sh a -> (ix, a)
argMaximum = ((ix, a), (ix, a)) -> (ix, a)
forall a b. (a, b) -> b
snd (((ix, a), (ix, a)) -> (ix, a))
-> (T shA sh a -> ((ix, a), (ix, a))) -> T shA sh a -> (ix, a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. T shA sh a -> ((ix, a), (ix, a))
forall shA sh ix a.
(C shA, InvIndexed sh, Index sh ~ ix, Real a) =>
T shA sh a -> ((ix, a), (ix, a))
argLimits

{- |
prop> forSliced real_ $ \xs -> isNonEmpty xs ==> VectorSlice.limits xs == Array.limits (VectorSlice.toVector xs)

In contrast to 'Array.limits'
this implementation is based on fast BLAS functions.
It should be faster than @Array.minimum@ and @Array.maximum@
although it is certainly not as fast as possible.
Boths limits share the precision of the limit with the larger absolute value.
This implies for example that you cannot rely on the property
that @raise (- minimum x) x@ has only non-negative elements.
-}
limits :: (Shape.C shA, Shape.C sh, Class.Real a) => T shA sh a -> (a,a)
limits :: forall shA sh a. (C shA, C sh, Real a) => T shA sh a -> (a, a)
limits T shA sh a
xs0 =
   let xs :: T shA (Deferred sh) a
xs = (sh -> Deferred sh) -> T shA sh a -> T shA (Deferred sh) a
forall slice0 slice1 sh a.
(slice0 -> slice1) -> T sh slice0 a -> T sh slice1 a
mapShape sh -> Deferred sh
forall sh. sh -> Deferred sh
Shape.Deferred T shA sh a
xs0
       x0 :: a
x0 = (DeferredIndex sh, a) -> a
forall a b. (a, b) -> b
snd ((DeferredIndex sh, a) -> a) -> (DeferredIndex sh, a) -> a
forall a b. (a -> b) -> a -> b
$ T shA (Deferred sh) a -> (Index (Deferred sh), a)
forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbs1Maximum T shA (Deferred sh) a
xs
       x1 :: a
x1 = T shA (Deferred sh) a
xs T shA (Deferred sh) a -> Index (Deferred sh) -> a
forall shA sh a.
(C shA, Indexed sh, Storable a) =>
T shA sh a -> Index sh -> a
! (DeferredIndex sh, a) -> DeferredIndex sh
forall a b. (a, b) -> a
fst (T (Deferred sh) (Deferred sh) a -> (Index (Deferred sh), a)
forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbs1Maximum (Vector (Deferred sh) a -> T (Deferred sh) (Deferred sh) a
forall sh a. C sh => Vector sh a -> T sh sh a
fromVector (a -> T shA (Deferred sh) a -> Vector (Deferred sh) a
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
raise (-a
x0) T shA (Deferred sh) a
xs)))
   in if a
x0a -> a -> Bool
forall a. Ord a => a -> a -> Bool
>=a
0 then (a
x1,a
x0) else (a
x0,a
x1)

{- |
For restrictions see 'limits'.
-}
argLimits ::
   (Shape.C shA, Shape.InvIndexed sh, Shape.Index sh ~ ix, Class.Real a) =>
   T shA sh a -> ((ix,a),(ix,a))
argLimits :: forall shA sh ix a.
(C shA, InvIndexed sh, Index sh ~ ix, Real a) =>
T shA sh a -> ((ix, a), (ix, a))
argLimits T shA sh a
xs =
   let p0 :: (Index sh, a)
p0@(Index sh
_i0,a
x0) = T shA sh a -> (Index sh, a)
forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbs1Maximum T shA sh a
xs
       p1 :: (ix, a)
p1 = (ix
i1,T shA sh a
xsT shA sh a -> Index sh -> a
forall shA sh a.
(C shA, Indexed sh, Storable a) =>
T shA sh a -> Index sh -> a
!ix
Index sh
i1); i1 :: ix
i1 = (ix, a) -> ix
forall a b. (a, b) -> a
fst ((ix, a) -> ix) -> (ix, a) -> ix
forall a b. (a -> b) -> a -> b
$ T sh sh a -> (Index sh, a)
forall sh a shA.
(InvIndexed sh, Floating a) =>
T shA sh a -> (Index sh, a)
argAbs1Maximum (T sh sh a -> (Index sh, a)) -> T sh sh a -> (Index sh, a)
forall a b. (a -> b) -> a -> b
$ Vector sh a -> T sh sh a
forall sh a. C sh => Vector sh a -> T sh sh a
fromVector (Vector sh a -> T sh sh a) -> Vector sh a -> T sh sh a
forall a b. (a -> b) -> a -> b
$ a -> T shA sh a -> Vector sh a
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
raise (-a
x0) T shA sh a
xs
   in if a
x0a -> a -> Bool
forall a. Ord a => a -> a -> Bool
>=a
0 then ((ix, a)
p1,(ix, a)
(Index sh, a)
p0) else ((ix, a)
(Index sh, a)
p0,(ix, a)
p1)


{- |
prop> forSliced number_ $ \xs -> VectorSlice.negate xs == VectorSlice.scale minusOne xs
prop> forSliced number_ $ \xs -> VectorSlice.scale 2 xs == VectorSlice.add xs xs
-}
scale, _scale ::
   (Shape.C sh, Class.Floating a) =>
   a -> T shA sh a -> Vector sh a

scale :: forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
scale a
alpha T shA sh a
x = sh -> (Int -> Ptr a -> IO ()) -> Array sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) ((Int -> Ptr a -> IO ()) -> Array sh a)
-> (Int -> Ptr a -> IO ()) -> Array sh a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
syPtr -> do
   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 a
alphaPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
alpha
      Ptr CInt
nPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      (Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x
      Ptr CInt
incyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 ()
Blas.copy Ptr CInt
nPtr Ptr a
sxPtr Ptr CInt
incxPtr Ptr a
syPtr Ptr CInt
incyPtr
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 a -> Ptr CInt -> IO ()
forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr a -> Ptr CInt -> IO ()
Blas.scal Ptr CInt
nPtr Ptr a
alphaPtr Ptr a
syPtr Ptr CInt
incyPtr

_scale :: forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
_scale a
a T shA sh a
b = sh -> (Int -> Ptr a -> IO ()) -> Array sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
b) ((Int -> Ptr a -> IO ()) -> Array sh a)
-> (Int -> Ptr a -> IO ()) -> Array sh a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
cPtr -> do
   let m :: Int
m = Int
1
   let k :: Int
k = 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
transaPtr <- Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'N'
      Ptr CChar
transbPtr <- Char -> FortranIO () (Ptr CChar)
forall r. Char -> FortranIO r (Ptr CChar)
Call.char Char
'N'
      Ptr CInt
mPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
m
      Ptr CInt
kPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
k
      Ptr CInt
nPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      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
Scalar.one
      Ptr a
aPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
a
      Ptr CInt
ldaPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
      (Ptr a
bPtr, Ptr CInt
ldbPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
b
      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
Scalar.zero
      Ptr CInt
ldcPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.leadingDim Int
m
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 CChar
-> 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 CChar
-> Ptr CInt
-> Ptr CInt
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr CInt
-> Ptr a
-> Ptr a
-> Ptr CInt
-> IO ()
Blas.gemm
            Ptr CChar
transaPtr Ptr CChar
transbPtr Ptr CInt
mPtr Ptr CInt
nPtr Ptr CInt
kPtr Ptr a
alphaPtr
            Ptr a
aPtr Ptr CInt
ldaPtr Ptr a
bPtr Ptr CInt
ldbPtr Ptr a
betaPtr Ptr a
cPtr Ptr CInt
ldcPtr


{- |
Complex implementation requires double number of multiplications
compared to 'Numeric.BLAS.Vector.scaleReal'.
-}
scaleReal ::
   (Shape.C sh, Class.Floating a) =>
   RealOf a -> T shA sh a -> Vector sh a
scaleReal :: forall sh a shA.
(C sh, Floating a) =>
RealOf a -> T shA sh a -> Vector sh a
scaleReal =
   ScaleReal (T shA sh) (Array sh) a
-> RealOf a -> T shA sh a -> Array sh a
forall (f :: * -> *) (g :: * -> *) a.
ScaleReal f g a -> RealOf a -> f a -> g a
getScaleReal (ScaleReal (T shA sh) (Array sh) a
 -> RealOf a -> T shA sh a -> Array sh a)
-> ScaleReal (T shA sh) (Array sh) a
-> RealOf a
-> T shA sh a
-> Array sh a
forall a b. (a -> b) -> a -> b
$
   ScaleReal (T shA sh) (Array sh) Float
-> ScaleReal (T shA sh) (Array sh) Double
-> ScaleReal (T shA sh) (Array sh) (Complex Float)
-> ScaleReal (T shA sh) (Array sh) (Complex Double)
-> ScaleReal (T shA sh) (Array sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((RealOf Float -> T shA sh Float -> Array sh Float)
-> ScaleReal (T shA sh) (Array sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(RealOf a -> f a -> g a) -> ScaleReal f g a
ScaleReal Float -> T shA sh Float -> Array sh Float
RealOf Float -> T shA sh Float -> Array sh Float
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
scale)
      ((RealOf Double -> T shA sh Double -> Array sh Double)
-> ScaleReal (T shA sh) (Array sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(RealOf a -> f a -> g a) -> ScaleReal f g a
ScaleReal Double -> T shA sh Double -> Array sh Double
RealOf Double -> T shA sh Double -> Array sh Double
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
scale)
      ((RealOf (Complex Float)
 -> T shA sh (Complex Float) -> Array sh (Complex Float))
-> ScaleReal (T shA sh) (Array sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(RealOf a -> f a -> g a) -> ScaleReal f g a
ScaleReal ((RealOf (Complex Float)
  -> T shA sh (Complex Float) -> Array sh (Complex Float))
 -> ScaleReal (T shA sh) (Array sh) (Complex Float))
-> (RealOf (Complex Float)
    -> T shA sh (Complex Float) -> Array sh (Complex Float))
-> ScaleReal (T shA sh) (Array sh) (Complex Float)
forall a b. (a -> b) -> a -> b
$ Complex Float
-> T shA sh (Complex Float) -> Array sh (Complex Float)
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
scale (Complex Float
 -> T shA sh (Complex Float) -> Array sh (Complex Float))
-> (Float -> Complex Float)
-> Float
-> T shA sh (Complex Float)
-> Array sh (Complex Float)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Float -> Complex Float
RealOf (Complex Float) -> Complex Float
forall a. Floating a => RealOf a -> a
Scalar.fromReal)
      ((RealOf (Complex Double)
 -> T shA sh (Complex Double) -> Array sh (Complex Double))
-> ScaleReal (T shA sh) (Array sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(RealOf a -> f a -> g a) -> ScaleReal f g a
ScaleReal ((RealOf (Complex Double)
  -> T shA sh (Complex Double) -> Array sh (Complex Double))
 -> ScaleReal (T shA sh) (Array sh) (Complex Double))
-> (RealOf (Complex Double)
    -> T shA sh (Complex Double) -> Array sh (Complex Double))
-> ScaleReal (T shA sh) (Array sh) (Complex Double)
forall a b. (a -> b) -> a -> b
$ Complex Double
-> T shA sh (Complex Double) -> Array sh (Complex Double)
forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
scale (Complex Double
 -> T shA sh (Complex Double) -> Array sh (Complex Double))
-> (Double -> Complex Double)
-> Double
-> T shA sh (Complex Double)
-> Array sh (Complex Double)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> Complex Double
RealOf (Complex Double) -> Complex Double
forall a. Floating a => RealOf a -> a
Scalar.fromReal)

newtype ScaleReal f g a = ScaleReal {forall (f :: * -> *) (g :: * -> *) a.
ScaleReal f g a -> RealOf a -> f a -> g a
getScaleReal :: RealOf a -> f a -> g a}



infixl 6 `add`, `sub`


{- |
prop> forSliced2 number_ $ \xs ys -> VectorSlice.add xs ys == VectorSlice.add ys xs
prop> forSliced2 number_ $ \xs ys -> VectorSlice.toVector xs == VectorSlice.sub xs ys |+| VectorSlice.toVector ys
-}
add, sub ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   T shA sh a -> T shB sh a -> Vector sh a
add :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> Vector sh a
add = a -> T shA sh a -> T shB sh a -> Vector sh a
forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
a -> T shA sh a -> T shB sh a -> Vector sh a
mac a
forall a. Floating a => a
Scalar.one
sub :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> Vector sh a
sub T shA sh a
x T shB sh a
y = a -> T shB sh a -> T shA sh a -> Vector sh a
forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
a -> T shA sh a -> T shB sh a -> Vector sh a
mac a
forall a. Floating a => a
Scalar.minusOne T shB sh a
y T shA sh a
x

mac ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   a -> T shA sh a -> T shB sh a -> Vector sh a
mac :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
a -> T shA sh a -> T shB sh a -> Vector sh a
mac a
alpha T shA sh a
x T shB sh a
y =
   sh -> (Int -> Ptr a -> IO ()) -> Array sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) ((Int -> Ptr a -> IO ()) -> Array sh a)
-> (Int -> Ptr a -> IO ()) -> Array sh a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
szPtr -> do
   String -> Bool -> IO ()
Call.assert String
"mac: shapes mismatch" (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== T shB sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh a
y)
   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 -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      Ptr a
saPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
alpha
      (Ptr a
sxPtr, Ptr CInt
incxPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x
      (Ptr a
syPtr, Ptr CInt
incyPtr) <- T shB sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shB sh a
y
      Ptr CInt
inczPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 ()
Blas.copy Ptr CInt
nPtr Ptr a
syPtr Ptr CInt
incyPtr Ptr a
szPtr Ptr CInt
inczPtr
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr a -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
Blas.axpy Ptr CInt
nPtr Ptr a
saPtr Ptr a
sxPtr Ptr CInt
incxPtr Ptr a
szPtr Ptr CInt
inczPtr


{- |
prop> forSliced number_ $ \xs -> VectorSlice.toVector xs == Vector.negate (VectorSlice.negate xs)
-}
negate :: (Shape.C sh, Class.Floating a) => T shA sh a -> Vector sh a
negate :: forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
negate =
   Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a
forall (f :: * -> *) (g :: * -> *) a. Conjugate f g a -> f a -> g a
getConjugate (Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a)
-> Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a
forall a b. (a -> b) -> a -> b
$
   Conjugate (T shA sh) (Array sh) Float
-> Conjugate (T shA sh) (Array sh) Double
-> Conjugate (T shA sh) (Array sh) (Complex Float)
-> Conjugate (T shA sh) (Array sh) (Complex Double)
-> Conjugate (T shA sh) (Array sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> Array sh Float)
-> Conjugate (T shA sh) (Array sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate ((T shA sh Float -> Array sh Float)
 -> Conjugate (T shA sh) (Array sh) Float)
-> (T shA sh Float -> Array sh Float)
-> Conjugate (T shA sh) (Array sh) Float
forall a b. (a -> b) -> a -> b
$ RealOf Float -> T shA sh Float -> Array sh Float
forall sh a shA.
(C sh, Floating a) =>
RealOf a -> T shA sh a -> Vector sh a
scaleReal Float
RealOf Float
forall a. Floating a => a
Scalar.minusOne)
      ((T shA sh Double -> Array sh Double)
-> Conjugate (T shA sh) (Array sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate ((T shA sh Double -> Array sh Double)
 -> Conjugate (T shA sh) (Array sh) Double)
-> (T shA sh Double -> Array sh Double)
-> Conjugate (T shA sh) (Array sh) Double
forall a b. (a -> b) -> a -> b
$ RealOf Double -> T shA sh Double -> Array sh Double
forall sh a shA.
(C sh, Floating a) =>
RealOf a -> T shA sh a -> Vector sh a
scaleReal Double
RealOf Double
forall a. Floating a => a
Scalar.minusOne)
      ((T shA sh (Complex Float) -> Array sh (Complex Float))
-> Conjugate (T shA sh) (Array sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate ((T shA sh (Complex Float) -> Array sh (Complex Float))
 -> Conjugate (T shA sh) (Array sh) (Complex Float))
-> (T shA sh (Complex Float) -> Array sh (Complex Float))
-> Conjugate (T shA sh) (Array sh) (Complex Float)
forall a b. (a -> b) -> a -> b
$ RealOf (Complex Float)
-> T shA sh (Complex Float) -> Array sh (Complex Float)
forall sh a shA.
(C sh, Floating a) =>
RealOf a -> T shA sh a -> Vector sh a
scaleReal Float
RealOf (Complex Float)
forall a. Floating a => a
Scalar.minusOne)
      ((T shA sh (Complex Double) -> Array sh (Complex Double))
-> Conjugate (T shA sh) (Array sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate ((T shA sh (Complex Double) -> Array sh (Complex Double))
 -> Conjugate (T shA sh) (Array sh) (Complex Double))
-> (T shA sh (Complex Double) -> Array sh (Complex Double))
-> Conjugate (T shA sh) (Array sh) (Complex Double)
forall a b. (a -> b) -> a -> b
$ RealOf (Complex Double)
-> T shA sh (Complex Double) -> Array sh (Complex Double)
forall sh a shA.
(C sh, Floating a) =>
RealOf a -> T shA sh a -> Vector sh a
scaleReal Double
RealOf (Complex Double)
forall a. Floating a => a
Scalar.minusOne)


{- |
prop> QC.forAll (genNumber maxElem) $ \d -> forSliced number_ $ \xs -> VectorSlice.toVector xs == Vector.raise (-d) (VectorSlice.raise d xs)
-}
raise :: (Shape.C sh, Class.Floating a) => a -> T shA sh a -> Vector sh a
raise :: forall sh a shA.
(C sh, Floating a) =>
a -> T shA sh a -> Vector sh a
raise a
c T shA sh a
x =
   sh -> (Int -> Ptr a -> IO ()) -> Array sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) ((Int -> Ptr a -> IO ()) -> Array sh a)
-> (Int -> Ptr a -> IO ()) -> Array 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 -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      Ptr a
cPtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
c
      Ptr a
onePtr <- a -> FortranIO () (Ptr a)
forall a r. Floating a => a -> FortranIO r (Ptr a)
Call.number a
forall a. Floating a => a
Scalar.one
      Ptr CInt
inccPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
0
      (Ptr a
xPtr, Ptr CInt
incxPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x
      Ptr CInt
inc1Ptr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
         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 ()
Blas.copy Ptr CInt
nPtr Ptr a
xPtr Ptr CInt
incxPtr Ptr a
yPtr Ptr CInt
inc1Ptr
         Ptr CInt
-> Ptr a -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
forall a.
Floating a =>
Ptr CInt
-> Ptr a -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
Blas.axpy Ptr CInt
nPtr Ptr a
onePtr Ptr a
cPtr Ptr CInt
inccPtr Ptr a
yPtr Ptr CInt
inc1Ptr


{- |
prop> forSliced2 number_ $ \xs ys -> VectorSlice.mul xs ys == VectorSlice.mul ys xs
-}
mul ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   T shA sh a -> T shB sh a -> Vector sh a
mul :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> Vector sh a
mul = Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
mulConjugation Conjugation
NonConjugated

{- |
prop> forSliced2 number_ $ \xs ys -> VectorSlice.mulConj xs ys == Vector.mul (VectorSlice.conjugate xs) (VectorSlice.toVector ys)
-}
mulConj ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   T shA sh a -> T shB sh a -> Vector sh a
mulConj :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
T shA sh a -> T shB sh a -> Vector sh a
mulConj = Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
mulConjugation Conjugation
Conjugated

mulConjugation ::
   (Shape.C sh, Eq sh, Class.Floating a) =>
   Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
mulConjugation :: forall sh a shA shB.
(C sh, Eq sh, Floating a) =>
Conjugation -> T shA sh a -> T shB sh a -> Vector sh a
mulConjugation Conjugation
conj T shA sh a
a T shB sh a
x =
      sh -> (Int -> Ptr a -> IO ()) -> Array sh a
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shB sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh a
x) ((Int -> Ptr a -> IO ()) -> Array sh a)
-> (Int -> Ptr a -> IO ()) -> Array sh a
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr a
yPtr -> do
   String -> Bool -> IO ()
Call.assert String
"mul: shapes mismatch" (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
a sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== T shB sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh a
x)
   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 a
aPtr <- T shA sh a -> FortranIO () (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shA sh a
a
      Ptr a
xPtr <- T shB sh a -> FortranIO () (Ptr a)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a)
startArg T shB sh a
x
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ Conjugation
-> Int -> Ptr a -> Int -> Ptr a -> Int -> Ptr a -> Int -> IO ()
forall a.
Floating a =>
Conjugation
-> Int -> Ptr a -> Int -> Ptr a -> Int -> Ptr a -> Int -> IO ()
Private.mul Conjugation
conj Int
n Ptr a
aPtr (T shA sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shA sh a
a) Ptr a
xPtr (T shB sh a -> Int
forall sh slice a. T sh slice a -> Int
increment T shB sh a
x) Ptr a
yPtr Int
1


newtype Conjugate f g a = Conjugate {forall (f :: * -> *) (g :: * -> *) a. Conjugate f g a -> f a -> g a
getConjugate :: f a -> g a}

conjugate ::
   (Shape.C sh, Class.Floating a) =>
   T shA sh a -> Vector sh a
conjugate :: forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
conjugate =
   Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a
forall (f :: * -> *) (g :: * -> *) a. Conjugate f g a -> f a -> g a
getConjugate (Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a)
-> Conjugate (T shA sh) (Array sh) a -> T shA sh a -> Array sh a
forall a b. (a -> b) -> a -> b
$
   Conjugate (T shA sh) (Array sh) Float
-> Conjugate (T shA sh) (Array sh) Double
-> Conjugate (T shA sh) (Array sh) (Complex Float)
-> Conjugate (T shA sh) (Array sh) (Complex Double)
-> Conjugate (T shA sh) (Array sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> Array sh Float)
-> Conjugate (T shA sh) (Array sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate T shA sh Float -> Array sh Float
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh Double -> Array sh Double)
-> Conjugate (T shA sh) (Array sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate T shA sh Double -> Array sh Double
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh (Complex Float) -> Array sh (Complex Float))
-> Conjugate (T shA sh) (Array sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate T shA sh (Complex Float) -> Array sh (Complex Float)
forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> Vector sh (Complex a)
complexConjugate)
      ((T shA sh (Complex Double) -> Array sh (Complex Double))
-> Conjugate (T shA sh) (Array sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g a) -> Conjugate f g a
Conjugate T shA sh (Complex Double) -> Array sh (Complex Double)
forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> Vector sh (Complex a)
complexConjugate)

complexConjugate ::
   (Shape.C sh, Class.Real a) =>
   T shA sh (Complex a) -> Vector sh (Complex a)
complexConjugate :: forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> Vector sh (Complex a)
complexConjugate T shA sh (Complex a)
x = sh -> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh (Complex a) -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh (Complex a)
x) ((Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a))
-> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr (Complex a)
syPtr ->
   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 -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      (Ptr (Complex a)
sxPtr, Ptr CInt
incxPtr) <- T shA sh (Complex a) -> FortranIO () (Ptr (Complex a), Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh (Complex a)
x
      Ptr CInt
incyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
1
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
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 (Complex a)
-> Ptr CInt
-> Ptr (Complex a)
-> Ptr CInt
-> IO ()
forall a.
Floating a =>
Ptr CInt -> Ptr a -> Ptr CInt -> Ptr a -> Ptr CInt -> IO ()
copyConjugate Ptr CInt
nPtr Ptr (Complex a)
sxPtr Ptr CInt
incxPtr Ptr (Complex a)
syPtr Ptr CInt
incyPtr


fromReal ::
   (Shape.C sh, Class.Floating a) => T shA sh (RealOf a) -> Vector sh a
fromReal :: forall sh a shA.
(C sh, Floating a) =>
T shA sh (RealOf a) -> Vector sh a
fromReal =
   FromReal (T shA sh) (Array sh) a
-> T shA sh (RealOf a) -> Array sh a
forall (f :: * -> *) (g :: * -> *) a.
FromReal f g a -> f (RealOf a) -> g a
getFromReal (FromReal (T shA sh) (Array sh) a
 -> T shA sh (RealOf a) -> Array sh a)
-> FromReal (T shA sh) (Array sh) a
-> T shA sh (RealOf a)
-> Array sh a
forall a b. (a -> b) -> a -> b
$
   FromReal (T shA sh) (Array sh) Float
-> FromReal (T shA sh) (Array sh) Double
-> FromReal (T shA sh) (Array sh) (Complex Float)
-> FromReal (T shA sh) (Array sh) (Complex Double)
-> FromReal (T shA sh) (Array sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh (RealOf Float) -> Array sh Float)
-> FromReal (T shA sh) (Array sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f (RealOf a) -> g a) -> FromReal f g a
FromReal T shA sh Float -> Array sh Float
T shA sh (RealOf Float) -> Array sh Float
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh (RealOf Double) -> Array sh Double)
-> FromReal (T shA sh) (Array sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f (RealOf a) -> g a) -> FromReal f g a
FromReal T shA sh Double -> Array sh Double
T shA sh (RealOf Double) -> Array sh Double
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh (RealOf (Complex Float)) -> Array sh (Complex Float))
-> FromReal (T shA sh) (Array sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f (RealOf a) -> g a) -> FromReal f g a
FromReal T shA sh Float -> Array sh (Complex Float)
T shA sh (RealOf (Complex Float)) -> Array sh (Complex Float)
forall sh a shA.
(C sh, Real a) =>
T shA sh a -> Vector sh (Complex a)
complexFromReal)
      ((T shA sh (RealOf (Complex Double)) -> Array sh (Complex Double))
-> FromReal (T shA sh) (Array sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f (RealOf a) -> g a) -> FromReal f g a
FromReal T shA sh Double -> Array sh (Complex Double)
T shA sh (RealOf (Complex Double)) -> Array sh (Complex Double)
forall sh a shA.
(C sh, Real a) =>
T shA sh a -> Vector sh (Complex a)
complexFromReal)

newtype FromReal f g a = FromReal {forall (f :: * -> *) (g :: * -> *) a.
FromReal f g a -> f (RealOf a) -> g a
getFromReal :: f (RealOf a) -> g a}

toComplex ::
   (Shape.C sh, Class.Floating a) => T shA sh a -> Vector sh (ComplexOf a)
toComplex :: forall sh a shA.
(C sh, Floating a) =>
T shA sh a -> Vector sh (ComplexOf a)
toComplex =
   ToComplex (T shA sh) (Array sh) a
-> T shA sh a -> Array sh (Complex (RealOf a))
forall (f :: * -> *) (g :: * -> *) a.
ToComplex f g a -> f a -> g (ComplexOf a)
getToComplex (ToComplex (T shA sh) (Array sh) a
 -> T shA sh a -> Array sh (Complex (RealOf a)))
-> ToComplex (T shA sh) (Array sh) a
-> T shA sh a
-> Array sh (Complex (RealOf a))
forall a b. (a -> b) -> a -> b
$
   ToComplex (T shA sh) (Array sh) Float
-> ToComplex (T shA sh) (Array sh) Double
-> ToComplex (T shA sh) (Array sh) (Complex Float)
-> ToComplex (T shA sh) (Array sh) (Complex Double)
-> ToComplex (T shA sh) (Array sh) a
forall a (f :: * -> *).
Floating a =>
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
forall (f :: * -> *).
f Float
-> f Double -> f (Complex Float) -> f (Complex Double) -> f a
Class.switchFloating
      ((T shA sh Float -> Array sh (ComplexOf Float))
-> ToComplex (T shA sh) (Array sh) Float
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g (ComplexOf a)) -> ToComplex f g a
ToComplex T shA sh Float -> Vector sh (Complex Float)
T shA sh Float -> Array sh (ComplexOf Float)
forall sh a shA.
(C sh, Real a) =>
T shA sh a -> Vector sh (Complex a)
complexFromReal)
      ((T shA sh Double -> Array sh (ComplexOf Double))
-> ToComplex (T shA sh) (Array sh) Double
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g (ComplexOf a)) -> ToComplex f g a
ToComplex T shA sh Double -> Vector sh (Complex Double)
T shA sh Double -> Array sh (ComplexOf Double)
forall sh a shA.
(C sh, Real a) =>
T shA sh a -> Vector sh (Complex a)
complexFromReal)
      ((T shA sh (Complex Float) -> Array sh (ComplexOf (Complex Float)))
-> ToComplex (T shA sh) (Array sh) (Complex Float)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g (ComplexOf a)) -> ToComplex f g a
ToComplex T shA sh (Complex Float) -> Vector sh (Complex Float)
T shA sh (Complex Float) -> Array sh (ComplexOf (Complex Float))
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)
      ((T shA sh (Complex Double)
 -> Array sh (ComplexOf (Complex Double)))
-> ToComplex (T shA sh) (Array sh) (Complex Double)
forall (f :: * -> *) (g :: * -> *) a.
(f a -> g (ComplexOf a)) -> ToComplex f g a
ToComplex T shA sh (Complex Double) -> Vector sh (Complex Double)
T shA sh (Complex Double) -> Array sh (ComplexOf (Complex Double))
forall slice a sh.
(C slice, Floating a) =>
T sh slice a -> Vector slice a
toVector)

newtype ToComplex f g a = ToComplex {forall (f :: * -> *) (g :: * -> *) a.
ToComplex f g a -> f a -> g (ComplexOf a)
getToComplex :: f a -> g (ComplexOf a)}

complexFromReal ::
   (Shape.C sh, Class.Real a) => T shA sh a -> Vector sh (Complex a)
complexFromReal :: forall sh a shA.
(C sh, Real a) =>
T shA sh a -> Vector sh (Complex a)
complexFromReal T shA sh a
x =
   sh -> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
x) ((Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a))
-> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr (Complex a)
yPtr ->
   case Ptr (Complex a) -> Ptr (RealOf (Complex a))
forall a. Ptr a -> Ptr (RealOf a)
realPtr Ptr (Complex a)
yPtr of
      Ptr (RealOf (Complex a))
yrPtr -> 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 -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
         (Ptr a
xPtr, Ptr CInt
incxPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
x
         Ptr CInt
incyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
2
         Ptr CInt
inczPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
0
         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
Scalar.zero
         IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
            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 ()
Blas.copy Ptr CInt
nPtr Ptr a
xPtr Ptr CInt
incxPtr Ptr a
Ptr (RealOf (Complex a))
yrPtr Ptr CInt
incyPtr
            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 ()
Blas.copy Ptr CInt
nPtr Ptr a
zPtr Ptr CInt
inczPtr (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
Ptr (RealOf (Complex a))
yrPtr Int
1) Ptr CInt
incyPtr


realFromComplexVector ::
   (Shape.C sh) =>
   Vector sh (Complex a) -> T (sh, ComplexShape) (sh, ComplexShape) a
realFromComplexVector :: forall sh a.
C sh =>
Vector sh (Complex a) -> T (sh, ComplexShape) (sh, ComplexShape) a
realFromComplexVector (Array sh
sh ForeignPtr (Complex a)
a) =
   let csh :: (sh, ComplexShape)
csh = (sh
sh, ComplexShape
forall sh. Static sh => sh
Shape.static) in
   T (sh, ComplexShape)
-> Vector (sh, ComplexShape) a
-> T (sh, ComplexShape) (sh, ComplexShape) a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons ((sh, ComplexShape) -> T (sh, ComplexShape)
forall sh. C sh => sh -> T sh
Slice.fromShape (sh, ComplexShape)
csh) ((sh, ComplexShape) -> ForeignPtr a -> Vector (sh, ComplexShape) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array (sh, ComplexShape)
csh (ForeignPtr (Complex a) -> ForeignPtr a
forall a b. ForeignPtr a -> ForeignPtr b
castForeignPtr ForeignPtr (Complex a)
a))


realPart ::
   (Shape.C sh, Class.Real a) =>
   T shA sh (Complex a) -> T (shA, ComplexShape) sh a
realPart :: forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> T (shA, ComplexShape) sh a
realPart (Cons (Slice.Cons Int
s Int
k sh
slc) (Array shA
sh ForeignPtr (Complex a)
a)) =
   T sh -> Vector (shA, ComplexShape) a -> T (shA, ComplexShape) sh a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons
      (Int -> Int -> sh -> T sh
forall sh. Int -> Int -> sh -> T sh
Slice.Cons (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
s) (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k) sh
slc)
      ((shA, ComplexShape) -> ForeignPtr a -> Vector (shA, ComplexShape) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array (shA
sh, ComplexShape
forall sh. Static sh => sh
Shape.static) (ForeignPtr (Complex a) -> ForeignPtr a
forall a b. ForeignPtr a -> ForeignPtr b
castForeignPtr ForeignPtr (Complex a)
a))

imaginaryPart ::
   (Shape.C sh, Class.Real a) =>
   T shA sh (Complex a) -> T (shA, ComplexShape) sh a
imaginaryPart :: forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> T (shA, ComplexShape) sh a
imaginaryPart (Cons (Slice.Cons Int
s Int
k sh
slc) (Array shA
sh ForeignPtr (Complex a)
a)) =
   T sh -> Vector (shA, ComplexShape) a -> T (shA, ComplexShape) sh a
forall sh slice a. T slice -> Vector sh a -> T sh slice a
Cons
      (Int -> Int -> sh -> T sh
forall sh. Int -> Int -> sh -> T sh
Slice.Cons (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k) sh
slc)
      ((shA, ComplexShape) -> ForeignPtr a -> Vector (shA, ComplexShape) a
forall sh a. sh -> ForeignPtr a -> Array sh a
Array (shA
sh, ComplexShape
forall sh. Static sh => sh
Shape.static) (ForeignPtr (Complex a) -> ForeignPtr a
forall a b. ForeignPtr a -> ForeignPtr b
castForeignPtr ForeignPtr (Complex a)
a))


zipComplex ::
   (Shape.C sh, Eq sh, Class.Real a) =>
   T shA sh a -> T shB sh a -> Vector sh (Complex a)
zipComplex :: forall sh a shA shB.
(C sh, Eq sh, Real a) =>
T shA sh a -> T shB sh a -> Vector sh (Complex a)
zipComplex T shA sh a
xr T shB sh a
xi =
   sh -> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall sh a.
(C sh, Storable a) =>
sh -> (Int -> Ptr a -> IO ()) -> Array sh a
Array.unsafeCreateWithSize (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
xr) ((Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a))
-> (Int -> Ptr (Complex a) -> IO ()) -> Array sh (Complex a)
forall a b. (a -> b) -> a -> b
$ \Int
n Ptr (Complex 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
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ String -> Bool -> IO ()
Call.assert String
"zipComplex: shapes mismatch" (T shA sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shA sh a
xr sh -> sh -> Bool
forall a. Eq a => a -> a -> Bool
== T shB sh a -> sh
forall sh slice a. T sh slice a -> slice
shape T shB sh a
xi)
      Ptr CInt
nPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
n
      (Ptr a
xrPtr, Ptr CInt
incxrPtr) <- T shA sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shA sh a
xr
      (Ptr a
xiPtr, Ptr CInt
incxiPtr) <- T shB sh a -> FortranIO () (Ptr a, Ptr CInt)
forall a sh slice r.
Storable a =>
T sh slice a -> FortranIO r (Ptr a, Ptr CInt)
sliceArg T shB sh a
xi
      let yrPtr :: Ptr (RealOf (Complex a))
yrPtr = Ptr (Complex a) -> Ptr (RealOf (Complex a))
forall a. Ptr a -> Ptr (RealOf a)
realPtr Ptr (Complex a)
yPtr
      Ptr CInt
incyPtr <- Int -> ContT () IO (Ptr CInt)
forall r. Int -> FortranIO r (Ptr CInt)
Call.cint Int
2
      IO () -> ContT () IO ()
forall a. IO a -> ContT () IO a
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> ContT () IO ()) -> IO () -> ContT () IO ()
forall a b. (a -> b) -> a -> b
$ do
         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 ()
Blas.copy Ptr CInt
nPtr Ptr a
xrPtr Ptr CInt
incxrPtr Ptr a
Ptr (RealOf (Complex a))
yrPtr Ptr CInt
incyPtr
         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 ()
Blas.copy Ptr CInt
nPtr Ptr a
xiPtr Ptr CInt
incxiPtr (Ptr a -> Int -> Ptr a
forall a. Storable a => Ptr a -> Int -> Ptr a
advancePtr Ptr a
Ptr (RealOf (Complex a))
yrPtr Int
1) Ptr CInt
incyPtr

{- |
prop> forSliced complex_ $ \xs -> approxReal 1e-2 (VectorSlice.norm2 xs) $ let (xrs,xis) = VectorSlice.unzipComplex xs in sqrt $ VectorSlice.norm2Squared xrs + VectorSlice.norm2Squared xis
-}
unzipComplex ::
   (Shape.C sh, Class.Real a) =>
   T shA sh (Complex a) ->
   (T (shA,ComplexShape) sh a, T (shA,ComplexShape) sh a)
unzipComplex :: forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a)
-> (T (shA, ComplexShape) sh a, T (shA, ComplexShape) sh a)
unzipComplex T shA sh (Complex a)
x = (T shA sh (Complex a) -> T (shA, ComplexShape) sh a
forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> T (shA, ComplexShape) sh a
realPart T shA sh (Complex a)
x, T shA sh (Complex a) -> T (shA, ComplexShape) sh a
forall sh a shA.
(C sh, Real a) =>
T shA sh (Complex a) -> T (shA, ComplexShape) sh a
imaginaryPart T shA sh (Complex a)
x)