From 5b6de561f131d75049fdb999e98a07939ec2e8e7 Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Sat, 24 May 2014 13:32:58 +0200 Subject: backward compatibility --- packages/base/src/Numeric/Container.hs | 273 +++------------------ packages/base/src/Numeric/HMatrix.hs | 158 ++++++++++++ packages/base/src/Numeric/LinearAlgebra.hs | 149 +---------- .../base/src/Numeric/LinearAlgebra/Algorithms.hs | 2 +- packages/base/src/Numeric/LinearAlgebra/Data.hs | 18 +- packages/base/src/Numeric/LinearAlgebra/Util.hs | 149 +++++++++-- packages/base/src/Numeric/LinearAlgebra/Util/CG.hs | 2 +- .../src/Numeric/LinearAlgebra/Util/Convolution.hs | 18 +- packages/base/src/Numeric/Sparse.hs | 6 +- 9 files changed, 351 insertions(+), 424 deletions(-) create mode 100644 packages/base/src/Numeric/HMatrix.hs (limited to 'packages/base/src/Numeric') diff --git a/packages/base/src/Numeric/Container.hs b/packages/base/src/Numeric/Container.hs index 067c5fa..6a841aa 100644 --- a/packages/base/src/Numeric/Container.hs +++ b/packages/base/src/Numeric/Container.hs @@ -1,258 +1,49 @@ -{-# LANGUAGE TypeFamilies #-} -{-# LANGUAGE FlexibleContexts #-} -{-# LANGUAGE FlexibleInstances #-} -{-# LANGUAGE MultiParamTypeClasses #-} -{-# LANGUAGE FunctionalDependencies #-} -{-# LANGUAGE UndecidableInstances #-} - ------------------------------------------------------------------------------ --- | --- Module : Numeric.Container --- Copyright : (c) Alberto Ruiz 2010-14 --- License : BSD3 --- Maintainer : Alberto Ruiz --- Stability : provisional --- --- Basic numeric operations on 'Vector' and 'Matrix', including conversion routines. --- --- The 'Container' class is used to define optimized generic functions which work --- on 'Vector' and 'Matrix' with real or complex elements. --- --- Some of these functions are also available in the instances of the standard --- numeric Haskell classes provided by "Numeric.LinearAlgebra". --- ------------------------------------------------------------------------------ {-# OPTIONS_HADDOCK hide #-} -module Numeric.Container ( - -- * Basic functions +module Numeric.Container( module Data.Packed, - konst, build, + constant, linspace, - diag, ident, + diag, + ident, ctrans, - -- * Generic operations - Container, - add, mul, sub, divide, equal, scaleRecip, addConstant, - scalar, conj, scale, arctan2, cmap, - atIndex, minIndex, maxIndex, minElement, maxElement, + Container(scaleRecip, addConstant,add, sub, mul, divide, equal), + scalar, + conj, + scale, + arctan2, + cmap, + Konst(..), + Build(..), + atIndex, + minIndex, maxIndex, minElement, maxElement, sumElements, prodElements, step, cond, find, assoc, accum, - Transposable(..), Linear(..), - -- * Matrix product - Product(..), udot, dot, (◇), - Mul(..), - Contraction(..),(<.>), + Element(..), + Product(..), optimiseMult, - mXm,mXv,vXm,LSDiv(..), + mXm, mXv, vXm, (<.>), + Mul(..), + LSDiv, (<\>), outer, kronecker, - -- * Random numbers RandDist(..), - randomVector, - gaussianSample, - uniformSample, - -- * Element conversion + randomVector, gaussianSample, uniformSample, + meanCov, Convert(..), - Complexable(), - RealElement(), - - RealOf, ComplexOf, SingleOf, DoubleOf, - - IndexOf, + Complexable, + RealElement, + RealOf, ComplexOf, SingleOf, DoubleOf, IndexOf, module Data.Complex, - -- * IO - module Data.Packed.IO, - -- * Misc - Testable(..) + dispf, disps, dispcf, vecdisp, latexFormat, format, + loadMatrix, saveMatrix, readMatrix ) where -import Data.Packed hiding (stepD, stepF, condD, condF, conjugateC, conjugateQ) -import Data.Packed.Internal.Numeric -import Data.Complex -import Numeric.LinearAlgebra.Algorithms(Field,linearSolveSVD) -import Data.Monoid(Monoid(mconcat)) -import Data.Packed.IO -import Numeric.LinearAlgebra.Random - ------------------------------------------------------------------- - -{- | Creates a real vector containing a range of values: - ->>> linspace 5 (-3,7::Double) -fromList [-3.0,-0.5,2.0,4.5,7.0]@ - ->>> linspace 5 (8,2+i) :: Vector (Complex Double) -fromList [8.0 :+ 0.0,6.5 :+ 0.25,5.0 :+ 0.5,3.5 :+ 0.75,2.0 :+ 1.0] - -Logarithmic spacing can be defined as follows: - -@logspace n (a,b) = 10 ** linspace n (a,b)@ --} -linspace :: (Container Vector e) => Int -> (e, e) -> Vector e -linspace 0 (a,b) = fromList[(a+b)/2] -linspace n (a,b) = addConstant a $ scale s $ fromList $ map fromIntegral [0 .. n-1] - where s = (b-a)/fromIntegral (n-1) - --------------------------------------------------------- - -{- | Matrix product, matrix - vector product, and dot product (equivalent to 'contraction') - -(This operator can also be written using the unicode symbol ◇ (25c7).) - -Examples: - ->>> let a = (3><4) [1..] :: Matrix Double ->>> let v = fromList [1,0,2,-1] :: Vector Double ->>> let u = fromList [1,2,3] :: Vector Double - ->>> a -(3><4) - [ 1.0, 2.0, 3.0, 4.0 - , 5.0, 6.0, 7.0, 8.0 - , 9.0, 10.0, 11.0, 12.0 ] - -matrix × matrix: - ->>> disp 2 (a <.> trans a) -3x3 - 30 70 110 - 70 174 278 -110 278 446 - -matrix × vector: - ->>> a <.> v -fromList [3.0,11.0,19.0] - -dot product: - ->>> u <.> fromList[3,2,1::Double] -10 - -For complex vectors the first argument is conjugated: - ->>> fromList [1,i] <.> fromList[2*i+1,3] -1.0 :+ (-1.0) - ->>> fromList [1,i,1-i] <.> complex a -fromList [10.0 :+ 4.0,12.0 :+ 4.0,14.0 :+ 4.0,16.0 :+ 4.0] --} -infixl 7 <.> -(<.>) :: Contraction a b c => a -> b -> c -(<.>) = contraction - -class Contraction a b c | a b -> c - where - -- | Matrix product, matrix - vector product, and dot product - contraction :: a -> b -> c - -instance (Product t, Container Vector t) => Contraction (Vector t) (Vector t) t where - u `contraction` v = conj u `udot` v - -instance Product t => Contraction (Matrix t) (Vector t) (Vector t) where - contraction = mXv - -instance (Container Vector t, Product t) => Contraction (Vector t) (Matrix t) (Vector t) where - contraction v m = (conj v) `vXm` m - -instance Product t => Contraction (Matrix t) (Matrix t) (Matrix t) where - contraction = mXm - - --------------------------------------------------------------------------------- - -class Mul a b c | a b -> c where - infixl 7 <> - -- | Matrix-matrix, matrix-vector, and vector-matrix products. - (<>) :: Product t => a t -> b t -> c t - -instance Mul Matrix Matrix Matrix where - (<>) = mXm - -instance Mul Matrix Vector Vector where - (<>) m v = flatten $ m <> asColumn v - -instance Mul Vector Matrix Vector where - (<>) v m = flatten $ asRow v <> m - --------------------------------------------------------------------------------- - -class LSDiv c where - infixl 7 <\> - -- | least squares solution of a linear system, similar to the \\ operator of Matlab\/Octave (based on linearSolveSVD) - (<\>) :: Field t => Matrix t -> c t -> c t - -instance LSDiv Vector where - m <\> v = flatten (linearSolveSVD m (reshape 1 v)) - -instance LSDiv Matrix where - (<\>) = linearSolveSVD - --------------------------------------------------------------------------------- - -class Konst e d c | d -> c, c -> d - where - -- | - -- >>> konst 7 3 :: Vector Float - -- fromList [7.0,7.0,7.0] - -- - -- >>> konst i (3::Int,4::Int) - -- (3><4) - -- [ 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0 - -- , 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0 - -- , 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0, 0.0 :+ 1.0 ] - -- - konst :: e -> d -> c e - -instance Container Vector e => Konst e Int Vector - where - konst = konst' - -instance Container Vector e => Konst e (Int,Int) Matrix - where - konst = konst' - --------------------------------------------------------------------------------- - -class Build d f c e | d -> c, c -> d, f -> e, f -> d, f -> c, c e -> f, d e -> f - where - -- | - -- >>> build 5 (**2) :: Vector Double - -- fromList [0.0,1.0,4.0,9.0,16.0] - -- - -- Hilbert matrix of order N: - -- - -- >>> let hilb n = build (n,n) (\i j -> 1/(i+j+1)) :: Matrix Double - -- >>> putStr . dispf 2 $ hilb 3 - -- 3x3 - -- 1.00 0.50 0.33 - -- 0.50 0.33 0.25 - -- 0.33 0.25 0.20 - -- - build :: d -> f -> c e - -instance Container Vector e => Build Int (e -> e) Vector e - where - build = build' - -instance Container Matrix e => Build (Int,Int) (e -> e -> e) Matrix e - where - build = build' - --------------------------------------------------------------------------------- - --- | alternative unicode symbol (25c7) for 'contraction' -(◇) :: Contraction a b c => a -> b -> c -infixl 7 ◇ -(◇) = contraction - --- | dot product: @cdot u v = 'udot' ('conj' u) v@ -dot :: (Container Vector t, Product t) => Vector t -> Vector t -> t -dot u v = udot (conj u) v - --------------------------------------------------------------------------------- +import Data.Packed.Numeric +import Data.Packed +import Data.Packed.Internal(constantD) +import Data.Complex -optimiseMult :: Monoid (Matrix t) => [Matrix t] -> Matrix t -optimiseMult = mconcat +constant :: Element a => a -> Int -> Vector a +constant = constantD diff --git a/packages/base/src/Numeric/HMatrix.hs b/packages/base/src/Numeric/HMatrix.hs new file mode 100644 index 0000000..a56c3d2 --- /dev/null +++ b/packages/base/src/Numeric/HMatrix.hs @@ -0,0 +1,158 @@ +----------------------------------------------------------------------------- +{- | +Module : Numeric.HMatrix +Copyright : (c) Alberto Ruiz 2006-14 +License : BSD3 +Maintainer : Alberto Ruiz +Stability : provisional + +-} +----------------------------------------------------------------------------- +module Numeric.HMatrix ( + + -- * Basic types and data processing + module Numeric.LinearAlgebra.Data, + + -- * Arithmetic and numeric classes + -- | + -- The standard numeric classes are defined elementwise: + -- + -- >>> fromList [1,2,3] * fromList [3,0,-2 :: Double] + -- fromList [3.0,0.0,-6.0] + -- + -- >>> (3><3) [1..9] * ident 3 :: Matrix Double + -- (3><3) + -- [ 1.0, 0.0, 0.0 + -- , 0.0, 5.0, 0.0 + -- , 0.0, 0.0, 9.0 ] + -- + -- In arithmetic operations single-element vectors and matrices + -- (created from numeric literals or using 'scalar') automatically + -- expand to match the dimensions of the other operand: + -- + -- >>> 5 + 2*ident 3 :: Matrix Double + -- (3><3) + -- [ 7.0, 5.0, 5.0 + -- , 5.0, 7.0, 5.0 + -- , 5.0, 5.0, 7.0 ] + -- + + -- * Matrix product + (<.>), + + -- | The overloaded multiplication operators may need type annotations to remove + -- ambiguity. In those cases we can also use the specific functions 'mXm', 'mXv', and 'dot'. + -- + -- The matrix x matrix product is also implemented in the "Data.Monoid" instance, where + -- single-element matrices (created from numeric literals or using 'scalar') + -- are used for scaling. + -- + -- >>> let m = (2><3)[1..] :: Matrix Double + -- >>> m <> 2 <> diagl[0.5,1,0] + -- (2><3) + -- [ 1.0, 4.0, 0.0 + -- , 4.0, 10.0, 0.0 ] + -- + -- 'mconcat' uses 'optimiseMult' to get the optimal association order. + + + -- * Other products + outer, kronecker, cross, + scale, + sumElements, prodElements, + + -- * Linear Systems + (<\>), + linearSolve, + linearSolveLS, + linearSolveSVD, + luSolve, + cholSolve, + cgSolve, + cgSolve', + + -- * Inverse and pseudoinverse + inv, pinv, pinvTol, + + -- * Determinant and rank + rcond, rank, ranksv, + det, invlndet, + + -- * Singular value decomposition + svd, + fullSVD, + thinSVD, + compactSVD, + singularValues, + leftSV, rightSV, + + -- * Eigensystems + eig, eigSH, eigSH', + eigenvalues, eigenvaluesSH, eigenvaluesSH', + geigSH', + + -- * QR + qr, rq, qrRaw, qrgr, + + -- * Cholesky + chol, cholSH, mbCholSH, + + -- * Hessenberg + hess, + + -- * Schur + schur, + + -- * LU + lu, luPacked, + + -- * Matrix functions + expm, + sqrtm, + matFunc, + + -- * Nullspace + nullspacePrec, + nullVector, + nullspaceSVD, + null1, null1sym, + + orth, + + -- * Norms + norm_0, norm_1, norm_2, norm_Inf, + mnorm_0, mnorm_1, mnorm_2, mnorm_Inf, + norm_Frob, norm_nuclear, + + -- * Correlation and convolution + corr, conv, corrMin, corr2, conv2, + + -- * Random arrays + + RandDist(..), randomVector, rand, randn, gaussianSample, uniformSample, + + -- * Misc + meanCov, peps, relativeError, haussholder, optimiseMult, dot, udot, mXm, mXv, smXv, (<>), (◇), Seed, checkT, + -- * Auxiliary classes + Element, Container, Product, Contraction, LSDiv, + Complexable, RealElement, + RealOf, ComplexOf, SingleOf, DoubleOf, + IndexOf, + Field, + Normed, + CGMat, Transposable, + ℕ,ℤ,ℝ,ℂ,ℝn,ℂn, 𝑖, i_C --ℍ +) where + +import Numeric.LinearAlgebra.Data + +import Numeric.Matrix() +import Numeric.Vector() +import Data.Packed.Numeric +import Numeric.LinearAlgebra.Algorithms +import Numeric.LinearAlgebra.Util +import Numeric.LinearAlgebra.Random +import Numeric.Sparse(smXv) +import Numeric.LinearAlgebra.Util.CG + + diff --git a/packages/base/src/Numeric/LinearAlgebra.hs b/packages/base/src/Numeric/LinearAlgebra.hs index 242122f..ad315e4 100644 --- a/packages/base/src/Numeric/LinearAlgebra.hs +++ b/packages/base/src/Numeric/LinearAlgebra.hs @@ -1,4 +1,4 @@ ------------------------------------------------------------------------------ +-------------------------------------------------------------------------------- {- | Module : Numeric.LinearAlgebra Copyright : (c) Alberto Ruiz 2006-14 @@ -7,149 +7,16 @@ Maintainer : Alberto Ruiz Stability : provisional -} ------------------------------------------------------------------------------ -module Numeric.LinearAlgebra ( - - -- * Basic types and data processing - module Numeric.LinearAlgebra.Data, - - -- * Arithmetic and numeric classes - -- | - -- The standard numeric classes are defined elementwise: - -- - -- >>> fromList [1,2,3] * fromList [3,0,-2 :: Double] - -- fromList [3.0,0.0,-6.0] - -- - -- >>> (3><3) [1..9] * ident 3 :: Matrix Double - -- (3><3) - -- [ 1.0, 0.0, 0.0 - -- , 0.0, 5.0, 0.0 - -- , 0.0, 0.0, 9.0 ] - -- - -- In arithmetic operations single-element vectors and matrices - -- (created from numeric literals or using 'scalar') automatically - -- expand to match the dimensions of the other operand: - -- - -- >>> 5 + 2*ident 3 :: Matrix Double - -- (3><3) - -- [ 7.0, 5.0, 5.0 - -- , 5.0, 7.0, 5.0 - -- , 5.0, 5.0, 7.0 ] - -- - - -- * Matrix product - (<.>), - - -- | The overloaded multiplication operators may need type annotations to remove - -- ambiguity. In those cases we can also use the specific functions 'mXm', 'mXv', and 'dot'. - -- - -- The matrix x matrix product is also implemented in the "Data.Monoid" instance, where - -- single-element matrices (created from numeric literals or using 'scalar') - -- are used for scaling. - -- - -- >>> let m = (2><3)[1..] :: Matrix Double - -- >>> m <> 2 <> diagl[0.5,1,0] - -- (2><3) - -- [ 1.0, 4.0, 0.0 - -- , 4.0, 10.0, 0.0 ] - -- - -- 'mconcat' uses 'optimiseMult' to get the optimal association order. - - - -- * Other products - outer, kronecker, cross, - scale, - sumElements, prodElements, absSum, - - -- * Linear Systems - (<\>), - linearSolve, - linearSolveLS, - linearSolveSVD, - luSolve, - cholSolve, - cgSolve, - cgSolve', - - -- * Inverse and pseudoinverse - inv, pinv, pinvTol, - - -- * Determinant and rank - rcond, rank, ranksv, - det, invlndet, - - -- * Singular value decomposition - svd, - fullSVD, - thinSVD, - compactSVD, - singularValues, - leftSV, rightSV, - - -- * Eigensystems - eig, eigSH, eigSH', - eigenvalues, eigenvaluesSH, eigenvaluesSH', - geigSH', - - -- * QR - qr, rq, qrRaw, qrgr, - - -- * Cholesky - chol, cholSH, mbCholSH, - - -- * Hessenberg - hess, - - -- * Schur - schur, - - -- * LU - lu, luPacked, - - -- * Matrix functions - expm, - sqrtm, - matFunc, +-------------------------------------------------------------------------------- +{-# OPTIONS_HADDOCK hide #-} - -- * Nullspace - nullspacePrec, - nullVector, - nullspaceSVD, - null1, null1sym, - - orth, - - -- * Norms - norm1, norm2, normInf, pnorm, NormType(..), - - -- * Correlation and convolution - corr, conv, corrMin, corr2, conv2, - - -- * Random arrays - - RandDist(..), randomVector, rand, randn, gaussianSample, uniformSample, - - -- * Misc - meanCov, peps, relativeError, haussholder, optimiseMult, dot, udot, mXm, mXv, smXv, (<>), (◇), Seed, checkT, - -- * Auxiliary classes - Element, Container, Product, Contraction, LSDiv, - Complexable(), RealElement(), - RealOf, ComplexOf, SingleOf, DoubleOf, - IndexOf, - Field, Normed, - CGMat, Transposable, - R,V +module Numeric.LinearAlgebra ( + module Numeric.Container, + module Numeric.LinearAlgebra.Algorithms ) where -import Numeric.LinearAlgebra.Data - -import Numeric.Matrix() -import Numeric.Vector() import Numeric.Container import Numeric.LinearAlgebra.Algorithms -import Numeric.LinearAlgebra.Util -import Numeric.LinearAlgebra.Random -import Numeric.Sparse(smXv) -import Numeric.LinearAlgebra.Util.CG - +import Numeric.Matrix() +import Numeric.Vector() diff --git a/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs b/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs index c7e7043..bbcc513 100644 --- a/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs +++ b/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs @@ -76,12 +76,12 @@ module Numeric.LinearAlgebra.Algorithms ( ) where -import Data.Packed.Development hiding ((//)) import Data.Packed import Numeric.LinearAlgebra.LAPACK as LAPACK import Data.List(foldl1') import Data.Array import Data.Packed.Internal.Numeric +import Data.Packed.Internal(shSize) {- | Generic linear algebra functions for double precision real and complex matrices. diff --git a/packages/base/src/Numeric/LinearAlgebra/Data.hs b/packages/base/src/Numeric/LinearAlgebra/Data.hs index e3cbe31..89bebbe 100644 --- a/packages/base/src/Numeric/LinearAlgebra/Data.hs +++ b/packages/base/src/Numeric/LinearAlgebra/Data.hs @@ -16,13 +16,19 @@ module Numeric.LinearAlgebra.Data( -- * Vector -- | 1D arrays are storable vectors from the vector package. - Vector, (|>), dim, (@>), - + vect, (|>), + -- * Matrix - Matrix, (><), size, (@@>), trans, ctrans, + + mat, (><), tr, + + -- * Indexing + + size, + Indexable(..), -- * Construction - scalar, konst, build, assoc, accum, linspace, -- ones, zeros, + scalar, Konst(..), Build(..), assoc, accum, linspace, -- ones, zeros, -- * Diagonal ident, diag, diagl, diagRect, takeDiag, @@ -62,11 +68,13 @@ module Numeric.LinearAlgebra.Data( module Data.Complex, + Vector, Matrix + ) where import Data.Packed.Vector import Data.Packed.Matrix -import Numeric.Container +import Data.Packed.Numeric import Numeric.LinearAlgebra.Util import Data.Complex import Numeric.Sparse diff --git a/packages/base/src/Numeric/LinearAlgebra/Util.hs b/packages/base/src/Numeric/LinearAlgebra/Util.hs index 2f91e18..a7d6946 100644 --- a/packages/base/src/Numeric/LinearAlgebra/Util.hs +++ b/packages/base/src/Numeric/LinearAlgebra/Util.hs @@ -1,4 +1,9 @@ {-# LANGUAGE FlexibleContexts #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE TypeFamilies #-} +{-# LANGUAGE MultiParamTypeClasses #-} +{-# LANGUAGE FunctionalDependencies #-} + ----------------------------------------------------------------------------- {- | Module : Numeric.LinearAlgebra.Util @@ -14,19 +19,24 @@ Stability : provisional module Numeric.LinearAlgebra.Util( -- * Convenience functions - size, disp, + vect, mat, + disp, zeros, ones, diagl, row, col, (&), (¦), (——), (#), (?), (¿), + Indexable(..), size, + rand, randn, cross, norm, + ℕ,ℤ,ℝ,ℂ,ℝn,ℂn,𝑖,i_C, --ℍ + norm_1, norm_2, norm_0, norm_Inf, norm_Frob, norm_nuclear, + mnorm_1, mnorm_2, mnorm_0, mnorm_Inf, unitary, mt, pairwiseD2, - meanCov, rowOuters, null1, null1sym, @@ -48,13 +58,49 @@ module Numeric.LinearAlgebra.Util( vtrans ) where -import Numeric.Container +import Data.Packed.Numeric import Numeric.LinearAlgebra.Algorithms hiding (i) import Numeric.Matrix() import Numeric.Vector() - +import Numeric.LinearAlgebra.Random import Numeric.LinearAlgebra.Util.Convolution +type ℝ = Double +type ℕ = Int +type ℤ = Int +type ℂ = Complex Double +type ℝn = Vector ℝ +type ℂn = Vector ℂ +--newtype ℍ m = H m + +i_C, 𝑖 :: ℂ +𝑖 = 0:+1 +i_C = 𝑖 + +{- | create a real vector + +>>> vect [1..5] +fromList [1.0,2.0,3.0,4.0,5.0] + +-} +vect :: [ℝ] -> ℝn +vect = fromList + +{- | create a real matrix + +>>> mat 5 [1..15] +(3><5) + [ 1.0, 2.0, 3.0, 4.0, 5.0 + , 6.0, 7.0, 8.0, 9.0, 10.0 + , 11.0, 12.0, 13.0, 14.0, 15.0 ] + +-} +mat + :: Int -- ^ columns + -> [ℝ] -- ^ elements + -> Matrix ℝ +mat c = reshape c . fromList + {- | print a real matrix with given number of digits after the decimal point >>> disp 5 $ ident 2 / 3 @@ -175,38 +221,97 @@ norm :: Vector Double -> Double -- ^ 2-norm of real vector norm = pnorm PNorm2 +norm_2 :: Normed Vector t => Vector t -> RealOf t +norm_2 = pnorm PNorm2 + +norm_1 :: Normed Vector t => Vector t -> RealOf t +norm_1 = pnorm PNorm1 + +norm_Inf :: Normed Vector t => Vector t -> RealOf t +norm_Inf = pnorm Infinity + +norm_0 :: Vector ℝ -> ℝ +norm_0 v = sumElements (step (abs v - scalar (eps*mx))) + where + mx = norm_Inf v + +norm_Frob :: Normed Matrix t => Matrix t -> RealOf t +norm_Frob = pnorm Frobenius + +norm_nuclear :: Field t => Matrix t -> ℝ +norm_nuclear = sumElements . singularValues + +mnorm_2 :: Normed Matrix t => Matrix t -> RealOf t +mnorm_2 = pnorm PNorm2 + +mnorm_1 :: Normed Matrix t => Matrix t -> RealOf t +mnorm_1 = pnorm PNorm1 + +mnorm_Inf :: Normed Matrix t => Matrix t -> RealOf t +mnorm_Inf = pnorm Infinity + +mnorm_0 :: Matrix ℝ -> ℝ +mnorm_0 = norm_0 . flatten -- | Obtains a vector in the same direction with 2-norm=1 unitary :: Vector Double -> Vector Double unitary v = v / scalar (norm v) --- | ('rows' &&& 'cols') -size :: Matrix t -> (Int, Int) -size m = (rows m, cols m) -- | trans . inv mt :: Matrix Double -> Matrix Double mt = trans . inv -------------------------------------------------------------------------------- +{- | + +>>> size $ fromList[1..10::Double] +10 +>>> size $ (2><5)[1..10::Double] +(2,5) + +-} +size :: Container c t => c t -> IndexOf c +size = size' -{- | Compute mean vector and covariance matrix of the rows of a matrix. +{- | + +>>> vect [1..10] ! 3 +4.0 + +>>> mat 5 [1..15] ! 1 +fromList [6.0,7.0,8.0,9.0,10.0] ->>> meanCov $ gaussianSample 666 1000 (fromList[4,5]) (diagl[2,3]) -(fromList [4.010341078059521,5.0197204699640405], -(2><2) - [ 1.9862461923890056, -1.0127225830525157e-2 - , -1.0127225830525157e-2, 3.0373954915729318 ]) +>>> mat 5 [1..15] ! 1 ! 3 +9.0 -} -meanCov :: Matrix Double -> (Vector Double, Matrix Double) -meanCov x = (med,cov) where - r = rows x - k = 1 / fromIntegral r - med = konst k r `vXm` x - meds = konst 1 r `outer` med - xc = x `sub` meds - cov = scale (recip (fromIntegral (r-1))) (trans xc `mXm` xc) +class Indexable c t | c -> t , t -> c + where + infixl 9 ! + (!) :: c -> Int -> t + +instance Indexable (Vector Double) Double + where + (!) = (@>) + +instance Indexable (Vector Float) Float + where + (!) = (@>) + +instance Indexable (Vector (Complex Double)) (Complex Double) + where + (!) = (@>) + +instance Indexable (Vector (Complex Float)) (Complex Float) + where + (!) = (@>) + +instance Element t => Indexable (Matrix t) (Vector t) + where + m!j = subVector (j*c) c (flatten m) + where + c = cols m -------------------------------------------------------------------------------- @@ -220,7 +325,7 @@ pairwiseD2 x y | ok = x2 `outer` oy + ox `outer` y2 - 2* x <> trans y ox = one (rows x) oy = one (rows y) oc = one (cols x) - one k = constant 1 k + one k = konst 1 k x2 = x * x <> oc y2 = y * y <> oc ok = cols x == cols y diff --git a/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs b/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs index d21602d..5e2ea84 100644 --- a/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs +++ b/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs @@ -6,7 +6,7 @@ module Numeric.LinearAlgebra.Util.CG( CGMat, CGState(..), R, V ) where -import Numeric.Container +import Data.Packed.Numeric import Numeric.Vector() {- diff --git a/packages/base/src/Numeric/LinearAlgebra/Util/Convolution.hs b/packages/base/src/Numeric/LinearAlgebra/Util/Convolution.hs index e4cba8f..c8c7536 100644 --- a/packages/base/src/Numeric/LinearAlgebra/Util/Convolution.hs +++ b/packages/base/src/Numeric/LinearAlgebra/Util/Convolution.hs @@ -16,16 +16,18 @@ module Numeric.LinearAlgebra.Util.Convolution( corr2, conv2, separable ) where -import Numeric.Container +import Data.Packed.Numeric vectSS :: Element t => Int -> Vector t -> Matrix t vectSS n v = fromRows [ subVector k n v | k <- [0 .. dim v - n] ] -corr :: Product t => Vector t -- ^ kernel - -> Vector t -- ^ source - -> Vector t +corr + :: (Container Vector t, Product t) + => Vector t -- ^ kernel + -> Vector t -- ^ source + -> Vector t {- ^ correlation >>> corr (fromList[1,2,3]) (fromList [1..10]) @@ -33,12 +35,12 @@ fromList [14.0,20.0,26.0,32.0,38.0,44.0,50.0,56.0] -} corr ker v - | dim ker == 0 = constant 0 (dim v) + | dim ker == 0 = konst 0 (dim v) | dim ker <= dim v = vectSS (dim ker) v <> ker | otherwise = error $ "corr: dim kernel ("++show (dim ker)++") > dim vector ("++show (dim v)++")" -conv :: (Product t, Num t) => Vector t -> Vector t -> Vector t +conv :: (Container Vector t, Product t, Num t) => Vector t -> Vector t -> Vector t {- ^ convolution ('corr' with reversed kernel and padded input, equivalent to polynomial product) >>> conv (fromList[1,1]) (fromList [-1,1]) @@ -46,12 +48,12 @@ fromList [-1.0,0.0,1.0] -} conv ker v - | dim ker == 0 = constant 0 (dim v) + | dim ker == 0 = konst 0 (dim v) | otherwise = corr ker' v' where ker' = (flatten.fliprl.asRow) ker v' = vjoin [z,v,z] - z = constant 0 (dim ker -1) + z = konst 0 (dim ker -1) corrMin :: (Container Vector t, RealElement t, Product t) => Vector t diff --git a/packages/base/src/Numeric/Sparse.hs b/packages/base/src/Numeric/Sparse.hs index 1957d3a..2df4578 100644 --- a/packages/base/src/Numeric/Sparse.hs +++ b/packages/base/src/Numeric/Sparse.hs @@ -10,7 +10,7 @@ module Numeric.Sparse( smXv )where -import Numeric.Container +import Data.Packed.Numeric import qualified Data.Vector.Storable as V import Data.Function(on) import Control.Arrow((***)) @@ -133,10 +133,6 @@ instance Transposable SMatrix tr (CSC vs rs cs n m) = CSR vs rs cs m n tr (Diag v n m) = Diag v m n -instance Transposable (Matrix Double) - where - tr = trans - instance CGMat SMatrix instance CGMat (Matrix Double) -- cgit v1.2.3