summaryrefslogtreecommitdiff
path: root/lib/Numeric/GSL
diff options
context:
space:
mode:
Diffstat (limited to 'lib/Numeric/GSL')
-rw-r--r--lib/Numeric/GSL/Vector.hs32
-rw-r--r--lib/Numeric/GSL/gsl-aux.c114
2 files changed, 142 insertions, 4 deletions
diff --git a/lib/Numeric/GSL/Vector.hs b/lib/Numeric/GSL/Vector.hs
index 28c3b06..d09323b 100644
--- a/lib/Numeric/GSL/Vector.hs
+++ b/lib/Numeric/GSL/Vector.hs
@@ -14,10 +14,10 @@
14----------------------------------------------------------------------------- 14-----------------------------------------------------------------------------
15 15
16module Numeric.GSL.Vector ( 16module Numeric.GSL.Vector (
17 FunCodeS(..), toScalarR, 17 FunCodeS(..), toScalarR, toScalarF,
18 FunCodeV(..), vectorMapR, vectorMapC, 18 FunCodeV(..), vectorMapR, vectorMapC, vectorMapF,
19 FunCodeSV(..), vectorMapValR, vectorMapValC, 19 FunCodeSV(..), vectorMapValR, vectorMapValC, vectorMapValF,
20 FunCodeVV(..), vectorZipR, vectorZipC, 20 FunCodeVV(..), vectorZipR, vectorZipC, vectorZipF,
21 RandDist(..), randomVector 21 RandDist(..), randomVector
22) where 22) where
23 23
@@ -106,6 +106,12 @@ toScalarR oper = toScalarAux c_toScalarR (fromei oper)
106 106
107foreign import ccall safe "gsl-aux.h toScalarR" c_toScalarR :: CInt -> TVV 107foreign import ccall safe "gsl-aux.h toScalarR" c_toScalarR :: CInt -> TVV
108 108
109-- | obtains different functions of a vector: norm1, norm2, max, min, posmax, posmin, etc.
110toScalarF :: FunCodeS -> Vector Float -> Float
111toScalarF oper = toScalarAux c_toScalarF (fromei oper)
112
113foreign import ccall safe "gsl-aux.h toScalarF" c_toScalarF :: CInt -> TFF
114
109------------------------------------------------------------------ 115------------------------------------------------------------------
110 116
111-- | map of real vectors with given function 117-- | map of real vectors with given function
@@ -120,6 +126,12 @@ vectorMapC oper = vectorMapAux c_vectorMapC (fromei oper)
120 126
121foreign import ccall safe "gsl-aux.h mapC" c_vectorMapC :: CInt -> TCVCV 127foreign import ccall safe "gsl-aux.h mapC" c_vectorMapC :: CInt -> TCVCV
122 128
129-- | map of real vectors with given function
130vectorMapF :: FunCodeV -> Vector Float -> Vector Float
131vectorMapF = vectorMapAux c_vectorMapF
132
133foreign import ccall safe "gsl-aux.h mapF" c_vectorMapF :: CInt -> TFF
134
123------------------------------------------------------------------- 135-------------------------------------------------------------------
124 136
125-- | map of real vectors with given function 137-- | map of real vectors with given function
@@ -134,6 +146,12 @@ vectorMapValC = vectorMapValAux c_vectorMapValC
134 146
135foreign import ccall safe "gsl-aux.h mapValC" c_vectorMapValC :: CInt -> Ptr (Complex Double) -> TCVCV 147foreign import ccall safe "gsl-aux.h mapValC" c_vectorMapValC :: CInt -> Ptr (Complex Double) -> TCVCV
136 148
149-- | map of real vectors with given function
150vectorMapValF :: FunCodeSV -> Float -> Vector Float -> Vector Float
151vectorMapValF oper = vectorMapValAux c_vectorMapValF (fromei oper)
152
153foreign import ccall safe "gsl-aux.h mapValF" c_vectorMapValF :: CInt -> Ptr Float -> TFF
154
137------------------------------------------------------------------- 155-------------------------------------------------------------------
138 156
139-- | elementwise operation on real vectors 157-- | elementwise operation on real vectors
@@ -148,6 +166,12 @@ vectorZipC = vectorZipAux c_vectorZipC
148 166
149foreign import ccall safe "gsl-aux.h zipC" c_vectorZipC :: CInt -> TCVCVCV 167foreign import ccall safe "gsl-aux.h zipC" c_vectorZipC :: CInt -> TCVCVCV
150 168
169-- | elementwise operation on real vectors
170vectorZipF :: FunCodeVV -> Vector Float -> Vector Float -> Vector Float
171vectorZipF = vectorZipAux c_vectorZipF
172
173foreign import ccall safe "gsl-aux.h zipF" c_vectorZipF :: CInt -> TFFF
174
151----------------------------------------------------------------------- 175-----------------------------------------------------------------------
152 176
153data RandDist = Uniform -- ^ uniform distribution in [0,1) 177data RandDist = Uniform -- ^ uniform distribution in [0,1)
diff --git a/lib/Numeric/GSL/gsl-aux.c b/lib/Numeric/GSL/gsl-aux.c
index f5a95f4..6bb16f0 100644
--- a/lib/Numeric/GSL/gsl-aux.c
+++ b/lib/Numeric/GSL/gsl-aux.c
@@ -10,6 +10,16 @@
10#define KCVEC(A) int A##n, const gsl_complex*A##p 10#define KCVEC(A) int A##n, const gsl_complex*A##p
11#define KCMAT(A) int A##r, int A##c, const gsl_complex* A##p 11#define KCMAT(A) int A##r, int A##c, const gsl_complex* A##p
12 12
13#define FVEC(A) int A##n, float*A##p
14#define FMAT(A) int A##r, int A##c, float* A##p
15#define KFVEC(A) int A##n, const float*A##p
16#define KFMAT(A) int A##r, int A##c, const float* A##p
17
18#define QVEC(A) int A##n, gsl_complex_float*A##p
19#define QMAT(A) int A##r, int A##c, gsl_complex_float* A##p
20#define KQVEC(A) int A##n, const gsl_complex_float*A##p
21#define KQMAT(A) int A##r, int A##c, const gsl_complex_float* A##p
22
13#include <gsl/gsl_blas.h> 23#include <gsl/gsl_blas.h>
14#include <gsl/gsl_math.h> 24#include <gsl/gsl_math.h>
15#include <gsl/gsl_errno.h> 25#include <gsl/gsl_errno.h>
@@ -64,12 +74,24 @@
64#define KCVVIEW(A) gsl_vector_complex_const_view A = gsl_vector_complex_const_view_array((double*)A##p,A##n) 74#define KCVVIEW(A) gsl_vector_complex_const_view A = gsl_vector_complex_const_view_array((double*)A##p,A##n)
65#define KCMVIEW(A) gsl_matrix_complex_const_view A = gsl_matrix_complex_const_view_array((double*)A##p,A##r,A##c) 75#define KCMVIEW(A) gsl_matrix_complex_const_view A = gsl_matrix_complex_const_view_array((double*)A##p,A##r,A##c)
66 76
77#define FVVIEW(A) gsl_vector_float_view A = gsl_vector_float_view_array(A##p,A##n)
78#define FMVIEW(A) gsl_matrix_float_view A = gsl_matrix_float_view_array(A##p,A##r,A##c)
79#define QVVIEW(A) gsl_vector_float_complex_view A = gsl_vector_float_complex_view_array((float*)A##p,A##n)
80#define QMVIEW(A) gsl_matrix_float_complex_view A = gsl_matrix_float_complex_view_array((float*)A##p,A##r,A##c)
81#define KFVVIEW(A) gsl_vector_float_const_view A = gsl_vector_float_const_view_array(A##p,A##n)
82#define KFMVIEW(A) gsl_matrix_float_const_view A = gsl_matrix_float_const_view_array(A##p,A##r,A##c)
83#define KQVVIEW(A) gsl_vector_float_complex_const_view A = gsl_vector_float_complex_const_view_array((float*)A##p,A##n)
84#define KQMVIEW(A) gsl_matrix_float_complex_const_view A = gsl_matrix_float_complex_const_view_array((float*)A##p,A##r,A##c)
85
67#define V(a) (&a.vector) 86#define V(a) (&a.vector)
68#define M(a) (&a.matrix) 87#define M(a) (&a.matrix)
69 88
70#define GCVEC(A) int A##n, gsl_complex*A##p 89#define GCVEC(A) int A##n, gsl_complex*A##p
71#define KGCVEC(A) int A##n, const gsl_complex*A##p 90#define KGCVEC(A) int A##n, const gsl_complex*A##p
72 91
92#define GQVEC(A) int A##n, gsl_complex_float*A##p
93#define KGQVEC(A) int A##n, const gsl_complex_float*A##p
94
73#define BAD_SIZE 2000 95#define BAD_SIZE 2000
74#define BAD_CODE 2001 96#define BAD_CODE 2001
75#define MEM 2002 97#define MEM 2002
@@ -99,6 +121,24 @@ int toScalarR(int code, KRVEC(x), RVEC(r)) {
99 OK 121 OK
100} 122}
101 123
124int toScalarF(int code, KFVEC(x), FVEC(r)) {
125 REQUIRES(rn==1,BAD_SIZE);
126 DEBUGMSG("toScalarF");
127 KFVVIEW(x);
128 float res;
129 switch(code) {
130 case 0: { res = gsl_blas_snrm2(V(x)); break; }
131 case 1: { res = gsl_blas_sasum(V(x)); break; }
132 case 2: { res = gsl_vector_float_max_index(V(x)); break; }
133 case 3: { res = gsl_vector_float_max(V(x)); break; }
134 case 4: { res = gsl_vector_float_min_index(V(x)); break; }
135 case 5: { res = gsl_vector_float_min(V(x)); break; }
136 default: ERROR(BAD_CODE);
137 }
138 rp[0] = res;
139 OK
140}
141
102 142
103inline double sign(double x) { 143inline double sign(double x) {
104 if(x>0) { 144 if(x>0) {
@@ -110,6 +150,16 @@ inline double sign(double x) {
110 } 150 }
111} 151}
112 152
153inline float float_sign(float x) {
154 if(x>0) {
155 return +1.0;
156 } else if (x<0) {
157 return -1.0;
158 } else {
159 return 0.0;
160 }
161}
162
113inline gsl_complex complex_abs(gsl_complex z) { 163inline gsl_complex complex_abs(gsl_complex z) {
114 gsl_complex r; 164 gsl_complex r;
115 r.dat[0] = gsl_complex_abs(z); 165 r.dat[0] = gsl_complex_abs(z);
@@ -159,6 +209,32 @@ int mapR(int code, KRVEC(x), RVEC(r)) {
159 } 209 }
160} 210}
161 211
212int mapF(int code, KFVEC(x), FVEC(r)) {
213 int k;
214 REQUIRES(xn == rn,BAD_SIZE);
215 DEBUGMSG("mapF");
216 switch (code) {
217 OP(0,sin)
218 OP(1,cos)
219 OP(2,tan)
220 OP(3,fabs)
221 OP(4,asin)
222 OP(5,acos)
223 OP(6,atan) /* atan2 mediante vectorZip */
224 OP(7,sinh)
225 OP(8,cosh)
226 OP(9,tanh)
227 OP(10,gsl_asinh)
228 OP(11,gsl_acosh)
229 OP(12,gsl_atanh)
230 OP(13,exp)
231 OP(14,log)
232 OP(15,sign)
233 OP(16,sqrt)
234 default: ERROR(BAD_CODE);
235 }
236}
237
162 238
163int mapCAux(int code, KGCVEC(x), GCVEC(r)) { 239int mapCAux(int code, KGCVEC(x), GCVEC(r)) {
164 int k; 240 int k;
@@ -210,6 +286,22 @@ int mapValR(int code, double* pval, KRVEC(x), RVEC(r)) {
210 } 286 }
211} 287}
212 288
289int mapValF(int code, float* pval, KFVEC(x), FVEC(r)) {
290 int k;
291 float val = *pval;
292 REQUIRES(xn == rn,BAD_SIZE);
293 DEBUGMSG("mapValF");
294 switch (code) {
295 OPV(0,val*xp[k])
296 OPV(1,val/xp[k])
297 OPV(2,val+xp[k])
298 OPV(3,val-xp[k])
299 OPV(4,pow(val,xp[k]))
300 OPV(5,pow(xp[k],val))
301 default: ERROR(BAD_CODE);
302 }
303}
304
213int mapValCAux(int code, gsl_complex* pval, KGCVEC(x), GCVEC(r)) { 305int mapValCAux(int code, gsl_complex* pval, KGCVEC(x), GCVEC(r)) {
214 int k; 306 int k;
215 gsl_complex val = *pval; 307 gsl_complex val = *pval;
@@ -255,6 +347,28 @@ int zipR(int code, KRVEC(a), KRVEC(b), RVEC(r)) {
255} 347}
256 348
257 349
350int zipF(int code, KFVEC(a), KFVEC(b), FVEC(r)) {
351 REQUIRES(an == bn && an == rn, BAD_SIZE);
352 int k;
353 switch(code) {
354 OPZE(4,"zipF Pow",pow)
355 OPZE(5,"zipF ATan2",atan2)
356 }
357 KFVVIEW(a);
358 KFVVIEW(b);
359 FVVIEW(r);
360 gsl_vector_float_memcpy(V(r),V(a));
361 int res;
362 switch(code) {
363 OPZV(0,"zipF Add",gsl_vector_float_add)
364 OPZV(1,"zipF Sub",gsl_vector_float_sub)
365 OPZV(2,"zipF Mul",gsl_vector_float_mul)
366 OPZV(3,"zipF Div",gsl_vector_float_div)
367 default: ERROR(BAD_CODE);
368 }
369}
370
371
258int zipCAux(int code, KGCVEC(a), KGCVEC(b), GCVEC(r)) { 372int zipCAux(int code, KGCVEC(a), KGCVEC(b), GCVEC(r)) {
259 REQUIRES(an == bn && an == rn, BAD_SIZE); 373 REQUIRES(an == bn && an == rn, BAD_SIZE);
260 int k; 374 int k;