[pure-lang-svn] SF.net SVN: pure-lang:[844] pure/trunk
Status: Beta
Brought to you by:
agraef
|
From: <ag...@us...> - 2008-09-23 23:25:32
|
Revision: 844
http://pure-lang.svn.sourceforge.net/pure-lang/?rev=844&view=rev
Author: agraef
Date: 2008-09-23 23:25:21 +0000 (Tue, 23 Sep 2008)
Log Message:
-----------
Overhaul of matrix-pointer operations.
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-23 22:03:18 UTC (rev 843)
+++ pure/trunk/lib/matrices.pure 2008-09-23 23:25:21 UTC (rev 844)
@@ -43,17 +43,36 @@
rowvectorp x = matrixp x && dim x!0==1;
colvectorp x = matrixp x && dim x!1==1;
-/* Convenience functions to create numeric matrices with the given dimensions
+/* Matrix comparisons. */
+
+x::matrix == y::matrix = x === y
+ if nmatrixp x && matrix_type x == matrix_type y;
+// mixed numeric cases
+ = cmatrix x === y if nmatrixp x && cmatrixp y;
+ = x === cmatrix y if cmatrixp x && nmatrixp y;
+ = dmatrix x === y if imatrixp x && dmatrixp y;
+ = x === dmatrix y if dmatrixp x && imatrixp y;
+// comparisons with symbolic matrices
+ = 0 if dim x != dim y;
+ = compare 0 with
+ compare i::int = 1 if i>=n;
+ = 0 if x!i != y!i;
+ = compare (i+1);
+ end when n::int = #x end;
+
+x::matrix != y::matrix = not x == y;
+
+/* Convenience functions to create zero matrices with the given dimensions
(either a pair denoting the number of rows and columns, or just the row
- size in order to create a row vector), and all elements zero. */
+ size in order to create a row vector). */
-dmatrix (n::int,m::int) = cdmatrix (pointer 0) (n,m);
-cmatrix (n::int,m::int) = ccmatrix (pointer 0) (n,m);
-imatrix (n::int,m::int) = cimatrix (pointer 0) (n,m);
+dmatrix (n::int,m::int) = dmatrix_from_array_dup (pointer 0) (n,m);
+cmatrix (n::int,m::int) = cmatrix_from_array_dup (pointer 0) (n,m);
+imatrix (n::int,m::int) = imatrix_from_array_dup (pointer 0) (n,m);
-dmatrix n::int = cdmatrix (pointer 0) n;
-cmatrix n::int = ccmatrix (pointer 0) n;
-imatrix n::int = cimatrix (pointer 0) n;
+dmatrix n::int = dmatrix (1,n);
+cmatrix n::int = cmatrix (1,n);
+imatrix n::int = imatrix (1,n);
/* Size of a matrix (#x) and its dimensions (dim x=n,m where n is the number
of rows, m the number of columns). Note that Pure supports empty matrices,
@@ -328,34 +347,34 @@
im x::matrix = matrix_im x if nmatrixp x;
/* 'Pack' a matrix. This creates a copy of the matrix which has the data in
- contiguous storage. If possible, it also frees up extra memory if the
- matrix was created as a slice from a bigger matrix. The 'packed' predicate
- can be used to verify whether a matrix is already packed. */
+ contiguous storage. It also frees up extra memory if the matrix was created
+ as a slice from a bigger matrix. The 'packed' predicate can be used to
+ verify whether a matrix is already packed. */
pack x::matrix = colcat [x,{}];
packed x::matrix = stride x==dim x!1;
-/* Convert a matrix to a pointer. This returns a pointer to the underlying C
- array, which allows you to change the data in-place by shelling out to C
- functions, so beware. You also need to be careful when passing such a
- pointer to C functions if the underlying data is non-contiguous; when in
- doubt use the 'pack' function from above to place the data in contiguous
- storage first. */
+/* Low-level operations for converting between matrices and raw pointers.
+ These are typically used to shovel around massive amounts of numeric data
+ between Pure and external C routines, when performance and throughput is an
+ important consideration (e.g., graphics, video and audio applications). The
+ usual caveats apply. */
+/* Convert a matrix to a pointer. Returns a pointer pointing 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. */
+
private pure_pointerval;
extern expr* pure_pointerval(expr*);
pointer x::matrix = pure_pointerval x;
-/* Create a numeric matrix from a pointer. The pointer to be converted must
- point to a malloc'ed, properly initialized memory area big enough to
- accommodate the requested dimensions. The pointer may also be NULL in which
- case a suitable pointer is allocated automatically. This memory is taken
- over by Pure and will be freed automatically when the matrix object is
- garbage-collected. CAVEAT: These are low-level operations and hence should
- be used with care. They are useful, in particular, if you need to handle
- massive amounts of matrix or vector data generated or processed by external
- C functions, such as routines dealing with graphics and sound data. */
+/* 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. */
private matrix_from_double_array;
private matrix_from_complex_array;
@@ -364,35 +383,41 @@
extern expr* matrix_from_complex_array(void* p, int n, int m);
extern expr* matrix_from_int_array(void* p, int n, int m);
-cdmatrix p::pointer (n::int,m::int)
+dmatrix_from_array p::pointer (n::int,m::int)
= matrix_from_double_array p n m;
-ccmatrix p::pointer (n::int,m::int)
+cmatrix_from_array p::pointer (n::int,m::int)
= matrix_from_complex_array p n m;
-cimatrix p::pointer (n::int,m::int)
+imatrix_from_array p::pointer (n::int,m::int)
= matrix_from_int_array p n m;
-cdmatrix p::pointer n::int
- = cdmatrix p (1,n);
-ccmatrix p::pointer n::int
- = ccmatrix p (1,n);
-cimatrix p::pointer n::int
- = cimatrix p (1,n);
+dmatrix_from_array p::pointer n::int
+ = dmatrix_from_array p (1,n);
+cmatrix_from_array p::pointer n::int
+ = cmatrix_from_array p (1,n);
+imatrix_from_array p::pointer n::int
+ = imatrix_from_array p (1,n);
-/* Matrix comparisons. */
+/* 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. */
-x::matrix == y::matrix = x === y
- if nmatrixp x && matrix_type x == matrix_type y;
-// mixed numeric cases
- = cmatrix x === y if nmatrixp x && cmatrixp y;
- = x === cmatrix y if cmatrixp x && nmatrixp y;
- = dmatrix x === y if imatrixp x && dmatrixp y;
- = x === dmatrix y if dmatrixp x && imatrixp y;
-// comparisons with symbolic matrices
- = 0 if dim x != dim y;
- = compare 0 with
- compare i::int = 1 if i>=n;
- = 0 if x!i != y!i;
- = compare (i+1);
- end when n::int = #x end;
+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);
-x::matrix != y::matrix = not x == y;
+dmatrix_from_array_dup p::pointer (n::int,m::int)
+ = matrix_from_double_array_dup p n m;
+cmatrix_from_array_dup p::pointer (n::int,m::int)
+ = matrix_from_complex_array_dup p n m;
+imatrix_from_array_dup p::pointer (n::int,m::int)
+ = matrix_from_int_array_dup p n m;
+
+dmatrix_from_array_dup p::pointer n::int
+ = dmatrix_from_array_dup p (1,n);
+cmatrix_from_array_dup p::pointer n::int
+ = cmatrix_from_array_dup p (1,n);
+imatrix_from_array_dup p::pointer n::int
+ = imatrix_from_array_dup p (1,n);
Modified: pure/trunk/runtime.cc
===================================================================
--- pure/trunk/runtime.cc 2008-09-23 22:03:18 UTC (rev 843)
+++ pure/trunk/runtime.cc 2008-09-23 23:25:21 UTC (rev 844)
@@ -388,19 +388,19 @@
#ifdef HAVE_GSL
case EXPR::DMATRIX: {
gsl_matrix *m = (gsl_matrix*)x->data.mat.p;
- m->owner = owner;
+ m->owner = owner && m->block;
gsl_matrix_free(m);
break;
}
case EXPR::CMATRIX: {
gsl_matrix_complex *m = (gsl_matrix_complex*)x->data.mat.p;
- m->owner = owner;
+ m->owner = owner && m->block;
gsl_matrix_complex_free(m);
break;
}
case EXPR::IMATRIX: {
gsl_matrix_int *m = (gsl_matrix_int*)x->data.mat.p;
- m->owner = owner;
+ m->owner = owner && m->block;
gsl_matrix_int_free(m);
break;
}
@@ -5016,14 +5016,94 @@
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
return pure_double_matrix(create_double_matrix(n1, n2));
- if (!p) p = calloc(n1*n2, sizeof(double));
if (!p) return 0;
gsl_matrix_view v = gsl_matrix_view_array((double*)p, n1, n2);
// take a copy of the view matrix
gsl_matrix *m = (gsl_matrix*)malloc(sizeof(gsl_matrix));
- gsl_block *b = (gsl_block*)malloc(sizeof(gsl_block));
assert(m && v.matrix.data);
*m = v.matrix;
+ pure_expr *x = new_expr();
+ x->tag = EXPR::DMATRIX;
+ x->data.mat.p = m;
+ x->data.mat.refc = new uint32_t;
+ *x->data.mat.refc = 1;
+ MEMDEBUG_NEW(x)
+ return x;
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+pure_expr *matrix_from_complex_array(void *p, uint32_t n1, uint32_t n2)
+{
+#ifdef HAVE_GSL
+ if (n1 == 0 || n2 == 0) // empty matrix
+ return pure_complex_matrix(create_complex_matrix(n1, n2));
+ if (!p) return 0;
+ gsl_matrix_complex_view v =
+ gsl_matrix_complex_view_array((double*)p, n1, n2);
+ // take a copy of the view matrix
+ gsl_matrix_complex *m =
+ (gsl_matrix_complex*)malloc(sizeof(gsl_matrix_complex));
+ assert(m && v.matrix.data);
+ *m = v.matrix;
+ pure_expr *x = new_expr();
+ x->tag = EXPR::CMATRIX;
+ x->data.mat.p = m;
+ x->data.mat.refc = new uint32_t;
+ *x->data.mat.refc = 1;
+ MEMDEBUG_NEW(x)
+ return x;
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+pure_expr *matrix_from_int_array(void *p, uint32_t n1, uint32_t n2)
+{
+#ifdef HAVE_GSL
+ if (n1 == 0 || n2 == 0) // empty matrix
+ return pure_int_matrix(create_int_matrix(n1, n2));
+ if (!p) return 0;
+ gsl_matrix_int_view v = gsl_matrix_int_view_array((int*)p, n1, n2);
+ // take a copy of the view matrix
+ gsl_matrix_int *m = (gsl_matrix_int*)malloc(sizeof(gsl_matrix_int));
+ assert(m && v.matrix.data);
+ *m = v.matrix;
+ pure_expr *x = new_expr();
+ x->tag = EXPR::IMATRIX;
+ x->data.mat.p = m;
+ x->data.mat.refc = new uint32_t;
+ *x->data.mat.refc = 1;
+ MEMDEBUG_NEW(x)
+ return x;
+#else
+ return 0;
+#endif
+}
+
+extern "C"
+pure_expr *matrix_from_double_array_dup(void *p, uint32_t n1, uint32_t n2)
+{
+#ifdef HAVE_GSL
+ if (n1 == 0 || n2 == 0) // empty matrix
+ return pure_double_matrix(create_double_matrix(n1, n2));
+ if (!p)
+ p = calloc(n1*n2, sizeof(double));
+ else {
+ void *q = malloc(n1*n2*sizeof(double));
+ memcpy(q, p, n1*n2*sizeof(double));
+ p = q;
+ }
+ if (!p) return 0;
+ gsl_matrix_view v = gsl_matrix_view_array((double*)p, n1, n2);
+ // take a copy of the view matrix
+ gsl_matrix *m = (gsl_matrix*)malloc(sizeof(gsl_matrix));
+ gsl_block *b = (gsl_block*)malloc(sizeof(gsl_block));
+ assert(m && b && v.matrix.data);
+ *m = v.matrix;
b->size = n1*n2;
b->data = m->data;
m->block = b;
@@ -5040,12 +5120,18 @@
}
extern "C"
-pure_expr *matrix_from_complex_array(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_complex_array_dup(void *p, uint32_t n1, uint32_t n2)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
return pure_complex_matrix(create_complex_matrix(n1, n2));
- if (!p) p = calloc(2*n1*n2, sizeof(double));
+ if (!p)
+ p = calloc(2*n1*n2, sizeof(double));
+ else {
+ void *q = malloc(2*n1*n2*sizeof(double));
+ memcpy(q, p, 2*n1*n2*sizeof(double));
+ p = q;
+ }
if (!p) return 0;
gsl_matrix_complex_view v =
gsl_matrix_complex_view_array((double*)p, n1, n2);
@@ -5053,7 +5139,7 @@
gsl_matrix_complex *m =
(gsl_matrix_complex*)malloc(sizeof(gsl_matrix_complex));
gsl_block_complex *b = (gsl_block_complex*)malloc(sizeof(gsl_block_complex));
- assert(m && v.matrix.data);
+ assert(m && b && v.matrix.data);
*m = v.matrix;
b->size = n1*n2;
b->data = m->data;
@@ -5071,18 +5157,24 @@
}
extern "C"
-pure_expr *matrix_from_int_array(void *p, uint32_t n1, uint32_t n2)
+pure_expr *matrix_from_int_array_dup(void *p, uint32_t n1, uint32_t n2)
{
#ifdef HAVE_GSL
if (n1 == 0 || n2 == 0) // empty matrix
return pure_int_matrix(create_int_matrix(n1, n2));
- if (!p) p = calloc(n1*n2, sizeof(int));
+ if (!p)
+ p = calloc(n1*n2, sizeof(int));
+ else {
+ void *q = malloc(n1*n2*sizeof(int));
+ memcpy(q, p, n1*n2*sizeof(int));
+ p = q;
+ }
if (!p) return 0;
gsl_matrix_int_view v = gsl_matrix_int_view_array((int*)p, n1, n2);
// take a copy of the view matrix
gsl_matrix_int *m = (gsl_matrix_int*)malloc(sizeof(gsl_matrix_int));
gsl_block_int *b = (gsl_block_int*)malloc(sizeof(gsl_block_int));
- assert(m && v.matrix.data);
+ assert(m && b && v.matrix.data);
*m = v.matrix;
b->size = n1*n2;
b->data = m->data;
Modified: pure/trunk/runtime.h
===================================================================
--- pure/trunk/runtime.h 2008-09-23 22:03:18 UTC (rev 843)
+++ pure/trunk/runtime.h 2008-09-23 23:25:21 UTC (rev 844)
@@ -734,17 +734,25 @@
pure_expr *matrix_im(pure_expr *x);
/* Create a matrix object of the given dimensions which uses the given pointer
- p as its underlying C array. There are no checks whatsoever, so the caller
- is responsible for making sure that the memory has the right size and is
- properly initialized. p must point to a malloc'ed memory area, which is
- taken over by Pure and will be freed automatically when the matrix is
- garbage-collected. p may also be NULL in which case a suitable pointer is
- allocated automatically. */
+ p as its underlying storage. There are no checks whatsoever and the data is
+ *not* copied, so the caller is responsible for making sure that the memory
+ has the right size, is properly initialized and is not freed while the
+ 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);
+/* 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.
+ 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);
+
/* 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.
|