From 29099e3bfb4eec87ac3d4d675d7cfc82234c20d6 Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Mon, 6 Sep 2010 07:37:18 +0000 Subject: working on conversion / linear --- lib/Data/Packed/Matrix.hs | 2 +- lib/Data/Packed/Random.hs | 3 +- lib/Graphics/Plot.hs | 1 + lib/Numeric/Container.hs | 51 ++++++++++------------- lib/Numeric/LinearAlgebra.hs | 6 ++- lib/Numeric/LinearAlgebra/Algorithms.hs | 2 + lib/Numeric/LinearAlgebra/Linear.hs | 72 +++++++++++++-------------------- lib/Numeric/Matrix.hs | 9 +++-- lib/Numeric/Vector.hs | 48 +++++++++++++++++++++- 9 files changed, 112 insertions(+), 82 deletions(-) diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs index 07258f8..b860cb1 100644 --- a/lib/Data/Packed/Matrix.hs +++ b/lib/Data/Packed/Matrix.hs @@ -19,7 +19,7 @@ ----------------------------------------------------------------------------- module Data.Packed.Matrix ( - Element(..), + Element, Matrix,rows,cols, (><), trans, diff --git a/lib/Data/Packed/Random.hs b/lib/Data/Packed/Random.hs index 579d13c..0bc5350 100644 --- a/lib/Data/Packed/Random.hs +++ b/lib/Data/Packed/Random.hs @@ -21,10 +21,11 @@ module Data.Packed.Random ( import Numeric.GSL.Vector import Data.Packed.Matrix -import Data.Packed.Vector +import Numeric.Vector import Numeric.LinearAlgebra.Linear import Numeric.LinearAlgebra.Algorithms + -- | Obtains a matrix whose rows are pseudorandom samples from a multivariate -- Gaussian distribution. gaussianSample :: Int -- ^ seed diff --git a/lib/Graphics/Plot.hs b/lib/Graphics/Plot.hs index 478b4ad..ae88a92 100644 --- a/lib/Graphics/Plot.hs +++ b/lib/Graphics/Plot.hs @@ -29,6 +29,7 @@ module Graphics.Plot( ) where import Data.Packed +import Numeric.Vector import Numeric.LinearAlgebra.Linear import Data.List(intersperse) import System.Process (system) diff --git a/lib/Numeric/Container.hs b/lib/Numeric/Container.hs index 0bec2e8..010235f 100644 --- a/lib/Numeric/Container.hs +++ b/lib/Numeric/Container.hs @@ -19,10 +19,14 @@ ----------------------------------------------------------------------------- module Numeric.Container ( - RealElement, AutoReal(..), - Container(..), Linear(..), - Convert(..), RealOf, ComplexOf, SingleOf, DoubleOf, ElementOf, IndexOf, - Precision(..), comp, + Container(..), RealElement, Precision, NumericContainer(..), comp, + Convert(..), AutoReal(..), + RealOf, ComplexOf, SingleOf, DoubleOf, + +-- ElementOf, + + IndexOf, + module Data.Complex ) where @@ -62,7 +66,7 @@ instance RealElement Double instance RealElement Float -- | Conversion utilities -class Container c where +class NumericContainer c where toComplex :: (RealElement e) => (c e, c e) -> c (Complex e) fromComplex :: (RealElement e) => c (Complex e) -> (c e, c e) complex' :: (RealElement e) => c e -> c (Complex e) @@ -71,9 +75,11 @@ class Container c where single' :: Precision a b => c b -> c a double' :: Precision a b => c a -> c b +-- | a synonym for "complex'" +comp :: (NumericContainer c, RealElement e) => c e -> c (Complex e) comp x = complex' x -instance Container Vector where +instance NumericContainer Vector where toComplex = toComplexV fromComplex = fromComplexV complex' v = toComplex (v,constantD 0 (dim v)) @@ -82,7 +88,7 @@ instance Container Vector where single' = double2FloatG double' = float2DoubleG -instance Container Matrix where +instance NumericContainer Matrix where toComplex = uncurry $ liftMatrix2 $ curry toComplex fromComplex z = (reshape c *** reshape c) . fromComplex . flatten $ z where c = cols z @@ -138,10 +144,10 @@ type instance IndexOf Matrix = (Int,Int) -- | generic conversion functions class Convert t where - real :: Container c => c (RealOf t) -> c t - complex :: Container c => c t -> c (ComplexOf t) - single :: Container c => c t -> c (SingleOf t) - double :: Container c => c t -> c (DoubleOf t) + real :: NumericContainer c => c (RealOf t) -> c t + complex :: NumericContainer c => c t -> c (ComplexOf t) + single :: NumericContainer c => c t -> c (SingleOf t) + double :: NumericContainer c => c t -> c (DoubleOf t) instance Convert Double where real = id @@ -171,8 +177,8 @@ instance Convert (Complex Float) where -- | to be replaced by Convert class Convert t => AutoReal t where - real'' :: Container c => c Double -> c t - complex'' :: Container c => c t -> c (Complex Double) + real'' :: NumericContainer c => c Double -> c t + complex'' :: NumericContainer c => c t -> c (Complex Double) instance AutoReal Double where real'' = real @@ -193,23 +199,7 @@ instance AutoReal (Complex Float) where ------------------------------------------------------------------- -- | Basic element-by-element functions. -class (Element e, Container c) => Linear c e where - -- | create a structure with a single element - scalar :: e -> c e - scale :: e -> c e -> c e - -- | scale the element by element reciprocal of the object: - -- - -- @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@ - scaleRecip :: e -> c e -> c e - addConstant :: e -> c e -> c e - add :: c e -> c e -> c e - sub :: c e -> c e -> c e - -- | element by element multiplication - mul :: c e -> c e -> c e - -- | element by element division - divide :: c e -> c e -> c e - equal :: c e -> c e -> Bool - -- +class (Element e) => Container c e where minIndex :: c e -> IndexOf c maxIndex :: c e -> IndexOf c minElement :: c e -> e @@ -217,3 +207,4 @@ class (Element e, Container c) => Linear c e where + diff --git a/lib/Numeric/LinearAlgebra.hs b/lib/Numeric/LinearAlgebra.hs index 3df9bd7..da4a4b5 100644 --- a/lib/Numeric/LinearAlgebra.hs +++ b/lib/Numeric/LinearAlgebra.hs @@ -15,9 +15,13 @@ This module reexports all normally required functions for Linear Algebra applica module Numeric.LinearAlgebra ( module Numeric.LinearAlgebra.Algorithms, module Numeric.LinearAlgebra.Interface, - module Numeric.LinearAlgebra.Linear + module Numeric.LinearAlgebra.Linear, + module Numeric.Matrix, + module Numeric.Vector ) where import Numeric.LinearAlgebra.Algorithms import Numeric.LinearAlgebra.Interface import Numeric.LinearAlgebra.Linear +import Numeric.Matrix +import Numeric.Vector diff --git a/lib/Numeric/LinearAlgebra/Algorithms.hs b/lib/Numeric/LinearAlgebra/Algorithms.hs index 14bf5d8..ac46847 100644 --- a/lib/Numeric/LinearAlgebra/Algorithms.hs +++ b/lib/Numeric/LinearAlgebra/Algorithms.hs @@ -81,6 +81,8 @@ import Numeric.LinearAlgebra.Linear import Numeric.LinearAlgebra.LAPACK as LAPACK import Data.List(foldl1') import Data.Array +import Numeric.Vector +import Numeric.Matrix() -- | Auxiliary typeclass used to define generic computations for both real and complex matrices. class (Product t, Linear Vector t, Linear Matrix t) => Field t where diff --git a/lib/Numeric/LinearAlgebra/Linear.hs b/lib/Numeric/LinearAlgebra/Linear.hs index 9048204..952661d 100644 --- a/lib/Numeric/LinearAlgebra/Linear.hs +++ b/lib/Numeric/LinearAlgebra/Linear.hs @@ -18,24 +18,25 @@ Basic optimized operations on vectors and matrices. module Numeric.LinearAlgebra.Linear ( -- * Linear Algebra Typeclasses - Vectors(..), + Vectors(..), Linear(..), -- * Products Product(..), mXm,mXv,vXm, outer, kronecker, -- * Modules - module Numeric.Vector, - module Numeric.Matrix, + --module Numeric.Vector, + --module Numeric.Matrix, module Numeric.Container ) where import Data.Packed.Internal.Common import Data.Packed.Matrix +import Data.Packed.Vector import Data.Complex import Numeric.Container -import Numeric.Vector -import Numeric.Matrix -import Numeric.GSL.Vector +--import Numeric.Vector +--import Numeric.Matrix +--import Numeric.GSL.Vector import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ) -- | basic Vector functions @@ -49,43 +50,6 @@ class Num e => Vectors a e where norm2 :: a e -> e normInf :: a e -> e - -instance Vectors Vector Float where - vectorSum = sumF - vectorProd = prodF - norm2 = toScalarF Norm2 - absSum = toScalarF AbsSum - dot = dotF - norm1 = toScalarF AbsSum - normInf = maxElement . vectorMapF Abs - -instance Vectors Vector Double where - vectorSum = sumR - vectorProd = prodR - norm2 = toScalarR Norm2 - absSum = toScalarR AbsSum - dot = dotR - norm1 = toScalarR AbsSum - normInf = maxElement . vectorMapR Abs - -instance Vectors Vector (Complex Float) where - vectorSum = sumQ - vectorProd = prodQ - norm2 = (:+ 0) . toScalarQ Norm2 - absSum = (:+ 0) . toScalarQ AbsSum - dot = dotQ - norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapQ Abs - normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapQ Abs - -instance Vectors Vector (Complex Double) where - vectorSum = sumC - vectorProd = prodC - norm2 = (:+ 0) . toScalarC Norm2 - absSum = (:+ 0) . toScalarC AbsSum - dot = dotC - norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapC Abs - normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapC Abs - ---------------------------------------------------- class Element t => Product t where @@ -162,4 +126,24 @@ kronecker a b = fromBlocks . toRows $ flatten a `outer` flatten b --------------------------------------------------- + +------------------------------------------------------------------- + + +-- | Basic element-by-element functions. +class (Element e, Container c e) => Linear c e where + -- | create a structure with a single element + scalar :: e -> c e + scale :: e -> c e -> c e + -- | scale the element by element reciprocal of the object: + -- + -- @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@ + scaleRecip :: e -> c e -> c e + addConstant :: e -> c e -> c e + add :: c e -> c e -> c e + sub :: c e -> c e -> c e + -- | element by element multiplication + mul :: c e -> c e -> c e + -- | element by element division + divide :: c e -> c e -> c e + equal :: c e -> c e -> Bool diff --git a/lib/Numeric/Matrix.hs b/lib/Numeric/Matrix.hs index 8d3764a..f240384 100644 --- a/lib/Numeric/Matrix.hs +++ b/lib/Numeric/Matrix.hs @@ -27,7 +27,8 @@ module Numeric.Matrix ( import Data.Packed.Vector import Data.Packed.Matrix import Numeric.Container -import Numeric.Vector() +import Numeric.LinearAlgebra.Linear +--import Numeric.Vector import Control.Monad(ap) @@ -74,7 +75,7 @@ instance (Linear Vector a, Floating (Vector a), Fractional (Matrix a)) => Floati --------------------------------------------------------------- -instance (Linear Vector a, Container Matrix) => (Linear Matrix a) where +instance (Linear Vector a, NumericContainer Matrix) => (Linear Matrix a) where scale x = liftMatrix (scale x) scaleRecip x = liftMatrix (scaleRecip x) addConstant x = liftMatrix (addConstant x) @@ -84,6 +85,9 @@ instance (Linear Vector a, Container Matrix) => (Linear Matrix a) where divide = liftMatrix2 divide equal a b = cols a == cols b && flatten a `equal` flatten b scalar x = (1><1) [x] + + +instance (Linear Vector a, NumericContainer Matrix) => (Container Matrix a) where minIndex m = let (r,c) = (rows m,cols m) i = 1 + (minIndex $ flatten m) in (i `div` r,i `mod` r) @@ -94,4 +98,3 @@ instance (Linear Vector a, Container Matrix) => (Linear Matrix a) where maxElement = ap (@@>) maxIndex ---------------------------------------------------- - diff --git a/lib/Numeric/Vector.hs b/lib/Numeric/Vector.hs index ced202f..d92a5e4 100644 --- a/lib/Numeric/Vector.hs +++ b/lib/Numeric/Vector.hs @@ -29,10 +29,11 @@ import Data.Complex import Control.Monad(ap) import Data.Packed.Vector -import Data.Packed.Matrix(Element(..)) +import Data.Packed.Internal.Matrix(Element(..)) import Numeric.GSL.Vector import Numeric.Container +import Numeric.LinearAlgebra.Linear ------------------------------------------------------------------- @@ -263,12 +264,13 @@ instance Linear Vector Float where divide = vectorZipF Div equal u v = dim u == dim v && maxElement (vectorMapF Abs (sub u v)) == 0.0 scalar x = fromList [x] + +instance Container Vector Float where minIndex = round . toScalarF MinIdx maxIndex = round . toScalarF MaxIdx minElement = toScalarF Min maxElement = toScalarF Max - instance Linear Vector Double where scale = vectorMapValR Scale scaleRecip = vectorMapValR Recip @@ -279,6 +281,8 @@ instance Linear Vector Double where divide = vectorZipR Div equal u v = dim u == dim v && maxElement (vectorMapR Abs (sub u v)) == 0.0 scalar x = fromList [x] + +instance Container Vector Double where minIndex = round . toScalarR MinIdx maxIndex = round . toScalarR MaxIdx minElement = toScalarR Min @@ -294,6 +298,8 @@ instance Linear Vector (Complex Double) where divide = vectorZipC Div equal u v = dim u == dim v && maxElement (mapVector magnitude (sub u v)) == 0.0 scalar x = fromList [x] + +instance Container Vector (Complex Double) where minIndex = minIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) maxIndex = maxIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) minElement = ap (@>) minIndex @@ -309,9 +315,47 @@ instance Linear Vector (Complex Float) where divide = vectorZipQ Div equal u v = dim u == dim v && maxElement (mapVector magnitude (sub u v)) == 0.0 scalar x = fromList [x] + +instance Container Vector (Complex Float) where minIndex = minIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) maxIndex = maxIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) minElement = ap (@>) minIndex maxElement = ap (@>) maxIndex --------------------------------------------------------------- + +instance Vectors Vector Float where + vectorSum = sumF + vectorProd = prodF + norm2 = toScalarF Norm2 + absSum = toScalarF AbsSum + dot = dotF + norm1 = toScalarF AbsSum + normInf = maxElement . vectorMapF Abs + +instance Vectors Vector Double where + vectorSum = sumR + vectorProd = prodR + norm2 = toScalarR Norm2 + absSum = toScalarR AbsSum + dot = dotR + norm1 = toScalarR AbsSum + normInf = maxElement . vectorMapR Abs + +instance Vectors Vector (Complex Float) where + vectorSum = sumQ + vectorProd = prodQ + norm2 = (:+ 0) . toScalarQ Norm2 + absSum = (:+ 0) . toScalarQ AbsSum + dot = dotQ + norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapQ Abs + normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapQ Abs + +instance Vectors Vector (Complex Double) where + vectorSum = sumC + vectorProd = prodC + norm2 = (:+ 0) . toScalarC Norm2 + absSum = (:+ 0) . toScalarC AbsSum + dot = dotC + norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapC Abs + normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapC Abs \ No newline at end of file -- cgit v1.2.3