[pure-lang-svn] SF.net SVN: pure-lang:[884] pure/trunk
Status: Beta
Brought to you by:
agraef
Menu
▾
▴
[pure-lang-svn] SF.net SVN: pure-lang:[884] pure/trunk
|
From: <ag...@us...> - 2008-09-27 16:03:53
|
Revision: 884
http://pure-lang.svn.sourceforge.net/pure-lang/?rev=884&view=rev
Author: agraef
Date: 2008-09-27 16:03:46 +0000 (Sat, 27 Sep 2008)
Log Message:
-----------
Added missing complex->double/int matrix conversions, overhaul of matrix<->pointer conversions.
Modified Paths:
--------------
pure/trunk/lib/matrices.pure
pure/trunk/runtime.cc
pure/trunk/runtime.h
Modified: pure/trunk/lib/matrices.pure
===================================================================
--- pure/trunk/lib/matrices.pure 2008-09-27 13:11:09 UTC (rev 883)
+++ pure/trunk/lib/matrices.pure 2008-09-27 16:03:46 UTC (rev 884)
@@ -217,9 +217,9 @@
(either a pair denoting the number of rows and columns, or just the row
size in order to create a row vector). */
-dmatrix (n::int,m::int) = dmatrix_dup (pointer 0,n,m);
-cmatrix (n::int,m::int) = cmatrix_dup (pointer 0,n,m);
-imatrix (n::int,m::int) = imatrix_dup (pointer 0,n,m);
+dmatrix (n::int,m::int) = double_matrix (n,m) (pointer 0);
+cmatrix (n::int,m::int) = complex_matrix (n,m) (pointer 0);
+imatrix (n::int,m::int) = int_matrix (n,m) (pointer 0);
dmatrix n::int = dmatrix (1,n);
cmatrix n::int = cmatrix (1,n);
@@ -233,8 +233,8 @@
extern expr* matrix_double(expr *x), expr* matrix_complex(expr *x),
expr* matrix_int(expr *x);
-dmatrix x::matrix = matrix_double x if imatrixp x || dmatrixp x;
-imatrix x::matrix = matrix_int x if imatrixp x || dmatrixp x;
+dmatrix x::matrix = matrix_double x if nmatrixp x;
+imatrix x::matrix = matrix_int x if nmatrixp x;
cmatrix x::matrix = matrix_complex x if nmatrixp x;
private matrix_re matrix_im;
@@ -374,94 +374,126 @@
important consideration (e.g., graphics, video and audio applications). The
usual caveats apply. */
-/* Convert a matrix to a pointer, which points directly to the underlying C
- array. You have to be careful when passing such a pointer to C functions if
- the underlying data is non-contiguous; when in doubt, first use the 'pack'
- function from above to place the data in contiguous storage. */
+/* Get a pointer to the underlying C array of a matrix. The data is *not*
+ copied. Hence you have to be careful when passing such a pointer to C
+ functions if the underlying data is non-contiguous; when in doubt, first
+ use the 'pack' function from above to place the data in contiguous storage,
+ or use one of the matrix-pointer conversion routines below. */
private pure_pointerval;
extern expr* pure_pointerval(expr*);
pointer x::matrix = pure_pointerval x;
-/* Create a numeric matrix from a pointer, without copying the data. The
- caller must ensure that the pointer points to properly initialized memory
- big enough to accommodate the requested dimensions, which persists for the
- entire lifetime of the matrix object. */
+/* Copy the contents of a matrix to a given pointer and return that pointer,
+ converting to the target data type on the fly if necessary. The given
+ pointer may also be NULL, in which case suitable memory is malloc'ed and
+ returned; otherwise the caller must ensure that the memory pointed to by p
+ is big enough for the contents of the given matrix. */
+private matrix_to_double_array;
+private matrix_to_float_array;
+private matrix_to_complex_array;
+private matrix_to_complex_float_array;
+private matrix_to_int_array;
+private matrix_to_short_array;
+private matrix_to_byte_array;
+extern void* matrix_to_double_array(void* p, expr* x);
+extern void* matrix_to_float_array(void* p, expr* x);
+extern void* matrix_to_complex_array(void* p, expr* x);
+extern void* matrix_to_complex_float_array(void* p, expr* x);
+extern void* matrix_to_int_array(void* p, expr* x);
+extern void* matrix_to_short_array(void* p, expr* x);
+extern void* matrix_to_byte_array(void* p, expr* x);
+
+double_pointer p::pointer x::matrix
+ = matrix_to_double_array p x if nmatrixp x;
+float_pointer p::pointer x::matrix
+ = matrix_to_float_array p x if nmatrixp x;
+complex_pointer p::pointer x::matrix
+ = matrix_to_complex_array p x if nmatrixp x;
+complex_float_pointer p::pointer x::matrix
+ = matrix_to_complex_float_array p x if nmatrixp x;
+int_pointer p::pointer x::matrix
+ = matrix_to_int_array p x if nmatrixp x;
+short_pointer p::pointer x::matrix
+ = matrix_to_short_array p x if nmatrixp x;
+byte_pointer p::pointer x::matrix
+ = matrix_to_byte_array p x if nmatrixp x;
+
+/* Create a numeric matrix from a pointer, copying the data and converting it
+ from the source type on the fly if necessary. The source pointer p may also
+ be NULL, in which case the new matrix is filled with zeros instead.
+ Otherwise the caller must ensure that the pointer points to properly
+ initialized memory big enough for the requested dimensions. */
+
private matrix_from_double_array;
+private matrix_from_float_array;
private matrix_from_complex_array;
+private matrix_from_complex_float_array;
private matrix_from_int_array;
-extern expr* matrix_from_double_array(void* p, int n, int m);
-extern expr* matrix_from_complex_array(void* p, int n, int m);
-extern expr* matrix_from_int_array(void* p, int n, int m);
+private matrix_from_short_array;
+private matrix_from_byte_array;
+extern expr* matrix_from_double_array(int n, int m, void* p);
+extern expr* matrix_from_float_array(int n, int m, void* p);
+extern expr* matrix_from_complex_array(int n, int m, void* p);
+extern expr* matrix_from_complex_float_array(int n, int m, void* p);
+extern expr* matrix_from_int_array(int n, int m, void* p);
+extern expr* matrix_from_short_array(int n, int m, void* p);
+extern expr* matrix_from_byte_array(int n, int m, void* p);
-dmatrix (p::pointer,n::int,m::int)
- = matrix_from_double_array p n m;
-cmatrix (p::pointer,n::int,m::int)
- = matrix_from_complex_array p n m;
-imatrix (p::pointer,n::int,m::int)
- = matrix_from_int_array p n m;
+double_matrix (n::int,m::int) p::pointer
+ = matrix_from_double_array n m p;
+float_matrix (n::int,m::int) p::pointer
+ = matrix_from_float_array n m p;
+complex_matrix (n::int,m::int) p::pointer
+ = matrix_from_complex_array n m p;
+complex_float_matrix (n::int,m::int) p::pointer
+ = matrix_from_complex_float_array n m p;
+int_matrix (n::int,m::int) p::pointer
+ = matrix_from_int_array n m p;
+short_matrix (n::int,m::int) p::pointer
+ = matrix_from_short_array n m p;
+byte_matrix (n::int,m::int) p::pointer
+ = matrix_from_byte_array n m p;
-dmatrix (p::pointer,n::int)
- = dmatrix (p,1,n);
-cmatrix (p::pointer,n::int)
- = cmatrix (p,1,n);
-imatrix (p::pointer,n::int)
- = imatrix (p,1,n);
+double_matrix n::int p::pointer
+ = double_matrix (1,n) p;
+float_matrix n::int p::pointer
+ = float_matrix (1,n) p;
+complex_matrix n::int p::pointer
+ = complex_matrix (1,n) p;
+complex_float_matrix n::int p::pointer
+ = complex_float_matrix (1,n) p;
+int_matrix n::int p::pointer
+ = int_matrix (1,n) p;
+short_matrix n::int p::pointer
+ = short_matrix (1,n) p;
+byte_matrix n::int p::pointer
+ = byte_matrix (1,n) p;
-/* Create a numeric matrix from a pointer, copying the data to fresh memory.
- The source pointer p may also be NULL, in which case the new matrix is
- filled with zeros instead. */
+/* Create a numeric matrix view of existing data, without copying the data.
+ The data must be double, complex or int, the pointer must not be NULL and
+ the caller must also ensure that the memory persists for the entire
+ lifetime of the matrix object. */
-private matrix_from_double_array_dup;
-private matrix_from_complex_array_dup;
-private matrix_from_int_array_dup;
-extern expr* matrix_from_double_array_dup(void* p, int n, int m);
-extern expr* matrix_from_complex_array_dup(void* p, int n, int m);
-extern expr* matrix_from_int_array_dup(void* p, int n, int m);
+private matrix_from_double_array_nodup;
+private matrix_from_complex_array_nodup;
+private matrix_from_int_array_nodup;
+extern expr* matrix_from_double_array_nodup(int n, int m, void* p);
+extern expr* matrix_from_complex_array_nodup(int n, int m, void* p);
+extern expr* matrix_from_int_array_nodup(int n, int m, void* p);
-dmatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_double_array_dup p n m;
-cmatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_complex_array_dup p n m;
-imatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_int_array_dup p n m;
+double_matrix_view (n::int,m::int) p::pointer
+ = matrix_from_double_array_nodup n m p;
+complex_matrix_view (n::int,m::int) p::pointer
+ = matrix_from_complex_array_nodup n m p;
+int_matrix_view (n::int,m::int) p::pointer
+ = matrix_from_int_array_nodup n m p;
-dmatrix_dup (p::pointer,n::int)
- = dmatrix_dup (p,1,n);
-cmatrix_dup (p::pointer,n::int)
- = cmatrix_dup (p,1,n);
-imatrix_dup (p::pointer,n::int)
- = imatrix_dup (p,1,n);
-
-/* Some additional functions for alternate base types. These work like the
- routines above, but initialize the data from float, complex float, short
- and byte arrays, respectively. */
-
-private matrix_from_float_array_dup;
-private matrix_from_complex_float_array_dup;
-private matrix_from_short_array_dup;
-private matrix_from_byte_array_dup;
-extern expr* matrix_from_float_array_dup(void* p, int n, int m);
-extern expr* matrix_from_complex_float_array_dup(void* p, int n, int m);
-extern expr* matrix_from_short_array_dup(void* p, int n, int m);
-extern expr* matrix_from_byte_array_dup(void* p, int n, int m);
-
-float_dmatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_float_array_dup p n m;
-float_cmatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_complex_float_array_dup p n m;
-short_imatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_short_array_dup p n m;
-byte_imatrix_dup (p::pointer,n::int,m::int)
- = matrix_from_byte_array_dup p n m;
-
-float_dmatrix_dup (p::pointer,n::int)
- = float_dmatrix_dup (p,1,n);
-float_cmatrix_dup (p::pointer,n::int)
- = float_cmatrix_dup (p,1,n);
-short_imatrix_dup (p::pointer,n::int)
- = short_imatrix_dup (p,1,n);
-byte_imatrix_dup (p::pointer,n::int)
- = byte_imatrix_dup (p,1,n);
+double_matrix_view n::int p::pointer
+ = double_matrix_view (1,n) p;
+complex_matrix_view n::int p::pointer
+ = complex_matrix_view (1,n) p;
+int_matrix_view n::int p::pointer
+ = int_matrix_view (1,n) p;
Modified: pure/trunk/runtime.cc
===================================================================
--- pure/trunk/runtime.cc 2008-09-27 13:11:09 UTC (rev 883)
+++ pure/trunk/runtime.cc 2008-09-27 16:03:46 UTC (rev 884)
@@ -4882,6 +4882,19 @@
m2->data[i*m2->tda+j] = (double)m1->data[i*m1->tda+j];
return pure_double_matrix(m2);
}
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ gsl_matrix *m2 = create_double_matrix(n, 2*m);
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t k = i*m2->tda+2*j;
+ size_t l = 2*(i*m1->tda+j);
+ m2->data[k] = m1->data[l];
+ m2->data[k+1] = m1->data[l+1];
+ }
+ return pure_double_matrix(m2);
+ }
default:
return 0;
}
@@ -4945,6 +4958,19 @@
}
case EXPR::IMATRIX:
return x;
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ gsl_matrix_int *m2 = create_int_matrix(n, 2*m);
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t k = i*m2->tda+2*j;
+ size_t l = 2*(i*m1->tda+j);
+ m2->data[k] = (int)m1->data[l];
+ m2->data[k+1] = (int)m1->data[l+1];
+ }
+ return pure_int_matrix(m2);
+ }
default:
return 0;
}
@@ -5019,7 +5045,7 @@
}
extern "C"
-pure_expr *matrix_from_double_array(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_double_array_nodup(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5043,7 +5069,7 @@
}
extern "C"
-pure_expr *matrix_from_complex_array(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_complex_array_nodup(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5069,7 +5095,7 @@
}
extern "C"
-pure_expr *matrix_from_int_array(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_int_array_nodup(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5093,7 +5119,7 @@
}
extern "C"
-pure_expr *matrix_from_double_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_double_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5128,7 +5154,7 @@
}
extern "C"
-pure_expr *matrix_from_complex_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_complex_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5165,7 +5191,7 @@
}
extern "C"
-pure_expr *matrix_from_int_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_int_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5200,9 +5226,178 @@
}
extern "C"
-pure_expr *matrix_from_float_array_dup(void *p, uint32_t n1, uint32_t n2)
+void *matrix_to_double_array(void *p, pure_expr *x)
{
#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = m1->data[l];
+ q[k++] = m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = m1->data[i*m1->tda+j];
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (double)m1->data[i*m1->tda+j];
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+void *matrix_to_complex_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = m1->data[l];
+ q[k++] = m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ q[k++] = m1->data[i*m1->tda+j];
+ q[k++] = 0.0;
+ }
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(double));
+ if (!p) return 0;
+ double *q = (double*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ q[k++] = (double)m1->data[i*m1->tda+j];
+ q[k++] = 0.0;
+ }
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+void *matrix_to_int_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(int));
+ if (!p) return 0;
+ int *q = (int*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = (int)m1->data[l];
+ q[k++] = (int)m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(int));
+ if (!p) return 0;
+ int *q = (int*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (int)m1->data[i*m1->tda+j];
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(int));
+ if (!p) return 0;
+ int *q = (int*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = m1->data[i*m1->tda+j];
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+pure_expr *matrix_from_float_array(uint32_t n1, uint32_t n2, void *p)
+{
+#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
return pure_double_matrix(create_double_matrix(n1, n2));
if (!p)
@@ -5236,7 +5431,7 @@
}
extern "C"
-pure_expr *matrix_from_complex_float_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_complex_float_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5274,7 +5469,7 @@
}
extern "C"
-pure_expr *matrix_from_short_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_short_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5310,7 +5505,7 @@
}
extern "C"
-pure_expr *matrix_from_byte_array_dup(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_byte_array(uint32_t n1, uint32_t n2, void *p)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
@@ -5345,6 +5540,230 @@
#endif
}
+extern "C"
+void *matrix_to_float_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = (float)m1->data[l];
+ q[k++] = (float)m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (float)m1->data[i*m1->tda+j];
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (float)m1->data[i*m1->tda+j];
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+void *matrix_to_complex_float_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = (float)m1->data[l];
+ q[k++] = (float)m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ q[k++] = (float)m1->data[i*m1->tda+j];
+ q[k++] = 0.0;
+ }
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(float));
+ if (!p) return 0;
+ float *q = (float*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ q[k++] = (float)m1->data[i*m1->tda+j];
+ q[k++] = 0.0;
+ }
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+void *matrix_to_short_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(short));
+ if (!p) return 0;
+ short *q = (short*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = (short)m1->data[l];
+ q[k++] = (short)m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(short));
+ if (!p) return 0;
+ short *q = (short*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (short)m1->data[i*m1->tda+j];
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(short));
+ if (!p) return 0;
+ short *q = (short*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (short)m1->data[i*m1->tda+j];
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+void *matrix_to_byte_array(void *p, pure_expr *x)
+{
+#ifdef HAVE_GSL
+ switch (x->tag) {
+ case EXPR::CMATRIX: {
+ gsl_matrix_complex *m1 = (gsl_matrix_complex*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(2*n*m*sizeof(int8_t));
+ if (!p) return 0;
+ int8_t *q = (int8_t*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++) {
+ size_t l = 2*(i*m1->tda+j);
+ q[k++] = (int8_t)m1->data[l];
+ q[k++] = (int8_t)m1->data[l+1];
+ }
+ return p;
+ }
+ case EXPR::DMATRIX: {
+ gsl_matrix *m1 = (gsl_matrix*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(int8_t));
+ if (!p) return 0;
+ int8_t *q = (int8_t*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (int8_t)m1->data[i*m1->tda+j];
+ return p;
+ }
+ case EXPR::IMATRIX: {
+ gsl_matrix_int *m1 = (gsl_matrix_int*)x->data.mat.p;
+ size_t n = m1->size1, m = m1->size2;
+ if (n==0 || m==0) return p;
+ if (!p) p = malloc(n*m*sizeof(int8_t));
+ if (!p) return 0;
+ int8_t *q = (int8_t*)p;
+ size_t k = 0;
+ for (size_t i = 0; i < n; i++)
+ for (size_t j = 0; j < m; j++)
+ q[k++] = (int8_t)m1->data[i*m1->tda+j];
+ return p;
+ }
+ default:
+ return 0;
+ }
+#else
+ return 0;
+#endif
+}
+
static uint32_t mpz_hash(const mpz_t z)
{
uint32_t h = 0;
Modified: pure/trunk/runtime.h
===================================================================
--- pure/trunk/runtime.h 2008-09-27 13:11:09 UTC (rev 883)
+++ pure/trunk/runtime.h 2008-09-27 16:03:46 UTC (rev 884)
@@ -715,10 +715,7 @@
pure_expr *matrix_transpose(pure_expr *x);
-/* Convert between different types of numeric matrices. Note that any numeric
- matrix can be converted to a complex matrix, but for the other conversions
- the input must be a a double or integer matrix (see matrix_re and matrix_im
- below to handle the complex->double case). */
+/* Convert between different types of numeric matrices. */
pure_expr *matrix_double(pure_expr *x);
pure_expr *matrix_complex(pure_expr *x);
@@ -740,25 +737,42 @@
matrix is still in use. The memory is *not* freed when the matrix is
garbage-collected but remains in the ownership of the caller. */
-pure_expr *matrix_from_double_array(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_complex_array(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_int_array(void *p, uint32_t n, uint32_t m);
+pure_expr *matrix_from_double_array_nodup(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_complex_array_nodup(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_int_array_nodup(uint32_t n, uint32_t m, void *p);
/* The following routines work like the above, but copy the data to newly
allocated memory, so the original data can be freed after the call.
- Additional routines are provided for various alternate data sizes.
Moreover, the source pointer p may also be NULL in which case the new
matrix is filled with zeros instead. */
-pure_expr *matrix_from_double_array_dup(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_complex_array_dup(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_int_array_dup(void *p, uint32_t n, uint32_t m);
+pure_expr *matrix_from_double_array(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_complex_array(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_int_array(uint32_t n, uint32_t m, void *p);
-pure_expr *matrix_from_float_array_dup(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_complex_float_array_dup(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_short_array_dup(void *p, uint32_t n, uint32_t m);
-pure_expr *matrix_from_byte_array_dup(void *p, uint32_t n, uint32_t m);
+/* Copy data from the given matrix to the given data pointer, which is then
+ returned. If p is NULL then memory of sufficient size is malloc'ed;
+ otherwise p must point to a memory area of sufficient size. */
+void *matrix_to_double_array(void *p, pure_expr *x);
+void *matrix_to_complex_array(void *p, pure_expr *x);
+void *matrix_to_int_array(void *p, pure_expr *x);
+
+/* Additional routines for alternative base types. These work like the
+ routines above but take data consisting of base types which are not
+ directly supported by Pure GSL matrices: float, complex float, short,
+ byte. */
+
+pure_expr *matrix_from_float_array(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_complex_float_array(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_short_array(uint32_t n, uint32_t m, void *p);
+pure_expr *matrix_from_byte_array(uint32_t n, uint32_t m, void *p);
+
+void *matrix_to_float_array(void *p, pure_expr *x);
+void *matrix_to_complex_float_array(void *p, pure_expr *x);
+void *matrix_to_short_array(void *p, pure_expr *x);
+void *matrix_to_byte_array(void *p, pure_expr *x);
+
/* Compute a 32 bit hash code of a Pure expression. This makes it possible to
use arbitary Pure values as keys in a hash table. */
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
