Actual source code: mathematica.c

  1: #include <petsc/private/viewerimpl.h>
  2: #include <petsc/private/pcimpl.h>
  3: #include <../src/mat/impls/aij/seq/aij.h>
  4: #include <mathematica.h>

  6: #if defined(PETSC_HAVE__SNPRINTF) && !defined(PETSC_HAVE_SNPRINTF)
  7:   #define snprintf _snprintf
  8: #endif

 10: PetscViewer  PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE = NULL;
 11: static void *mathematicaEnv                         = NULL;

 13: static PetscBool PetscViewerMathematicaPackageInitialized = PETSC_FALSE;
 14: /*@C
 15:   PetscViewerMathematicaFinalizePackage - This function destroys everything in the Petsc interface to Mathematica. It is
 16:   called from PetscFinalize().

 18:   Level: developer

 20: .seealso: `PetscFinalize()`
 21: @*/
 22: PetscErrorCode PetscViewerMathematicaFinalizePackage(void)
 23: {
 24:   PetscFunctionBegin;
 25:   if (mathematicaEnv) MLDeinitialize((MLEnvironment)mathematicaEnv);
 26:   PetscViewerMathematicaPackageInitialized = PETSC_TRUE;
 27:   PetscFunctionReturn(PETSC_SUCCESS);
 28: }

 30: /*@C
 31:   PetscViewerMathematicaInitializePackage - This function initializes everything in the Petsc interface to Mathematica. It is
 32:   called from `PetscViewerInitializePackage()`.

 34:   Level: developer

 36: .seealso: `PetscSysInitializePackage()`, `PetscInitialize()`
 37: @*/
 38: PetscErrorCode PetscViewerMathematicaInitializePackage(void)
 39: {
 40:   PetscError ierr;

 42:   PetscFunctionBegin;
 43:   if (PetscViewerMathematicaPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
 44:   PetscViewerMathematicaPackageInitialized = PETSC_TRUE;

 46:   mathematicaEnv = (void *)MLInitialize(0);

 48:   PetscCall(PetscRegisterFinalize(PetscViewerMathematicaFinalizePackage));
 49:   PetscFunctionReturn(PETSC_SUCCESS);
 50: }

 52: static PetscErrorCode PetscViewerInitializeMathematicaWorld_Private()
 53: {
 54:   PetscFunctionBegin;
 55:   if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) PetscFunctionReturn(PETSC_SUCCESS);
 56:   PetscCall(PetscViewerMathematicaOpen(PETSC_COMM_WORLD, PETSC_DECIDE, NULL, NULL, &PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE));
 57:   PetscFunctionReturn(PETSC_SUCCESS);
 58: }

 60: static PetscErrorCode PetscViewerDestroy_Mathematica(PetscViewer viewer)
 61: {
 62:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;

 64:   PetscFunctionBegin;
 65:   MLClose(vmath->link);
 66:   PetscCall(PetscFree(vmath->linkname));
 67:   PetscCall(PetscFree(vmath->linkhost));
 68:   PetscCall(PetscFree(vmath));
 69:   PetscFunctionReturn(PETSC_SUCCESS);
 70: }

 72: static PetscErrorCode PetscViewerDestroyMathematica_Private(void)
 73: {
 74:   PetscFunctionBegin;
 75:   if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) PetscCall(PetscViewerDestroy(PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE));
 76:   PetscFunctionReturn(PETSC_SUCCESS);
 77: }

 79: static PetscErrorCode PetscViewerMathematicaSetupConnection_Private(PetscViewer v)
 80: {
 81:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
 82: #if defined(MATHEMATICA_3_0)
 83:   int   argc = 6;
 84:   char *argv[6];
 85: #else
 86:   int   argc = 5;
 87:   char *argv[5];
 88: #endif
 89:   char hostname[256];
 90:   long lerr;

 92:   PetscFunctionBegin;
 93:   /* Link name */
 94:   argv[0] = "-linkname";
 95:   if (!vmath->linkname) argv[1] = "math -mathlink";
 96:   else argv[1] = vmath->linkname;

 98:   /* Link host */
 99:   argv[2] = "-linkhost";
100:   if (!vmath->linkhost) {
101:     PetscCall(PetscGetHostName(hostname, sizeof(hostname)));
102:     argv[3] = hostname;
103:   } else argv[3] = vmath->linkhost;

105:     /* Link mode */
106: #if defined(MATHEMATICA_3_0)
107:   argv[4] = "-linkmode";
108:   switch (vmath->linkmode) {
109:   case MATHEMATICA_LINK_CREATE:
110:     argv[5] = "Create";
111:     break;
112:   case MATHEMATICA_LINK_CONNECT:
113:     argv[5] = "Connect";
114:     break;
115:   case MATHEMATICA_LINK_LAUNCH:
116:     argv[5] = "Launch";
117:     break;
118:   }
119: #else
120:   switch (vmath->linkmode) {
121:   case MATHEMATICA_LINK_CREATE:
122:     argv[4] = "-linkcreate";
123:     break;
124:   case MATHEMATICA_LINK_CONNECT:
125:     argv[4] = "-linkconnect";
126:     break;
127:   case MATHEMATICA_LINK_LAUNCH:
128:     argv[4] = "-linklaunch";
129:     break;
130:   }
131: #endif
132:   vmath->link = MLOpenInEnv(mathematicaEnv, argc, argv, &lerr);
133: #endif
134:   PetscFunctionReturn(PETSC_SUCCESS);
135: }

