GridSolve: GridSolve/src/matlab_client/matlab_gs

Go to the documentation of this file.
00001 
00011 /* $Id: matlab_gs_utils.c,v 1.20 2008/05/01 13:53:05 yarkhan Exp $ */
00012 /* $UTK_Copyright: $ */
00013 
00014 #include <stdio.h>
00015 #include <string.h>
00016 
00017 #include <matrix.h>
00018 #include <mex.h>
00019 
00020 #include "grpc.h"
00021 #include "comm_data.h"
00022 #include "matlab_gs.h"
00023 
00024 #define MAX_FILES_IN_PACKEDFILE 1024
00025 
00026 typedef enum {DO_NOT_ALLOC, MAYBE_ALLOC} gs_alloc_t;
00027 
00028 #define COPY_A_TO_B(TYPEA, a, aoffset, TYPEB, b, boffset, size)  do {int ss = 0; \
00029 TYPEA* ss_a = (TYPEA*) a; \
00030 TYPEB* ss_b = (TYPEB*) b; \
00031 ss_a+=aoffset; ss_b+=boffset; \
00032 for (ss = 0; ss < size; ss++) \
00033 *((ss_b)++) = (TYPEB) *((ss_a)++);\
00034 } while (0)
00035 
00036 
00037 #define COPY_MXARRAY_TO_GSCOMPLEX(MX_DATA_T, mxr, mxi, mxoffset, GS_COMPLEX_T, gs_cplx, gsoffset, GS_DATA_T, size)  \
00038 do {int ss = 0; \
00039 MX_DATA_T* ss_mxr = (MX_DATA_T*) mxr; \
00040 MX_DATA_T* ss_mxi = (MX_DATA_T*) mxi; \
00041 GS_COMPLEX_T* ss_gs_cplx = (GS_COMPLEX_T*) gs_cplx; \
00042 ss_mxr+=mxoffset; ss_mxi+=mxoffset; ss_gs_cplx+=gsoffset; \
00043 for (ss = 0; ss < size; ss++) { \
00044   ss_gs_cplx[ss].r = (GS_DATA_T) *((ss_mxr)++); \
00045   if (mxi == NULL) \
00046     ss_gs_cplx[ss].i = (GS_DATA_T) 0.0; \
00047   else \
00048     ss_gs_cplx[ss].i = (GS_DATA_T) *((ss_mxi)++); \
00049   } \
00050 } while (0)
00051 
00052 
00064 void *matlab_gs_calloc(size_t N, size_t S)
00065 {
00066   void *ptr;
00067   /* ptr = mxCalloc(N, S); */
00068   ptr = calloc(N, S);
00069   if (ptr == NULL) {
00070     ERRPRINTF("Matlab is out of memory in GridSolve client\n");
00071     return NULL;
00072   }
00073   /* else mexMakeMemoryPersistent(ptr); */
00074   return ptr;
00075 }
00076 
00084 void *matlab_gs_free(void *ptr)
00085 {
00086   /* if ( ptr ) mxFree(ptr); */
00087   if ( ptr ) free(ptr);  
00088   return NULL;
00089 }
00090 
00100 int matlab_gs_arg_is_sparse_index_or_pointer(gs_argument_t *argp, gs_problem_t *problem)
00101 {
00102   int arg_is_sparse_attr = 0;
00103   gs_argument_t* arg_attr = NULL;
00104   
00105   for (arg_attr = problem->arglist; arg_attr != NULL; arg_attr = arg_attr->next)
00106     if ((arg_attr->objecttype == GS_SPARSEMATRIX)  && 
00107     (!strcmp(argp->name, arg_attr->sparse_attr.indices) || 
00108      !strcmp(argp->name, arg_attr->sparse_attr.pointer))) {
00109       arg_is_sparse_attr = 1;
00110       break;
00111     }
00112   return arg_is_sparse_attr;
00113 }
00114 
00115 
00127 static mxArray *matlab_gs_create_matrix(
00128   int objecttype, int M, int N, int nnz, int datatype)
00129 {
00130   mxComplexity realorcomplex = mxREAL;
00131   mxClassID classid = mxDOUBLE_CLASS;
00132   mxArray *ptr;
00133 
00134   switch (datatype) {
00135   case GS_INT: 
00136   case GS_FLOAT:
00137   case GS_DOUBLE:
00138     classid = mxDOUBLE_CLASS;
00139     realorcomplex = mxREAL;
00140     break;
00141   case GS_SCOMPLEX:
00142   case GS_DCOMPLEX:
00143     classid = mxDOUBLE_CLASS;
00144     realorcomplex = mxCOMPLEX;
00145     break;
00146   case GS_CHAR:
00147     break;
00148   case GS_BAD_DTYPE:
00149     return NULL;
00150     break;
00151   default:
00152     return NULL;
00153     break;
00154   }
00155 
00156   if (objecttype == GS_SPARSEMATRIX) {
00157     ptr = mxCreateSparse(M, N ,nnz, realorcomplex);
00158   } else {
00159     ptr = mxCreateNumericMatrix(M, N, classid, realorcomplex); 
00160   }
00161   return ptr;
00162 }
00163 
00164 
00178 static int transfer_mx_to_gs_int(mxClassID mxclassid, void *pr, int poffset, void** datapp, int doffset, int nelem, gs_alloc_t allocp)
00179 {
00180   /* Allocate space, and copy data over to it */
00181   if (allocp != DO_NOT_ALLOC) 
00182     if ((*datapp = (int*)matlab_gs_calloc(nelem, sizeof(int)))==NULL) return -1;
00183   
00184   if (mxclassid == mxDOUBLE_CLASS) 
00185     COPY_A_TO_B(double, pr, poffset, int, *datapp, doffset, nelem);
00186   else if (mxclassid == mxSINGLE_CLASS)
00187     COPY_A_TO_B(float, pr, poffset, int, *datapp, doffset, nelem);
00188   else if (mxclassid == mxINT8_CLASS)
00189     COPY_A_TO_B(int8_t, pr, poffset, int, *datapp, doffset, nelem);
00190   else if (mxclassid == mxUINT8_CLASS)
00191     COPY_A_TO_B(uint8_t, pr, poffset, int, *datapp, doffset, nelem);
00192   else if (mxclassid == mxINT16_CLASS)
00193     COPY_A_TO_B(int16_t, pr, poffset, int, *datapp, doffset, nelem);
00194   else if (mxclassid == mxUINT16_CLASS)
00195     COPY_A_TO_B(uint16_t, pr, poffset, int, *datapp, doffset, nelem);
00196   else if (mxclassid == mxINT32_CLASS)
00197     COPY_A_TO_B(int32_t, pr, poffset, int, *datapp, doffset, nelem);
00198   else if (mxclassid == mxUINT32_CLASS)
00199     COPY_A_TO_B(uint32_t, pr, poffset, int, *datapp, doffset, nelem);
00200   else if (mxclassid == mxCHAR_CLASS)
00201     COPY_A_TO_B(char, pr, poffset, int, *datapp, doffset, nelem);
00202 /*   else if (mxclassid == mxINT64_CLASS) */
00203 /*     COPY_A_TO_B(int64_t, pr, poffset, int, *datapp, doffset, nelem); */
00204 /*   else if (mxclassid == mxUINT64_CLASS) */
00205 /*     COPY_A_TO_B(int64_t, pr, poffset, int, *datapp, doffset, nelem); */
00206   else {
00207     ERRPRINTF("Could not convert matlab data to GridSolve data type (int)\n");
00208     return -1;
00209   }
00210   return 1;
00211 }
00212 
00226 static int transfer_mx_to_gs_double(mxClassID mxclassid, void *pr, int poffset, void** datapp, int doffset, int nelem, gs_alloc_t allocp)
00227 {
00228   /* Allocate space if needed, and copy data over to it */
00229   if (allocp != DO_NOT_ALLOC) 
00230     if ((*datapp = (double*)matlab_gs_calloc(nelem, sizeof(double)))==NULL) return -1;
00231 
00232   if (mxclassid == mxDOUBLE_CLASS) 
00233     COPY_A_TO_B(double, pr, poffset, double, *datapp, doffset, nelem);
00234   else if (mxclassid == mxSINGLE_CLASS)
00235     COPY_A_TO_B(float, pr, poffset, double, *datapp, doffset, nelem);
00236   else if (mxclassid == mxINT8_CLASS)
00237     COPY_A_TO_B(int8_t, pr, poffset, double, *datapp, doffset, nelem);
00238   else if (mxclassid == mxUINT8_CLASS)
00239     COPY_A_TO_B(uint8_t, pr, poffset, double, *datapp, doffset, nelem);
00240   else if (mxclassid == mxINT16_CLASS)
00241     COPY_A_TO_B(int16_t, pr, poffset, double, *datapp, doffset, nelem);
00242   else if (mxclassid == mxUINT16_CLASS)
00243     COPY_A_TO_B(uint16_t, pr, poffset, double, *datapp, doffset, nelem);
00244   else if (mxclassid == mxINT32_CLASS)
00245     COPY_A_TO_B(int32_t, pr, poffset, double, *datapp, doffset, nelem);
00246   else if (mxclassid == mxUINT32_CLASS)
00247     COPY_A_TO_B(uint32_t, pr, poffset, double, *datapp, doffset, nelem);
00248 /*   else if (mxclassid == mxINT64_CLASS) */
00249 /*     COPY_A_TO_B(int64_t, pr, poffset, double, *datapp, doffset, nelem); */
00250 /*   else if (mxclassid == mxUINT64_CLASS) */
00251 /*     COPY_A_TO_B(uint64_t, pr, poffset, double, *datapp, doffset, nelem); */
00252   else {
00253     ERRPRINTF("Could not convert matlab data to GridSolve data type (double)\n");
00254     return -1;
00255   }
00256 
00257   return 1;
00258 }
00259 
00273 static int transfer_mx_to_gs_float(mxClassID mxclassid, void *pr, int poffset, void** datapp, int doffset, int nelem, gs_alloc_t allocp)
00274 {
00275   /* Allocate space if needed, and copy data over to it */
00276   if (allocp != DO_NOT_ALLOC) 
00277     if ((*datapp = (double*)matlab_gs_calloc(nelem, sizeof(float)))==NULL) return -1;
00278 
00279   if (mxclassid == mxDOUBLE_CLASS) 
00280     COPY_A_TO_B(double, pr, poffset, float, *datapp, doffset, nelem);
00281   else if (mxclassid == mxSINGLE_CLASS)
00282     COPY_A_TO_B(float, pr, poffset, float, *datapp, doffset, nelem);
00283   else if (mxclassid == mxINT8_CLASS)
00284     COPY_A_TO_B(int8_t, pr, poffset, float, *datapp, doffset, nelem);
00285   else if (mxclassid == mxUINT8_CLASS)
00286     COPY_A_TO_B(uint8_t, pr, poffset, float, *datapp, doffset, nelem);
00287   else if (mxclassid == mxINT16_CLASS)
00288     COPY_A_TO_B(int16_t, pr, poffset, float, *datapp, doffset, nelem);
00289   else if (mxclassid == mxUINT16_CLASS)
00290     COPY_A_TO_B(uint16_t, pr, poffset, float, *datapp, doffset, nelem);
00291   else if (mxclassid == mxINT32_CLASS)
00292     COPY_A_TO_B(int32_t, pr, poffset, float, *datapp, doffset, nelem);
00293   else if (mxclassid == mxUINT32_CLASS)
00294     COPY_A_TO_B(uint32_t, pr, poffset, float, *datapp, doffset, nelem);
00295 /*   else if (mxclassid == mxINT64_CLASS) */
00296 /*     COPY_A_TO_B(int64_t, pr, poffset, float, *datapp, doffset, nelem); */
00297 /*   else if (mxclassid == mxUINT64_CLASS) */
00298 /*     COPY_A_TO_B(uint64_t, pr, poffset, float, *datapp, doffset, nelem); */
00299   else {
00300     ERRPRINTF("Could not convert matlab data to GridSolve data type (float)\n");
00301     return -1;
00302   }
00303   return 1;
00304 }
00305 
00306 
00321 static int transfer_mx_to_gs_scomplex(mxClassID mxclassid, void *pr, void *pi, int poffset, void** datapp, int doffset, int nelem, gs_alloc_t allocp)
00322 {
00323   /* Allocate space, and copy data over to it */
00324   if (allocp != DO_NOT_ALLOC) 
00325     if ((*datapp = (gs_scomplex*)matlab_gs_calloc(nelem, sizeof(gs_scomplex)))==NULL) return -1; 
00326   
00327   if (mxclassid == mxDOUBLE_CLASS) 
00328     COPY_MXARRAY_TO_GSCOMPLEX(double, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00329   else if (mxclassid == mxSINGLE_CLASS)
00330     COPY_MXARRAY_TO_GSCOMPLEX(float, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00331   else if (mxclassid == mxINT8_CLASS)
00332     COPY_MXARRAY_TO_GSCOMPLEX(int8_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00333   else if (mxclassid == mxUINT8_CLASS)
00334     COPY_MXARRAY_TO_GSCOMPLEX(uint8_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00335   else if (mxclassid == mxINT16_CLASS)
00336     COPY_MXARRAY_TO_GSCOMPLEX(int16_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00337   else if (mxclassid == mxUINT16_CLASS)
00338     COPY_MXARRAY_TO_GSCOMPLEX(uint16_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00339   else if (mxclassid == mxINT32_CLASS)
00340     COPY_MXARRAY_TO_GSCOMPLEX(int32_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00341   else if (mxclassid == mxUINT32_CLASS)
00342     COPY_MXARRAY_TO_GSCOMPLEX(uint32_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem);
00343 /*   else if (mxclassid == mxINT64_CLASS) */
00344 /*     COPY_MXARRAY_TO_GSCOMPLEX(int64_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem); */
00345 /*   else if (mxclassid == mxUINT64_CLASS) */
00346 /*     COPY_MXARRAY_TO_GSCOMPLEX(uint64_t, pr, pi, poffset, gs_scomplex, *datapp, doffset, float, nelem); */
00347   else {
00348     ERRPRINTF("Could not convert matlab data to GridSolve data type (single precision complex)\n");
00349     return -1;
00350   }
00351   return 1;
00352 }
00353 
00354 
00369 static int transfer_mx_to_gs_dcomplex(mxClassID mxclassid, void *pr, void *pi, int poffset, void** datapp, int doffset, int nelem, gs_alloc_t allocp)
00370 {
00371   /* Allocate space, and copy data over to it */
00372   if (allocp != DO_NOT_ALLOC) 
00373     if ((*datapp = (gs_dcomplex*)matlab_gs_calloc(nelem, sizeof(gs_dcomplex)))==NULL) return -1; 
00374   
00375   if (mxclassid == mxDOUBLE_CLASS) 
00376     COPY_MXARRAY_TO_GSCOMPLEX(double, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00377   else if (mxclassid == mxSINGLE_CLASS)
00378     COPY_MXARRAY_TO_GSCOMPLEX(float, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00379   else if (mxclassid == mxINT8_CLASS)
00380     COPY_MXARRAY_TO_GSCOMPLEX(int8_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00381   else if (mxclassid == mxUINT8_CLASS)
00382     COPY_MXARRAY_TO_GSCOMPLEX(uint8_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00383   else if (mxclassid == mxINT16_CLASS)
00384     COPY_MXARRAY_TO_GSCOMPLEX(int16_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00385   else if (mxclassid == mxUINT16_CLASS)
00386     COPY_MXARRAY_TO_GSCOMPLEX(uint16_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00387   else if (mxclassid == mxINT32_CLASS)
00388     COPY_MXARRAY_TO_GSCOMPLEX(int32_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00389   else if (mxclassid == mxUINT32_CLASS)
00390     COPY_MXARRAY_TO_GSCOMPLEX(uint32_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem);
00391 /*   else if (mxclassid == mxINT64_CLASS) */
00392 /*     COPY_MXARRAY_TO_GSCOMPLEX(int64_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem); */
00393 /*   else if (mxclassid == mxUINT64_CLASS) */
00394 /*     COPY_MXARRAY_TO_GSCOMPLEX(uint64_t, pr, pi, poffset, gs_dcomplex, *datapp, doffset, double, nelem); */
00395   else {
00396     ERRPRINTF("Could not convert matlab data to GridSolve data type (double precision complex)\n");
00397     return -1;
00398   }
00399   return 1;
00400 }
00401 
00402 
00413 static mxArray *transfer_gs_to_mx(gs_argument_t *argp)
00414 {
00415   int nelements, N, M, j;
00416   mxArray *mxarray;
00417   void *datap;
00418 
00419   if (argp->inout==GS_VAROUT) 
00420     datap = *((void **)argp->data);
00421   else 
00422     datap = argp->data;
00423 
00424   nelements = argp->rows*argp->cols; 
00425   M = argp->rows;
00426   N = argp->cols;
00427 
00428   switch (argp->datatype) {
00429   case GS_DOUBLE:
00430     if ((mxarray = matlab_gs_create_matrix(argp->objecttype, M, N, nelements, argp->datatype)) == NULL) goto error_output;
00431     for (j=0; j<nelements; j++) 
00432       mxGetPr(mxarray)[j] = ((double*)((double *)datap))[j];
00433     break;
00434   case GS_INT:
00435     if ((mxarray = matlab_gs_create_matrix(argp->objecttype, M, N, nelements, argp->datatype)) == NULL) goto error_output;
00436     for (j=0; j<nelements; j++) 
00437       mxGetPr(mxarray)[j] = (double)(((int*)(datap))[j]);
00438     break;
00439   case GS_FLOAT:
00440     if ((mxarray = matlab_gs_create_matrix(argp->objecttype, M, N, nelements, argp->datatype)) == NULL) goto error_output;
00441     for (j=0; j<nelements; j++) 
00442       mxGetPr(mxarray)[j] = (double)(((float*)(datap))[j]);
00443     break;
00444   case GS_SCOMPLEX:
00445     if ((mxarray = matlab_gs_create_matrix(argp->objecttype, M, N, nelements, argp->datatype)) == NULL) goto error_output;
00446     for (j=0; j<nelements; j++) {
00447       mxGetPr(mxarray)[j] = (double)(((gs_scomplex*)(datap))[j].r);
00448       mxGetPi(mxarray)[j] = (double)(((gs_scomplex*)(datap))[j].i);
00449     }
00450     break;
00451   case GS_DCOMPLEX:
00452     if ((mxarray = matlab_gs_create_matrix(argp->objecttype, M, N, nelements, argp->datatype)) == NULL) goto error_output;
00453     for (j=0; j<nelements; j++) {
00454       mxGetPr(mxarray)[j] = (double)(((gs_dcomplex*)(datap))[j].r);
00455       mxGetPi(mxarray)[j] = (double)(((gs_dcomplex*)(datap))[j].i);
00456     }
00457     break;
00458   case GS_CHAR:
00459   {
00460     mxChar *charData;
00461     int dims[2] = {argp->rows, argp->cols};
00462     if ((mxarray = mxCreateCharArray(2, dims)) == NULL) goto error_output;
00463     charData = (mxChar *)mxGetData(mxarray);
00464     for(j=0; j<argp->rows*argp->cols; j++) 
00465       charData[j] = (mxChar)((char*)datap)[j]; 
00466     break;
00467   }
00468   default:
00469     ERRPRINTF("Unknown datatype %d\n", argp->datatype);
00470     goto error_output;
00471     break;
00472   }
00473   return mxarray;
00474 
00475  error_output:
00476   ERRPRINTF("GridSolve: Error transfering  GridSolve output (arg %s type %s) to Matlab", argp->name, gs_c_datatype[argp->datatype]);
00477   return NULL;
00478 
00479 
00480 }
00481 
00482 
00492 static int arg_data_same_as_mx_dat(const mxArray *mxarray, gs_argument_t *argp) 
00493 {
00494   int gsdatatype = argp->datatype;
00495   mxClassID mxclassid = mxGetClassID(mxarray);
00496 
00497   if (argp->inout != GS_IN) return 0;
00498 
00499   if (mxclassid==mxDOUBLE_CLASS && gsdatatype==GS_DOUBLE) 
00500     return 1;
00501   else if (mxclassid==mxSINGLE_CLASS && gsdatatype==GS_FLOAT) 
00502     return 1;
00503   else if (mxclassid==mxINT32_CLASS && gsdatatype==GS_INT) 
00504     return 1;
00505   else 
00506     return 0;
00507 }
00508 
00509 
00535 static int matlab_gs_sparse_matcsc_to_csr(
00536   int m, int n, int nnz, 
00537   int mxclassid, void *ar, void *ai, int *colind, int *rowptr,
00538   int gsdatatype, void *gsat, int *rowind, int *colptr)
00539 {
00540   register int i, j, col, relpos;
00541   int *marker;
00542   
00543   marker = (int *) matlab_gs_calloc(n, sizeof(int));
00544   if (marker == NULL) return -1;
00545  
00546   /* Get counts of each column of A, and set up column pointers */
00547   for (i = 0; i < m; ++i)
00548     for (j = rowptr[i]; j < rowptr[i+1]; ++j) 
00549       ++marker[colind[j]];
00550   
00551   colptr[0] = 0;
00552 
00553   for (j = 0; j < n; ++j) {
00554     colptr[j+1] = colptr[j] + marker[j];
00555     marker[j] = colptr[j];
00556   }
00557  
00558   /* Transfer the matrix into the compressed column storage. */
00559   for (i = 0; i < m; ++i) {
00560     for (j = rowptr[i]; j < rowptr[i+1]; ++j) {
00561       col = colind[j];
00562       relpos = marker[col];
00563       rowind[relpos] = i;
00564 
00565       /* Copy one data item at a time */
00566       switch(gsdatatype) {
00567       case GS_INT:
00568     /* Transfer 1 item from ar (starting at j) to gsat (starting at relpos) */
00569     if (transfer_mx_to_gs_int(mxclassid, ar, j, &gsat, relpos, 1, DO_NOT_ALLOC) == -1) goto error;
00570     break;
00571       case GS_FLOAT:
00572     if (transfer_mx_to_gs_float(mxclassid, ar, j, &gsat, relpos, 1, DO_NOT_ALLOC) == -1) goto error;
00573     break;
00574       case GS_DOUBLE:
00575     if (transfer_mx_to_gs_double(mxclassid, ar, j, &gsat, relpos, 1, DO_NOT_ALLOC) == -1) goto error;
00576     break;
00577       case GS_SCOMPLEX:
00578     if (transfer_mx_to_gs_scomplex(mxclassid, ar, ai, j, &gsat, relpos, 1, DO_NOT_ALLOC) == -1) goto error;
00579     break;
00580       case GS_DCOMPLEX:
00581     if (transfer_mx_to_gs_dcomplex(mxclassid, ar, ai, j, &gsat, relpos, 1, DO_NOT_ALLOC) == -1) goto error;
00582     break;
00583       default:
00584     ERRPRINTF("GridSolve sparse of datatype (%s) not supported \n", gs_c_datatype[gsdatatype]);
00585     /* Return an error  */
00586     return -1;
00587     break;
00588       }
00589 
00590       ++marker[col];
00591     }
00592   }
00593   marker = matlab_gs_free(marker);
00594   return 0;
00595 
00596  error:
00597   ERRPRINTF("Could not convert matlab data to GridSolve sparse CSR format (%s)\n", gs_c_datatype[gsdatatype]);
00598   return -1;
00599 }
00600 
00626 static int matlab_gs_sparse_csr_to_matcsc(
00627   int m, int n, int nnz, 
00628   int datatype, void *a, int *colind, int *rowptr,
00629   double *at, double *ati, int *rowind, int *colptr)
00630 {
00631   register int i, j, col, relpos;
00632   int *marker;
00633 
00634   marker = (int *)matlab_gs_calloc(n, sizeof(int));
00635   if (marker == NULL) return -1;
00636  
00637   /* Get counts of each column of A, and set up column pointers */
00638   for (i = 0; i < m; ++i)
00639     for (j = rowptr[i]; j < rowptr[i+1]; ++j) ++marker[colind[j]];
00640 
00641   colptr[0] = 0;
00642 
00643   for (j = 0; j < n; ++j) {
00644     colptr[j+1] = colptr[j] + marker[j];
00645     marker[j] = colptr[j];
00646   }
00647  
00648   /* Transfer the matrix into the compressed column storage. */
00649   for (i = 0; i < m; ++i) {
00650     for (j = rowptr[i]; j < rowptr[i+1]; ++j) {
00651       col = colind[j];
00652       relpos = marker[col];
00653       rowind[relpos] = i;
00654       
00655       switch(datatype) {
00656       case GS_INT:
00657     ((double *)at)[relpos] = (double)((int *)a)[j];
00658     break;
00659       case GS_FLOAT:
00660     ((double *)at)[relpos] = (double)((float *)a)[j];
00661     break;
00662       case GS_DOUBLE:
00663     ((double *)at)[relpos] = ((double *)a)[j];
00664     break;
00665       case GS_SCOMPLEX:
00666     ((double *)at)[relpos] = (double)(((gs_scomplex *)a)[j].r);
00667     ((double *)ati)[relpos] = (double)(((gs_scomplex *)a)[j].i);
00668     break;
00669       case GS_DCOMPLEX:
00670     ((double *)at)[relpos] = (double)((gs_dcomplex *)a)[j].r;
00671     ((double *)ati)[relpos] = (double)((gs_dcomplex *)a)[j].i;
00672     break;
00673       default:
00674     ERRPRINTF("GridSolve sparse of datatype (%s) not supported \n", gs_c_datatype[datatype]);
00675     /* Return an error  */
00676     return -1;
00677     break;
00678       }
00679       ++marker[col];
00680     }
00681   }
00682   marker = matlab_gs_free(marker);
00683   return 0;
00684 }
00685 
00686 
00702 static int matlab_gs_setup_scalar_arguments(
00703   gs_problem_t* problem, int nlhs, mxArray *plhs[],
00704   int nrhs, const mxArray *prhs[])
00705 {
00706   gs_argument_t* argp = NULL;
00707   double *pr = NULL;
00708   double *pi = NULL;
00709   int irhs = 0;
00710 
00711   for (argp=problem->arglist,irhs=2; argp!=NULL && irhs<=nrhs; argp=argp->next) {
00712     
00713     /* Skip sparse indices and pointers, scalar nnz values must be provided */
00714     if (matlab_gs_arg_is_sparse_index_or_pointer(argp, problem)) {
00715       DBGPRINTF("Skipping sparse attribute: %s \n", argp->name);
00716       continue;
00717     }
00718 
00719     /* Grab scalar input values so that we can calulate sizes of non scalars */
00720     if (argp->inout==GS_IN || argp->inout==GS_INOUT ) {
00721       switch (argp->objecttype) {
00722       case GS_SCALAR:
00723     DBGPRINTF("Transfer scalar Matlab arg %d (%s %d) to GridSolve arg (%s %s)\n", irhs-2, mxGetClassName(prhs[irhs]), mxGetNumberOfElements(prhs[irhs]),argp->name, gs_c_datatype[argp->datatype]);
00724     if (mxGetNumberOfElements(prhs[irhs])!=1) { DBGPRINTF("numelements mismatch\n"); goto error; }
00725     pr = mxGetPr(prhs[irhs]); 
00726     pi = mxGetPi(prhs[irhs]); 
00727     argp->data = &(argp->scalar_val);
00728     switch (argp->datatype) {
00729     case GS_DOUBLE:
00730       if (transfer_mx_to_gs_double(mxGetClassID(prhs[irhs]), pr, 0, (void **)&(argp->data), 0, 1, DO_NOT_ALLOC)<0) goto error;
00731       break;
00732     case GS_INT:
00733       if (transfer_mx_to_gs_int(mxGetClassID(prhs[irhs]), pr, 0, (void **)&(argp->data), 0, 1, DO_NOT_ALLOC)<0) goto error;
00734       break;
00735     case GS_FLOAT:
00736       if (transfer_mx_to_gs_float(mxGetClassID(prhs[irhs]), pr, 0, (void **)&(argp->data), 0, 1, DO_NOT_ALLOC)<0) goto error;
00737       break;
00738     case GS_SCOMPLEX:
00739       if (transfer_mx_to_gs_scomplex(mxGetClassID(prhs[irhs]), pr, pi, 0, (void **)&(argp->data), 0, 1, DO_NOT_ALLOC)<0) goto error;
00740       break;
00741     case GS_DCOMPLEX:
00742       if (transfer_mx_to_gs_dcomplex(mxGetClassID(prhs[irhs]), pr, pi, 0, (void **)&(argp->data), 0, 1, DO_NOT_ALLOC)<0) goto error;
00743       break;
00744     case GS_CHAR:
00745       argp->scalar_val.char_val = (char)mxGetScalar(prhs[irhs]);
00746       break;
00747     default:
00748       argp->data = NULL;
00749       ERRPRINTF("GridSolve datatype (%s) is not handled yet\n", gs_c_datatype[argp->datatype]);
00750       goto error;
00751       break;
00752     }
00753     break;
00754       default:
00755     DBGPRINTF("Vector, matrix or other non-scalar input (rhs arg %d name %s )\n", irhs-2, argp->name);
00756     break;
00757       }
00758       irhs++;
00759       
00760     } else if (argp->inout==GS_OUT) {
00761       /* For scalar output, the argp->data will point to arg->scalar_val */
00762       switch (argp->objecttype) {
00763       case GS_SCALAR:
00764     DBGPRINTF("Scalar output (name %s)\n", argp->name);
00765     argp->data = &(argp->scalar_val);
00766     break;
00767       default:
00768     DBGPRINTF("Non scalar output (arg %d name %s )\n", irhs-2, argp->name);
00769     break;
00770       }
00771       
00772     } else if (argp->inout==GS_VAROUT) {
00773       /* For scalar varout, the argp->data will point to arg->scalar_val */
00774       switch (argp->objecttype) {
00775       case GS_SCALAR:
00776     ERRPRINTF("GridSolve scalars arguments cannot be marked as VAROUT in the interface (matlab input arg %d)\n", irhs-2);
00777     goto error;
00778     break;
00779       default:
00780     DBGPRINTF("Non scalar varout (arg %d name %s )\n", irhs-2, argp->name);
00781     break;
00782       }
00783 
00784     } else if (argp->inout==GS_WORKSPACE) {
00785       /* Do nothing, just move to next argument */
00786       argp->data = NULL;
00787       
00788     } else {
00789       ERRPRINTF("GridSolve interface in/out type is unknown (not IN, OUT, INOUT, or VAROUT or WORKSPACE) (matlab input argument %d) \n", irhs-2);
00790       goto error;
00791     }
00792   }
00793   return 0;
00794 
00795  error:
00796   ERRPRINTF("GridSolve: Scalar argument conversion problem: Matlab arg num %d (type %s size %d) could not be converted to GridSolve arg (name %s type %s)\n", irhs-2, mxGetClassName(prhs[irhs]), mxGetNumberOfElements(prhs[irhs]),argp->name, gs_c_datatype[argp->datatype]);
00797   return -1;
00798 }
00799 
00800 
00818 int
00819 matlab_gs_setup_args(gs_problem_t* problem, int nlhs, mxArray *plhs[],
00820              int nrhs, const mxArray *prhs[])
00821 {
00822   gs_argument_t* argp = NULL;
00823   double *pr = NULL;
00824   double *pi = NULL;
00825   mxChar *mxcharp = NULL;
00826   int nelements;
00827   int i,j,m,n,irhs=0;
00828   char buffer[256];
00829   FILE *fp;
00830   char *tmpname;
00831   char **filenames = NULL;
00832   gs_argument_t* arg_indices = NULL;
00833   gs_argument_t* arg_pointer = NULL;
00834   int *ir = NULL;
00835   int *jc = NULL;
00836   int nnz;
00837   
00838   if (problem==NULL) goto error_unknown;
00839 
00840   /* TODO Check that we have correct number of input and output args */
00841 
00842   /* Setup scalar arguments, copying from Matlab as needed */
00843   if (matlab_gs_setup_scalar_arguments(problem, nlhs, plhs, nrhs, prhs) < 0) goto error_scalar;
00844 
00845   /* Compute sizes for all other arguments using scalar arguments */
00846   if (gs_receiver_compute_arg_sizes(problem, GS_IN) < 0) {
00847     ERRPRINTF("GridSolve: Could not compute argument sizes");  
00848     return -1;
00849   }
00850   
00851   /* For each non-scalar argument, check datatype, and attach data to the argp->data */
00852   DBGPRINTF("For each argument, check datatype, attach data for non-scalar args\n");
00853   for (argp=problem->arglist,irhs=2; argp!=NULL && irhs<=nrhs; argp=argp->next) {
00854     
00855     /* skip this if it's a sparse matrix attribute, since the sparse
00856        attributes (for sparse arguments) (IN, INOUT or OUT) are not
00857        sent from matlab, they will be handled when the actual sparse
00858        argument is handled  */
00859     if (matlab_gs_arg_is_sparse_index_or_pointer(argp, problem)) {
00860       DBGPRINTF("Skipping sparse attribute: %s \n", argp->name);
00861       continue;
00862     }
00863 
00864     /* ***************************** */
00865     if (argp->inout==GS_IN || argp->inout==GS_INOUT) {
00866     /* ***************************** */
00867 
00868       m = mxGetM(prhs[irhs]); 
00869       n = mxGetN(prhs[irhs]); 
00870       pr = mxGetPr(prhs[irhs]); 
00871       pi = mxGetPi(prhs[irhs]); 
00872 
00873       /* Setup input non scalars, copying args from Matlab to GridSolve, allocating space for arguments */
00874       switch (argp->objecttype) {
00875 
00876       case GS_SCALAR:
00877     /* Nothing to do for scalars, they are already handled */
00878     break;
00879 
00880       case GS_SPARSEMATRIX:
00881     DBGPRINTF("Transfer mx (rhs %d) to gs (%s [%d,%d] %s %s %s) \n", irhs-2, argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype]); fflush(0);
00882     nnz = mxGetNzmax(prhs[irhs]);
00883     /* Lookup which gridsolve arg should hold indices and pointers */
00884     arg_indices = gs_arg_lookup_by_name(argp->sparse_attr.indices, problem->arglist);
00885     arg_pointer = gs_arg_lookup_by_name(argp->sparse_attr.pointer, problem->arglist);
00886     if (arg_indices==NULL || arg_pointer==NULL) return -1;
00887     /* Allocate space for indices, pointers, data */
00888     if ((arg_indices->data = matlab_gs_calloc(nnz, sizeof(int))) == NULL) return -1;
00889     if ((arg_pointer->data = matlab_gs_calloc(n+1, sizeof(int))) == NULL) return -1;
00890     if ((argp->data = matlab_gs_calloc(nnz, gs_datatype_sizeof(argp->datatype))) == NULL) return -1;
00891     /* Copy CSC to CSR */
00892     ir = mxGetIr(prhs[irhs]);
00893     jc = mxGetJc(prhs[irhs]);
00894     if (matlab_gs_sparse_matcsc_to_csr(m, n, nnz, mxGetClassID(prhs[irhs]), pr, pi, ir, jc, argp->datatype, argp->data, arg_indices->data, arg_pointer->data) == -1) goto error_input;
00895     break;
00896     
00897       case GS_MATRIX:
00898       case GS_VECTOR:
00899     DBGPRINTF("Transfer mx (rhs %d) to gs (%s [%d,%d] %s %s %s) \n", irhs-2, argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype]); fflush(0);
00900     nelements = argp->rows*argp->cols;
00901 
00902     if (arg_data_same_as_mx_dat(prhs[irhs], argp)) {
00903       /* If argument datatypes match then we can simply point to
00904        * the input provided by Matlab */
00905       argp->data = pr;
00906     } else {
00907       /* Copy data from Matlab to GridSolve */
00908       switch (argp->datatype) {
00909       case GS_DOUBLE:
00910         if (transfer_mx_to_gs_double(mxGetClassID(prhs[irhs]), pr, 0, &argp->data, 0, nelements, MAYBE_ALLOC) < 0)  goto error_input;
00911         break;
00912       case GS_INT:
00913         if (transfer_mx_to_gs_int(mxGetClassID(prhs[irhs]), pr, 0, &argp->data, 0, nelements, MAYBE_ALLOC) < 0) goto error_input;
00914         break;
00915       case GS_FLOAT:
00916         if (transfer_mx_to_gs_float(mxGetClassID(prhs[irhs]), pr, 0, &argp->data, 0, nelements, MAYBE_ALLOC) < 0) goto error_input;
00917         break;
00918       case GS_SCOMPLEX:
00919         if (transfer_mx_to_gs_scomplex(mxGetClassID(prhs[irhs]), pr, pi, 0, &argp->data, 0, nelements, MAYBE_ALLOC) < 0) goto error_input;
00920         break;
00921       case GS_DCOMPLEX:
00922         if (transfer_mx_to_gs_dcomplex(mxGetClassID(prhs[irhs]), pr, pi, 0, &argp->data, 0, nelements, MAYBE_ALLOC) < 0) goto error_input;
00923         break;
00924       case GS_CHAR:
00925         if ((argp->data = (void*)matlab_gs_calloc(nelements, gs_datatype_sizeof(argp->datatype)))==NULL) goto error_input;
00926         for (j=0; j<nelements; j++) ((char*)argp->data)[j] = (char)((mxChar*)mxGetData(prhs[irhs]))[j];
00927         break;
00928       default:
00929         DBGPRINTF("Datatype not handled\n");
00930         break;          
00931       } /* switch(argp->datatype) */
00932     } /* else arg_data_same_as_mx_dat */
00933     break;
00934 
00935       case GS_FILE:
00936     DBGPRINTF("Transfer mx (rhs %d) to gs (%s [%d,%d] %s %s %s) \n", irhs-2, argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype]); fflush(0);
00937     if (!mxIsChar(prhs[irhs])) {
00938       ERRPRINTF("Arg %s is of class %s but should be a string\n", argp->name, mxGetClassName(prhs[irhs]));
00939       goto error_input;
00940     }
00941     if ((argp->data = (char *)matlab_gs_calloc(m*n+1, sizeof(char))) == NULL) return -1;
00942     mxGetString(prhs[irhs], argp->data, n*m+1);
00943     break;
00944     
00945       case GS_PACKEDFILE:
00946     DBGPRINTF("Transfer mx (rhs %d) to gs (%s [%d,%d] %s %s %s) \n", irhs-2, argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype]); fflush(0);
00947     if (!mxIsChar(prhs[irhs])) {
00948       ERRPRINTF("Arg %s is of class %s but should be a matrix of strings\n", argp->name, mxGetClassName(prhs[irhs]));
00949       goto error_input;
00950     }
00951     mxcharp = (mxChar*)mxGetPr(prhs[irhs]);
00952     if ((filenames = (char **)matlab_gs_calloc(m, sizeof(char*))) == NULL) return -1;
00953     for (i=0; i<m; i++) {
00954       if ((filenames[i] = (char *)matlab_gs_calloc(n+1, sizeof(char))) == NULL) return -1;
00955       /* Copy over full name, going from column major to row major */
00956       for (j=0; j<n; j++) filenames[i][j] = (char)((char)mxcharp[i+j*m]);
00957       /* Make last char blank */
00958       filenames[i][n] = ' ';
00959       /* Backing up from end, get rid of blanks */
00960       for (j=n; j>0; j++) 
00961         if (filenames[i][j] == ' ') 
00962           filenames[i][j] = '\0';
00963         else break;
00964       DBGPRINTF("Packed input/inout file %d **%s** \n", i, filenames[i]);
00965     }
00966     argp->data = filenames;
00967     break;
00968     
00969       default:
00970     ERRPRINTF("Unknown objectype\n");   
00971     goto error_input;
00972     break;
00973     
00974       } /* argp->objectype */
00975       
00976       /* increment irhs index */
00977       irhs++;           
00978       /* increment irhs index */
00979       
00980       
00981       /* ***************************** */
00982     } else if (argp->inout==GS_WORKSPACE) {
00983       /* ***************************** */
00984 
00985       /* Do nothing, just move to next argument */
00986       argp->data = NULL;
00987 
00988       /* ***************************** */
00989     } else if (argp->inout==GS_OUT) {
00990       /* ***************************** */
00991       
00992       DBGPRINTF("Handling GS_OUT: %s %s %s %s \n", argp->name, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype]); fflush(0);
00993       
00994       /* Setup output non scalars, allocating space in GridSolve as needed */
00995       switch (argp->objecttype) {
00996     
00997       case GS_SCALAR:
00998     break;
00999     
01000       case GS_SPARSEMATRIX:
01001     DBGPRINTF("Attaching output sparse attributes \n");
01002     /* Lookup which gridsolve arg should hold indices and pointers */
01003     arg_indices = gs_arg_lookup_by_name(argp->sparse_attr.indices, problem->arglist);
01004     arg_pointer = gs_arg_lookup_by_name(argp->sparse_attr.pointer, problem->arglist);
01005     if (arg_indices==NULL || arg_pointer==NULL) goto error_output;
01006     /* Allocate space for indices, pointer and data */
01007     if ((arg_indices->data = matlab_gs_calloc(argp->sparse_attr.nnz, sizeof(int))) == NULL) goto error_memory;
01008     if ((arg_pointer->data = matlab_gs_calloc(argp->sparse_attr.cols_val_saved+1, sizeof(int))) == NULL) goto error_memory;
01009     if ((argp->data = matlab_gs_calloc(argp->sparse_attr.nnz, gs_datatype_sizeof(argp->datatype))) == NULL) goto error_memory;
01010     break;
01011     
01012       case GS_VECTOR:
01013       case GS_MATRIX:
01014     /* Compute number of elements to be allocated and allocate space for output */
01015     nelements = argp->rows*argp->cols;
01016     if ((argp->data = (double*)matlab_gs_calloc(nelements, gs_datatype_sizeof(argp->datatype)))==NULL) goto error_memory;
01017     break;
01018     
01019       case GS_FILE: 
01020     /* In matlab, the client does not pass in the GS_FILE name for
01021        an output argument, so use problem_arg_XXXXXX as the
01022        filename */
01023     snprintf(buffer, 256, "./%s_%s_XXXXXX", problem->name, argp->name);
01024     tmpname = mktemp(buffer);
01025     if (tmpname == NULL) { ERRPRINTF("Could not make tmp file\n"); goto error_output; }
01026     fp = fopen(tmpname, "w");
01027     if (fp == NULL) { ERRPRINTF("Could not make tmp file\n"); goto error_output; }
01028     fclose(fp);
01029     argp->data = strdup(buffer);
01030     if (argp->data==NULL) goto error_memory;
01031     DBGPRINTF("output file name %s\n", (char*)argp->data);
01032     break;
01033     
01034       case GS_PACKEDFILE:
01035     /* In matlab, the client does not pass in the GS_FILE name for
01036        an output argument, so use problem_arg_%d_XXXXXX as the
01037        filename */
01038     nelements = argp->rows*argp->cols;
01039     filenames = (char **)matlab_gs_calloc(MAX_FILES_IN_PACKEDFILE, sizeof(char*));
01040     /* TODO Free memory */
01041     if (filenames==NULL) goto error_memory;
01042     for (i=0; i<nelements; i++) {
01043       snprintf(buffer, 256, "./%s_%s_%d_XXXXXX", problem->name, argp->name, i+1);
01044       tmpname = mktemp(buffer);
01045       if (tmpname == NULL) { ERRPRINTF("Could not make tmp file\n"); goto error_output; }
01046       fp = fopen(tmpname, "w");
01047       if (fp == NULL) { ERRPRINTF("Could not make tmp file\n"); goto error_output; }
01048       fclose(fp);
01049       filenames[i] = strdup(buffer);
01050       if (filenames[i]==NULL) goto error_memory;
01051       DBGPRINTF("packed output filename %d %s\n", i, filenames[i]);
01052     }
01053     argp->data = filenames;
01054     break;
01055     
01056       default:
01057     DBGPRINTF("Object type not handled %d\n", argp->objecttype);
01058     break;
01059       } /* switch(argp->objectype) */
01060       
01061 
01062       /* ***************************** */
01063     } else if (argp->inout==GS_VAROUT) {
01064       /* ***************************** */
01065       
01066       switch (argp->objecttype) {
01067       case GS_SPARSEMATRIX: 
01068     ERRPRINTF("VAROUT SPARSEMATIX is not supported \n");
01069     return -1;
01070     break;
01071       default:
01072     /* If VAROUT, create a (char **) and use it as the arguments data space  */
01073     argp->data = (char **)matlab_gs_calloc(1,sizeof(char*));
01074     if (!argp->data) goto error_memory;
01075     break;
01076       } /* switch objectype */
01077       
01078       /* ***************************** */
01079     } else {
01080       /* ***************************** */
01081 
01082       ERRPRINTF("Unknown inout mode\n");
01083       return -1;
01084     }
01085     
01086   }
01087   DBGPRINTF("Done\n");
01088   return 0;
01089 
01090 
01091  error_scalar:
01092   ERRPRINTF("GridSolve: Could not transfer scalar arguments");
01093   return -1;
01094 
01095  error_input:
01096   ERRPRINTF("GridSolve: Data conversion problem between matlab input data (arg %d type %s) and GridSolve (arg %s type %s)", irhs-2, mxGetClassName(prhs[irhs]), argp->name, gs_c_datatype[argp->datatype]);
01097   return -1;
01098 
01099  error_output:
01100   ERRPRINTF("GridSolve: Data preparation problem for GridSolve output (arg %s type %s)", argp->name, gs_c_datatype[argp->datatype]);
01101   return -1;
01102 
01103  error_memory:
01104   ERRPRINTF("GridSolve: Memory allocation problem for GridSolve argument (arg %s type %s)", argp->name, gs_c_datatype[argp->datatype]);
01105   return -1;
01106 
01107  error_unknown:
01108   ERRPRINTF("GridSolve: Unknown error occured");
01109   return -1;
01110 
01111 }
01112 
01113 
01126 int 
01127 matlab_gs_get_output(gs_problem_t* problem, int nlhs, mxArray *plhs[],
01128               int nrhs, const mxArray *prhs[])
01129 {
01130   gs_argument_t *argp = NULL;
01131   gs_argument_t *arg_pointer = NULL;  
01132   gs_argument_t *arg_indices = NULL;  
01133   gs_argument_t *row_arg, *col_arg;
01134   int M, N, ilhs;
01135   
01136   DBGPRINTF("Transferring output arguments to matlab\n");
01137   for (argp=problem->arglist,ilhs=0; argp!=NULL && ilhs<nlhs; argp=argp->next) {
01138     
01139     if (argp->inout==GS_OUT || argp->inout==GS_INOUT || argp->inout==GS_VAROUT) {
01140 
01141       /* skip this if it's a sparse matrix attribute, since the sparse
01142      attributes are not sent as matlab parameter, they will be
01143      handled when the actual sparse argument is handled  */
01144       if (matlab_gs_arg_is_sparse_index_or_pointer(argp, problem)) {
01145     DBGPRINTF("Skipping sparse attribute: %s \n", argp->name);
01146     continue;
01147       }
01148 
01149       switch (argp->objecttype) {
01150     
01151       case GS_FILE:
01152     DBGPRINTF("Transfer GS_FILE %s using filename %s\n", argp->name, (char *)argp->data);
01153     plhs[ilhs] = mxCreateString((char*)argp->data);
01154     break;
01155     
01156       case GS_PACKEDFILE:
01157     DBGPRINTF("Transfer GS_PACKEDFILE\n");
01158     plhs[ilhs] = mxCreateCharMatrixFromStrings(argp->rows, (const char**)argp->data);
01159     break;
01160 
01161       case GS_SPARSEMATRIX:
01162     DBGPRINTF("Transfer gs (%s [%d,%d] %s %s %s) to mx (ilhs %d) \n", argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype], ilhs); fflush(0);
01163     if (argp->inout==GS_VAROUT) { ERRPRINTF("VAROUT SPARSEMATRIX not supported\n");  goto error_output; }
01164     /* Lookup which gridsolve arg should hold indices and pointers */
01165     arg_indices = gs_arg_lookup_by_name(argp->sparse_attr.indices, problem->arglist);
01166     arg_pointer = gs_arg_lookup_by_name(argp->sparse_attr.pointer, problem->arglist);
01167     M = argp->sparse_attr.rows_val_saved;
01168     N = argp->sparse_attr.cols_val_saved;
01169     /* Copy CSR to CSC */
01170     if ((plhs[ilhs] = matlab_gs_create_matrix(argp->objecttype, M, N, argp->sparse_attr.nnz, argp->datatype)) == NULL) goto error_output;
01171     if (matlab_gs_sparse_csr_to_matcsc(M, N, argp->sparse_attr.nnz, argp->datatype, argp->data, arg_indices->data, arg_pointer->data, mxGetPr(plhs[ilhs]), mxGetPi(plhs[ilhs]), mxGetIr(plhs[ilhs]), mxGetJc(plhs[ilhs])) == -1) goto error_output;
01172     break;
01173     
01174       case GS_SCALAR:
01175     DBGPRINTF("Transfer gs (%s [%d,%d] %s %s %s) to mx (ilhs %d) \n", argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype], ilhs); fflush(0);
01176     plhs[ilhs] = transfer_gs_to_mx(argp);
01177     if (plhs[ilhs] == NULL) return -1;
01178     break;
01179       
01180       case GS_VECTOR:
01181       case GS_MATRIX:
01182     DBGPRINTF("Transfer gs (%s [%d,%d] %s %s %s) to mx (ilhs %d) \n", argp->name, argp->rows, argp->cols, gs_inout[argp->inout], gs_objecttype[argp->objecttype], gs_c_datatype[argp->datatype], ilhs); fflush(0);
01183     if (argp->inout==GS_VAROUT) { 
01184       /* If VAROUT, find row and col value holders, get the value
01185        * and transfer those values to the argument */
01186       argp->rows = 1;
01187       argp->cols = 1;
01188       row_arg = gs_arg_lookup_by_name(argp->rowexp, problem->arglist);
01189       if (row_arg==NULL || row_arg->inout!=GS_OUT || row_arg->datatype!=GS_INT 
01190           || row_arg->objecttype!=GS_SCALAR) goto error_output;
01191       argp->rows = row_arg->scalar_val.int_val;
01192       col_arg = gs_arg_lookup_by_name(argp->colexp, problem->arglist);
01193       if (col_arg==NULL) 
01194         argp->cols = 1;
01195       else {        /* if col_arg != NULL */
01196         if ( col_arg->inout!=GS_OUT || col_arg->datatype!=GS_INT
01197          || col_arg->objecttype!=GS_SCALAR) goto error_output;
01198         argp->cols = col_arg->scalar_val.int_val;
01199       }
01200     }
01201     plhs[ilhs] = transfer_gs_to_mx(argp);
01202     if (plhs[ilhs] == NULL) return -1;
01203     break;
01204     
01205       default:
01206     DBGPRINTF("Unknown object type name %s type %d\n", argp->name, argp->objecttype);
01207     goto error_output;
01208     break;
01209       }
01210 
01211       /* Increment index for output arguments */
01212       ilhs++;           
01213 
01214     }
01215   }
01216   
01217   return 0;
01218 
01219  error_output:
01220   ERRPRINTF("GridSolve: Could not transfer output to Matlab (output arg %d) for GridSolve output (arg %s type %s)", ilhs, argp->name, gs_c_datatype[argp->datatype]);
01221   return -1;
01222 
01223 }
01224 
01225 
01238 int 
01239 matlab_gs_free_input_args(gs_problem_t* problem, int nlhs, mxArray *plhs[],int nrhs, const mxArray *prhs[])
01240 {
01241   int irhs, i;
01242   gs_argument_t* argp = NULL;
01243   if (problem==NULL) return 0;
01244 
01245   DBGPRINTF("Freeing input args\n");
01246 
01247   for (argp=problem->arglist,irhs=2; argp!=NULL && irhs<=nrhs; argp=argp->next) {
01248     if (argp->inout==GS_IN) {
01249       switch(argp->objecttype) {
01250       case GS_SCALAR:
01251     break;
01252       case GS_VECTOR:
01253       case GS_MATRIX:
01254     if (arg_data_same_as_mx_dat(prhs[irhs], argp)) 
01255       argp->data = NULL;
01256     else 
01257       argp->data = matlab_gs_free(argp->data);
01258     break;
01259       case GS_SPARSEMATRIX:
01260     argp->data = matlab_gs_free(argp->data);
01261     break;
01262       case GS_FILE:
01263     argp->data = matlab_gs_free(argp->data);    
01264     break;
01265       case GS_PACKEDFILE:
01266     for (i=0; i<argp->rows*argp->cols; i++) 
01267       matlab_gs_free(((char **)argp->data)[i]);
01268     argp->data = matlab_gs_free(argp->data);
01269     break;
01270       default:
01271     break;
01272       }
01273 
01274       /* increment irhs index if not sparse related, since sparse
01275        * index/pointer are not passed in in Matlab */
01276       if (!matlab_gs_arg_is_sparse_index_or_pointer(argp, problem))
01277     irhs++;         
01278     }
01279   }
01280 
01281   DBGPRINTF("Finished\n");
01282   return 0;
01283 
01284 }
01285 
01286 
01298 int 
01299 matlab_gs_free_args(gs_problem_t* problem, int nlhs, mxArray *plhs[],
01300             int nrhs, const mxArray *prhs[])
01301 {
01302   int ilhs = 0;
01303   int i;
01304   gs_argument_t* argp = NULL;
01305 
01306   if (problem==NULL) return 0;
01307 
01308   DBGPRINTF("Freeing all output args\n");
01309 
01310   for (argp=problem->arglist,ilhs=0; argp!=NULL && ilhs<nlhs; argp=argp->next) {
01311     if (argp->inout==GS_INOUT || argp->inout==GS_OUT || argp->inout==GS_VAROUT ) {
01312       switch(argp->objecttype) {
01313       case GS_SCALAR:
01314     break;
01315       case GS_VECTOR:
01316       case GS_MATRIX:
01317     argp->data = matlab_gs_free(argp->data);
01318     break;
01319       case GS_SPARSEMATRIX:
01320     argp->data = matlab_gs_free(argp->data);
01321     break;
01322       case GS_FILE:
01323     argp->data = matlab_gs_free(argp->data);
01324     break;
01325       case GS_PACKEDFILE:
01326     for (i=0; i<argp->rows*argp->cols; i++) 
01327       matlab_gs_free(((char **)argp->data)[i]);
01328     argp->data = matlab_gs_free(argp->data);
01329     break;
01330       default:
01331     break;
01332     
01333       }
01334     }
01335   }
01336   
01337   DBGPRINTF("Finished\n");
01338   return 0;
01339 
01340 }
01341