Actual source code: mathematica.c
1: #include <petsc/private/viewerimpl.h>
2: #include <mathematica.h>
4: #if defined(PETSC_HAVE__SNPRINTF) && !defined(PETSC_HAVE_SNPRINTF)
5: #define snprintf _snprintf
6: #endif
8: PetscViewer PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE = NULL;
9: static void *mathematicaEnv = NULL;
11: static PetscBool PetscViewerMathematicaPackageInitialized = PETSC_FALSE;
12: /*@C
13: PetscViewerMathematicaFinalizePackage - This function destroys everything in the PETSc interface to Mathematica. It is
14: called from PetscFinalize().
16: Level: developer
18: .seealso: `PetscFinalize()`
19: @*/
20: PetscErrorCode PetscViewerMathematicaFinalizePackage(void)
21: {
22: PetscFunctionBegin;
23: if (mathematicaEnv) MLDeinitialize((MLEnvironment)mathematicaEnv);
24: PetscViewerMathematicaPackageInitialized = PETSC_TRUE;
25: PetscFunctionReturn(PETSC_SUCCESS);
26: }
28: /*@C
29: PetscViewerMathematicaInitializePackage - This function initializes everything in the PETSc interface to Mathematica. It is
30: called from `PetscViewerInitializePackage()`.
32: Level: developer
34: .seealso: `PetscSysInitializePackage()`, `PetscInitialize()`
35: @*/
36: PetscErrorCode PetscViewerMathematicaInitializePackage(void)
37: {
38: PetscError ierr;
40: PetscFunctionBegin;
41: if (PetscViewerMathematicaPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
42: PetscViewerMathematicaPackageInitialized = PETSC_TRUE;
44: mathematicaEnv = (void *)MLInitialize(0);
46: PetscCall(PetscRegisterFinalize(PetscViewerMathematicaFinalizePackage));
47: PetscFunctionReturn(PETSC_SUCCESS);
48: }
50: static PetscErrorCode PetscViewerInitializeMathematicaWorld_Private()
51: {
52: PetscFunctionBegin;
53: if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) PetscFunctionReturn(PETSC_SUCCESS);
54: PetscCall(PetscViewerMathematicaOpen(PETSC_COMM_WORLD, PETSC_DECIDE, NULL, NULL, &PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE));
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }
58: static PetscErrorCode PetscViewerDestroy_Mathematica(PetscViewer viewer)
59: {
60: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
62: PetscFunctionBegin;
63: MLClose(vmath->link);
64: PetscCall(PetscFree(vmath->linkname));
65: PetscCall(PetscFree(vmath->linkhost));
66: PetscCall(PetscFree(vmath));
67: PetscFunctionReturn(PETSC_SUCCESS);
68: }
70: static PetscErrorCode PetscViewerDestroyMathematica_Private(void)
71: {
72: PetscFunctionBegin;
73: if (PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE) PetscCall(PetscViewerDestroy(PETSC_VIEWER_MATHEMATICA_WORLD_PRIVATE));
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
77: static PetscErrorCode PetscViewerMathematicaSetupConnection_Private(PetscViewer v)
78: {
79: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
80: #if defined(MATHEMATICA_3_0)
81: int argc = 6;
82: char *argv[6];
83: #else
84: int argc = 5;
85: char *argv[5];
86: #endif
87: char hostname[256];
88: long lerr;
90: PetscFunctionBegin;
91: /* Link name */
92: argv[0] = "-linkname";
93: if (!vmath->linkname) argv[1] = "math -mathlink";
94: else argv[1] = vmath->linkname;
96: /* Link host */
97: argv[2] = "-linkhost";
98: if (!vmath->linkhost) {
99: PetscCall(PetscGetHostName(hostname, sizeof(hostname)));
100: argv[3] = hostname;
101: } else argv[3] = vmath->linkhost;
103: /* Link mode */
104: #if defined(MATHEMATICA_3_0)
105: argv[4] = "-linkmode";
106: switch (vmath->linkmode) {
107: case MATHEMATICA_LINK_CREATE:
108: argv[5] = "Create";
109: break;
110: case MATHEMATICA_LINK_CONNECT:
111: argv[5] = "Connect";
112: break;
113: case MATHEMATICA_LINK_LAUNCH:
114: argv[5] = "Launch";
115: break;
116: }
117: #else
118: switch (vmath->linkmode) {
119: case MATHEMATICA_LINK_CREATE:
120: argv[4] = "-linkcreate";
121: break;
122: case MATHEMATICA_LINK_CONNECT:
123: argv[4] = "-linkconnect";
124: break;
125: case MATHEMATICA_LINK_LAUNCH:
126: argv[4] = "-linklaunch";
127: break;
128: }
129: #endif
130: vmath->link = MLOpenInEnv(mathematicaEnv, argc, argv, &lerr);
131: #endif
132: PetscFunctionReturn(PETSC_SUCCESS);
133: }
135: PETSC_EXTERN PetscErrorCode PetscViewerCreate_Mathematica(PetscViewer v)
136: {
137: PetscViewer_Mathematica *vmath;
139: PetscFunctionBegin;
140: PetscCall(PetscViewerMathematicaInitializePackage());
142: PetscCall(PetscNew(&vmath));
143: v->data = (void *)vmath;
144: v->ops->destroy = PetscViewerDestroy_Mathematica;
145: v->ops->flush = 0;
146: PetscCall(PetscStrallocpy(PETSC_VIEWER_MATHEMATICA, &((PetscObject)v)->type_name));
148: vmath->linkname = NULL;
149: vmath->linkhost = NULL;
150: vmath->linkmode = MATHEMATICA_LINK_CONNECT;
151: vmath->graphicsType = GRAPHICS_MOTIF;
152: vmath->plotType = MATHEMATICA_TRIANGULATION_PLOT;
153: vmath->objName = NULL;
155: PetscCall(PetscViewerMathematicaSetFromOptions(v));
156: PetscCall(PetscViewerMathematicaSetupConnection_Private(v));
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: static PetscErrorCode PetscViewerMathematicaParseLinkMode(char *modename, LinkMode *mode)
161: {
162: PetscBool isCreate, isConnect, isLaunch;
164: PetscFunctionBegin;
165: PetscCall(PetscStrcasecmp(modename, "Create", &isCreate));
166: PetscCall(PetscStrcasecmp(modename, "Connect", &isConnect));
167: PetscCall(PetscStrcasecmp(modename, "Launch", &isLaunch));
168: if (isCreate) *mode = MATHEMATICA_LINK_CREATE;
169: else if (isConnect) *mode = MATHEMATICA_LINK_CONNECT;
170: else if (isLaunch) *mode = MATHEMATICA_LINK_LAUNCH;
171: else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid Mathematica link mode: %s", modename);
172: PetscFunctionReturn(PETSC_SUCCESS);
173: }
175: PetscErrorCode PetscViewerMathematicaSetFromOptions(PetscViewer v)
176: {
177: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
178: char linkname[256];
179: char modename[256];
180: char hostname[256];
181: char type[256];
182: PetscInt numPorts;
183: PetscInt *ports;
184: PetscInt numHosts;
185: int h;
186: char **hosts;
187: PetscMPIInt size, rank;
188: PetscBool opt;
190: PetscFunctionBegin;
191: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)v), &size));
192: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)v), &rank));
194: /* Get link name */
195: PetscCall(PetscOptionsGetString("viewer_", "-math_linkname", linkname, sizeof(linkname), &opt));
196: if (opt) PetscCall(PetscViewerMathematicaSetLinkName(v, linkname));
197: /* Get link port */
198: numPorts = size;
199: PetscCall(PetscMalloc1(size, &ports));
200: PetscCall(PetscOptionsGetIntArray("viewer_", "-math_linkport", ports, &numPorts, &opt));
201: if (opt) {
202: if (numPorts > rank) snprintf(linkname, sizeof(linkname), "%6d", ports[rank]);
203: else snprintf(linkname, sizeof(linkname), "%6d", ports[0]);
204: PetscCall(PetscViewerMathematicaSetLinkName(v, linkname));
205: }
206: PetscCall(PetscFree(ports));
207: /* Get link host */
208: numHosts = size;
209: PetscCall(PetscMalloc1(size, &hosts));
210: PetscCall(PetscOptionsGetStringArray("viewer_", "-math_linkhost", hosts, &numHosts, &opt));
211: if (opt) {
212: if (numHosts > rank) {
213: PetscCall(PetscStrncpy(hostname, hosts[rank], sizeof(hostname)));
214: } else {
215: PetscCall(PetscStrncpy(hostname, hosts[0], sizeof(hostname)));
216: }
217: PetscCall(PetscViewerMathematicaSetLinkHost(v, hostname));
218: }
219: for (h = 0; h < numHosts; h++) PetscCall(PetscFree(hosts[h]));
220: PetscCall(PetscFree(hosts));
221: /* Get link mode */
222: PetscCall(PetscOptionsGetString("viewer_", "-math_linkmode", modename, sizeof(modename), &opt));
223: if (opt) {
224: LinkMode mode;
226: PetscCall(PetscViewerMathematicaParseLinkMode(modename, &mode));
227: PetscCall(PetscViewerMathematicaSetLinkMode(v, mode));
228: }
229: /* Get graphics type */
230: PetscCall(PetscOptionsGetString("viewer_", "-math_graphics", type, sizeof(type), &opt));
231: if (opt) {
232: PetscBool isMotif, isPS, isPSFile;
234: PetscCall(PetscStrcasecmp(type, "Motif", &isMotif));
235: PetscCall(PetscStrcasecmp(type, "PS", &isPS));
236: PetscCall(PetscStrcasecmp(type, "PSFile", &isPSFile));
237: if (isMotif) vmath->graphicsType = GRAPHICS_MOTIF;
238: else if (isPS) vmath->graphicsType = GRAPHICS_PS_STDOUT;
239: else if (isPSFile) vmath->graphicsType = GRAPHICS_PS_FILE;
240: }
241: /* Get plot type */
242: PetscCall(PetscOptionsGetString("viewer_", "-math_type", type, sizeof(type), &opt));
243: if (opt) {
244: PetscBool isTri, isVecTri, isVec, isSurface;
246: PetscCall(PetscStrcasecmp(type, "Triangulation", &isTri));
247: PetscCall(PetscStrcasecmp(type, "VectorTriangulation", &isVecTri));
248: PetscCall(PetscStrcasecmp(type, "Vector", &isVec));
249: PetscCall(PetscStrcasecmp(type, "Surface", &isSurface));
250: if (isTri) vmath->plotType = MATHEMATICA_TRIANGULATION_PLOT;
251: else if (isVecTri) vmath->plotType = MATHEMATICA_VECTOR_TRIANGULATION_PLOT;
252: else if (isVec) vmath->plotType = MATHEMATICA_VECTOR_PLOT;
253: else if (isSurface) vmath->plotType = MATHEMATICA_SURFACE_PLOT;
254: }
255: PetscFunctionReturn(PETSC_SUCCESS);
256: }
258: PetscErrorCode PetscViewerMathematicaSetLinkName(PetscViewer v, const char *name)
259: {
260: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
262: PetscFunctionBegin;
264: PetscAssertPointer(name, 2);
265: PetscCall(PetscStrallocpy(name, &vmath->linkname));
266: PetscFunctionReturn(PETSC_SUCCESS);
267: }
269: PetscErrorCode PetscViewerMathematicaSetLinkPort(PetscViewer v, int port)
270: {
271: char name[16];
273: PetscFunctionBegin;
274: snprintf(name, 16, "%6d", port);
275: PetscCall(PetscViewerMathematicaSetLinkName(v, name));
276: PetscFunctionReturn(PETSC_SUCCESS);
277: }
279: PetscErrorCode PetscViewerMathematicaSetLinkHost(PetscViewer v, const char *host)
280: {
281: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
283: PetscFunctionBegin;
285: PetscAssertPointer(host, 2);
286: PetscCall(PetscStrallocpy(host, &vmath->linkhost));
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
290: PetscErrorCode PetscViewerMathematicaSetLinkMode(PetscViewer v, LinkMode mode)
291: {
292: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)v->data;
294: PetscFunctionBegin;
295: vmath->linkmode = mode;
296: PetscFunctionReturn(PETSC_SUCCESS);
297: }
299: /*----------------------------------------- Public Functions --------------------------------------------------------*/
300: /*@
301: PetscViewerMathematicaOpen - Communicates with Mathemtica using MathLink.
303: Collective
305: Input Parameters:
306: + comm - The MPI communicator
307: . port - [optional] The port to connect on, or PETSC_DECIDE
308: . machine - [optional] The machine to run Mathematica on, or NULL
309: - mode - [optional] The connection mode, or NULL
311: Output Parameter:
312: . v - The Mathematica viewer
314: Options Database Keys:
315: + -viewer_math_linkhost <machine> - The host machine for the kernel
316: . -viewer_math_linkname <name> - The full link name for the connection
317: . -viewer_math_linkport <port> - The port for the connection
318: . -viewer_math_mode <mode> - The mode, e.g. Launch, Connect
319: . -viewer_math_type <type> - The plot type, e.g. Triangulation, Vector
320: - -viewer_math_graphics <output> - The output type, e.g. Motif, PS, PSFile
322: Level: intermediate
324: Note:
325: Most users should employ the following commands to access the
326: Mathematica viewers
327: .vb
328: PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
329: MatView(Mat matrix, PetscViewer viewer)
331: or
333: PetscViewerMathematicaOpen(MPI_Comm comm, int port, char *machine, char *mode, PetscViewer &viewer)
334: VecView(Vec vector, PetscViewer viewer)
335: .ve
337: .seealso: `PETSCVIEWERMATHEMATICA`, `MatView()`, `VecView()`
338: @*/
339: PetscErrorCode PetscViewerMathematicaOpen(MPI_Comm comm, int port, const char machine[], const char mode[], PetscViewer *v)
340: {
341: PetscFunctionBegin;
342: PetscCall(PetscViewerCreate(comm, v));
343: #if 0
344: LinkMode linkmode;
345: PetscCall(PetscViewerMathematicaSetLinkPort(*v, port));
346: PetscCall(PetscViewerMathematicaSetLinkHost(*v, machine));
347: PetscCall(PetscViewerMathematicaParseLinkMode(mode, &linkmode));
348: PetscCall(PetscViewerMathematicaSetLinkMode(*v, linkmode));
349: #endif
350: PetscCall(PetscViewerSetType(*v, PETSC_VIEWER_MATHEMATICA));
351: PetscFunctionReturn(PETSC_SUCCESS);
352: }
354: /*
355: PetscViewerMathematicaGetLink - Returns the link to Mathematica from a `PETSCVIEWERMATHEMATICA`
357: Input Parameters:
358: + viewer - The Mathematica viewer
359: - link - The link to Mathematica
361: Level: intermediate
363: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaOpen()`
364: */
365: static PetscErrorCode PetscViewerMathematicaGetLink(PetscViewer viewer, MLINK *link)
366: {
367: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
369: PetscFunctionBegin;
371: *link = vmath->link;
372: PetscFunctionReturn(PETSC_SUCCESS);
373: }
375: /*@C
376: PetscViewerMathematicaSkipPackets - Discard packets sent by Mathematica until a certain packet type is received
378: Input Parameters:
379: + viewer - The Mathematica viewer
380: - type - The packet type to search for, e.g RETURNPKT
382: Level: advanced
384: .seealso: `PetscViewerMathematicaSetName()`, `PetscViewerMathematicaGetVector()`
385: @*/
386: PetscErrorCode PetscViewerMathematicaSkipPackets(PetscViewer viewer, int type)
387: {
388: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
389: MLINK link = vmath->link; /* The link to Mathematica */
390: int pkt; /* The packet type */
392: PetscFunctionBegin;
393: while ((pkt = MLNextPacket(link)) && (pkt != type)) MLNewPacket(link);
394: if (!pkt) {
395: MLClearError(link);
396: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, (char *)MLErrorMessage(link));
397: }
398: PetscFunctionReturn(PETSC_SUCCESS);
399: }
401: /*@C
402: PetscViewerMathematicaGetName - Retrieve the default name for objects communicated to Mathematica via `PETSCVIEWERMATHEMATICA`
404: Input Parameter:
405: . viewer - The Mathematica viewer
407: Output Parameter:
408: . name - The name for new objects created in Mathematica
410: Level: intermediate
412: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaSetName()`, `PetscViewerMathematicaClearName()`
413: @*/
414: PetscErrorCode PetscViewerMathematicaGetName(PetscViewer viewer, const char **name)
415: {
416: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
418: PetscFunctionBegin;
420: PetscAssertPointer(name, 2);
421: *name = vmath->objName;
422: PetscFunctionReturn(PETSC_SUCCESS);
423: }
425: /*@C
426: PetscViewerMathematicaSetName - Override the default name for objects communicated to Mathematica via `PETSCVIEWERMATHEMATICA`
428: Input Parameters:
429: + viewer - The Mathematica viewer
430: - name - The name for new objects created in Mathematica
432: Level: intermediate
434: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaClearName()`
435: @*/
436: PetscErrorCode PetscViewerMathematicaSetName(PetscViewer viewer, const char name[])
437: {
438: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
440: PetscFunctionBegin;
442: PetscAssertPointer(name, 2);
443: vmath->objName = name;
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: /*@
448: PetscViewerMathematicaClearName - Use the default name for objects communicated to Mathematica
450: Input Parameter:
451: . viewer - The Mathematica viewer
453: Level: intermediate
455: .seealso: `PETSCVIEWERMATHEMATICA`, `PetscViewerMathematicaGetName()`, `PetscViewerMathematicaSetName()`
456: @*/
457: PetscErrorCode PetscViewerMathematicaClearName(PetscViewer viewer)
458: {
459: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
461: PetscFunctionBegin;
463: vmath->objName = NULL;
464: PetscFunctionReturn(PETSC_SUCCESS);
465: }
467: PetscErrorCode PetscViewerMathematicaPutMatrix(PetscViewer viewer, int m, int n, PetscReal *a)
468: {
469: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
470: MLINK link = vmath->link; /* The link to Mathematica */
471: char *name;
473: PetscFunctionBegin;
474: /* Determine the object name */
475: if (!vmath->objName) name = "mat";
476: else name = (char *)vmath->objName;
478: /* Send the dense matrix object */
479: MLPutFunction(link, "EvaluatePacket", 1);
480: MLPutFunction(link, "Set", 2);
481: MLPutSymbol(link, name);
482: MLPutFunction(link, "Transpose", 1);
483: MLPutFunction(link, "Partition", 2);
484: MLPutRealList(link, a, m * n);
485: MLPutInteger(link, m);
486: MLEndPacket(link);
487: /* Skip packets until ReturnPacket */
488: PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
489: /* Skip ReturnPacket */
490: MLNewPacket(link);
491: PetscFunctionReturn(PETSC_SUCCESS);
492: }
494: PetscErrorCode PetscViewerMathematicaPutCSRMatrix(PetscViewer viewer, int m, int n, int *i, int *j, PetscReal *a)
495: {
496: PetscViewer_Mathematica *vmath = (PetscViewer_Mathematica *)viewer->data;
497: MLINK link = vmath->link; /* The link to Mathematica */
498: const char *symbol;
499: char *name;
500: PetscBool match;
502: PetscFunctionBegin;
503: /* Determine the object name */
504: if (!vmath->objName) name = "mat";
505: else name = (char *)vmath->objName;
507: /* Make sure Mathematica recognizes sparse matrices */
508: MLPutFunction(link, "EvaluatePacket", 1);
509: MLPutFunction(link, "Needs", 1);
510: MLPutString(link, "LinearAlgebra`CSRMatrix`");
511: MLEndPacket(link);
512: /* Skip packets until ReturnPacket */
513: PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
514: /* Skip ReturnPacket */
515: MLNewPacket(link);
517: /* Send the CSRMatrix object */
518: MLPutFunction(link, "EvaluatePacket", 1);
519: MLPutFunction(link, "Set", 2);
520: MLPutSymbol(link, name);
521: MLPutFunction(link, "CSRMatrix", 5);
522: MLPutInteger(link, m);
523: MLPutInteger(link, n);
524: MLPutFunction(link, "Plus", 2);
525: MLPutIntegerList(link, i, m + 1);
526: MLPutInteger(link, 1);
527: MLPutFunction(link, "Plus", 2);
528: MLPutIntegerList(link, j, i[m]);
529: MLPutInteger(link, 1);
530: MLPutRealList(link, a, i[m]);
531: MLEndPacket(link);
532: /* Skip packets until ReturnPacket */
533: PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
534: /* Skip ReturnPacket */
535: MLNewPacket(link);
537: /* Check that matrix is valid */
538: MLPutFunction(link, "EvaluatePacket", 1);
539: MLPutFunction(link, "ValidQ", 1);
540: MLPutSymbol(link, name);
541: MLEndPacket(link);
542: PetscCall(PetscViewerMathematicaSkipPackets(viewer, RETURNPKT));
543: MLGetSymbol(link, &symbol);
544: PetscCall(PetscStrcmp("True", (char *)symbol, &match));
545: if (!match) {
546: MLDisownSymbol(link, symbol);
547: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid CSR matrix in Mathematica");
548: }
549: MLDisownSymbol(link, symbol);
550: /* Skip ReturnPacket */
551: MLNewPacket(link);
552: PetscFunctionReturn(PETSC_SUCCESS);
553: }