137: PETSC_EXTERN PetscErrorCode PetscViewerCreate_Mathematica(PetscViewer v)
138: {
139:   PetscViewer_Mathematica *vmath;

141:   PetscFunctionBegin;
142:   PetscCall(PetscViewerMathematicaInitializePackage());

144:   PetscCall(PetscNew(&vmath));
145:   v->data         = (void *)vmath;
146:   v->ops->destroy = PetscViewerDestroy_Mathematica;
147:   v->ops->flush   = 0;
148:   PetscCall(PetscStrallocpy(PETSC_VIEWER_MATHEMATICA, &((PetscObject)v)->type_name));

150:   vmath->linkname     = NULL;
151:   vmath->linkhost     = NULL;
152:   vmath->linkmode     = MATHEMATICA_LINK_CONNECT;
153:   vmath->graphicsType = GRAPHICS_MOTIF;
154:   vmath->plotType     = MATHEMATICA_TRIANGULATION_PLOT;
155:   vmath->objName      = NULL;

157:   PetscCall(PetscViewerMathematicaSetFromOptions(v));
158:   PetscCall(PetscViewerMathematicaSetupConnection_Private(v));
159:   PetscFunctionReturn(PETSC_SUCCESS);
160: }

162: static PetscErrorCode PetscViewerMathematicaParseLinkMode(char *modename, LinkMode *mode)
163: {
164:   PetscBool isCreate, isConnect, isLaunch;

166:   PetscFunctionBegin;
167:   PetscCall(PetscStrcasecmp(modename, "Create", &isCreate));
168:   PetscCall(PetscStrcasecmp(modename, "Connect", &isConnect));
169:   PetscCall(PetscStrcasecmp(modename, "Launch", &isLaunch));
170:   if (isCreate) *mode = MATHEMATICA_LINK_CREATE;
171:   else if (isConnect) *mode = MATHEMATICA_LINK_CONNECT;
172:   else if (isLaunch) *mode = MATHEMATICA_LINK_LAUNCH;
173:   else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid Mathematica link mode: %s", modename);
174:   PetscFunctionReturn(PETSC_SUCCESS);
175: }

177: PetscErrorCode PetscViewerMathematicaSetFromOptions(PetscViewer v)
178: {
179:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
180:   char                     linkname[256];
181:   char                     modename[256];
182:   char                     hostname[256];
183:   char                     type[256];
184:   PetscInt                 numPorts;
185:   PetscInt                *ports;
186:   PetscInt                 numHosts;
187:   int                      h;
188:   char                   **hosts;
189:   PetscMPIInt              size, rank;
190:   PetscBool                opt;

192:   PetscFunctionBegin;
193:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)v), &size));
194:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)v), &rank));

