diff --git a/CHANGELOG b/CHANGELOG --- a/CHANGELOG +++ b/CHANGELOG @@ -1,3 +1,10 @@ +0.16.1.0 +-------- + + * Added (|||) and (===) for "besides" and "above" + + * rowOuters + 0.16.0.0 -------- diff --git a/THANKS.md b/THANKS.md --- a/THANKS.md +++ b/THANKS.md @@ -103,7 +103,7 @@ deprecation warnings, used more explicit imports, and updated to ghc-7.4. - Tom Nielsen discovered a problem in Config.hs, exposed by link problems - in Ubuntu 11.10 beta. + in Ubuntu 11.10 beta, and fixed the link options on freebsd. - Daniel Fischer reported some Haddock markup errors. @@ -170,4 +170,10 @@ - Ian Ross reported the max/minIndex bug. - Niklas Hambüchen improved the documentation. + +- "erdeszt" optimized "conv" using a direct vector reverse. + +- John Shahbazian added support for openBLAS. + +- "yongqli" reported the bug in randomVector (rand() is not thread-safe). diff --git a/hmatrix.cabal b/hmatrix.cabal --- a/hmatrix.cabal +++ b/hmatrix.cabal @@ -1,5 +1,5 @@ Name: hmatrix -Version: 0.16.0.6 +Version: 0.16.1.0 License: BSD3 License-file: LICENSE Author: Alberto Ruiz @@ -32,6 +32,10 @@ extra-source-files: src/C/lapack-aux.h +flag openblas + description: Link with OpenBLAS (https://github.com/xianyi/OpenBLAS) optimized libraries. + default: False + library Build-Depends: base >= 4 && < 5, @@ -100,7 +104,11 @@ cpp-options: -DBINARY - extra-libraries: blas lapack + if flag(openblas) + extra-lib-dirs: /usr/lib/openblas/lib + extra-libraries: openblas + else + extra-libraries: blas lapack if os(OSX) extra-lib-dirs: /opt/local/lib/ @@ -114,7 +122,7 @@ if os(freebsd) extra-lib-dirs: /usr/local/lib include-dirs: /usr/local/include - extra-libraries: blas lapack + extra-libraries: blas lapack gfortran if os(windows) extra-libraries: blas lapack diff --git a/src/C/vector-aux.c b/src/C/vector-aux.c --- a/src/C/vector-aux.c +++ b/src/C/vector-aux.c @@ -700,17 +700,24 @@ //////////////////////////////////////////////////////////////////////////////// +inline double urandom(struct drand48_data * buffer) { + double res; + drand48_r(buffer,&res); + return res; +} + + // http://c-faq.com/lib/gaussian.html -double gaussrand() +double gaussrand(struct drand48_data *buffer, + int *phase, double *pV1, double *pV2, double *pS) { - static double V1, V2, S; - static int phase = 0; + double V1=*pV1, V2=*pV2, S=*pS; double X; - if(phase == 0) { + if(*phase == 0) { do { - double U1 = (double)rand() / RAND_MAX; - double U2 = (double)rand() / RAND_MAX; + double U1 = urandom(buffer); + double U2 = urandom(buffer); V1 = 2 * U1 - 1; V2 = 2 * U2 - 1; @@ -721,24 +728,29 @@ } else X = V2 * sqrt(-2 * log(S) / S); - phase = 1 - phase; + *phase = 1 - *phase; + *pV1=V1; *pV2=V2; *pS=S; return X; } -int random_vector(int seed, int code, DVEC(r)) { - srand(seed); +int random_vector(unsigned int seed, int code, DVEC(r)) { + struct drand48_data buffer; + srand48_r(seed,&buffer); + int phase = 0; + double V1,V2,S; + int k; switch (code) { case 0: { // uniform for (k=0; k>> mapMatrixWithIndex (\\(i,j) v -> 100*v + 10*fromIntegral i + fromIntegral j) (ident 3:: Matrix Double) +>>> mapMatrixWithIndex (\(i,j) v -> 100*v + 10*fromIntegral i + fromIntegral j) (ident 3:: Matrix Double) (3><3) [ 100.0, 1.0, 2.0 , 10.0, 111.0, 12.0 diff --git a/src/Numeric/LinearAlgebra/Algorithms.hs b/src/Numeric/LinearAlgebra/Algorithms.hs --- a/src/Numeric/LinearAlgebra/Algorithms.hs +++ b/src/Numeric/LinearAlgebra/Algorithms.hs @@ -809,7 +809,7 @@ -------------------------------------------------------------- {- | Matrix square root. Currently it uses a simple iterative algorithm described in Wikipedia. -It only works with invertible matrices that have a real solution. For diagonalizable matrices you can try @matFunc sqrt@. +It only works with invertible matrices that have a real solution. @m = (2><2) [4,9 ,0,4] :: Matrix Double@ @@ -818,6 +818,13 @@ (2><2) [ 2.0, 2.25 , 0.0, 2.0 ] + +For diagonalizable matrices you can try 'matFunc' @sqrt@: + +>>> matFunc sqrt ((2><2) [1,0,0,-1]) +(2><2) + [ 1.0 :+ 0.0, 0.0 :+ 0.0 + , 0.0 :+ 0.0, 0.0 :+ 1.0 ] -} sqrtm :: Field t => Matrix t -> Matrix t diff --git a/src/Numeric/LinearAlgebra/Data.hs b/src/Numeric/LinearAlgebra/Data.hs --- a/src/Numeric/LinearAlgebra/Data.hs +++ b/src/Numeric/LinearAlgebra/Data.hs @@ -40,7 +40,7 @@ takeRows, dropRows, takeColumns, dropColumns, subMatrix, (?), (¿), fliprl, flipud, -- * Block matrix - fromBlocks, (¦), (——), diagBlock, repmat, toBlocks, toBlocksEvery, + fromBlocks, (|||), (===), diagBlock, repmat, toBlocks, toBlocksEvery, -- * Mapping functions conj, cmap, step, cond, @@ -66,7 +66,7 @@ arctan2, rows, cols, separable, - + (¦),(——), module Data.Complex, Vector, Matrix, GMatrix, nRows, nCols diff --git a/src/Numeric/LinearAlgebra/HMatrix.hs b/src/Numeric/LinearAlgebra/HMatrix.hs --- a/src/Numeric/LinearAlgebra/HMatrix.hs +++ b/src/Numeric/LinearAlgebra/HMatrix.hs @@ -134,7 +134,7 @@ Seed, RandDist(..), randomVector, rand, randn, gaussianSample, uniformSample, -- * Misc - meanCov, peps, relativeError, haussholder, optimiseMult, udot, nullspaceSVD, orthSVD, ranksv, + meanCov, rowOuters, peps, relativeError, haussholder, optimiseMult, udot, nullspaceSVD, orthSVD, ranksv, ℝ,ℂ,iC, -- * Auxiliary classes Element, Container, Product, Numeric, LSDiv, diff --git a/src/Numeric/LinearAlgebra/Util.hs b/src/Numeric/LinearAlgebra/Util.hs --- a/src/Numeric/LinearAlgebra/Util.hs +++ b/src/Numeric/LinearAlgebra/Util.hs @@ -33,7 +33,7 @@ diagl, row, col, - (&), (¦), (——), (#), + (&), (¦), (|||), (——), (===), (#), (?), (¿), Indexable(..), size, Numeric, @@ -157,26 +157,35 @@ {- | horizontal concatenation of real matrices - (unicode 0x00a6, broken bar) - ->>> ident 3 ¦ konst 7 (3,4) +>>> ident 3 ||| konst 7 (3,4) (3><7) [ 1.0, 0.0, 0.0, 7.0, 7.0, 7.0, 7.0 , 0.0, 1.0, 0.0, 7.0, 7.0, 7.0, 7.0 , 0.0, 0.0, 1.0, 7.0, 7.0, 7.0, 7.0 ] -} +infixl 3 ||| +(|||) :: Matrix Double -> Matrix Double -> Matrix Double +a ||| b = fromBlocks [[a,b]] + +-- | a synonym for ('|||') (unicode 0x00a6, broken bar) infixl 3 ¦ (¦) :: Matrix Double -> Matrix Double -> Matrix Double -a ¦ b = fromBlocks [[a,b]] +(¦) = (|||) + -- | vertical concatenation of real matrices -- --- (unicode 0x2014, em dash) +(===) :: Matrix Double -> Matrix Double -> Matrix Double +infixl 2 === +a === b = fromBlocks [[a],[b]] + +-- | a synonym for ('===') (unicode 0x2014, em dash) (——) :: Matrix Double -> Matrix Double -> Matrix Double infixl 2 —— -a —— b = fromBlocks [[a],[b]] +(——) = (===) + (#) :: Matrix Double -> Matrix Double -> Matrix Double infixl 2 # a # b = fromBlocks [[a],[b]] @@ -356,7 +365,26 @@ -------------------------------------------------------------------------------- --- | outer products of rows +{- | outer products of rows + +>>> a +(3><2) + [ 1.0, 2.0 + , 10.0, 20.0 + , 100.0, 200.0 ] +>>> b +(3><3) + [ 1.0, 2.0, 3.0 + , 4.0, 5.0, 6.0 + , 7.0, 8.0, 9.0 ] + +>>> rowOuters a (b ||| 1) +(3><8) + [ 1.0, 2.0, 3.0, 1.0, 2.0, 4.0, 6.0, 2.0 + , 40.0, 50.0, 60.0, 10.0, 80.0, 100.0, 120.0, 20.0 + , 700.0, 800.0, 900.0, 100.0, 1400.0, 1600.0, 1800.0, 200.0 ] + +-} rowOuters :: Matrix Double -> Matrix Double -> Matrix Double rowOuters a b = a' * b' where diff --git a/src/Numeric/LinearAlgebra/Util/Convolution.hs b/src/Numeric/LinearAlgebra/Util/Convolution.hs --- a/src/Numeric/LinearAlgebra/Util/Convolution.hs +++ b/src/Numeric/LinearAlgebra/Util/Convolution.hs @@ -16,6 +16,7 @@ corr2, conv2, separable ) where +import qualified Data.Vector.Storable as SV import Data.Packed.Numeric @@ -51,7 +52,7 @@ | dim ker == 0 = konst 0 (dim v) | otherwise = corr ker' v' where - ker' = (flatten.fliprl.asRow) ker + ker' = SV.reverse ker v' = vjoin [z,v,z] z = konst 0 (dim ker -1)