From 0ff13d993b880739295de343bca62f06fac5ca0c Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Wed, 12 Sep 2007 19:09:47 +0000 Subject: documentation --- lib/Data/Packed/Internal.hs | 5 +- lib/Data/Packed/Internal/Common.hs | 24 +++++--- lib/Data/Packed/Internal/Matrix.hs | 116 +++++++++++++++++++++++++------------ lib/Data/Packed/Internal/Vector.hs | 62 ++++++++++++++++---- lib/Data/Packed/Internal/aux.c | 8 +-- lib/Data/Packed/Matrix.hs | 62 ++++++++++---------- lib/Data/Packed/Plot.hs | 1 + lib/Data/Packed/Vector.hs | 35 ++--------- 8 files changed, 190 insertions(+), 123 deletions(-) (limited to 'lib/Data') diff --git a/lib/Data/Packed/Internal.hs b/lib/Data/Packed/Internal.hs index 822bb1b..e5028d3 100644 --- a/lib/Data/Packed/Internal.hs +++ b/lib/Data/Packed/Internal.hs @@ -11,15 +11,14 @@ -- Reexports all internal modules -- ----------------------------------------------------------------------------- +-- #hide module Data.Packed.Internal ( module Data.Packed.Internal.Common, module Data.Packed.Internal.Vector, module Data.Packed.Internal.Matrix, --- module Data.Packed.Internal.Tensor ) where import Data.Packed.Internal.Common import Data.Packed.Internal.Vector -import Data.Packed.Internal.Matrix ---import Data.Packed.Internal.Tensor \ No newline at end of file +import Data.Packed.Internal.Matrix \ No newline at end of file diff --git a/lib/Data/Packed/Internal/Common.hs b/lib/Data/Packed/Internal/Common.hs index 5548285..218fb6b 100644 --- a/lib/Data/Packed/Internal/Common.hs +++ b/lib/Data/Packed/Internal/Common.hs @@ -9,9 +9,10 @@ -- Stability : provisional -- Portability : portable (uses FFI) -- --- Common tools +-- Development utilities. -- ----------------------------------------------------------------------------- +-- #hide module Data.Packed.Internal.Common where @@ -33,11 +34,15 @@ instance (Storable a, RealFloat a) => Storable (Complex a) where -- poke p (a :+ b) = pokeArray (castPtr p) [a,b] -- ---------------------------------------------------------------------- +-- | @debug x = trace (show x) x@ +debug :: (Show a) => a -> a debug x = trace (show x) x +-- | useful for expressions like @sortBy (compare \`on\` length)@ on :: (a -> a -> b) -> (t -> a) -> t -> t -> b on f g = \x y -> f (g x) (g y) +-- | @partit 3 [1..9] == [[1,2,3],[4,5,6],[7,8,9]]@ partit :: Int -> [a] -> [[a]] partit _ [] = [] partit n l = take n l : partit n (drop n l) @@ -50,19 +55,20 @@ common f = commonval . map f where commonval [a] = Just a commonval (a:b:xs) = if a==b then commonval (b:xs) else Nothing +-- | postfix function application (@flip ($)@) (//) :: x -> (x -> y) -> y infixl 0 // (//) = flip ($) --- our codes should start from 1024 - +-- GSL error codes are <= 1024 +-- | error codes for the auxiliary functions required by the wrappers errorCode :: Int -> String -errorCode 1000 = "bad size" -errorCode 1001 = "bad function code" -errorCode 1002 = "memory problem" -errorCode 1003 = "bad file" -errorCode 1004 = "singular" -errorCode 1005 = "didn't converge" +errorCode 2000 = "bad size" +errorCode 2001 = "bad function code" +errorCode 2002 = "memory problem" +errorCode 2003 = "bad file" +errorCode 2004 = "singular" +errorCode 2005 = "didn't converge" errorCode n = "code "++show n {- | conversion of Haskell functions into function pointers that can be used in the C side diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs index 2db4838..f9dd9a9 100644 --- a/lib/Data/Packed/Internal/Matrix.hs +++ b/lib/Data/Packed/Internal/Matrix.hs @@ -12,7 +12,7 @@ -- Internal matrix representation -- ----------------------------------------------------------------------------- --- --#hide +-- #hide module Data.Packed.Internal.Matrix where @@ -57,10 +57,14 @@ import Data.Maybe(fromJust) data MatrixOrder = RowMajor | ColumnMajor deriving (Show,Eq) +-- | Matrix representation suitable for GSL and LAPACK computations. data Matrix t = MC { rows :: Int, cols :: Int, cdat :: Vector t, fdat :: Vector t } | MF { rows :: Int, cols :: Int, fdat :: Vector t, cdat :: Vector t } --- transposition just changes the data order +-- MC: preferred by C, fdat may require a transposition +-- MF: preferred by LAPACK, cdat may require a transposition + +-- | matrix transpose trans :: Matrix t -> Matrix t trans MC {rows = r, cols = c, cdat = d, fdat = dt } = MF {rows = c, cols = r, fdat = d, cdat = dt } trans MF {rows = r, cols = c, fdat = d, cdat = dt } = MC {rows = c, cols = r, cdat = d, fdat = dt } @@ -166,32 +170,29 @@ compat m1 m2 = rows m1 == rows m2 && cols m1 == cols m2 ---------------------------------------------------------------- --- | element types for which optimized matrix computations are provided +-- | Optimized matrix computations are provided for elements in the Field class. class Storable a => Field a where - -- | @constant val n@ creates a vector with @n@ elements, all equal to @val@. - constant :: a -> Int -> Vector a + constantD :: a -> Int -> Vector a transdata :: Int -> Vector a -> Int -> Vector a - multiplyD :: MatrixOrder -> Matrix a -> Matrix a -> Matrix a - -- | extracts a submatrix froma a matrix - subMatrix :: (Int,Int) -- ^ (r0,c0) starting position - -> (Int,Int) -- ^ (rt,ct) dimensions of submatrix - -> Matrix a -> Matrix a - -- | creates a square matrix with the given diagonal - diag :: Vector a -> Matrix a + multiplyD :: Matrix a -> Matrix a -> Matrix a + subMatrixD :: (Int,Int) -- ^ (r0,c0) starting position + -> (Int,Int) -- ^ (rt,ct) dimensions of submatrix + -> Matrix a -> Matrix a + diagD :: Vector a -> Matrix a instance Field Double where - constant = constantR + constantD = constantR transdata = transdataR - multiplyD = multiplyR - subMatrix = subMatrixR - diag = diagR + multiplyD = multiplyR + subMatrixD = subMatrixR + diagD = diagR instance Field (Complex Double) where - constant = constantC - transdata = transdataC - multiplyD = multiplyC - subMatrix = subMatrixC - diag = diagC + constantD = constantC + transdata = transdataC + multiplyD = multiplyC + subMatrixD = subMatrixC + diagD = diagC ------------------------------------------------------------------ @@ -209,6 +210,15 @@ dsp as = (++" ]") . (" ["++) . init . drop 2 . unlines . map (" , "++) . map unw ------------------------------------------------------------------ +(>|<) :: (Field a) => Int -> Int -> [a] -> Matrix a +r >|< c = f where + f l | dim v == r*c = matrixFromVector ColumnMajor c v + | otherwise = error $ "inconsistent list size = " + ++show (dim v) ++" in ("++show r++"><"++show c++")" + where v = fromList l + +------------------------------------------------------------------- + transdataR :: Int -> Vector Double -> Int -> Vector Double transdataR = transdataAux ctransR @@ -237,10 +247,10 @@ foreign import ccall safe "aux.h transC" gmatC MF {rows = r, cols = c, fdat = d} f = f 1 c r (ptr d) gmatC MC {rows = r, cols = c, cdat = d} f = f 0 r c (ptr d) -multiplyAux fun order a b = unsafePerformIO $ do +multiplyAux fun a b = unsafePerformIO $ do when (cols a /= rows b) $ error $ "inconsistent dimensions in contraction "++ show (rows a,cols a) ++ " x " ++ show (rows b, cols b) - r <- createMatrix order (rows a) (cols b) + r <- createMatrix RowMajor (rows a) (cols b) fun // gmatC a // gmatC b // mat dat r // check "multiplyAux" [dat a, dat b] return r @@ -258,33 +268,41 @@ foreign import ccall safe "aux.h multiplyC" -> Int -> Int -> Ptr (Complex Double) -> IO Int -multiply :: (Field a) => MatrixOrder -> Matrix a -> Matrix a -> Matrix a -multiply RowMajor a b = multiplyD RowMajor a b -multiply ColumnMajor a b = MF {rows = c, cols = r, fdat = d, cdat = dt } - where MC {rows = r, cols = c, cdat = d, fdat = dt } = multiplyD RowMajor (trans b) (trans a) --- FIXME using MatrixFromVector +multiply' :: (Field a) => MatrixOrder -> Matrix a -> Matrix a -> Matrix a +multiply' RowMajor a b = multiplyD a b +multiply' ColumnMajor a b = trans $ multiplyD (trans b) (trans a) + + +-- | matrix product +multiply :: (Field a) => Matrix a -> Matrix a -> Matrix a +multiply = multiplyD ---------------------------------------------------------------------- --- | extraction of a submatrix of a real matrix -subMatrixR :: (Int,Int) -- ^ (r0,c0) starting position - -> (Int,Int) -- ^ (rt,ct) dimensions of submatrix - -> Matrix Double -> Matrix Double +-- | extraction of a submatrix from a real matrix +subMatrixR :: (Int,Int) -> (Int,Int) -> Matrix Double -> Matrix Double subMatrixR (r0,c0) (rt,ct) x = unsafePerformIO $ do r <- createMatrix RowMajor rt ct c_submatrixR r0 (r0+rt-1) c0 (c0+ct-1) // mat cdat x // mat dat r // check "subMatrixR" [dat r] return r foreign import ccall "aux.h submatrixR" c_submatrixR :: Int -> Int -> Int -> Int -> TMM --- | extraction of a submatrix of a complex matrix -subMatrixC :: (Int,Int) -- ^ (r0,c0) starting position - -> (Int,Int) -- ^ (rt,ct) dimensions of submatrix - -> Matrix (Complex Double) -> Matrix (Complex Double) +-- | extraction of a submatrix from a complex matrix +subMatrixC :: (Int,Int) -> (Int,Int) -> Matrix (Complex Double) -> Matrix (Complex Double) subMatrixC (r0,c0) (rt,ct) x = reshape ct . asComplex . cdat . subMatrixR (r0,2*c0) (rt,2*ct) . reshape (2*cols x) . asReal . cdat $ x +-- | Extracts a submatrix from a matrix. +subMatrix :: Field a + => (Int,Int) -- ^ (r0,c0) starting position + -> (Int,Int) -- ^ (rt,ct) dimensions of submatrix + -> Matrix a -- ^ input matrix + -> Matrix a -- ^ result +subMatrix = subMatrixD + + --------------------------------------------------------------------- diagAux fun msg (v@V {dim = n}) = unsafePerformIO $ do @@ -302,6 +320,10 @@ diagC :: Vector (Complex Double) -> Matrix (Complex Double) diagC = diagAux c_diagC "diagC" foreign import ccall "aux.h diagC" c_diagC :: TCVCM +-- | creates a square matrix with the given diagonal +diag :: Field a => Vector a -> Matrix a +diag = diagD + ------------------------------------------------------------------------ constantAux fun x n = unsafePerformIO $ do @@ -321,6 +343,28 @@ constantC = constantAux cconstantC foreign import ccall safe "aux.h constantC" cconstantC :: Ptr (Complex Double) -> TCV -- Complex Double :> IO Int +{- | creates a vector with a given number of equal components: + +@> constant 2 7 +7 |> [2.0,2.0,2.0,2.0,2.0,2.0,2.0]@ +-} +constant :: Field a => a -> Int -> Vector a +constant = constantD + +-------------------------------------------------------------------------- + +-- | obtains the complex conjugate of a complex vector +conj :: Vector (Complex Double) -> Vector (Complex Double) +conj v = asComplex $ cdat $ reshape 2 (asReal v) `multiply` diag (fromList [1,-1]) + +-- | creates a complex vector from vectors with real and imaginary parts +toComplex :: (Vector Double, Vector Double) -> Vector (Complex Double) +toComplex (r,i) = asComplex $ cdat $ fromColumns [r,i] + +-- | converts a real vector into a complex representation (with zero imaginary parts) +comp :: Vector Double -> Vector (Complex Double) +comp v = toComplex (v,constant 0 (dim v)) + ------------------------------------------------------------------------- -- Generic definitions diff --git a/lib/Data/Packed/Internal/Vector.hs b/lib/Data/Packed/Internal/Vector.hs index 0d9dc70..f0ef8b6 100644 --- a/lib/Data/Packed/Internal/Vector.hs +++ b/lib/Data/Packed/Internal/Vector.hs @@ -12,6 +12,7 @@ -- Vector implementation -- ----------------------------------------------------------------------------- +-- #hide module Data.Packed.Internal.Vector where @@ -24,34 +25,39 @@ import Debug.Trace(trace) import Foreign.C.String(peekCString) import Foreign.C.Types - -data Vector t = V { dim :: Int - , fptr :: ForeignPtr t - , ptr :: Ptr t +-- | A one-dimensional array of objects stored in a contiguous memory block. +data Vector t = V { dim :: Int -- ^ number of elements + , fptr :: ForeignPtr t -- ^ foreign pointer to the memory block + , ptr :: Ptr t -- ^ ordinary pointer to the actual starting address (usually the same) } +-- | check the error code and touch foreign ptr of vector arguments (if any) check :: String -> [Vector a] -> IO Int -> IO () check msg ls f = do err <- f - when (err/=0) $ if err > 999 -- FIXME, it should be 1024 - then (error (msg++": "++errorCode err)) - else do + when (err/=0) $ if err > 1024 + then (error (msg++": "++errorCode err)) -- our errors + else do -- GSL errors ps <- gsl_strerror err s <- peekCString ps error (msg++": "++s) mapM_ (touchForeignPtr . fptr) ls return () +-- | description of GSL error codes foreign import ccall "aux.h gsl_strerror" gsl_strerror :: Int -> IO (Ptr CChar) +-- | signature of foreign functions admitting C-style vectors type Vc t s = Int -> Ptr t -> s -- not yet admitted by my haddock version -- infixr 5 :> -- type t :> s = Vc t s +-- | adaptation of our vectors to be admitted by foreign functions: @f \/\/ vec v@ vec :: Vector t -> (Vc t s) -> s vec v f = f (dim v) (ptr v) +-- | allocates memory for a new vector createVector :: Storable a => Int -> IO (Vector a) createVector n = do when (n <= 0) $ error ("trying to createVector of dim "++show n) @@ -60,6 +66,12 @@ createVector n = do --putStrLn ("\n---------> V"++show n) return $ V n fp p +{- | creates a Vector from a list: + +@> fromList [2,3,5,7] +4 |> [2.0,3.0,5.0,7.0]@ + +-} fromList :: Storable a => [a] -> Vector a fromList l = unsafePerformIO $ do v <- createVector (length l) @@ -67,14 +79,25 @@ fromList l = unsafePerformIO $ do f // vec v // check "fromList" [] return v +{- | extracts the Vector elements to a list + +@> toList (linspace 5 (1,10)) +[1.0,3.25,5.5,7.75,10.0]@ + +-} toList :: Storable a => Vector a -> [a] toList v = unsafePerformIO $ peekArray (dim v) (ptr v) +-- | an alternative to 'fromList' with explicit dimension, used also in the instances for Show (Vector a). +(|>) :: (Storable a) => Int -> [a] -> Vector a +infixl 9 |> n |> l = if length l == n then fromList l else error "|> with wrong size" +-- | access to Vector elements without range checking at' :: Storable a => Vector a -> Int -> a at' v n = unsafePerformIO $ peekElemOff (ptr v) n +-- | access to Vector elements with range checking. at :: Storable a => Vector a -> Int -> a at v n | n >= 0 && n < dim v = at' v n | otherwise = error "vector index out of range" @@ -82,9 +105,14 @@ at v n | n >= 0 && n < dim v = at' v n instance (Show a, Storable a) => (Show (Vector a)) where show v = (show (dim v))++" |> " ++ show (toList v) --- | creates a Vector taking a number of consecutive toList from another Vector +{- | takes a number of consecutive elements from a Vector + +@> subVector 2 3 (fromList [1..10]) +3 |> [3.0,4.0,5.0]@ + +-} subVector :: Storable t => Int -- ^ index of the starting element - -> Int -- ^ number of toList to extract + -> Int -- ^ number of elements to extract -> Vector t -- ^ source -> Vector t -- ^ result subVector k l (v@V {dim=n, ptr=p, fptr=fp}) @@ -100,13 +128,23 @@ subVector' k l (v@V {dim=n, ptr=p, fptr=fp}) | otherwise = v {dim=l, ptr=advancePtr p k} --- | Reads a vector position. +{- | Reads a vector position: + +@> fromList [0..9] \@\> 7 +7.0@ + +-} (@>) :: Storable t => Vector t -> Int -> t infixl 9 @> (@>) = at --- | creates a new Vector by joining a list of Vectors +{- | creates a new Vector by joining a list of Vectors + +@> join [fromList [1..5], constant 1 3] +8 |> [1.0,2.0,3.0,4.0,5.0,1.0,1.0,1.0]@ + +-} join :: Storable t => [Vector t] -> Vector t join [] = error "joining zero vectors" join as = unsafePerformIO $ do @@ -131,9 +169,11 @@ asComplex v = V { dim = dim v `div` 2, fptr = castForeignPtr (fptr v), ptr = ca ---------------------------------------------------------------- +-- | map on Vectors liftVector :: (Storable a, Storable b) => (a-> b) -> Vector a -> Vector b liftVector f = fromList . map f . toList +-- | zipWith for Vectors liftVector2 :: (Storable a, Storable b, Storable c) => (a-> b -> c) -> Vector a -> Vector b -> Vector c liftVector2 f u v = fromList $ zipWith f (toList u) (toList v) diff --git a/lib/Data/Packed/Internal/aux.c b/lib/Data/Packed/Internal/aux.c index fe611e2..3db3535 100644 --- a/lib/Data/Packed/Internal/aux.c +++ b/lib/Data/Packed/Internal/aux.c @@ -58,10 +58,10 @@ #define GCVEC(A) int A##n, gsl_complex*A##p #define KGCVEC(A) int A##n, const gsl_complex*A##p -#define BAD_SIZE 1000 -#define BAD_CODE 1001 -#define MEM 1002 -#define BAD_FILE 1003 +#define BAD_SIZE 2000 +#define BAD_CODE 2001 +#define MEM 2002 +#define BAD_FILE 2003 int transR(KRMAT(x),RMAT(t)) { REQUIRES(xr==tc && xc==tr,BAD_SIZE); diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs index 45aaaba..01e8133 100644 --- a/lib/Data/Packed/Matrix.hs +++ b/lib/Data/Packed/Matrix.hs @@ -13,19 +13,19 @@ ----------------------------------------------------------------------------- module Data.Packed.Matrix ( - Matrix(rows,cols), - fromLists, toLists, (><), (>|<), (@@>), + Field, + Matrix,rows,cols, + (><), reshape, flatten, + fromLists, toLists, + (@@>), trans, conjTrans, - reshape, flatten, asRow, asColumn, - fromRows, toRows, fromColumns, toColumns, fromBlocks, - joinVert, joinHoriz, + asRow, asColumn, + fromRows, toRows, fromColumns, toColumns, + fromBlocks, joinVert, joinHoriz, flipud, fliprl, + subMatrix, takeRows, dropRows, takeColumns, dropColumns, + diag, takeDiag, diagRect, ident, liftMatrix, liftMatrix2, - multiply, - outer, - subMatrix, - takeRows, dropRows, takeColumns, dropColumns, - diag, takeDiag, diagRect, ident ) where import Data.Packed.Internal @@ -83,6 +83,18 @@ takeDiag m = fromList [cdat m `at` (k*cols m+k) | k <- [0 .. min (rows m) (cols ident :: (Num t, Field t) => Int -> Matrix t ident n = diag (constant 1 n) +------------------------------------------------------------ + +{- | An easy way to create a matrix: + +@\> (2><3)[1..6] +(2><3) + [ 1.0, 2.0, 3.0 + , 4.0, 5.0, 6.0 ]@ + +This is the format produced by the instances of Show (Matrix a), which +can also be used for input. +-} (><) :: (Field a) => Int -> Int -> [a] -> Matrix a r >< c = f where f l | dim v == r*c = matrixFromVector RowMajor c v @@ -90,13 +102,6 @@ r >< c = f where ++show (dim v) ++" in ("++show r++"><"++show c++")" where v = fromList l -(>|<) :: (Field a) => Int -> Int -> [a] -> Matrix a -r >|< c = f where - f l | dim v == r*c = matrixFromVector ColumnMajor c v - | otherwise = error $ "inconsistent list size = " - ++show (dim v) ++" in ("++show r++"><"++show c++")" - where v = fromList l - ---------------------------------------------------------------- -- | Creates a matrix with the first n rows of another matrix @@ -117,12 +122,19 @@ dropColumns n mat = subMatrix (0,n) (rows mat, cols mat - n) mat {- | Creates a vector by concatenation of rows @\> flatten ('ident' 3) -9 # [1.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,1.0]@ +9 |> [1.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,1.0]@ -} flatten :: Field t => Matrix t -> Vector t flatten = cdat --- | Creates a 'Matrix' from a list of lists (considered as rows). +{- | Creates a 'Matrix' from a list of lists (considered as rows). + +@\> fromLists [[1,2],[3,4],[5,6]] +(3><2) + [ 1.0, 2.0 + , 3.0, 4.0 + , 5.0, 6.0 ]@ +-} fromLists :: Field t => [[t]] -> Matrix t fromLists = fromRows . map fromList @@ -137,15 +149,3 @@ asColumn v = reshape 1 v ------------------------------------------------ -{- | Outer product of two vectors. - -@\> 'fromList' [1,2,3] \`outer\` 'fromList' [5,2,3] -(3><3) - [ 5.0, 2.0, 3.0 - , 10.0, 4.0, 6.0 - , 15.0, 6.0, 9.0 ]@ --} -outer :: (Num t, Field t) => Vector t -> Vector t -> Matrix t -outer u v = multiply RowMajor r c - where r = matrixFromVector RowMajor 1 u - c = matrixFromVector RowMajor (dim v) v diff --git a/lib/Data/Packed/Plot.hs b/lib/Data/Packed/Plot.hs index a0a4aae..7d63426 100644 --- a/lib/Data/Packed/Plot.hs +++ b/lib/Data/Packed/Plot.hs @@ -26,6 +26,7 @@ module Data.Packed.Plot( import Data.Packed.Vector import Data.Packed.Matrix +import LinearAlgebra.Linear(outer) import GSL.Vector(FunCodeS(Max,Min),toScalarR) import Data.List(intersperse) import System diff --git a/lib/Data/Packed/Vector.hs b/lib/Data/Packed/Vector.hs index 867b77b..9818d58 100644 --- a/lib/Data/Packed/Vector.hs +++ b/lib/Data/Packed/Vector.hs @@ -13,33 +13,17 @@ ----------------------------------------------------------------------------- module Data.Packed.Vector ( - Vector(dim), Field, - fromList, toList, - (@>), + Vector, + fromList, (|>), toList, + dim, (@>), subVector, join, - constant, - toComplex, comp, - conj, - dot, - linspace, + constant, linspace, + liftVector, liftVector2 ) where import Data.Packed.Internal import Complex ---import GSL.Vector - --- | creates a complex vector from vectors with real and imaginary parts -toComplex :: (Vector Double, Vector Double) -> Vector (Complex Double) -toComplex (r,i) = asComplex $ cdat $ fromColumns [r,i] - --- | obtains the complex conjugate of a complex vector -conj :: Vector (Complex Double) -> Vector (Complex Double) -conj v = asComplex $ cdat $ reshape 2 (asReal v) `mulC` diag (fromList [1,-1]) - where mulC = multiply RowMajor - -comp :: Vector Double -> Vector (Complex Double) -comp v = toComplex (v,constant 0 (dim v)) {- | Creates a real vector containing a range of values: @@ -47,12 +31,5 @@ comp v = toComplex (v,constant 0 (dim v)) 5 |> [-3.0,-0.5,2.0,4.5,7.0]@ -} linspace :: Int -> (Double, Double) -> Vector Double -linspace n (a,b) = fromList [a::Double,a+delta .. b] +linspace n (a,b) = fromList [a, a+delta .. b] where delta = (b-a)/(fromIntegral n -1) - - -dot :: (Field t) => Vector t -> Vector t -> t -dot u v = dat (multiply RowMajor r c) `at` 0 - where r = matrixFromVector RowMajor (dim u) u - c = matrixFromVector RowMajor 1 v - -- cgit v1.2.3