196:   /* Get link name */
197:   PetscCall(PetscOptionsGetString("viewer_", "-math_linkname", linkname, sizeof(linkname), &opt));
198:   if (opt) PetscCall(PetscViewerMathematicaSetLinkName(v, linkname));
199:   /* Get link port */
200:   numPorts = size;
201:   PetscCall(PetscMalloc1(size, &ports));
202:   PetscCall(PetscOptionsGetIntArray("viewer_", "-math_linkport", ports, &numPorts, &opt));
203:   if (opt) {
204:     if (numPorts > rank) snprintf(linkname, sizeof(linkname), "%6d", ports[rank]);
205:     else snprintf(linkname, sizeof(linkname), "%6d", ports[0]);
206:     PetscCall(PetscViewerMathematicaSetLinkName(v, linkname));
207:   }
208:   PetscCall(PetscFree(ports));
209:   /* Get link host */
210:   numHosts = size;
211:   PetscCall(PetscMalloc1(size, &hosts));
212:   PetscCall(PetscOptionsGetStringArray("viewer_", "-math_linkhost", hosts, &numHosts, &opt));
213:   if (opt) {
214:     if (numHosts > rank) {
215:       PetscCall(PetscStrncpy(hostname, hosts[rank], sizeof(hostname)));
216:     } else {
217:       PetscCall(PetscStrncpy(hostname, hosts[0], sizeof(hostname)));
218:     }
219:     PetscCall(PetscViewerMathematicaSetLinkHost(v, hostname));
220:   }
221:   for (h = 0; h < numHosts; h++) PetscCall(PetscFree(hosts[h]));
222:   PetscCall(PetscFree(hosts));
223:   /* Get link mode */
224:   PetscCall(PetscOptionsGetString("viewer_", "-math_linkmode", modename, sizeof(modename), &opt));
225:   if (opt) {
226:     LinkMode mode;

228:     PetscCall(PetscViewerMathematicaParseLinkMode(modename, &mode));
229:     PetscCall(PetscViewerMathematicaSetLinkMode(v, mode));
230:   }
231:   /* Get graphics type */
232:   PetscCall(PetscOptionsGetString("viewer_", "-math_graphics", type, sizeof(type), &opt));
233:   if (opt) {
234:     PetscBool isMotif, isPS, isPSFile;

236:     PetscCall(PetscStrcasecmp(type, "Motif", &isMotif));
237:     PetscCall(PetscStrcasecmp(type, "PS", &isPS));
238:     PetscCall(PetscStrcasecmp(type, "PSFile", &isPSFile));
239:     if (isMotif) vmath->graphicsType = GRAPHICS_MOTIF;
240:     else if (isPS) vmath->graphicsType = GRAPHICS_PS_STDOUT;
241:     else if (isPSFile) vmath->graphicsType = GRAPHICS_PS_FILE;
242:   }
243:   /* Get plot type */
244:   PetscCall(PetscOptionsGetString("viewer_", "-math_type", type, sizeof(type), &opt));
245:   if (opt) {
246:     PetscBool isTri, isVecTri, isVec, isSurface;

248:     PetscCall(PetscStrcasecmp(type, "Triangulation", &isTri));
249:     PetscCall(PetscStrcasecmp(type, "VectorTriangulation", &isVecTri));
250:     PetscCall(PetscStrcasecmp(type, "Vector", &isVec));
251:     PetscCall(PetscStrcasecmp(type, "Surface", &isSurface));
252:     if (isTri) vmath->plotType = MATHEMATICA_TRIANGULATION_PLOT;
253:     else if (isVecTri) vmath->plotType = MATHEMATICA_VECTOR_TRIANGULATION_PLOT;
254:     else if (isVec) vmath->plotType = MATHEMATICA_VECTOR_PLOT;
255:     else if (isSurface) vmath->plotType = MATHEMATICA_SURFACE_PLOT;
256:   }
257:   PetscFunctionReturn(PETSC_SUCCESS);
258: }

260: PetscErrorCode PetscViewerMathematicaSetLinkName(PetscViewer v, const char *name)
261: {
262:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;

264:   PetscFunctionBegin;
266:   PetscAssertPointer(name, 2);
267:   PetscCall(PetscStrallocpy(name, &vmath->linkname));
268:   PetscFunctionReturn(PETSC_SUCCESS);
269: }

271: PetscErrorCode PetscViewerMathematicaSetLinkPort(PetscViewer v, int port)
272: {
273:   char name[16];

275:   PetscFunctionBegin;
276:   snprintf(name, 16, "%6d", port);
277:   PetscCall(PetscViewerMathematicaSetLinkName(v, name));
278:   PetscFunctionReturn(PETSC_SUCCESS);
279: }

281: PetscErrorCode PetscViewerMathematicaSetLinkHost(PetscViewer v, const char *host)
282: {
283:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;

285:   PetscFunctionBegin;
287:   PetscAssertPointer(host, 2);
288:   PetscCall(PetscStrallocpy(host, &vmath->linkhost));
289:   PetscFunctionReturn(PETSC_SUCCESS);
290: }

