[pure-lang-svn] SF.net SVN: pure-lang:[782] pure/trunk
Status: Beta
Brought to you by:
agraef
From: <ag...@us...> - 2008-09-18 11:44:19
|
Revision: 782 http://pure-lang.svn.sourceforge.net/pure-lang/?rev=782&view=rev Author: agraef Date: 2008-09-18 11:44:30 +0000 (Thu, 18 Sep 2008) Log Message: ----------- Implement the remaining matrix functions in the public runtime API. This completes basic GSL matrix support, but note that not all functionality has been fully tested yet, and operations to inspect matrix values are still missing. Modified Paths: -------------- pure/trunk/runtime.cc pure/trunk/runtime.h Modified: pure/trunk/runtime.cc =================================================================== --- pure/trunk/runtime.cc 2008-09-18 11:41:11 UTC (rev 781) +++ pure/trunk/runtime.cc 2008-09-18 11:44:30 UTC (rev 782) @@ -1128,19 +1128,97 @@ pure_expr *pure_matrix_rowsl(uint32_t n, ...) { #ifdef HAVE_GSL - // XXXTODO - return 0; + va_list ap; + pure_expr **xs = (pure_expr**)alloca(n*sizeof(pure_expr*)); + va_start(ap, n); + for (size_t i = 0; i < n; i++) + xs[i] = va_arg(ap, pure_expr*); + va_end(ap); + return pure_matrix_rowsv(n, xs); #else return 0; #endif } extern "C" -pure_expr *pure_matrix_rowsv(uint32_t n, pure_expr **elems) +pure_expr *pure_matrix_rowsv(uint32_t n, pure_expr **xs) { #ifdef HAVE_GSL - // XXXTODO - return 0; + int k = -1; + size_t nrows = 0, ncols = 0; + int32_t target = EXPR::IMATRIX; + bool have_matrix = false; + for (size_t i = 0; i < n; i++) { + pure_expr *x = xs[i]; + switch (x->tag) { + case EXPR::DBL: + if (target == EXPR::IMATRIX) target = EXPR::MATRIX; + case EXPR::INT: + case EXPR::BIGINT: + if (k >= 0 && k != 1) return 0; + nrows++; k = 1; + break; + case EXPR::APP: { + double a, b; + if (!get_complex(x, a, b)) return 0; + if (k >= 0 && k != 1) return 0; + target = EXPR::CMATRIX; + nrows++; k = 1; + break; + } + case EXPR::MATRIX: { + gsl_matrix *mp = (gsl_matrix*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size2 != (size_t)k) + return 0; + nrows += mp->size1; k = mp->size2; + } else if (k>0) + return 0; + if (target == EXPR::IMATRIX) target = EXPR::MATRIX; + have_matrix = true; + break; + } + case EXPR::CMATRIX: { + gsl_matrix_complex *mp = (gsl_matrix_complex*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size2 != (size_t)k) + return 0; + nrows += mp->size1; k = mp->size2; + } else if (k>0) + return 0; + target = EXPR::CMATRIX; + have_matrix = true; + break; + } + case EXPR::IMATRIX: { + gsl_matrix_int *mp = (gsl_matrix_int*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size2 != (size_t)k) + return 0; + nrows += mp->size1; k = mp->size2; + } else if (k>0) + return 0; + have_matrix = true; + break; + } + default: + return 0; + } + } + if (n == 1 && have_matrix) return xs[0]; + if (k < 0) k = 0; + ncols = k; + if (nrows == 0 && ncols == 0) target = EXPR::MATRIX; + switch (target) { + case EXPR::MATRIX: + return double_matrix_rows(nrows, ncols, n, xs); + case EXPR::CMATRIX: + return complex_matrix_rows(nrows, ncols, n, xs); + case EXPR::IMATRIX: + return int_matrix_rows(nrows, ncols, n, xs); + default: + return 0; + } #else return 0; #endif @@ -1150,19 +1228,97 @@ pure_expr *pure_matrix_columnsl(uint32_t n, ...) { #ifdef HAVE_GSL - // XXXTODO - return 0; + va_list ap; + pure_expr **xs = (pure_expr**)alloca(n*sizeof(pure_expr*)); + va_start(ap, n); + for (size_t i = 0; i < n; i++) + xs[i] = va_arg(ap, pure_expr*); + va_end(ap); + return pure_matrix_columnsv(n, xs); #else return 0; #endif } extern "C" -pure_expr *pure_matrix_columnsv(uint32_t n, pure_expr **elems) +pure_expr *pure_matrix_columnsv(uint32_t n, pure_expr **xs) { #ifdef HAVE_GSL - // XXXTODO - return 0; + int k = -1; + size_t nrows = 0, ncols = 0; + int32_t target = EXPR::IMATRIX; + bool have_matrix = false; + for (size_t i = 0; i < n; i++) { + pure_expr *x = xs[i]; + switch (x->tag) { + case EXPR::DBL: + if (target == EXPR::IMATRIX) target = EXPR::MATRIX; + case EXPR::INT: + case EXPR::BIGINT: + if (k >= 0 && k != 1) return 0; + ncols++; k = 1; + break; + case EXPR::APP: { + double a, b; + if (!get_complex(x, a, b)) return 0; + if (k >= 0 && k != 1) return 0; + target = EXPR::CMATRIX; + ncols++; k = 1; + break; + } + case EXPR::MATRIX: { + gsl_matrix *mp = (gsl_matrix*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size1 != (size_t)k) + return 0; + ncols += mp->size2; k = mp->size1; + } else if (k>0) + return 0; + if (target == EXPR::IMATRIX) target = EXPR::MATRIX; + have_matrix = true; + break; + } + case EXPR::CMATRIX: { + gsl_matrix_complex *mp = (gsl_matrix_complex*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size1 != (size_t)k) + return 0; + ncols += mp->size2; k = mp->size1; + } else if (k>0) + return 0; + target = EXPR::CMATRIX; + have_matrix = true; + break; + } + case EXPR::IMATRIX: { + gsl_matrix_int *mp = (gsl_matrix_int*)x->data.mat.p; + if (mp) { + if (k >= 0 && mp->size1 != (size_t)k) + return 0; + ncols += mp->size2; k = mp->size1; + } else if (k>0) + return 0; + have_matrix = true; + break; + } + default: + return 0; + } + } + if (n == 1 && have_matrix) return xs[0]; + if (k < 0) k = 0; + nrows = k; + if (nrows == 0 && ncols == 0) target = EXPR::MATRIX; + switch (target) { + case EXPR::MATRIX: + return double_matrix_columns(nrows, ncols, n, xs); + case EXPR::CMATRIX: + return complex_matrix_columns(nrows, ncols, n, xs); + case EXPR::IMATRIX: + return int_matrix_columns(nrows, ncols, n, xs); + default: + return 0; + } #else return 0; #endif Modified: pure/trunk/runtime.h =================================================================== --- pure/trunk/runtime.h 2008-09-18 11:41:11 UTC (rev 781) +++ pure/trunk/runtime.h 2008-09-18 11:44:30 UTC (rev 782) @@ -149,10 +149,10 @@ The pure_matrix_rows functions arrange the elements vertically, while the pure_matrix_columns functions arrange them horizontally, given that the other dimensions match. A null expression is returned in case of an error - (no matrix support, dimension mismatch, or invalid element type). Otherwise - a new matrix expression is returned. Temporary element expressions are - taken over by the callee and will be garbage-collected, but the elems - vectors are owned by the caller and won't be freed. */ + (no matrix support, dimension mismatch, or invalid element type), leaving + the input expressions untouched. Otherwise a new matrix expression is + returned and temporary element expressions are garbage-collected. In any + case, the elems vectors are owned by the caller and won't be freed. */ pure_expr *pure_matrix_rowsl(uint32_t n, ...); pure_expr *pure_matrix_rowsv(uint32_t n, pure_expr **elems); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |