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: }