summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--lib/Data/Packed/Internal/Matrix.hs33
-rw-r--r--lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c22
-rw-r--r--lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h4
3 files changed, 57 insertions, 2 deletions
diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs
index 7a17ef0..090826d 100644
--- a/lib/Data/Packed/Internal/Matrix.hs
+++ b/lib/Data/Packed/Internal/Matrix.hs
@@ -252,10 +252,11 @@ class (Storable a, Floating a) => Element a where
252 -> Matrix a -> Matrix a 252 -> Matrix a -> Matrix a
253 subMatrixD = subMatrix' 253 subMatrixD = subMatrix'
254 transdata :: Int -> Vector a -> Int -> Vector a 254 transdata :: Int -> Vector a -> Int -> Vector a
255 transdata = transdata' 255 transdata = transdataP -- transdata'
256 constantD :: a -> Int -> Vector a 256 constantD :: a -> Int -> Vector a
257 constantD = constant' 257 constantD = constantP -- constant'
258 conjugateD :: Vector a -> Vector a 258 conjugateD :: Vector a -> Vector a
259 conjugateD = id
259 260
260instance Element Float where 261instance Element Float where
261 transdata = transdataAux ctransF 262 transdata = transdataAux ctransF
@@ -320,10 +321,27 @@ transdataAux fun c1 d c2 =
320 r2 = dim d `div` c2 321 r2 = dim d `div` c2
321 noneed = r1 == 1 || c1 == 1 322 noneed = r1 == 1 || c1 == 1
322 323
324transdataP :: Storable a => Int -> Vector a -> Int -> Vector a
325transdataP c1 d c2 =
326 if noneed
327 then d
328 else unsafePerformIO $ do
329 v <- createVector (dim d)
330 unsafeWith d $ \pd ->
331 unsafeWith v $ \pv ->
332 ctransP (fi r1) (fi c1) (castPtr pd) (fi sz) (fi r2) (fi c2) (castPtr pv) (fi sz) // check "transdataStorable"
333 return v
334 where r1 = dim d `div` c1
335 r2 = dim d `div` c2
336 sz = sizeOf (d @> 0)
337 noneed = r1 == 1 || c1 == 1
338
323foreign import ccall "transF" ctransF :: TFMFM 339foreign import ccall "transF" ctransF :: TFMFM
324foreign import ccall "transR" ctransR :: TMM 340foreign import ccall "transR" ctransR :: TMM
325foreign import ccall "transQ" ctransQ :: TQMQM 341foreign import ccall "transQ" ctransQ :: TQMQM
326foreign import ccall "transC" ctransC :: TCMCM 342foreign import ccall "transC" ctransC :: TCMCM
343foreign import ccall "transP" ctransP :: CInt -> CInt -> Ptr () -> CInt -> CInt -> CInt -> Ptr () -> CInt -> IO CInt
344
327---------------------------------------------------------------------- 345----------------------------------------------------------------------
328 346
329constant' v n = unsafePerformIO $ do 347constant' v n = unsafePerformIO $ do
@@ -359,6 +377,17 @@ constantC :: Complex Double -> Int -> Vector (Complex Double)
359constantC = constantAux cconstantC 377constantC = constantAux cconstantC
360foreign import ccall "constantC" cconstantC :: Ptr (Complex Double) -> TCV 378foreign import ccall "constantC" cconstantC :: Ptr (Complex Double) -> TCV
361 379
380constantP :: Storable a => a -> Int -> Vector a
381constantP a n = unsafePerformIO $ do
382 let sz = sizeOf a
383 v <- createVector n
384 unsafeWith v $ \p -> do
385 alloca $ \k -> do
386 poke k a
387 cconstantP (castPtr k) (fi n) (castPtr p) (fi sz) // check "constantP"
388 return v
389foreign import ccall "constantP" cconstantP :: Ptr () -> CInt -> Ptr () -> CInt -> IO CInt
390
362--------------------------------------- 391---------------------------------------
363 392
364conjugateAux fun x = unsafePerformIO $ do 393conjugateAux fun x = unsafePerformIO $ do
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
index 09979cd..e8bbbdb 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
+++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c
@@ -1139,6 +1139,19 @@ int transC(KCMAT(x),CMAT(t)) {
1139 OK 1139 OK
1140} 1140}
1141 1141
1142int transP(KPMAT(x), PMAT(t)) {
1143 REQUIRES(xr==tc && xc==tr,BAD_SIZE);
1144 REQUIRES(xs==ts,NOCONVER);
1145 DEBUGMSG("transP");
1146 int i,j;
1147 for (i=0; i<tr; i++) {
1148 for (j=0; j<tc; j++) {
1149 memcpy(tp+(i*tc+j)*xs,xp +(j*xc+i)*xs,xs);
1150 }
1151 }
1152 OK
1153}
1154
1142//////////////////// constant ///////////////////////// 1155//////////////////// constant /////////////////////////
1143 1156
1144int constantF(float * pval, FVEC(r)) { 1157int constantF(float * pval, FVEC(r)) {
@@ -1181,6 +1194,15 @@ int constantC(doublecomplex* pval, CVEC(r)) {
1181 OK 1194 OK
1182} 1195}
1183 1196
1197int constantP(void* pval, PVEC(r)) {
1198 DEBUGMSG("constantP")
1199 int k;
1200 for(k=0;k<rn;k++) {
1201 memcpy(rp+k*rs,pval,rs);
1202 }
1203 OK
1204}
1205
1184//////////////////// float-double conversion ///////////////////////// 1206//////////////////// float-double conversion /////////////////////////
1185 1207
1186int float2double(FVEC(x),DVEC(y)) { 1208int float2double(FVEC(x),DVEC(y)) {
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
index 3de5b19..12dcde4 100644
--- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
+++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.h
@@ -44,19 +44,23 @@ typedef short ftnlen;
44#define DVEC(A) int A##n, double*A##p 44#define DVEC(A) int A##n, double*A##p
45#define QVEC(A) int A##n, complex*A##p 45#define QVEC(A) int A##n, complex*A##p
46#define CVEC(A) int A##n, doublecomplex*A##p 46#define CVEC(A) int A##n, doublecomplex*A##p
47#define PVEC(A) int A##n, void* A##p, int A##s
47#define FMAT(A) int A##r, int A##c, float* A##p 48#define FMAT(A) int A##r, int A##c, float* A##p
48#define DMAT(A) int A##r, int A##c, double* A##p 49#define DMAT(A) int A##r, int A##c, double* A##p
49#define QMAT(A) int A##r, int A##c, complex* A##p 50#define QMAT(A) int A##r, int A##c, complex* A##p
50#define CMAT(A) int A##r, int A##c, doublecomplex* A##p 51#define CMAT(A) int A##r, int A##c, doublecomplex* A##p
52#define PMAT(A) int A##r, int A##c, void* A##p, int A##s
51 53
52#define KFVEC(A) int A##n, const float*A##p 54#define KFVEC(A) int A##n, const float*A##p
53#define KDVEC(A) int A##n, const double*A##p 55#define KDVEC(A) int A##n, const double*A##p
54#define KQVEC(A) int A##n, const complex*A##p 56#define KQVEC(A) int A##n, const complex*A##p
55#define KCVEC(A) int A##n, const doublecomplex*A##p 57#define KCVEC(A) int A##n, const doublecomplex*A##p
58#define KPVEC(A) int A##n, const void* A##p, int A##s
56#define KFMAT(A) int A##r, int A##c, const float* A##p 59#define KFMAT(A) int A##r, int A##c, const float* A##p
57#define KDMAT(A) int A##r, int A##c, const double* A##p 60#define KDMAT(A) int A##r, int A##c, const double* A##p
58#define KQMAT(A) int A##r, int A##c, const complex* A##p 61#define KQMAT(A) int A##r, int A##c, const complex* A##p
59#define KCMAT(A) int A##r, int A##c, const doublecomplex* A##p 62#define KCMAT(A) int A##r, int A##c, const doublecomplex* A##p
63#define KPMAT(A) int A##r, int A##c, const void* A##p, int A##s
60 64
61/********************************************************/ 65/********************************************************/
62 66