Actual source code: mathematica.c

      1: #define PETSC_DLL
      2: /*
      3:         Written by Matt Knepley, knepley@cs.purdue.edu 7/23/97
      4:         Major overhall for interactivity               11/14/97
      5:         Reorganized                                    11/8/98
      6: */
     7:  #include src/sys/src/viewer/viewerimpl.h
     8:  #include src/ksp/pc/pcimpl.h
     9:  #include src/mat/impls/aij/seq/aij.h
     10:  #include mathematica.h
     11: #include "petscfix.h"

     13: #if defined (PETSC_HAVE__SNPRINTF) && !defined(PETSC_HAVE_SNPRINTF)
     14: #define snprintf _snprintf
     15: #endif

     17: PetscViewer  PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE = PETSC_NULL;
     18: static void *mathematicaEnv                   = PETSC_NULL;

     22: /*@C
     23:   PetscViewerMathematicaInitializePackage - This function initializes everything in the Petsc interface to Mathematica. It is
     24:   called from PetscDLLibraryRegister() when using dynamic libraries, and on the call to PetscInitialize()
     25:   when using static libraries.

     27:   Input Parameter:
     28:   path - The dynamic library path, or PETSC_NULL

     30:   Level: developer

     32: .keywords: Petsc, initialize, package, PLAPACK
     33: .seealso: PetscInitializePackage(), PetscInitialize()
     34: @*/
     35: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaInitializePackage(char *path)
     36: {
     37:   static PetscTruth initialized = PETSC_FALSE;

     40:   if (initialized) return(0);
     41:   initialized = PETSC_TRUE;
     42:   mathematicaEnv = (void*) MLInitialize(0);
     43:   return(0);
     44: }

     48: /*@C
     49:   PetscViewerMathematicaDestroyPackage - This function destroys everything in the Petsc interface to Mathematica. It is
     50:   called from PetscFinalize().

     52:   Level: developer

     54: .keywords: Petsc, destroy, package, mathematica
     55: .seealso: PetscFinalize()
     56: @*/
     57: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaFinalizePackage(void)
     58: {
     60:   if (mathematicaEnv) MLDeinitialize((MLEnvironment) mathematicaEnv);
     61:   return(0);
     62: }

     66: PetscErrorCode PetscViewerInitializeMathematicaWorld_Private()
     67: {

     71:   if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) return(0);
     72:   PetscViewerMathematicaOpen(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_NULL, PETSC_NULL, &PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE);
     73:   return(0);
     74: }

     78: static PetscErrorCode PetscViewerDestroy_Mathematica(PetscViewer viewer)
     79: {
     80:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
     81:   PetscErrorCode          ierr;

     84:   MLClose(vmath->link);
     85:   if (vmath->linkname) {
     86:     PetscFree(vmath->linkname);
     87:   }
     88:   if (vmath->linkhost) {
     89:     PetscFree(vmath->linkhost);
     90:   }
     91:   PetscFree(vmath);
     92:   return(0);
     93: }

     97: PetscErrorCode PetscViewerDestroyMathematica_Private(void)
     98: {

    102:   if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) {
    103:     PetscViewerDestroy(PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE);
    104:   }
    105:   return(0);
    106: }

    110: PetscErrorCode PetscViewerMathematicaSetupConnection_Private(PetscViewer v)
    111: {
    112:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) v->data;
    113: #ifdef MATHEMATICA_3_0
    114:   int                     argc = 6;
    115:   char                    *argv[6];
    116: #else
    117:   int                     argc = 5;
    118:   char                    *argv[5];
    119: #endif
    120:   char                    hostname[256];
    121:   long                    lerr;
    122:   PetscErrorCode          ierr;

    125:   /* Link name */
    126:   argv[0] = "-linkname";
    127:   if (!vmath->linkname) {
    128:     argv[1] = "math -mathlink";
    129:   } else {
    130:     argv[1] = vmath->linkname;
    131:   }

    133:   /* Link host */
    134:   argv[2] = "-linkhost";
    135:   if (!vmath->linkhost) {
    136:     PetscGetHostName(hostname, 255);
    137:     argv[3] = hostname;
    138:   } else {
    139:     argv[3] = vmath->linkhost;
    140:   }

    142:   /* Link mode */
    143: #ifdef MATHEMATICA_3_0
    144:   argv[4] = "-linkmode";
    145:   switch(vmath->linkmode) {
    146:   case MATHEMATICA_LINK_CREATE:
    147:     argv[5] = "Create";
    148:     break;
    149:   case MATHEMATICA_LINK_CONNECT:
    150:     argv[5] = "Connect";
    151:     break;
    152:   case MATHEMATICA_LINK_LAUNCH:
    153:     argv[5] = "Launch";
    154:     break;
    155:   }
    156: #else
    157:   switch(vmath->linkmode) {
    158:   case MATHEMATICA_LINK_CREATE:
    159:     argv[4] = "-linkcreate";
    160:     break;
    161:   case MATHEMATICA_LINK_CONNECT:
    162:     argv[4] = "-linkconnect";
    163:     break;
    164:   case MATHEMATICA_LINK_LAUNCH:
    165:     argv[4] = "-linklaunch";
    166:     break;
    167:   }
    168: #endif
    169:   vmath->link = MLOpenInEnv(mathematicaEnv, argc, argv, &lerr);
    170: #endif
    171:   return(0);
    172: }

    177: PetscErrorCode PETSC_DLLEXPORT PetscViewerCreate_Mathematica(PetscViewer v)
    178: {
    179:   PetscViewer_Mathematica *vmath;
    180:   PetscErrorCode          ierr;

    184:   PetscNew(PetscViewer_Mathematica, &vmath);
    185:   v->data         = (void*) vmath;
    186:   v->ops->destroy = PetscViewerDestroy_Mathematica;
    187:   v->ops->flush   = 0;
    188:   PetscStrallocpy(PETSC_VIEWER_MATHEMATICA, &v->type_name);

    190:   vmath->linkname         = PETSC_NULL;
    191:   vmath->linkhost         = PETSC_NULL;
    192:   vmath->linkmode         = MATHEMATICA_LINK_CONNECT;
    193:   vmath->graphicsType     = GRAPHICS_MOTIF;
    194:   vmath->plotType         = MATHEMATICA_TRIANGULATION_PLOT;
    195:   vmath->objName          = PETSC_NULL;

    197:   PetscViewerMathematicaSetFromOptions(v);
    198:   PetscViewerMathematicaSetupConnection_Private(v);
    199:   return(0);
    200: }

    205: PetscErrorCode PetscViewerMathematicaParseLinkMode_Private(char *modename, LinkMode *mode) {
    206:   PetscTruth     isCreate, isConnect, isLaunch;

    210:   PetscStrcasecmp(modename, "Create",  &isCreate);
    211:   PetscStrcasecmp(modename, "Connect", &isConnect);
    212:   PetscStrcasecmp(modename, "Launch",  &isLaunch);
    213:   if (isCreate) {
    214:     *mode = MATHEMATICA_LINK_CREATE;
    215:   } else if (isConnect) {
    216:     *mode = MATHEMATICA_LINK_CONNECT;
    217:   } else if (isLaunch) {
    218:     *mode = MATHEMATICA_LINK_LAUNCH;
    219:   } else {
    220:     SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid Mathematica link mode: %s", modename);
    221:   }
    222:   return(0);
    223: }

    227: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetFromOptions(PetscViewer v)
    228: {
    229:   PetscViewer_Mathematica  *vmath = (PetscViewer_Mathematica *) v->data;
    230:   char                     linkname[256];
    231:   char                     modename[256];
    232:   char                     hostname[256];
    233:   char                     type[256];
    234:   PetscInt                 numPorts;
    235:   PetscInt                 *ports;
    236:   PetscInt                 numHosts;
    237:   int                      h;
    238:   char                     **hosts;
    239:   PetscMPIInt              size, rank;
    240:   PetscTruth               opt;
    241:   PetscErrorCode           ierr;

    244:   MPI_Comm_size(v->comm, &size);
    245:   MPI_Comm_rank(v->comm, &rank);

    247:   /* Get link name */
    248:   PetscOptionsGetString("viewer_", "-math_linkname", linkname, 255, &opt);
    249:   if (opt) {
    250:     PetscViewerMathematicaSetLinkName(v, linkname);
    251:   }
    252:   /* Get link port */
    253:   numPorts = size;
    254:   PetscMalloc(size * sizeof(int), &ports);
    255:   PetscOptionsGetIntArray("viewer_", "-math_linkport", ports, &numPorts, &opt);
    256:   if (opt) {
    257:     if (numPorts > rank) {
    258:       snprintf(linkname, 255, "%6d", ports[rank]);
    259:     } else {
    260:       snprintf(linkname, 255, "%6d", ports[0]);
    261:     }
    262:     PetscViewerMathematicaSetLinkName(v, linkname);
    263:   }
    264:   PetscFree(ports);
    265:   /* Get link host */
    266:   numHosts = size;
    267:   PetscMalloc(size * sizeof(char *), &hosts);
    268:   PetscOptionsGetStringArray("viewer_", "-math_linkhost", hosts, &numHosts, &opt);
    269:   if (opt) {
    270:     if (numHosts > rank) {
    271:       PetscStrncpy(hostname, hosts[rank], 255);
    272:     } else {
    273:       PetscStrncpy(hostname, hosts[0], 255);
    274:     }
    275:     PetscViewerMathematicaSetLinkHost(v, hostname);
    276:   }
    277:   for(h = 0; h < numHosts; h++) {
    278:     PetscFree(hosts[h]);
    279:   }
    280:   PetscFree(hosts);
    281:   /* Get link mode */
    282:   PetscOptionsGetString("viewer_", "-math_linkmode", modename, 255, &opt);
    283:   if (opt) {
    284:     LinkMode mode;

    286:     PetscViewerMathematicaParseLinkMode_Private(modename, &mode);
    287:     PetscViewerMathematicaSetLinkMode(v, mode);
    288:   }
    289:   /* Get graphics type */
    290:   PetscOptionsGetString("viewer_", "-math_graphics", type, 255, &opt);
    291:   if (opt) {
    292:     PetscTruth isMotif, isPS, isPSFile;

    294:     PetscStrcasecmp(type, "Motif",  &isMotif);
    295:     PetscStrcasecmp(type, "PS",     &isPS);
    296:     PetscStrcasecmp(type, "PSFile", &isPSFile);
    297:     if (isMotif) {
    298:       vmath->graphicsType = GRAPHICS_MOTIF;
    299:     } else if (isPS) {
    300:       vmath->graphicsType = GRAPHICS_PS_STDOUT;
    301:     } else if (isPSFile) {
    302:       vmath->graphicsType = GRAPHICS_PS_FILE;
    303:     }
    304:   }
    305:   /* Get plot type */
    306:   PetscOptionsGetString("viewer_", "-math_type", type, 255, &opt);
    307:   if (opt) {
    308:     PetscTruth isTri, isVecTri, isVec, isSurface;

    310:     PetscStrcasecmp(type, "Triangulation",       &isTri);
    311:     PetscStrcasecmp(type, "VectorTriangulation", &isVecTri);
    312:     PetscStrcasecmp(type, "Vector",              &isVec);
    313:     PetscStrcasecmp(type, "Surface",             &isSurface);
    314:     if (isTri) {
    315:       vmath->plotType     = MATHEMATICA_TRIANGULATION_PLOT;
    316:     } else if (isVecTri) {
    317:       vmath->plotType     = MATHEMATICA_VECTOR_TRIANGULATION_PLOT;
    318:     } else if (isVec) {
    319:       vmath->plotType     = MATHEMATICA_VECTOR_PLOT;
    320:     } else if (isSurface) {
    321:       vmath->plotType     = MATHEMATICA_SURFACE_PLOT;
    322:     }
    323:   }
    324:   return(0);
    325: }

    329: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetLinkName(PetscViewer v, const char *name) {
    330:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) v->data;
    331:   PetscErrorCode          ierr;

    336:   PetscStrallocpy(name, &vmath->linkname);
    337:   return(0);
    338: }

    342: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetLinkPort(PetscViewer v, int port) {
    343:   char           name[16];

    347:   snprintf(name, 16, "%6d", port);
    348:   PetscViewerMathematicaSetLinkName(v, name);
    349:   return(0);
    350: }

    354: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetLinkHost(PetscViewer v, const char *host) {
    355:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) v->data;
    356:   PetscErrorCode          ierr;

    361:   PetscStrallocpy(host, &vmath->linkhost);
    362:   return(0);
    363: }

    367: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetLinkMode(PetscViewer v, LinkMode mode) {
    368:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) v->data;

    371:   vmath->linkmode = mode;
    372:   return(0);
    373: }

    375: /*----------------------------------------- Public Functions --------------------------------------------------------*/
    378: /*@C
    379:   PetscViewerMathematicaOpen - Communicates with Mathemtica using MathLink.

    381:   Collective on comm

    383:   Input Parameters:
    384: + comm    - The MPI communicator
    385: . port    - [optional] The port to connect on, or PETSC_DECIDE
    386: . machine - [optional] The machine to run Mathematica on, or PETSC_NULL
    387: - mode    - [optional] The connection mode, or PETSC_NULL

    389:   Output Parameter:
    390: . viewer  - The Mathematica viewer

    392:   Level: intermediate

    394:   Notes:
    395:   Most users should employ the following commands to access the
    396:   Mathematica viewers
    397: $
    398: $    PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
    399: $    MatView(Mat matrix, PetscViewer viewer)
    400: $
    401: $                or
    402: $
    403: $    PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
    404: $    VecView(Vec vector, PetscViewer viewer)

    406:    Options Database Keys:
    407: $    -viewer_math_linkhost <machine> - The host machine for the kernel
    408: $    -viewer_math_linkname <name>    - The full link name for the connection
    409: $    -viewer_math_linkport <port>    - The port for the connection
    410: $    -viewer_math_mode <mode>        - The mode, e.g. Launch, Connect
    411: $    -viewer_math_type <type>        - The plot type, e.g. Triangulation, Vector
    412: $    -viewer_math_graphics <output>  - The output type, e.g. Motif, PS, PSFile

    414: .keywords: PetscViewer, Mathematica, open

    416: .seealso: MatView(), VecView()
    417: @*/
    418: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaOpen(MPI_Comm comm, int port, const char machine[], const char mode[], PetscViewer *v)
    419: {

    423:   PetscViewerCreate(comm, v);
    424: #if 0
    425:   LinkMode linkmode;
    426:   PetscViewerMathematicaSetLinkPort(*v, port);
    427:   PetscViewerMathematicaSetLinkHost(*v, machine);
    428:   PetscViewerMathematicaParseLinkMode_Private(mode, &linkmode);
    429:   PetscViewerMathematicaSetLinkMode(*v, linkmode);
    430: #endif
    431:   PetscViewerSetType(*v, PETSC_VIEWER_MATHEMATICA);
    432:   return(0);
    433: }

    437: /*@C
    438:   PetscViewerMathematicaGetLink - Returns the link to Mathematica

    440:   Input Parameters:
    441: . viewer - The Mathematica viewer
    442: . link   - The link to Mathematica

    444:   Level: intermediate

    446: .keywords PetscViewer, Mathematica, link
    447: .seealso PetscViewerMathematicaOpen()
    448: @*/
    449: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaGetLink(PetscViewer viewer, MLINK *link)
    450: {
    451:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;

    455:   *link = vmath->link;
    456:   return(0);
    457: }

    461: /*@C
    462:   PetscViewerMathematicaSkipPackets - Discard packets sent by Mathematica until a certain packet type is received

    464:   Input Parameters:
    465: . viewer - The Mathematica viewer
    466: . type   - The packet type to search for, e.g RETURNPKT

    468:   Level: advanced

    470: .keywords PetscViewer, Mathematica, packets
    471: .seealso PetscViewerMathematicaSetName(), PetscViewerMathematicaGetVector()
    472: @*/
    473: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSkipPackets(PetscViewer viewer, int type)
    474: {
    475:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
    476:   MLINK                   link  = vmath->link; /* The link to Mathematica */
    477:   int                     pkt;                 /* The packet type */

    480:   while((pkt = MLNextPacket(link)) && (pkt != type))
    481:     MLNewPacket(link);
    482:   if (!pkt) {
    483:     MLClearError(link);
    484:     SETERRQ(PETSC_ERR_LIB, (char *) MLErrorMessage(link));
    485:   }
    486:   return(0);
    487: }

    491: /*@C
    492:   PetscViewerMathematicaGetName - Retrieve the default name for objects communicated to Mathematica

    494:   Input Parameter:
    495: . viewer - The Mathematica viewer

    497:   Output Parameter:
    498: . name   - The name for new objects created in Mathematica

    500:   Level: intermediate

    502: .keywords PetscViewer, Mathematica, name
    503: .seealso PetscViewerMathematicaSetName(), PetscViewerMathematicaClearName()
    504: @*/
    505: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaGetName(PetscViewer viewer, const char **name)
    506: {
    507:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;

    512:   *name = vmath->objName;
    513:   return(0);
    514: }

    518: /*@C
    519:   PetscViewerMathematicaSetName - Override the default name for objects communicated to Mathematica

    521:   Input Parameters:
    522: . viewer - The Mathematica viewer
    523: . name   - The name for new objects created in Mathematica

    525:   Level: intermediate

    527: .keywords PetscViewer, Mathematica, name
    528: .seealso PetscViewerMathematicaSetName(), PetscViewerMathematicaClearName()
    529: @*/
    530: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaSetName(PetscViewer viewer, const char name[])
    531: {
    532:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;

    537:   vmath->objName = name;
    538:   return(0);
    539: }

    543: /*@C
    544:   PetscViewerMathematicaClearName - Use the default name for objects communicated to Mathematica

    546:   Input Parameter:
    547: . viewer - The Mathematica viewer

    549:   Level: intermediate

    551: .keywords PetscViewer, Mathematica, name
    552: .seealso PetscViewerMathematicaGetName(), PetscViewerMathematicaSetName()
    553: @*/
    554: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaClearName(PetscViewer viewer)
    555: {
    556:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;

    560:   vmath->objName = PETSC_NULL;
    561:   return(0);
    562: }

    566: /*@C
    567:   PetscViewerMathematicaGetVector - Retrieve a vector from Mathematica

    569:   Input Parameter:
    570: . viewer - The Mathematica viewer

    572:   Output Parameter:
    573: . v      - The vector

    575:   Level: intermediate

    577: .keywords PetscViewer, Mathematica, vector
    578: .seealso VecView(), PetscViewerMathematicaPutVector()
    579: @*/
    580: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaGetVector(PetscViewer viewer, Vec v) {
    581:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
    582:   MLINK                   link;   /* The link to Mathematica */
    583:   char                    *name;
    584:   PetscScalar             *mArray,*array;
    585:   long                    mSize;
    586:   int                     n;
    587:   PetscErrorCode          ierr;

    593:   /* Determine the object name */
    594:   if (!vmath->objName) {
    595:     name = "vec";
    596:   } else {
    597:     name = (char *) vmath->objName;
    598:   }

    600:   link = vmath->link;
    601:   VecGetLocalSize(v, &n);
    602:   VecGetArray(v, &array);
    603:   MLPutFunction(link, "EvaluatePacket", 1);
    604:     MLPutSymbol(link, name);
    605:   MLEndPacket(link);
    606:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    607:   MLGetRealList(link, &mArray, &mSize);
    608:   if (n != mSize) SETERRQ2(PETSC_ERR_ARG_WRONG, "Incompatible vector sizes %d %d",n,mSize);
    609:   PetscMemcpy(array, mArray, mSize * sizeof(double));
    610:   MLDisownRealList(link, mArray, mSize);
    611:   VecRestoreArray(v, &array);

    613:   return(0);
    614: }

    618: /*@C
    619:   PetscViewerMathematicaPutVector - Send a vector to Mathematica

    621:   Input Parameters:
    622: + viewer - The Mathematica viewer
    623: - v      - The vector

    625:   Level: intermediate

    627: .keywords PetscViewer, Mathematica, vector
    628: .seealso VecView(), PetscViewerMathematicaGetVector()
    629: @*/
    630: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaPutVector(PetscViewer viewer, Vec v)
    631: {
    632:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
    633:   MLINK                   link  = vmath->link; /* The link to Mathematica */
    634:   char                    *name;
    635:   PetscScalar             *array;
    636:   int                     n;
    637:   PetscErrorCode          ierr;

    640:   /* Determine the object name */
    641:   if (!vmath->objName) {
    642:     name = "vec";
    643:   } else {
    644:     name = (char *) vmath->objName;
    645:   }
    646:   VecGetLocalSize(v, &n);
    647:   VecGetArray(v, &array);

    649:   /* Send the Vector object */
    650:   MLPutFunction(link, "EvaluatePacket", 1);
    651:     MLPutFunction(link, "Set", 2);
    652:       MLPutSymbol(link, name);
    653:       MLPutRealList(link, array, n);
    654:   MLEndPacket(link);
    655:   /* Skip packets until ReturnPacket */
    656:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    657:   /* Skip ReturnPacket */
    658:   MLNewPacket(link);

    660:   VecRestoreArray(v, &array);
    661:   return(0);
    662: }

    664: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaPutMatrix(PetscViewer viewer, int m, int n, PetscReal *a)
    665: {
    666:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
    667:   MLINK                   link  = vmath->link; /* The link to Mathematica */
    668:   char                    *name;
    669:   PetscErrorCode          ierr;

    672:   /* Determine the object name */
    673:   if (!vmath->objName) {
    674:     name = "mat";
    675:   } else {
    676:     name = (char *) vmath->objName;
    677:   }

    679:   /* Send the dense matrix object */
    680:   MLPutFunction(link, "EvaluatePacket", 1);
    681:     MLPutFunction(link, "Set", 2);
    682:       MLPutSymbol(link, name);
    683:       MLPutFunction(link, "Transpose", 1);
    684:         MLPutFunction(link, "Partition", 2);
    685:           MLPutRealList(link, a, m*n);
    686:           MLPutInteger(link, m);
    687:   MLEndPacket(link);
    688:   /* Skip packets until ReturnPacket */
    689:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    690:   /* Skip ReturnPacket */
    691:   MLNewPacket(link);

    693:   return(0);
    694: }

    696: PetscErrorCode PETSC_DLLEXPORT PetscViewerMathematicaPutCSRMatrix(PetscViewer viewer, int m, int n, int *i, int *j, PetscReal *a)
    697: {
    698:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *) viewer->data;
    699:   MLINK                   link  = vmath->link; /* The link to Mathematica */
    700:   const char              *symbol;
    701:   char                    *name;
    702:   PetscTruth              match;
    703:   PetscErrorCode          ierr;

    706:   /* Determine the object name */
    707:   if (!vmath->objName) {
    708:     name = "mat";
    709:   } else {
    710:     name = (char *) vmath->objName;
    711:   }

    713:   /* Make sure Mathematica recognizes sparse matrices */
    714:   MLPutFunction(link, "EvaluatePacket", 1);
    715:     MLPutFunction(link, "Needs", 1);
    716:       MLPutString(link, "LinearAlgebra`CSRMatrix`");
    717:   MLEndPacket(link);
    718:   /* Skip packets until ReturnPacket */
    719:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    720:   /* Skip ReturnPacket */
    721:   MLNewPacket(link);

    723:   /* Send the CSRMatrix object */
    724:   MLPutFunction(link, "EvaluatePacket", 1);
    725:     MLPutFunction(link, "Set", 2);
    726:       MLPutSymbol(link, name);
    727:       MLPutFunction(link, "CSRMatrix", 5);
    728:         MLPutInteger(link, m);
    729:         MLPutInteger(link, n);
    730:         MLPutFunction(link, "Plus", 2);
    731:           MLPutIntegerList(link, i, m+1);
    732:           MLPutInteger(link, 1);
    733:         MLPutFunction(link, "Plus", 2);
    734:           MLPutIntegerList(link, j, i[m]);
    735:           MLPutInteger(link, 1);
    736:         MLPutRealList(link, a, i[m]);
    737:   MLEndPacket(link);
    738:   /* Skip packets until ReturnPacket */
    739:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    740:   /* Skip ReturnPacket */
    741:   MLNewPacket(link);

    743:   /* Check that matrix is valid */
    744:   MLPutFunction(link, "EvaluatePacket", 1);
    745:     MLPutFunction(link, "ValidQ", 1);
    746:       MLPutSymbol(link, name);
    747:   MLEndPacket(link);
    748:   PetscViewerMathematicaSkipPackets(viewer, RETURNPKT);
    749:   MLGetSymbol(link, &symbol);
    750:   PetscStrcmp("True", (char *) symbol, &match);
    751:   if (!match) {
    752:     MLDisownSymbol(link, symbol);
    753:     SETERRQ(PETSC_ERR_PLIB, "Invalid CSR matrix in Mathematica");
    754:   }
    755:   MLDisownSymbol(link, symbol);
    756:   /* Skip ReturnPacket */
    757:   MLNewPacket(link);

    759:   return(0);
    760: }