292: PetscErrorCode PetscViewerMathematicaSetLinkMode(PetscViewer v, LinkMode mode)
293: {
294:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;

296:   PetscFunctionBegin;
297:   vmath->linkmode = mode;
298:   PetscFunctionReturn(PETSC_SUCCESS);
299: }

301: /*----------------------------------------- Public Functions --------------------------------------------------------*/
302: /*@C
303:   PetscViewerMathematicaOpen - Communicates with Mathemtica using MathLink.

305:   Collective

307:   Input Parameters:
308: + comm    - The MPI communicator
309: . port    - [optional] The port to connect on, or PETSC_DECIDE
310: . machine - [optional] The machine to run Mathematica on, or NULL
311: - mode    - [optional] The connection mode, or NULL

313:   Output Parameter:
314: . v - The Mathematica viewer

316:   Options Database Keys:
317: + -viewer_math_linkhost <machine> - The host machine for the kernel
318: . -viewer_math_linkname <name>    - The full link name for the connection
319: . -viewer_math_linkport <port>    - The port for the connection
320: . -viewer_math_mode <mode>        - The mode, e.g. Launch, Connect
321: . -viewer_math_type <type>        - The plot type, e.g. Triangulation, Vector
322: - -viewer_math_graphics <output>  - The output type, e.g. Motif, PS, PSFile

324:   Level: intermediate

326:   Note:
327:   Most users should employ the following commands to access the
328:   Mathematica viewers
329: .vb
330:     PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
331:     MatView(Mat matrix, PetscViewer viewer)

333:                 or

335:     PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
336:     VecView(Vec vector, PetscViewer viewer)
337: .ve

339: .seealso: `PETSCVIEWERMATHEMATICA`, `MatView()`, `VecView()`
340: @*/
341: PetscErrorCode PetscViewerMathematicaOpen(MPI_Comm comm, int port, const char machine[], const char mode[], PetscViewer *v)
342: {
343:   PetscFunctionBegin;
344:   PetscCall(PetscViewerCreate(comm, v));
345: #if 0
346:   LinkMode linkmode;
347:   PetscCall(PetscViewerMathematicaSetLinkPort(*v, port));
348:   PetscCall(PetscViewerMathematicaSetLinkHost(*v, machine));
349:   PetscCall(PetscViewerMathematicaParseLinkMode(mode, &linkmode));
350:   PetscCall(PetscViewerMathematicaSetLinkMode(*v, linkmode));
351: #endif
352:   PetscCall(PetscViewerSetType(*v, PETSC_VIEWER_MATHEMATICA));
353:   PetscFunctionReturn(PETSC_SUCCESS);
354: }

356: /*
357:   PetscViewerMathematicaGetLink - Returns the link to Mathematica from a `PETSCVIEWERMATHEMATICA`

359:   Input Parameters:
360: + viewer - The Mathematica viewer
361: - link   - The link to Mathematica

363:   Level: intermediate

365: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaOpen()`
366: */
367: static PetscErrorCode PetscViewerMathematicaGetLink(PetscViewer viewer, MLINK *link)
368: {
369:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;

371:   PetscFunctionBegin;
373:   *link = vmath->link;
374:   PetscFunctionReturn(PETSC_SUCCESS);
375: }

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

380:   Input Parameters:
381: + viewer - The Mathematica viewer
382: - type   - The packet type to search for, e.g RETURNPKT

384:   Level: advanced

386: .seealso: `PetscViewerMathematicaSetName()`, `PetscViewerMathematicaGetVector()`
387: @*/
388: PetscErrorCode PetscViewerMathematicaSkipPackets(PetscViewer viewer, int type)
389: {
390:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
391:   MLINK                    link  = vmath->link; /* The link to Mathematica */
392:   int                      pkt;                 /* The packet type */

394:   PetscFunctionBegin;
395:   while ((pkt = MLNextPacket(link)) && (pkt != type)) MLNewPacket(link);
396:   if (!pkt) {
397:     MLClearError(link);
398:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, (char *)MLErrorMessage(link));
399:   }
400:   PetscFunctionReturn(PETSC_SUCCESS);
401: }

