From 5db2ed78986bbc737b82e428392ee63999c8abfd Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Wed, 17 Jun 2009 20:52:14 +0000 Subject: vector fread/fwrite, fscanf/fprinf --- lib/Data/Packed/Internal/Matrix.hs | 9 +------ lib/Data/Packed/Internal/Vector.hs | 47 ++++++++++++++++++++++++++++++++- lib/Data/Packed/Vector.hs | 6 +++-- lib/Graphics/Plot.hs | 6 ++--- lib/Numeric/GSL/gsl-aux.c | 54 +++++++++++++++++++++++++++++++++----- lib/Numeric/GSL/gsl-aux.h | 40 ---------------------------- lib/Numeric/LinearAlgebra.hs | 2 +- 7 files changed, 103 insertions(+), 61 deletions(-) delete mode 100644 lib/Numeric/GSL/gsl-aux.h diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs index ccc652a..01d2ccf 100644 --- a/lib/Data/Packed/Internal/Matrix.hs +++ b/lib/Data/Packed/Internal/Matrix.hs @@ -22,7 +22,6 @@ import Data.Packed.Internal.Vector import Foreign hiding (xor) import Complex -import Foreign.C.String import Foreign.C.Types ----------------------------------------------------------------- @@ -353,10 +352,4 @@ fromComplex z = (r,i) where -- | loads a matrix from an ASCII file (the number of rows and columns must be known in advance). fromFile :: FilePath -> (Int,Int) -> IO (Matrix Double) -fromFile filename (r,c) = do - charname <- newCString filename - res <- createMatrix RowMajor r c - app1 (c_gslReadMatrix charname) mat res "gslReadMatrix" - free charname - return res -foreign import ccall "matrix_fscanf" c_gslReadMatrix:: Ptr CChar -> TM +fromFile filename (r,c) = reshape c `fmap` fscanfVector filename (r*c) diff --git a/lib/Data/Packed/Internal/Vector.hs b/lib/Data/Packed/Internal/Vector.hs index 1b572a5..5784861 100644 --- a/lib/Data/Packed/Internal/Vector.hs +++ b/lib/Data/Packed/Internal/Vector.hs @@ -18,7 +18,8 @@ module Data.Packed.Internal.Vector where import Data.Packed.Internal.Common import Foreign -import Foreign.C.Types(CInt) +import Foreign.C.String +import Foreign.C.Types(CInt,CChar) import Complex import Control.Monad(when) @@ -255,3 +256,47 @@ foldVectorG f s0 v = foldLoop g s0 (dim v) where g !k !s = f k (at' v) s {-# INLINE g #-} -- Thanks to Ryan Ingram (http://permalink.gmane.org/gmane.comp.lang.haskell.cafe/46479) {-# INLINE foldVectorG #-} + +------------------------------------------------------------------- + +-- | Loads a vector from an ASCII file (the number of elements must be known in advance). +fscanfVector :: FilePath -> Int -> IO (Vector Double) +fscanfVector filename n = do + charname <- newCString filename + res <- createVector n + app1 (gsl_vector_fscanf charname) vec res "gsl_vector_fscanf" + free charname + return res + +foreign import ccall "vector_fscanf" gsl_vector_fscanf:: Ptr CChar -> TV + +-- | Saves the elements of a vector, with a given format (%f, %e, %g), to an ASCII file. +fprintfVector :: FilePath -> String -> Vector Double -> IO () +fprintfVector filename fmt v = do + charname <- newCString filename + charfmt <- newCString fmt + app1 (gsl_vector_fprintf charname charfmt) vec v "gsl_vector_fprintf" + free charname + free charfmt + +foreign import ccall "vector_fprintf" gsl_vector_fprintf :: Ptr CChar -> Ptr CChar -> TV + +-- | Loads a vector from a binary file (the number of elements must be known in advance). +freadVector :: FilePath -> Int -> IO (Vector Double) +freadVector filename n = do + charname <- newCString filename + res <- createVector n + app1 (gsl_vector_fread charname) vec res "gsl_vector_fread" + free charname + return res + +foreign import ccall "vector_fread" gsl_vector_fread:: Ptr CChar -> TV + +-- | Saves the elements of a vector to a binary file. +fwriteVector :: FilePath -> Vector Double -> IO () +fwriteVector filename v = do + charname <- newCString filename + app1 (gsl_vector_fwrite charname) vec v "gsl_vector_fwrite" + free charname + +foreign import ccall "vector_fwrite" gsl_vector_fwrite :: Ptr CChar -> TV diff --git a/lib/Data/Packed/Vector.hs b/lib/Data/Packed/Vector.hs index e53d455..c657e10 100644 --- a/lib/Data/Packed/Vector.hs +++ b/lib/Data/Packed/Vector.hs @@ -8,7 +8,7 @@ -- Stability : provisional -- Portability : portable -- --- A representation of 1D arrays suitable for numeric computations using external libraries. +-- 1D arrays suitable for numeric computations using external libraries. -- ----------------------------------------------------------------------------- @@ -19,8 +19,10 @@ module Data.Packed.Vector ( subVector, join, constant, linspace, vectorMax, vectorMin, vectorMaxIndex, vectorMinIndex, + mapVector, zipVector, + fscanfVector, fprintfVector, freadVector, fwriteVector, liftVector, liftVector2, - foldLoop, foldVector, foldVectorG, mapVector, zipVector + foldLoop, foldVector, foldVectorG ) where import Data.Packed.Internal diff --git a/lib/Graphics/Plot.hs b/lib/Graphics/Plot.hs index 9b64a1d..e6a9098 100644 --- a/lib/Graphics/Plot.hs +++ b/lib/Graphics/Plot.hs @@ -40,8 +40,8 @@ size = dim --fromFile filename = readFile filename >>= return . readMatrix read -- | Saves a real matrix to a formatted ascii text file -toFile :: FilePath -> Matrix Double -> IO () -toFile filename matrix = writeFile filename (unlines . map unwords. map (map show) . toLists $ matrix) +toFile' :: FilePath -> Matrix Double -> IO () +toFile' filename matrix = writeFile filename (unlines . map unwords. map (map show) . toLists $ matrix) ------------------------------------------------------------------------ @@ -74,7 +74,7 @@ mesh m = gnuplotX (command++dat) where mesh' :: Matrix Double -> IO () mesh' m = do writeFile "splot-gnu-command" "splot \"splot-tmp.txt\" matrix with lines; pause -1"; - toFile "splot-tmp.txt" m + toFile' "splot-tmp.txt" m putStr "Press [Return] to close the graphic and continue... " system "gnuplot -persist splot-gnu-command" system "rm splot-tmp.txt splot-gnu-command" diff --git a/lib/Numeric/GSL/gsl-aux.c b/lib/Numeric/GSL/gsl-aux.c index 0904a67..b3b524d 100644 --- a/lib/Numeric/GSL/gsl-aux.c +++ b/lib/Numeric/GSL/gsl-aux.c @@ -1,4 +1,15 @@ -#include "gsl-aux.h" +#include + +#define RVEC(A) int A##n, double*A##p +#define RMAT(A) int A##r, int A##c, double* A##p +#define KRVEC(A) int A##n, const double*A##p +#define KRMAT(A) int A##r, int A##c, const double* A##p + +#define CVEC(A) int A##n, gsl_complex*A##p +#define CMAT(A) int A##r, int A##c, gsl_complex* A##p +#define KCVEC(A) int A##n, const gsl_complex*A##p +#define KCMAT(A) int A##r, int A##c, const gsl_complex* A##p + #include #include #include @@ -8,7 +19,6 @@ #include #include #include -#include #include #include #include @@ -340,13 +350,45 @@ int polySolve(KRVEC(a), CVEC(z)) { OK; } -int matrix_fscanf(char*filename, RMAT(a)) { +int vector_fscanf(char*filename, RVEC(a)) { DEBUGMSG("gsl_matrix_fscanf"); - //printf(filename); printf("\n"); - DMVIEW(a); + DVVIEW(a); FILE * f = fopen(filename,"r"); CHECK(!f,BAD_FILE); - int res = gsl_matrix_fscanf(f, M(a)); + int res = gsl_vector_fscanf(f,V(a)); + CHECK(res,res); + fclose (f); + OK +} + +int vector_fprintf(char*filename, char*fmt, RVEC(a)) { + DEBUGMSG("gsl_vector_fprintf"); + DVVIEW(a); + FILE * f = fopen(filename,"w"); + CHECK(!f,BAD_FILE); + int res = gsl_vector_fprintf(f,V(a),fmt); + CHECK(res,res); + fclose (f); + OK +} + +int vector_fread(char*filename, RVEC(a)) { + DEBUGMSG("gsl_matrix_fscanf"); + DVVIEW(a); + FILE * f = fopen(filename,"r"); + CHECK(!f,BAD_FILE); + int res = gsl_vector_fread(f,V(a)); + CHECK(res,res); + fclose (f); + OK +} + +int vector_fwrite(char*filename, RVEC(a)) { + DEBUGMSG("gsl_vector_fprintf"); + DVVIEW(a); + FILE * f = fopen(filename,"w"); + CHECK(!f,BAD_FILE); + int res = gsl_vector_fwrite(f,V(a)); CHECK(res,res); fclose (f); OK diff --git a/lib/Numeric/GSL/gsl-aux.h b/lib/Numeric/GSL/gsl-aux.h deleted file mode 100644 index 881d0d0..0000000 --- a/lib/Numeric/GSL/gsl-aux.h +++ /dev/null @@ -1,40 +0,0 @@ -#include - -#define RVEC(A) int A##n, double*A##p -#define RMAT(A) int A##r, int A##c, double* A##p -#define KRVEC(A) int A##n, const double*A##p -#define KRMAT(A) int A##r, int A##c, const double* A##p - -#define CVEC(A) int A##n, gsl_complex*A##p -#define CMAT(A) int A##r, int A##c, gsl_complex* A##p -#define KCVEC(A) int A##n, const gsl_complex*A##p -#define KCMAT(A) int A##r, int A##c, const gsl_complex* A##p - -void no_abort_on_error(); - -int toScalarR(int code, KRVEC(x), RVEC(r)); -/* norm2, absdif, maximum, posmax, etc. */ - -int mapR(int code, KRVEC(x), RVEC(r)); -int mapC(int code, KCVEC(x), CVEC(r)); -/* sin cos tan etc. */ - -int mapValR(int code, double*, KRVEC(x), RVEC(r)); -int mapValC(int code, gsl_complex*, KCVEC(x), CVEC(r)); - -int zipR(int code, KRVEC(a), KRVEC(b), RVEC(r)); -int zipC(int code, KCVEC(a), KCVEC(b), CVEC(r)); - - -int fft(int code, KCVEC(a), CVEC(b)); - -int deriv(int code, double f(double, void*), double x, double h, double * result, double * abserr); - -int integrate_qng(double f(double, void*), double a, double b, double prec, - double *result, double*error); - -int integrate_qags(double f(double,void*), double a, double b, double prec, int w, - double *result, double* error); - -int polySolve(KRVEC(a), CVEC(z)); - diff --git a/lib/Numeric/LinearAlgebra.hs b/lib/Numeric/LinearAlgebra.hs index b8fd01e..f92c40d 100644 --- a/lib/Numeric/LinearAlgebra.hs +++ b/lib/Numeric/LinearAlgebra.hs @@ -10,7 +10,7 @@ Portability : uses ffi Basic matrix computations implemented by BLAS, LAPACK and GSL. -This is module reexports the most comon functions (including "Numeric.LinearAlgebra.Instances"). +This module reexports the most comon functions (including "Numeric.LinearAlgebra.Instances"). -} ----------------------------------------------------------------------------- -- cgit v1.2.3