packages feed

lapack 0.5.0.3 → 0.5.1

raw patch · 35 files changed

+131/−939 lines, 35 filesdep +comfort-blasdep ~comfort-arraydep ~deepseqdep ~text

Dependencies added: comfort-blas

Dependency ranges changed: comfort-array, deepseq, text

Files

lapack.cabal view
@@ -1,6 +1,6 @@ Cabal-Version:    2.2 Name:             lapack-Version:          0.5.0.3+Version:          0.5.1 License:          BSD-3-Clause License-File:     LICENSE Author:           Henning Thielemann <haskell@henning-thielemann.de>@@ -53,7 +53,7 @@   test-module.list  Source-Repository this-  Tag:         0.5.0.3+  Tag:         0.5.1   Type:        darcs   Location:    https://hub.darcs.net/thielema/lapack/ @@ -67,11 +67,12 @@  Library   Build-Depends:+    comfort-blas >=0.0 && <0.1,     lapack-ffi >=0.0.3 && <0.1,     blas-ffi >=0.0 && <0.2,     netlib-ffi >=0.1.1 && <0.2,     comfort-array-shape >=0.0 && <0.1,-    comfort-array >=0.5 && <0.6,+    comfort-array >=0.5.2 && <0.6,     guarded-allocation >=0.0 && <0.1,     hyper >=0.1 && <0.3,     blaze-html >=0.7 && <0.10,@@ -186,7 +187,6 @@     Numeric.LAPACK.Matrix.Array.Private     Numeric.LAPACK.Matrix.Lazy.UpperTriangular     Numeric.LAPACK.Matrix.ModifierTyped-    Numeric.LAPACK.Matrix.Modifier     Numeric.LAPACK.Matrix.RowMajor     Numeric.LAPACK.Matrix.Private     Numeric.LAPACK.Vector.Private
src/Numeric/LAPACK/Linear/LowerUpper.hs view
@@ -43,7 +43,7 @@ import qualified Numeric.LAPACK.Matrix as Matrix import qualified Numeric.LAPACK.Shape as ExtShape import Numeric.LAPACK.Matrix.Array.Private (Full)-import Numeric.LAPACK.Matrix.Modifier (Transposition, Conjugation, Inversion)+import Numeric.BLAS.Matrix.Modifier (Transposition, Conjugation, Inversion)  import qualified Numeric.Netlib.Class as Class
src/Numeric/LAPACK/Linear/Plain.hs view
@@ -57,7 +57,7 @@ import Numeric.LAPACK.Matrix.Triangular.Basic (Lower, Upper) import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor, ColumnMajor), Triangle(Triangle))-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed, Transposed),           Conjugation(NonConjugated, Conjugated),           Inversion(NonInverted, Inverted))
src/Numeric/LAPACK/Matrix.hs view
@@ -136,7 +136,7 @@ import qualified Numeric.LAPACK.Matrix.Array.Private as ArrMatrix import qualified Numeric.LAPACK.Matrix.Type.Private as Matrix import qualified Numeric.LAPACK.Matrix.Plain as Plain-import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Matrix.Divide as Divide import qualified Numeric.LAPACK.Matrix.Multiply as Multiply import qualified Numeric.LAPACK.Matrix.Shape.Omni as Omni
src/Numeric/LAPACK/Matrix/Array/Multiply.hs view
@@ -32,7 +32,7 @@ import Numeric.LAPACK.Matrix.Array.Banded (Banded) import Numeric.LAPACK.Matrix.Array.Private (ArrayMatrix, diagTag) import Numeric.LAPACK.Matrix.Type.Private (Matrix)-import Numeric.LAPACK.Matrix.Modifier (Transposition(NonTransposed,Transposed))+import Numeric.BLAS.Matrix.Modifier (Transposition(NonTransposed,Transposed)) import Numeric.LAPACK.Matrix.Shape.Omni (StripSingleton(StripBands, StripFilled)) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Shape.Private (Unchecked(Unchecked, deconsUnchecked))
src/Numeric/LAPACK/Matrix/Banded/Basic.hs view
@@ -46,7 +46,7 @@ import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor,ColumnMajor), transposeFromOrder, swapOnRowMajor,           UnaryProxy, addOffDiagonals)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated)) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (zero, one) import Numeric.LAPACK.Matrix.Extent.Private (Extent)
src/Numeric/LAPACK/Matrix/BandedHermitian.hs view
@@ -47,7 +47,7 @@ import Numeric.LAPACK.Matrix.Array.Mosaic (FlexHermitian) import Numeric.LAPACK.Matrix.Array.Private (Full) import Numeric.LAPACK.Matrix.Layout.Private (Order, UnaryProxy, natFromProxy)-import Numeric.LAPACK.Matrix.Modifier (Transposition(NonTransposed, Transposed))+import Numeric.BLAS.Matrix.Modifier (Transposition(NonTransposed, Transposed)) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (RealOf)
src/Numeric/LAPACK/Matrix/BandedHermitian/Basic.hs view
@@ -37,7 +37,7 @@ import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor,ColumnMajor), uploFromOrder,           UnaryProxy, natFromProxy)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed, Transposed), transposeOrder,           Conjugation(NonConjugated), conjugatedOnRowMajor) import Numeric.LAPACK.Vector (Vector)@@ -434,10 +434,13 @@          case transposed of             Transposed -> (Banded.transpose upper, Banded.transpose lower)             NonTransposed -> (lower,upper)-   in VectorPriv.mac one-         (Banded.multiplyFull+   in VectorPriv.recheck $+      Vector.mac one+         (VectorPriv.uncheck $+          Banded.multiplyFull             (Banded.mapExtent ExtentPriv.fromSquare lowerT) b)-         (Banded.multiplyFull+         (VectorPriv.uncheck $+          Banded.multiplyFull             (Banded.mapExtent ExtentPriv.fromSquare upperT) b)  takeUpper ::
src/Numeric/LAPACK/Matrix/BandedHermitian/Eigen.hs view
@@ -12,7 +12,7 @@ import qualified Numeric.LAPACK.Shape as ExtShape import Numeric.LAPACK.Matrix.Hermitian.Private (TakeDiagonal(..)) import Numeric.LAPACK.Matrix.Layout.Private (uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (conjugatedOnRowMajor)+import Numeric.BLAS.Matrix.Modifier (conjugatedOnRowMajor) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (RealOf) import Numeric.LAPACK.Private
src/Numeric/LAPACK/Matrix/BandedHermitianPositiveDefinite/Linear.hs view
@@ -16,7 +16,7 @@ import Numeric.LAPACK.Matrix.Hermitian.Private (Determinant(..)) import Numeric.LAPACK.Matrix.Mosaic.Private (copyTriangleToTemp) import Numeric.LAPACK.Matrix.Layout.Private (uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(Conjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(Conjugated)) import Numeric.LAPACK.Matrix.Private (Full) import Numeric.LAPACK.Scalar (RealOf, realPart) import Numeric.LAPACK.Private (copyBlock, withInfo, rankMsg, definiteMsg)
src/Numeric/LAPACK/Matrix/Basic.hs view
@@ -9,7 +9,7 @@ import qualified Numeric.LAPACK.Private as Private import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor, ColumnMajor), transposeFromOrder, flipOrder)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated)) import Numeric.LAPACK.Matrix.Private          (Full, Tall, Wide, Square, General, fromFull, ShapeInt, revealOrder) import Numeric.LAPACK.Vector (Vector)@@ -417,7 +417,7 @@       (RowMajor.recomplex .        RowMajor.scaleColumns          (RowMajor.tensorProduct (Left NonConjugated) x-            (Vector.one Shape.Enumeration)) .+            (Vector.one Shape.static)) .        RowMajor.decomplex)  scaleColumnsComplex ::
src/Numeric/LAPACK/Matrix/Divide.hs view
@@ -20,7 +20,7 @@ import qualified Numeric.LAPACK.Vector as Vector import Numeric.LAPACK.Matrix.Array.Private (Full) import Numeric.LAPACK.Matrix.Type.Private (scaleWithCheck)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed,Transposed),           Inversion(Inverted)) import Numeric.LAPACK.Vector (Vector)
src/Numeric/LAPACK/Matrix/Hermitian.hs view
@@ -87,7 +87,7 @@           SymmetricP) import Numeric.LAPACK.Matrix.Array.Private (Full, General, Square, packTag) import Numeric.LAPACK.Matrix.Layout.Private (Order)-import Numeric.LAPACK.Matrix.Modifier (Transposition(NonTransposed, Transposed))+import Numeric.BLAS.Matrix.Modifier (Transposition(NonTransposed, Transposed)) import Numeric.LAPACK.Matrix.Private (ShapeInt) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (RealOf, one)
src/Numeric/LAPACK/Matrix/Hermitian/Basic.hs view
@@ -22,7 +22,7 @@           withPacking, noLabel, applyFuncPair, triArg) import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor,ColumnMajor), uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed, Transposed),           Conjugation(Conjugated), conjugatedOnRowMajor) import Numeric.LAPACK.Vector (Vector)
src/Numeric/LAPACK/Matrix/Hermitian/Eigen.hs view
@@ -12,7 +12,7 @@ import Numeric.LAPACK.Matrix.Hermitian.Basic (Hermitian, HermitianP) import Numeric.LAPACK.Matrix.Square.Basic (Square) import Numeric.LAPACK.Matrix.Layout.Private (Order(ColumnMajor), uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (conjugatedOnRowMajor)+import Numeric.BLAS.Matrix.Modifier (conjugatedOnRowMajor) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (RealOf) import Numeric.LAPACK.Private
src/Numeric/LAPACK/Matrix/HermitianPositiveDefinite/Linear.hs view
@@ -18,7 +18,7 @@          (withPackingLinear, label, applyFuncPair, triArg, copyTriangleToTemp) import Numeric.LAPACK.Matrix.Mosaic.Basic (takeDiagonal) import Numeric.LAPACK.Matrix.Layout.Private (uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(Conjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(Conjugated)) import Numeric.LAPACK.Matrix.Private (Full) import Numeric.LAPACK.Linear.Private (solver) import Numeric.LAPACK.Scalar (RealOf, realPart, zero)
src/Numeric/LAPACK/Matrix/Layout/Private.hs view
@@ -3,10 +3,14 @@ {-# LANGUAGE EmptyDataDecls #-} {-# LANGUAGE StandaloneDeriving #-} {-# LANGUAGE GADTs #-}-module Numeric.LAPACK.Matrix.Layout.Private where+module Numeric.LAPACK.Matrix.Layout.Private (+   module Numeric.LAPACK.Matrix.Layout.Private,+   module Numeric.BLAS.Matrix.Layout,+   ) where  import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent import Numeric.LAPACK.Matrix.Extent.Private (Extent)+import Numeric.BLAS.Matrix.Layout (Order(..), flipOrder, transposeFromOrder)  import qualified Type.Data.Num.Unary.Literal as TypeNum import qualified Type.Data.Num.Unary.Proof as Proof@@ -25,21 +29,6 @@ import Data.Tuple.HT (mapSnd, swap) import Data.Bool.HT (if') --data Order = RowMajor | ColumnMajor-   deriving (Eq, Show)--instance NFData Order where-   rnf RowMajor = ()-   rnf ColumnMajor = ()--flipOrder :: Order -> Order-flipOrder RowMajor = ColumnMajor-flipOrder ColumnMajor = RowMajor--transposeFromOrder :: Order -> Char-transposeFromOrder RowMajor = 'T'-transposeFromOrder ColumnMajor = 'N'  swapOnRowMajor :: Order -> (a,a) -> (a,a) swapOnRowMajor order =
− src/Numeric/LAPACK/Matrix/Modifier.hs 
@@ -1,50 +0,0 @@-module Numeric.LAPACK.Matrix.Modifier where--import Numeric.LAPACK.Matrix.Layout.Private-         (Order(RowMajor,ColumnMajor), flipOrder)--import Data.Monoid (Monoid, mempty, mappend)-import Data.Semigroup (Semigroup, (<>))----data Transposition = NonTransposed | Transposed-   deriving (Eq, Show, Enum, Bounded)--instance Semigroup Transposition where-   x<>y = if x==y then NonTransposed else Transposed--instance Monoid Transposition where-   mempty = NonTransposed-   mappend = (<>)--transposeOrder :: Transposition -> Order -> Order-transposeOrder NonTransposed = id-transposeOrder Transposed = flipOrder----data Conjugation = NonConjugated | Conjugated-   deriving (Eq, Show, Enum, Bounded)--instance Semigroup Conjugation where-   x<>y = if x==y then NonConjugated else Conjugated--instance Monoid Conjugation where-   mempty = NonConjugated-   mappend = (<>)--conjugatedOnRowMajor :: Order -> Conjugation-conjugatedOnRowMajor RowMajor = Conjugated-conjugatedOnRowMajor ColumnMajor = NonConjugated---data Inversion = NonInverted | Inverted-   deriving (Eq, Show, Enum, Bounded)--instance Semigroup Inversion where-   x<>y = if x==y then NonInverted else Inverted--instance Monoid Inversion where-   mempty = NonInverted-   mappend = (<>)
src/Numeric/LAPACK/Matrix/ModifierTyped.hs view
@@ -1,7 +1,7 @@ {-# LANGUAGE GADTs #-} module Numeric.LAPACK.Matrix.ModifierTyped where -import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent import Numeric.LAPACK.Matrix.Extent.Private (Extent, Small)
src/Numeric/LAPACK/Matrix/Mosaic/Private.hs view
@@ -9,7 +9,7 @@ import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor,ColumnMajor), flipOrder, uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated)) import Numeric.LAPACK.Matrix.Private (Full) import Numeric.LAPACK.Scalar (zero) import Numeric.LAPACK.Shape.Private (Unchecked(Unchecked))
src/Numeric/LAPACK/Matrix/Multiply.hs view
@@ -14,7 +14,7 @@ import qualified Numeric.LAPACK.Matrix.Class as MatrixClass import qualified Numeric.LAPACK.Matrix.Type.Private as Matrix import qualified Numeric.LAPACK.Matrix.Basic as FullBasic-import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Matrix.Shape.Omni as Omni import qualified Numeric.LAPACK.Matrix.Shape as MatrixShape import qualified Numeric.LAPACK.Matrix.Layout.Private as Layout@@ -24,7 +24,7 @@ import qualified Numeric.LAPACK.Vector as Vector import Numeric.LAPACK.Matrix.Array.Private (Full) import Numeric.LAPACK.Matrix.Type.Private (Matrix, scaleWithCheck)-import Numeric.LAPACK.Matrix.Modifier (Transposition(NonTransposed,Transposed))+import Numeric.BLAS.Matrix.Modifier (Transposition(NonTransposed,Transposed)) import Numeric.LAPACK.Vector (Vector)  import qualified Numeric.Netlib.Class as Class
src/Numeric/LAPACK/Matrix/Permutation.hs view
@@ -20,7 +20,7 @@ import qualified Numeric.LAPACK.Matrix.Shape.Omni as Omni import qualified Numeric.LAPACK.Matrix.Layout.Private as Layout import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent-import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Permutation.Private as Plain import qualified Numeric.LAPACK.Permutation as Perm import Numeric.LAPACK.Matrix.Type.Private (Matrix(Permutation))
src/Numeric/LAPACK/Matrix/Plain.hs view
@@ -66,7 +66,7 @@ import Numeric.LAPACK.Matrix.Private          (Full, Tall, Wide, General, ShapeInt, shapeInt,           mapExtent, fromFull, generalizeTall, generalizeWide)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated,Conjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated,Conjugated)) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (zero, one) import Numeric.LAPACK.Private (pointerSeq, fill, copySubMatrix, copyBlock)
src/Numeric/LAPACK/Matrix/RowMajor.hs view
@@ -1,23 +1,33 @@ {-# LANGUAGE TypeFamilies #-} {-# LANGUAGE TypeOperators #-}-module Numeric.LAPACK.Matrix.RowMajor where+module Numeric.LAPACK.Matrix.RowMajor (+   Matrix,+   Matrix.takeRow,+   Matrix.takeColumn,+   Matrix.fromRows,+   Matrix.tensorProduct,+   Matrix.decomplex,+   Matrix.recomplex,+   Matrix.scaleRows,+   Matrix.scaleColumns,+   kronecker,+   ) where  import qualified Numeric.LAPACK.Matrix.Layout.Private as Layout import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent-import qualified Numeric.LAPACK.Private as Private-import Numeric.LAPACK.Matrix.Layout.Private (Order(RowMajor, ColumnMajor))-import Numeric.LAPACK.Matrix.Private (Full, ShapeInt, shapeInt)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated,Conjugated))-import Numeric.LAPACK.Scalar (zero, one)-import Numeric.LAPACK.Private (ComplexPart, pointerSeq)+import Numeric.LAPACK.Matrix.Private (Full) +import qualified Numeric.BLAS.Matrix.RowMajor as Matrix+import Numeric.BLAS.Matrix.RowMajor (Matrix)+import Numeric.BLAS.Matrix.Layout (Order(RowMajor, ColumnMajor))+import Numeric.BLAS.Scalar (zero, one)+ import qualified Numeric.BLAS.FFI.Generic as BlasGen import qualified Numeric.Netlib.Utility as Call import qualified Numeric.Netlib.Class as Class -import Foreign.Marshal.Array (copyArray, advancePtr)-import Foreign.ForeignPtr (withForeignPtr, castForeignPtr)-import Foreign.Storable (Storable)+import Foreign.Marshal.Array (advancePtr)+import Foreign.ForeignPtr (withForeignPtr)  import Control.Monad.Trans.Cont (ContT(ContT), evalContT) import Control.Monad.IO.Class (liftIO)@@ -27,150 +37,7 @@ import qualified Data.Array.Comfort.Shape as Shape import Data.Array.Comfort.Storable.Unchecked (Array(Array)) -import Data.Complex (Complex) import Data.Foldable (forM_)---type Matrix height width = Array (height,width)-type Vector = Array--takeRow ::-   (Shape.Indexed height, Shape.C width, Shape.Index height ~ ix,-    Storable a) =>-   ix -> Matrix height width a -> Vector width a-takeRow ix (Array (height,width) x) =-   Array.unsafeCreateWithSize width $ \n yPtr ->-   withForeignPtr x $ \xPtr ->-      copyArray yPtr (advancePtr xPtr (n * Shape.offset height ix)) n--takeColumn ::-   (Shape.C height, Shape.Indexed width, Shape.Index width ~ ix,-    Class.Floating a) =>-   ix -> Matrix height width a -> Vector height a-takeColumn ix (Array (height,width) x) =-   Array.unsafeCreateWithSize height $ \n yPtr -> evalContT $ do-      let offset = Shape.offset width ix-      nPtr <- Call.cint n-      xPtr <- ContT $ withForeignPtr x-      incxPtr <- Call.cint $ Shape.size width-      incyPtr <- Call.cint 1-      liftIO $ BlasGen.copy nPtr (advancePtr xPtr offset) incxPtr yPtr incyPtr---fromRows ::-   (Shape.C width, Eq width, Storable a) =>-   width -> [Vector width a] -> Matrix ShapeInt width a-fromRows width rows =-   Array.unsafeCreate (shapeInt $ length rows, width) $ \dstPtr ->-   let widthSize = Shape.size width-   in forM_ (zip (pointerSeq widthSize dstPtr) rows) $-         \(dstRowPtr, Array.Array rowWidth srcFPtr) ->-         withForeignPtr srcFPtr $ \srcPtr -> do-            Call.assert-               "Matrix.fromRows: non-matching vector size"-               (width == rowWidth)-            copyArray dstRowPtr srcPtr widthSize---tensorProduct ::-   (Shape.C height, Shape.C width, Class.Floating a) =>-   Either Conjugation Conjugation ->-   Vector height a -> Vector width a -> Matrix height width a-tensorProduct side (Array height x) (Array width y) =-   Array.unsafeCreate (height,width) $ \cPtr -> do-   let m = Shape.size width-   let n = Shape.size height-   let trans conjugated =-         case conjugated of NonConjugated -> 'T'; Conjugated -> 'C'-   let ((transa,transb),(lda,ldb)) =-         case side of-            Left c -> ((trans c, 'N'),(1,1))-            Right c -> (('N', trans c),(m,n))-   evalContT $ do-      transaPtr <- Call.char transa-      transbPtr <- Call.char transb-      mPtr <- Call.cint m-      nPtr <- Call.cint n-      kPtr <- Call.cint 1-      alphaPtr <- Call.number one-      aPtr <- ContT $ withForeignPtr y-      ldaPtr <- Call.leadingDim lda-      bPtr <- ContT $ withForeignPtr x-      ldbPtr <- Call.leadingDim ldb-      betaPtr <- Call.number zero-      ldcPtr <- Call.leadingDim m-      liftIO $-         BlasGen.gemm-            transaPtr transbPtr mPtr nPtr kPtr alphaPtr-            aPtr ldaPtr bPtr ldbPtr betaPtr cPtr ldcPtr---decomplex ::-   (Class.Real a) =>-   Matrix height width (Complex a) ->-   Matrix height (width, Shape.Enumeration ComplexPart) a-decomplex (Array (height,width) a) =-   Array (height, (width, Shape.Enumeration)) (castForeignPtr a)--recomplex ::-   (Class.Real a) =>-   Matrix height (width, Shape.Enumeration ComplexPart) a ->-   Matrix height width (Complex a)-recomplex (Array (height, (width, Shape.Enumeration)) a) =-   Array (height,width) (castForeignPtr a)---scaleRows ::-   (Shape.C height, Eq height, Shape.C width, Class.Floating a) =>-   Vector height a -> Matrix height width a -> Matrix height width a-scaleRows (Array heightX x) (Array shape@(height,width) a) =-      Array.unsafeCreate shape $ \bPtr -> do-   Call.assert "scaleRows: sizes mismatch" (heightX == height)-   evalContT $ do-      let m = Shape.size height-      let n = Shape.size width-      nPtr <- Call.cint n-      xPtr <- ContT $ withForeignPtr x-      aPtr <- ContT $ withForeignPtr a-      incaPtr <- Call.cint 1-      incbPtr <- Call.cint 1-      liftIO $ sequence_ $ take m $-         zipWith3-            (\xkPtr akPtr bkPtr -> do-               BlasGen.copy nPtr akPtr incaPtr bkPtr incbPtr-               BlasGen.scal nPtr xkPtr bkPtr incbPtr)-            (pointerSeq 1 xPtr)-            (pointerSeq n aPtr)-            (pointerSeq n bPtr)--scaleColumns ::-   (Shape.C height, Shape.C width, Eq width, Class.Floating a) =>-   Vector width a -> Matrix height width a -> Matrix height width a-scaleColumns (Array widthX x) (Array shape@(height,width) a) =-      Array.unsafeCreate shape $ \bPtr -> do-   Call.assert "scaleColumns: sizes mismatch" (widthX == width)-   evalContT $ do-      let m = Shape.size height-      let n = Shape.size width-      transPtr <- Call.char 'N'-      nPtr <- Call.cint n-      klPtr <- Call.cint 0-      kuPtr <- Call.cint 0-      alphaPtr <- Call.number one-      xPtr <- ContT $ withForeignPtr x-      ldxPtr <- Call.leadingDim 1-      aPtr <- ContT $ withForeignPtr a-      incaPtr <- Call.cint 1-      betaPtr <- Call.number zero-      incbPtr <- Call.cint 1-      liftIO $ sequence_ $ take m $-         zipWith-            (\akPtr bkPtr ->-               Private.gbmv transPtr-                  nPtr nPtr klPtr kuPtr alphaPtr xPtr ldxPtr-                  akPtr incaPtr betaPtr bkPtr incbPtr)-            (pointerSeq n aPtr)-            (pointerSeq n bPtr)   kronecker ::
src/Numeric/LAPACK/Matrix/Symmetric/Basic.hs view
@@ -18,7 +18,7 @@          (withPacking, noLabel, applyFuncPair, triArg) import Numeric.LAPACK.Matrix.Layout.Private          (MirrorSingleton(SimpleMirror), Order, uploFromOrder)-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated)) import Numeric.LAPACK.Vector (Vector) import Numeric.LAPACK.Scalar (zero) import Numeric.LAPACK.Private (fill)
src/Numeric/LAPACK/Matrix/Symmetric/Unified.hs view
@@ -21,7 +21,7 @@           Labelled(Labelled)) import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor,ColumnMajor), uploFromOrder, uploOrder)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Conjugation(NonConjugated, Conjugated), conjugatedOnRowMajor,           Transposition(NonTransposed, Transposed), transposeOrder) import Numeric.LAPACK.Matrix.Private (Full, General, Square)
src/Numeric/LAPACK/Orthogonal/Householder.hs view
@@ -31,12 +31,12 @@ import qualified Numeric.LAPACK.Matrix.Array.Private as ArrMatrix import qualified Numeric.LAPACK.Matrix.Extent.Strict as ExtentStrict import qualified Numeric.LAPACK.Matrix.Extent as Extent-import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Shape as ExtShape import Numeric.LAPACK.Orthogonal.Plain (Householder) import Numeric.LAPACK.Matrix.Array.Mosaic (Upper) import Numeric.LAPACK.Matrix.Array.Private (Full, Square)-import Numeric.LAPACK.Matrix.Modifier (Transposition, Conjugation)+import Numeric.BLAS.Matrix.Modifier (Transposition, Conjugation)  import qualified Numeric.Netlib.Class as Class
src/Numeric/LAPACK/Orthogonal/Plain.hs view
@@ -29,7 +29,7 @@ import Numeric.LAPACK.Matrix.Layout.Private          (Order(RowMajor, ColumnMajor), sideSwapFromOrder) import Numeric.LAPACK.Matrix.Extent.Private (Extent)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed, Transposed),           Conjugation(NonConjugated, Conjugated)) import Numeric.LAPACK.Matrix.Private (Full)
src/Numeric/LAPACK/Permutation.hs view
@@ -18,7 +18,7 @@  import qualified Numeric.LAPACK.Permutation.Private as Plain import qualified Numeric.LAPACK.Matrix.Array.Private as ArrMatrix-import qualified Numeric.LAPACK.Matrix.Modifier as Mod+import qualified Numeric.BLAS.Matrix.Modifier as Mod import qualified Numeric.LAPACK.Matrix.Extent.Private as Extent import Numeric.LAPACK.Permutation.Private (Permutation) import Numeric.LAPACK.Matrix.Array.Private (Full, Square)
src/Numeric/LAPACK/Permutation/Private.hs view
@@ -7,7 +7,7 @@ import qualified Numeric.LAPACK.Output as Output import Numeric.LAPACK.Output (Output, formatAligned) import Numeric.LAPACK.Matrix.Layout.Private (Order(RowMajor, ColumnMajor))-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition(NonTransposed,Transposed),           Inversion(NonInverted,Inverted)) import Numeric.LAPACK.Matrix.Private (Full, Square)
src/Numeric/LAPACK/Private.hs view
@@ -11,7 +11,7 @@ import qualified Numeric.BLAS.FFI.Generic as BlasGen import qualified Numeric.Netlib.Utility as Call import qualified Numeric.Netlib.Class as Class-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated, Conjugated))+import Numeric.BLAS.Matrix.Modifier (Conjugation(NonConjugated, Conjugated)) import Numeric.LAPACK.Scalar (RealOf, zero, one, isZero)  import qualified Foreign.Marshal.Array.Guarded as ForeignArray@@ -36,7 +36,7 @@ import Data.Array.Comfort.Storable (Array)  import qualified Data.Complex as Complex-import Data.Complex (Complex)+import Data.Complex (Complex((:+))) import Data.Tuple.HT (swap)  import Prelude hiding (sum)@@ -514,5 +514,17 @@ getConstFunc _ = getConst  -data ComplexPart = RealPart | ImaginaryPart-   deriving (Eq, Ord, Show, Enum, Bounded)+{-+ToDo:++type ComplexShape =+         Shape.NestedTuple Shape.TupleAccessor (Complex Shape.Element)++This would allow the use of Complex.realPart as accessor,+but it requires GHC>7.6.3 or so, where realPart has no RealFloat constraint.+-}+type ComplexShape = Shape.NestedTuple Shape.TupleIndex (Complex Shape.Element)++ixReal, ixImaginary :: Shape.ElementIndex (Complex Shape.Element)+ixReal :+ ixImaginary =+   Shape.indexTupleFromShape (Shape.static :: ComplexShape)
src/Numeric/LAPACK/Scalar.hs view
@@ -1,6 +1,3 @@-{-# LANGUAGE TypeFamilies #-}-{-# LANGUAGE TypeOperators #-}-{-# LANGUAGE GADTs #-} module Numeric.LAPACK.Scalar (    RealOf,    ComplexOf,@@ -31,128 +28,4 @@    conjugate,    ) where -import qualified Numeric.Netlib.Class as Class--import Data.Functor.Identity (Identity(Identity, runIdentity))--import qualified Data.Complex as Complex-import Data.Complex (Complex((:+)))---type family RealOf x--type instance RealOf Float = Float-type instance RealOf Double = Double-type instance RealOf (Complex.Complex a) = a---type ComplexOf x = Complex.Complex (RealOf x)---data ComplexSingleton a where-   Real :: (Class.Real a, RealOf a ~ a) => ComplexSingleton a-   Complex :: (Class.Real a) => ComplexSingleton (Complex.Complex a)--complexSingleton :: (Class.Floating a) => ComplexSingleton a-complexSingleton = Class.switchFloating Real Real Complex Complex--complexSingletonOf :: (Class.Floating a) => a -> ComplexSingleton a-complexSingletonOf = const complexSingleton--complexSingletonOfFunctor :: (Class.Floating a) => f a -> ComplexSingleton a-complexSingletonOfFunctor = const complexSingleton--withComplexSingleton :: (Class.Floating a) => (ComplexSingleton a -> a) -> a-withComplexSingleton f = f complexSingleton---data PrecisionSingleton a where-   Float :: PrecisionSingleton Float-   Double :: PrecisionSingleton Double--precisionSingleton :: (Class.Real a) => PrecisionSingleton a-precisionSingleton = Class.switchReal Float Double--precisionOf :: (Class.Real a) => a -> PrecisionSingleton a-precisionOf _ = precisionSingleton--precisionOfFunctor :: (Class.Real a) => f a -> PrecisionSingleton a-precisionOfFunctor _ = precisionSingleton----- move to netlib-carray:Utility or netlib-ffi:Class-zero, one, minusOne :: Class.Floating a => a-zero = selectFloating 0 0 0 0-one = selectFloating 1 1 1 1-minusOne = selectFloating (-1) (-1) (-1) (-1)--selectReal :: (Class.Real a) => Float -> Double -> a-selectReal rf rd = runIdentity $ Class.switchReal (Identity rf) (Identity rd)--selectFloating ::-   (Class.Floating a) =>-   Float -> Double -> Complex.Complex Float -> Complex.Complex Double -> a-selectFloating rf rd cf cd =-   withComplexSingleton $ \sw ->-      case sw of-         Real -> selectReal rf rd-         Complex -> Class.switchReal cf cd----equal :: (Class.Floating a) => a -> a -> Bool-equal a b =-   case complexSingletonOf a of-      Real -> a==b-      Complex -> a==b---isZero :: (Class.Floating a) => a -> Bool-isZero = equal zero---fromReal :: (Class.Floating a) => RealOf a -> a-fromReal a =-   withComplexSingleton $ \sw ->-      case sw of-         Real -> a-         Complex -> a:+0--toComplex :: (Class.Floating a) => a -> ComplexOf a-toComplex a =-   case complexSingletonOf a of-      Real -> a:+0-      Complex -> a--realPart :: (Class.Floating a) => a -> RealOf a-realPart a =-   case complexSingletonOf a of-      Real -> a-      Complex -> Complex.realPart a--absolute :: (Class.Floating a) => a -> RealOf a-absolute a =-   case complexSingletonOf a of-      Real -> abs a-      Complex -> Complex.magnitude a---norm1 :: (Class.Floating a) => a -> RealOf a-norm1 a =-   case complexSingletonOf a of-      Real -> abs a-      Complex -> case a of r:+i -> abs r + abs i---absoluteSquared :: (Class.Floating a) => a -> RealOf a-absoluteSquared a =-   case complexSingletonOf a of-      Real -> a*a-      Complex -> case a of r:+i -> r*r+i*i---conjugate :: (Class.Floating a) => a -> a-conjugate a =-   case complexSingletonOf a of-      Real -> a-      Complex -> Complex.conjugate a+import Numeric.BLAS.Scalar
src/Numeric/LAPACK/Split.hs view
@@ -14,7 +14,7 @@           swapOnRowMajor, sideSwapFromOrder,           Triangle, uploFromOrder, flipOrder) import Numeric.LAPACK.Matrix.Extent.Private (Extent)-import Numeric.LAPACK.Matrix.Modifier+import Numeric.BLAS.Matrix.Modifier          (Transposition, transposeOrder,           Conjugation(NonConjugated, Conjugated)) import Numeric.LAPACK.Matrix.Private (Full)
src/Numeric/LAPACK/Vector.hs view
@@ -4,245 +4,101 @@    Vector,    RealOf,    ComplexOf,-   toList,-   fromList,-   autoFromList,-   CheckedArray.append, (+++),+   Vector.toList,+   Vector.fromList,+   Vector.autoFromList,+   CheckedArray.append, (Vector.+++),    CheckedArray.take, CheckedArray.drop,    CheckedArray.takeLeft, CheckedArray.takeRight,-   swap,+   Vector.swap,    CheckedArray.singleton,-   constant,-   zero,-   one,-   unit,-   dot, inner, (-*|),-   sum,-   absSum,+   Vector.constant,+   Vector.zero,+   Vector.one,+   Vector.unit,+   Vector.dot, Vector.inner, (Vector.-*|),+   Vector.sum,+   Vector.absSum,    norm1,-   norm2,-   norm2Squared,+   Vector.norm2,+   Vector.norm2Squared,    normInf,-   normInf1,+   Vector.normInf1,    argAbsMaximum,-   argAbs1Maximum,-   product,-   scale, scaleReal, (.*|),-   add, sub, (|+|), (|-|),-   negate, raise,-   mac,-   mul, mulConj,+   Vector.argAbs1Maximum,+   Vector.product,+   Vector.scale, Vector.scaleReal, (Vector..*|),+   Vector.add, Vector.sub, (Vector.|+|), (Vector.|-|),+   Vector.negate, Vector.raise,+   Vector.mac,+   Vector.mul, Vector.mulConj,    divide, recip,-   minimum, argMinimum,-   maximum, argMaximum,-   limits, argLimits,+   Vector.minimum, Vector.argMinimum,+   Vector.maximum, Vector.argMaximum,+   Vector.limits, Vector.argLimits,    CheckedArray.foldl,    CheckedArray.foldl1,    CheckedArray.foldMap,     conjugate,-   fromReal,-   toComplex,-   realPart,-   imaginaryPart,-   zipComplex,-   unzipComplex,+   Vector.fromReal,+   Vector.toComplex,+   Vector.realPart,+   Vector.imaginaryPart,+   Vector.zipComplex,+   Vector.unzipComplex,     random, RandomDistribution(..),    ) where -import qualified Numeric.LAPACK.Matrix.RowMajor as RowMajor-import qualified Numeric.LAPACK.Vector.Private as Vector+import qualified Numeric.LAPACK.Vector.Private as VectorPriv import qualified Numeric.LAPACK.Scalar as Scalar import qualified Numeric.LAPACK.Private as Private-import Numeric.LAPACK.Matrix.Hermitian.Private-         (Determinant(Determinant, getDeterminant))-import Numeric.LAPACK.Matrix.Modifier (Conjugation(NonConjugated, Conjugated))+import qualified Numeric.BLAS.Vector as Vector import Numeric.LAPACK.Linear.Private (withInfo)-import Numeric.LAPACK.Scalar (ComplexOf, RealOf, minusOne, absolute)-import Numeric.LAPACK.Private-         (ComplexPart(RealPart, ImaginaryPart), fill, copyConjugate, realPtr)+import Numeric.LAPACK.Scalar (ComplexOf, RealOf, absolute)+import Numeric.LAPACK.Private (copyConjugate)  import qualified Numeric.LAPACK.FFI.Generic as LapackGen import qualified Numeric.LAPACK.FFI.Complex as LapackComplex-import qualified Numeric.BLAS.FFI.Generic as BlasGen-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 Foreign.Marshal.Array (advancePtr)-import Foreign.ForeignPtr (withForeignPtr, castForeignPtr)+import Foreign.ForeignPtr (withForeignPtr) import Foreign.Ptr (Ptr)-import Foreign.Storable (Storable, peek, peekElemOff, pokeElemOff)+import Foreign.Storable (Storable, peekElemOff, pokeElemOff) 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.ST (runST) import Control.Monad (fmap, return, (=<<)) import Control.Applicative (liftA3, (<$>)) -import qualified Data.Array.Comfort.Storable.Mutable.Unchecked as UMutArray-import qualified Data.Array.Comfort.Storable.Mutable as MutArray import qualified Data.Array.Comfort.Storable.Unchecked as Array import qualified Data.Array.Comfort.Storable as CheckedArray import qualified Data.Array.Comfort.Shape as Shape-import Data.Array.Comfort.Storable.Unchecked (Array(Array), append, (!))-import Data.Array.Comfort.Shape ((::+))+import Data.Array.Comfort.Storable.Unchecked (Array(Array)) -import Data.Function (id, flip, ($), (.))+import Data.Function (id, ($), (.)) import Data.Complex (Complex) import Data.Maybe (Maybe(Nothing,Just), maybe) import Data.Tuple.HT (mapFst, uncurry3)-import Data.Tuple (fst, snd)+import Data.Tuple (snd) import Data.Word (Word64) import Data.Bits (shiftR, (.&.))-import Data.Ord (Ord, (>=))+import Data.Ord (Ord) import Data.Eq (Eq, (==)) import Data.Bool (Bool(False,True)) -import Prelude (Int, fromIntegral, (+), (-), (*), Char, Show, Enum, error, IO)+import Prelude (Int, fromIntegral, (+), (-), (*), Show, Enum, error, IO)   type Vector = Array  -toList :: (Shape.C sh, Storable a) => Vector sh a -> [a]-toList = Array.toList--fromList :: (Shape.C sh, Storable a) => sh -> [a] -> Vector sh a-fromList = CheckedArray.fromList--autoFromList :: (Storable a) => [a] -> Vector (Shape.ZeroBased Int) a-autoFromList = Array.vectorFromList---{- |-> constant () = singleton--However, singleton does not need 'Class.Floating' constraint.--}-constant :: (Shape.C sh, Class.Floating a) => sh -> a -> Vector sh a-constant sh a = Array.unsafeCreateWithSize sh $ fill a--zero :: (Shape.C sh, Class.Floating a) => sh -> Vector sh a-zero = flip constant Scalar.zero--one :: (Shape.C sh, Class.Floating a) => sh -> Vector sh a-one = flip constant Scalar.one--unit ::-   (Shape.Indexed sh, Class.Floating a) =>-   sh -> Shape.Index sh -> Vector sh a-unit sh ix = Array.unsafeCreateWithSize sh $ \n xPtr -> do-   fill Scalar.zero n xPtr-   pokeElemOff xPtr (Shape.offset sh ix) Scalar.one---{- |-Precedence and associativity (right) of (List.++).-This also matches '(::+)'.--}-infixr 5 +++--(+++) ::-   (Shape.C shx, Shape.C shy, Storable a) =>-   Vector shx a -> Vector shy a -> Vector (shx::+shy) a-(+++) = append---swap ::-   (Shape.Indexed sh, Storable a) =>-   Shape.Index sh -> Shape.Index sh -> Vector sh a -> Vector sh a-swap i j x =-   runST (do-      y <- MutArray.thaw x-      xi <- MutArray.read y i-      xj <- MutArray.read y j-      MutArray.write y i xj-      MutArray.write y j xi-      UMutArray.unsafeFreeze y)---infixl 7 -*|, .*|--newtype Dot f a = Dot {runDot :: f a -> f a -> a}--{- |-> dot x y = Matrix.toScalar (singleRow x <#> singleColumn y)--}-dot, (-*|) ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Vector sh a -> Vector sh a -> a-(-*|) = dot-dot =-   runDot $-   Class.switchFloating-      (Dot dotReal)-      (Dot dotReal)-      (Dot $ dotComplex 'T')-      (Dot $ dotComplex 'T')--{- |-> inner x y = dot (conjugate x) y--}-inner ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Vector sh a -> Vector sh a -> a-inner =-   runDot $-   Class.switchFloating-      (Dot dotReal)-      (Dot dotReal)-      (Dot $ dotComplex 'C')-      (Dot $ dotComplex 'C')--dotReal ::-   (Shape.C sh, Eq sh, Class.Real a) =>-   Vector sh a -> Vector sh a -> a-dotReal arrX@(Array shX _x) (Array shY y) = unsafePerformIO $ do-   Call.assert "dot: shapes mismatch" (shX == shY)-   evalContT $ do-      (nPtr, sxPtr, incxPtr) <- vectorArgs arrX-      syPtr <- ContT $ withForeignPtr y-      incyPtr <- Call.cint 1-      liftIO $ BlasReal.dot nPtr sxPtr incxPtr syPtr 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) =>-   Char -> Vector sh (Complex a) -> Vector sh (Complex a) -> Complex a-dotComplex trans (Array shX x) (Array shY y) = unsafePerformIO $ do-   Call.assert "dot: shapes mismatch" (shX == shY)-   evalContT $ do-      let m = Shape.size shX-      transPtr <- Call.char trans-      mPtr <- Call.cint m-      nPtr <- Call.cint 1-      alphaPtr <- Call.number Scalar.one-      xPtr <- ContT $ withForeignPtr x-      ldxPtr <- Call.leadingDim m-      yPtr <- ContT $ withForeignPtr y-      incyPtr <- Call.cint 1-      betaPtr <- Call.number Scalar.zero-      zPtr <- Call.alloca-      inczPtr <- Call.cint 1-      liftIO $-         Private.gemv-            transPtr mPtr nPtr alphaPtr xPtr ldxPtr-            yPtr incyPtr betaPtr zPtr inczPtr-      liftIO $ peek zPtr--sum :: (Shape.C sh, Class.Floating a) => Vector sh a -> a-sum (Array sh x) = unsafePerformIO $-   withForeignPtr x $ \xPtr -> Private.sum (Shape.size sh) xPtr 1- norm1 :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a norm1 arr = unsafePerformIO $    evalContT $ liftIO . uncurry3 csum1 =<< vectorArgs arr@@ -256,78 +112,17 @@       (Norm LapackComplex.sum1)       (Norm LapackComplex.sum1) --{- |-Sum of the absolute values of real numbers or components of complex numbers.-For real numbers it is equivalent to 'norm1'.--}-absSum :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a-absSum arr = unsafePerformIO $-   evalContT $ liftIO . uncurry3 asum =<< vectorArgs arr--asum :: Class.Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)-asum =-   getNorm $-   Class.switchFloating-      (Norm BlasReal.asum) (Norm BlasReal.asum)-      (Norm BlasComplex.casum) (Norm BlasComplex.casum)---{- |-Euclidean norm of a vector or Frobenius norm of a matrix.--}-norm2 :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a-norm2 arr = unsafePerformIO $-   evalContT $ liftIO . uncurry3 nrm2 =<< vectorArgs arr--nrm2 :: Class.Floating a => Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)-nrm2 =-   getNorm $-   Class.switchFloating-      (Norm BlasReal.nrm2) (Norm BlasReal.nrm2)-      (Norm BlasComplex.cnrm2) (Norm BlasComplex.cnrm2)- newtype Norm a =    Norm {getNorm :: Ptr CInt -> Ptr a -> Ptr CInt -> IO (RealOf a)}  -norm2Squared :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a-norm2Squared =-   getDeterminant $-   Class.switchFloating-      (Determinant norm2SquaredReal)-      (Determinant norm2SquaredReal)-      (Determinant $ norm2SquaredReal . decomplex)-      (Determinant $ norm2SquaredReal . decomplex)--norm2SquaredReal :: (Shape.C sh, Class.Real a) => Vector sh a -> a-norm2SquaredReal arr =-   unsafePerformIO $ evalContT $ do-      (nPtr, sxPtr, incxPtr) <- vectorArgs arr-      liftIO $ BlasReal.dot nPtr sxPtr incxPtr sxPtr incxPtr-- normInf :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a normInf arr = unsafePerformIO $    evalContT $ do       (nPtr, sxPtr, incxPtr) <- vectorArgs arr       liftIO $          fmap (absolute . maybe Scalar.zero snd) $-         peekElemOff1 sxPtr =<< Vector.absMax nPtr sxPtr incxPtr--{- |-Computes (almost) the infinity norm of the vector.-For complex numbers every element is replaced-by the sum of the absolute component values first.--}-normInf1 :: (Shape.C sh, Class.Floating a) => Vector sh a -> RealOf a-normInf1 arr = unsafePerformIO $-   evalContT $ do-      (nPtr, sxPtr, incxPtr) <- vectorArgs arr-      liftIO $-         fmap (Scalar.norm1 . maybe Scalar.zero snd) $-         peekElemOff1 sxPtr =<< BlasGen.iamax nPtr sxPtr incxPtr-+         peekElemOff1 sxPtr =<< VectorPriv.absMax nPtr sxPtr incxPtr  {- | Returns the index and value of the element with the maximal absolute value.@@ -344,28 +139,9 @@             (maybe                (error "Vector.argAbsMaximum: empty vector")                (mapFst (Shape.uncheckedIndexFromOffset $ Array.shape arr))) $-         peekElemOff1 sxPtr =<< Vector.absMax nPtr sxPtr incxPtr+         peekElemOff1 sxPtr =<< VectorPriv.absMax nPtr sxPtr incxPtr  -{- |-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!--}-argAbs1Maximum ::-   (Shape.InvIndexed sh, Class.Floating a) =>-   Vector sh a -> (Shape.Index sh, a)-argAbs1Maximum arr = unsafePerformIO $-   evalContT $ do-      (nPtr, sxPtr, incxPtr) <- vectorArgs arr-      liftIO $-         fmap-            (maybe-               (error "Vector.argAbs1Maximum: empty vector")-               (mapFst (Shape.uncheckedIndexFromOffset $ Array.shape arr))) $-         peekElemOff1 sxPtr =<< BlasGen.iamax nPtr sxPtr incxPtr- vectorArgs ::    (Shape.C sh) => Array sh a -> ContT r IO (Ptr CInt, Ptr a, Ptr CInt) vectorArgs (Array sh x) =@@ -383,184 +159,6 @@          else Just . (,) ki <$> peekElemOff ptr ki  -product :: (Shape.C sh, Class.Floating a) => Vector sh a -> a-product (Array sh x) = unsafePerformIO $-   withForeignPtr x $ \xPtr -> Private.product (Shape.size sh) xPtr 1---{- |-For restrictions see 'limits'.--}-minimum, maximum :: (Shape.C sh, Class.Real a) => Vector sh a -> a-minimum = fst . limits-maximum = snd . limits--{- |-For restrictions see 'limits'.--}-argMinimum, argMaximum ::-   (Shape.InvIndexed sh, Shape.Index sh ~ ix, Class.Real a) =>-   Vector sh a -> (ix,a)-argMinimum = fst . argLimits-argMaximum = snd . argLimits--{- |-It should hold @limits x = Array.limits x@.-The function 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.-It is less precise if minimum and maximum differ considerably in magnitude-and there are several minimum or maximum candidates of similar value.-E.g. you cannot rely on the property-that @raise (- minimum x) x@ has only non-negative elements.--}-limits :: (Shape.C sh, Class.Real a) => Vector sh a -> (a,a)-limits xs0 =-   let xs = Array.mapShape Shape.Deferred xs0-       x0 = snd $ argAbsMaximum xs-       x1 = xs ! fst (argAbsMaximum (raise (-x0) xs))-   in if x0>=0 then (x1,x0) else (x0,x1)--argLimits ::-   (Shape.InvIndexed sh, Shape.Index sh ~ ix, Class.Real a) =>-   Vector sh a -> ((ix,a),(ix,a))-argLimits xs =-   let p0@(_i0,x0) = argAbsMaximum xs-       p1 = (i1,xs!i1); i1 = fst $ argAbsMaximum $ raise (-x0) xs-   in if x0>=0 then (p1,p0) else (p0,p1)---scale, _scale, (.*|) ::-   (Shape.C sh, Class.Floating a) =>-   a -> Vector sh a -> Vector sh a-(.*|) = scale--scale alpha (Array sh x) = Array.unsafeCreateWithSize sh $ \n syPtr -> do-   evalContT $ do-      alphaPtr <- Call.number alpha-      nPtr <- Call.cint n-      sxPtr <- ContT $ withForeignPtr x-      incxPtr <- Call.cint 1-      incyPtr <- Call.cint 1-      liftIO $ BlasGen.copy nPtr sxPtr incxPtr syPtr incyPtr-      liftIO $ BlasGen.scal nPtr alphaPtr syPtr incyPtr--_scale a (Array sh b) = Array.unsafeCreateWithSize sh $ \n cPtr -> do-   let m = 1-   let k = 1-   evalContT $ do-      transaPtr <- Call.char 'N'-      transbPtr <- Call.char 'N'-      mPtr <- Call.cint m-      kPtr <- Call.cint k-      nPtr <- Call.cint n-      alphaPtr <- Call.number Scalar.one-      aPtr <- Call.number a-      ldaPtr <- Call.leadingDim m-      bPtr <- ContT $ withForeignPtr b-      ldbPtr <- Call.leadingDim k-      betaPtr <- Call.number Scalar.zero-      ldcPtr <- Call.leadingDim m-      liftIO $-         BlasGen.gemm-            transaPtr transbPtr mPtr nPtr kPtr alphaPtr-            aPtr ldaPtr bPtr ldbPtr betaPtr cPtr ldcPtr---scaleReal ::-   (Shape.C sh, Class.Floating a) =>-   RealOf a -> Vector sh a -> Vector sh a-scaleReal =-   getScaleReal $-   Class.switchFloating-      (ScaleReal scale)-      (ScaleReal scale)-      (ScaleReal $ \x -> recomplex . scale x . decomplex)-      (ScaleReal $ \x -> recomplex . scale x . decomplex)--newtype ScaleReal f a = ScaleReal {getScaleReal :: RealOf a -> f a -> f a}---decomplex ::-   (Class.Real a) =>-   Vector sh (Complex a) -> Vector (sh, Shape.Enumeration ComplexPart) a-decomplex (Array sh a) = Array (sh, Shape.Enumeration) (castForeignPtr a)--recomplex ::-   (Class.Real a) =>-   Vector (sh, Shape.Enumeration ComplexPart) a -> Vector sh (Complex a)-recomplex (Array (sh, Shape.Enumeration) a) = Array sh (castForeignPtr a)----infixl 6 |+|, |-|, `add`, `sub`---add, sub, (|+|), (|-|) ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Vector sh a -> Vector sh a -> Vector sh a-add = mac Scalar.one-sub x y = mac minusOne y x--(|+|) = add-(|-|) = sub--mac ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   a -> Vector sh a -> Vector sh a -> Vector sh a-mac alpha x y =-   if Array.shape x == Array.shape y-      then Vector.mac alpha x y-      else error "mac: shapes mismatch"--negate :: (Shape.C sh, Class.Floating a) => Vector sh a -> Vector sh a-negate =-   getConjugate $-   Class.switchFloating-      (Conjugate $ scaleReal Scalar.minusOne)-      (Conjugate $ scaleReal Scalar.minusOne)-      (Conjugate $ scaleReal Scalar.minusOne)-      (Conjugate $ scaleReal Scalar.minusOne)---raise :: (Shape.C sh, Class.Floating a) => a -> Array sh a -> Array sh a-raise c (Array shX x) =-   Array.unsafeCreateWithSize shX $ \n yPtr -> evalContT $ do-      nPtr <- Call.cint n-      cPtr <- Call.number c-      onePtr <- Call.number Scalar.one-      inccPtr <- Call.cint 0-      xPtr <- ContT $ withForeignPtr x-      inc1Ptr <- Call.cint 1-      liftIO $ do-         BlasGen.copy nPtr xPtr inc1Ptr yPtr inc1Ptr-         BlasGen.axpy nPtr onePtr cPtr inccPtr yPtr inc1Ptr---mul ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Vector sh a -> Vector sh a -> Vector sh a-mul = mulConjugation NonConjugated--{- |-prop> \xs ys -> mulConj xs ys == mul (conjugate xs) ys--}-mulConj ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Vector sh a -> Vector sh a -> Vector sh a-mulConj = mulConjugation Conjugated--mulConjugation ::-   (Shape.C sh, Eq sh, Class.Floating a) =>-   Conjugation -> Vector sh a -> Vector sh a -> Vector sh a-mulConjugation conj (Array shA a) (Array shX x) =-      Array.unsafeCreateWithSize shX $ \n yPtr -> do-   Call.assert "mul: shapes mismatch" (shA == shX)-   evalContT $ do-      aPtr <- ContT $ withForeignPtr a-      xPtr <- ContT $ withForeignPtr x-      liftIO $ Private.mul conj n aPtr 1 xPtr 1 yPtr 1- divide ::    (Shape.C sh, Eq sh, Class.Floating a) =>    Vector sh a -> Vector sh a -> Vector sh a@@ -585,8 +183,10 @@  recip :: (Shape.C sh, Class.Floating a) => Vector sh a -> Vector sh a recip x =-   Vector.recheck $-   divide (Vector.uncheck $ one $ Array.shape x) (Vector.uncheck x)+   VectorPriv.recheck $+   divide+      (VectorPriv.uncheck $ Vector.one $ Array.shape x)+      (VectorPriv.uncheck x)   newtype Conjugate f a = Conjugate {getConjugate :: f a -> f a}@@ -612,86 +212,6 @@       incxPtr <- Call.cint 1       incyPtr <- Call.cint 1       liftIO $ copyConjugate nPtr sxPtr incxPtr syPtr incyPtr---fromReal ::-   (Shape.C sh, Class.Floating a) => Vector sh (RealOf a) -> Vector sh a-fromReal =-   getFromReal $-   Class.switchFloating-      (FromReal id)-      (FromReal id)-      (FromReal complexFromReal)-      (FromReal complexFromReal)--newtype FromReal f a = FromReal {getFromReal :: f (RealOf a) -> f a}--toComplex ::-   (Shape.C sh, Class.Floating a) => Vector sh a -> Vector sh (ComplexOf a)-toComplex =-   getToComplex $-   Class.switchFloating-      (ToComplex complexFromReal)-      (ToComplex complexFromReal)-      (ToComplex id)-      (ToComplex id)--newtype ToComplex f a = ToComplex {getToComplex :: f a -> f (ComplexOf a)}--complexFromReal ::-   (Shape.C sh, Class.Real a) => Vector sh a -> Vector sh (Complex a)-complexFromReal (Array sh x) =-   Array.unsafeCreateWithSize sh $ \n yPtr ->-   case realPtr yPtr of-      yrPtr -> evalContT $ do-         nPtr <- Call.cint n-         xPtr <- ContT $ withForeignPtr x-         incxPtr <- Call.cint 1-         incyPtr <- Call.cint 2-         inczPtr <- Call.cint 0-         zPtr <- Call.number Scalar.zero-         liftIO $ do-            BlasGen.copy nPtr xPtr incxPtr yrPtr incyPtr-            BlasGen.copy nPtr zPtr inczPtr (advancePtr yrPtr 1) incyPtr---realPart ::-   (Shape.C sh, Class.Floating a) => Vector sh a -> Vector sh (RealOf a)-realPart =-   getToReal $-   Class.switchFloating-      (ToReal id)-      (ToReal id)-      (ToReal $ RowMajor.takeColumn RealPart . decomplex)-      (ToReal $ RowMajor.takeColumn RealPart . decomplex)--newtype ToReal f a = ToReal {getToReal :: f a -> f (RealOf a)}--imaginaryPart ::-   (Shape.C sh, Class.Real a) => Vector sh (Complex a) -> Vector sh a-imaginaryPart = RowMajor.takeColumn ImaginaryPart . decomplex---zipComplex ::-   (Shape.C sh, Eq sh, Class.Real a) =>-   Vector sh a -> Vector sh a -> Vector sh (Complex a)-zipComplex (Array shr xr) (Array shi xi) =-   Array.unsafeCreateWithSize shr $ \n yPtr -> evalContT $ do-      liftIO $ Call.assert "zipComplex: shapes mismatch" (shr==shi)-      nPtr <- Call.cint n-      xrPtr <- ContT $ withForeignPtr xr-      xiPtr <- ContT $ withForeignPtr xi-      let yrPtr = realPtr yPtr-      incxPtr <- Call.cint 1-      incyPtr <- Call.cint 2-      liftIO $ do-         BlasGen.copy nPtr xrPtr incxPtr yrPtr incyPtr-         BlasGen.copy nPtr xiPtr incxPtr (advancePtr yrPtr 1) incyPtr--unzipComplex ::-   (Shape.C sh, Class.Real a) =>-   Vector sh (Complex a) -> (Vector sh a, Vector sh a)-unzipComplex x = (realPart x, imaginaryPart x)   data RandomDistribution =
src/Numeric/LAPACK/Vector/Private.hs view
@@ -4,35 +4,13 @@  import qualified Numeric.LAPACK.FFI.Complex as LapackComplex import qualified Numeric.BLAS.FFI.Generic as BlasGen-import qualified Numeric.Netlib.Utility as Call import qualified Numeric.Netlib.Class as Class -import Foreign.ForeignPtr (withForeignPtr) import Foreign.Ptr (Ptr) import Foreign.C.Types (CInt) -import Control.Monad.Trans.Cont (ContT(ContT), evalContT)-import Control.Monad.IO.Class (liftIO)- import qualified Data.Array.Comfort.Storable.Unchecked as Array-import qualified Data.Array.Comfort.Shape as Shape-import Data.Array.Comfort.Storable.Unchecked (Array(Array))---mac ::-   (Shape.C sh, Class.Floating a) =>-   a -> Array sh a -> Array sh a -> Array sh a-mac alpha (Array shX x) (Array _shY y) =-   Array.unsafeCreateWithSize shX $ \n szPtr -> evalContT $ do-      nPtr <- Call.cint n-      saPtr <- Call.number alpha-      sxPtr <- ContT $ withForeignPtr x-      incxPtr <- Call.cint 1-      syPtr <- ContT $ withForeignPtr y-      incyPtr <- Call.cint 1-      inczPtr <- Call.cint 1-      liftIO $ BlasGen.copy nPtr syPtr incyPtr szPtr inczPtr-      liftIO $ BlasGen.axpy nPtr saPtr sxPtr incxPtr szPtr inczPtr+import Data.Array.Comfort.Storable.Unchecked (Array)   newtype ArgMaximum a =