403: /*@C
404:   PetscViewerMathematicaGetName - Retrieve the default name for objects communicated to Mathematica via `PETSCVIEWERMATHEMATICA`

406:   Input Parameter:
407: . viewer - The Mathematica viewer

409:   Output Parameter:
410: . name - The name for new objects created in Mathematica

412:   Level: intermediate

414: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaSetName()`, `PetscViewerMathematicaClearName()`
415: @*/
416: PetscErrorCode PetscViewerMathematicaGetName(PetscViewer viewer, const char **name)
417: {
418:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;

420:   PetscFunctionBegin;
422:   PetscAssertPointer(name, 2);
423:   *name = vmath->objName;
424:   PetscFunctionReturn(PETSC_SUCCESS);
425: }

427: /*@C
428:   PetscViewerMathematicaSetName - Override the default name for objects communicated to Mathematica via `PETSCVIEWERMATHEMATICA`

430:   Input Parameters:
431: + viewer - The Mathematica viewer
432: - name   - The name for new objects created in Mathematica

434:   Level: intermediate

436: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaClearName()`
437: @*/
438: PetscErrorCode PetscViewerMathematicaSetName(PetscViewer viewer, const char name[])
439: {
440:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;

442:   PetscFunctionBegin;
444:   PetscAssertPointer(name, 2);
445:   vmath->objName = name;
446:   PetscFunctionReturn(PETSC_SUCCESS);
447: }

449: /*@
450:   PetscViewerMathematicaClearName - Use the default name for objects communicated to Mathematica

452:   Input Parameter:
453: . viewer - The Mathematica viewer

455:   Level: intermediate

457: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaGetName()`, `PetscViewerMathematicaSetName()`
458: @*/
459: PetscErrorCode PetscViewerMathematicaClearName(PetscViewer viewer)
460: {
461:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;

463:   PetscFunctionBegin;
465:   vmath->objName = NULL;
466:   PetscFunctionReturn(PETSC_SUCCESS);
467: }

469: /*@
470:   PetscViewerMathematicaGetVector - Retrieve a vector from Mathematica via a `PETSCVIEWERMATHEMATICA`

472:   Input Parameter:
473: . viewer - The Mathematica viewer

475:   Output Parameter:
476: . v - The vector

478:   Level: intermediate

480: .seealso: `PETSCVIEWERMATHEMATICA`, `VecView()`, `PetscViewerMathematicaPutVector()`
481: @*/
482: PetscErrorCode PetscViewerMathematicaGetVector(PetscViewer viewer, Vec v)
483: {
484:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
485:   MLINK                    link; /* The link to Mathematica */
486:   char                    *name;
487:   PetscScalar             *mArray, *array;
488:   long                     mSize;
489:   int                      n;

491:   PetscFunctionBegin;

495:   /* Determine the object name */
496:   if (!vmath->objName) name = "vec";
497:   else name = (char *)vmath->objName;

499:   link = vmath->link;
500:   PetscCall(VecGetLocalSize(v, &n));
501:   PetscCall(VecGetArray(v, &array));
502:   MLPutFunction(link, "EvaluatePacket", 1);
503:   MLPutSymbol(link, name);
504:   MLEndPacket(link);
505:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
506:   MLGetRealList(link, &mArray, &mSize);
507:   PetscCheck(n == mSize, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Incompatible vector sizes %d %d", n, mSize);
508:   PetscCall(PetscArraycpy(array, mArray, mSize));
509:   MLDisownRealList(link, mArray, mSize);
510:   PetscCall(VecRestoreArray(v, &array));
511:   PetscFunctionReturn(PETSC_SUCCESS);
512: }

514: /*@
515:   PetscViewerMathematicaPutVector - Send a vector to Mathematica via a `PETSCVIEWERMATHEMATICA` `PetscViewer`

517:   Input Parameters:
518: + viewer - The Mathematica viewer
519: - v      - The vector

521:   Level: intermediate

523: .seealso: `PETSCVIEWERMATHEMATICA`, `VecView()`, `PetscViewerMathematicaGetVector()`
524: @*/
525: PetscErrorCode PetscViewerMathematicaPutVector(PetscViewer viewer, Vec v)
526: {
527:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
528:   MLINK                    link  = vmath->link; /* The link to Mathematica */
529:   char                    *name;
530:   PetscScalar             *array;
531:   int                      n;

533:   PetscFunctionBegin;
534:   /* Determine the object name */
535:   if (!vmath->objName) name = "vec";
536:   else name = (char *)vmath->objName;

538:   PetscCall(VecGetLocalSize(v, &n));
539:   PetscCall(VecGetArray(v, &array));

541:   /* Send the Vector object */
542:   MLPutFunction(link, "EvaluatePacket", 1);
543:   MLPutFunction(link, "Set", 2);
544:   MLPutSymbol(link, name);
545:   MLPutRealList(link, array, n);
546:   MLEndPacket(link);
547:   /* Skip packets until ReturnPacket */
548:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
549:   /* Skip ReturnPacket */
550:   MLNewPacket(link);

552:   PetscCall(VecRestoreArray(v, &array));
553:   PetscFunctionReturn(PETSC_SUCCESS);
554: }

556: PetscErrorCode PetscViewerMathematicaPutMatrix(PetscViewer viewer, int m, int n, PetscReal *a)
557: {
558:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
559:   MLINK                    link  = vmath->link; /* The link to Mathematica */
560:   char                    *name;

562:   PetscFunctionBegin;
563:   /* Determine the object name */
564:   if (!vmath->objName) name = "mat";
565:   else name = (char *)vmath->objName;

567:   /* Send the dense matrix object */
568:   MLPutFunction(link, "EvaluatePacket", 1);
569:   MLPutFunction(link, "Set", 2);
570:   MLPutSymbol(link, name);
571:   MLPutFunction(link, "Transpose", 1);
572:   MLPutFunction(link, "Partition", 2);
573:   MLPutRealList(link, a, m * n);
574:   MLPutInteger(link, m);
575:   MLEndPacket(link);
576:   /* Skip packets until ReturnPacket */
577:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
578:   /* Skip ReturnPacket */
579:   MLNewPacket(link);
580:   PetscFunctionReturn(PETSC_SUCCESS);
581: }

583: PetscErrorCode PetscViewerMathematicaPutCSRMatrix(PetscViewer viewer, int m, int n, int *i, int *j, PetscReal *a)
584: {
585:   PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
586:   MLINK                    link  = vmath->link; /* The link to Mathematica */
587:   const char              *symbol;
588:   char                    *name;
589:   PetscBool                match;

591:   PetscFunctionBegin;
592:   /* Determine the object name */
593:   if (!vmath->objName) name = "mat";
594:   else name = (char *)vmath->objName;

596:   /* Make sure Mathematica recognizes sparse matrices */
597:   MLPutFunction(link, "EvaluatePacket", 1);
598:   MLPutFunction(link, "Needs", 1);
599:   MLPutString(link, "LinearAlgebra`CSRMatrix`");
600:   MLEndPacket(link);
601:   /* Skip packets until ReturnPacket */
602:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
603:   /* Skip ReturnPacket */
604:   MLNewPacket(link);

606:   /* Send the CSRMatrix object */
607:   MLPutFunction(link, "EvaluatePacket", 1);
608:   MLPutFunction(link, "Set", 2);
609:   MLPutSymbol(link, name);
610:   MLPutFunction(link, "CSRMatrix", 5);
611:   MLPutInteger(link, m);
612:   MLPutInteger(link, n);
613:   MLPutFunction(link, "Plus", 2);
614:   MLPutIntegerList(link, i, m + 1);
615:   MLPutInteger(link, 1);
616:   MLPutFunction(link, "Plus", 2);
617:   MLPutIntegerList(link, j, i[m]);
618:   MLPutInteger(link, 1);
619:   MLPutRealList(link, a, i[m]);
620:   MLEndPacket(link);
621:   /* Skip packets until ReturnPacket */
622:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
623:   /* Skip ReturnPacket */
624:   MLNewPacket(link);

626:   /* Check that matrix is valid */
627:   MLPutFunction(link, "EvaluatePacket", 1);
628:   MLPutFunction(link, "ValidQ", 1);
629:   MLPutSymbol(link, name);
630:   MLEndPacket(link);
631:   PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
632:   MLGetSymbol(link, &symbol);
633:   PetscCall(PetscStrcmp("True", (char *)symbol, &match));
634:   if (!match) {
635:     MLDisownSymbol(link, symbol);
636:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid CSR matrix in Mathematica");
637:   }
638:   MLDisownSymbol(link, symbol);
639:   /* Skip ReturnPacket */
640:   MLNewPacket(link);
641:   PetscFunctionReturn(PETSC_SUCCESS);
642: }