Actual source code: vmatlab.c

  1: #include <petsc/private/viewerimpl.h>
  2: #include <mat.h>                      /*I    "petscmat.h"      I*/

  4: typedef struct {
  5:   MATFile      *ep;
  6:   PetscMPIInt   rank;
  7:   PetscFileMode btype;
  8: } PetscViewer_Matlab;

 10: /*@
 11:   PetscViewerMatlabPutArray - Puts an array into the `PETSCVIEWERMATLAB` viewer.

 13:   Not Collective, only processor zero saves `array`

 15:   Input Parameters:
 16: + mfile - the viewer
 17: . m     - the first dimensions of `array`
 18: . n     - the second dimensions of `array`
 19: . array - the array (represented in one dimension)
 20: - name  - the MATLAB name of `array`

 22:   Level: advanced

 24:   Note:
 25:   Only writes `array` values on processor 0.

 27: .seealso: `PETSCVIEWERMATLAB`, `PetscViewerMatlabGetArray()`
 28: @*/
 29: PetscErrorCode PetscViewerMatlabPutArray(PetscViewer mfile, int m, int n, const PetscScalar *array, const char *name)
 30: {
 31:   PetscViewer_Matlab *ml;
 32:   mxArray            *mat;

 34:   PetscFunctionBegin;
 35:   PetscCheck(mfile, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Null argument: probably PETSC_VIEWER_MATLAB_() failed");
 36:   ml = (PetscViewer_Matlab *)mfile->data;
 37:   if (!ml->rank) {
 38:     PetscCall(PetscInfo(mfile, "Putting MATLAB array %s\n", name));
 39: #if !defined(PETSC_USE_COMPLEX)
 40:     mat = mxCreateDoubleMatrix(m, n, mxREAL);
 41: #else
 42:     mat = mxCreateDoubleMatrix(m, n, mxCOMPLEX);
 43: #endif
 44:     PetscCall(PetscArraycpy(mxGetPr(mat), array, m * n));
 45:     matPutVariable(ml->ep, name, mat);

 47:     PetscCall(PetscInfo(mfile, "Put MATLAB array %s\n", name));
 48:   }
 49:   PetscFunctionReturn(PETSC_SUCCESS);
 50: }

 52: PetscErrorCode PetscViewerMatlabPutVariable(PetscViewer viewer, const char *name, void *mat)
 53: {
 54:   PetscViewer_Matlab *ml = (PetscViewer_Matlab *)viewer->data;

 56:   PetscFunctionBegin;
 57:   matPutVariable(ml->ep, name, (mxArray *)mat);
 58:   PetscFunctionReturn(PETSC_SUCCESS);
 59: }

 61: /*@
 62:   PetscViewerMatlabGetArray - Gets a variable from a `PETSCVIEWERMATLAB` viewer into an array

 64:   Not Collective; only processor zero reads in the array

 66:   Input Parameters:
 67: + mfile - the MATLAB file viewer
 68: . m     - the first dimensions of `array`
 69: . n     - the second dimensions of `array`
 70: . array - the array (represented in one dimension), must of be length `m` * `n`
 71: - name  - the MATLAB name of `array`

 73:   Level: advanced

 75:   Note:
 76:   Only reads in `array` values on processor 0.

 78: .seealso: `PETSCVIEWERMATLAB`, `PetscViewerMatlabPutArray()`
 79: @*/
 80: PetscErrorCode PetscViewerMatlabGetArray(PetscViewer mfile, int m, int n, PetscScalar array[], const char *name)
 81: {
 82:   PetscViewer_Matlab *ml;
 83:   mxArray            *mat;

 85:   PetscFunctionBegin;
 86:   PetscCheck(mfile, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Null argument: probably PETSC_VIEWER_MATLAB_() failed");
 87:   ml = (PetscViewer_Matlab *)mfile->data;
 88:   if (!ml->rank) {
 89:     PetscCall(PetscInfo(mfile, "Getting MATLAB array %s\n", name));
 90:     mat = matGetVariable(ml->ep, name);
 91:     PetscCheck(mat, PETSC_COMM_SELF, PETSC_ERR_LIB, "Unable to get array %s from matlab", name);
 92:     PetscCall(PetscArraycpy(array, mxGetPr(mat), m * n));
 93:     PetscCall(PetscInfo(mfile, "Got MATLAB array %s\n", name));
 94:   }
 95:   PetscFunctionReturn(PETSC_SUCCESS);
 96: }

 98: static PetscErrorCode PetscViewerFileSetMode_Matlab(PetscViewer viewer, PetscFileMode type)
 99: {
100:   PetscViewer_Matlab *vmatlab = (PetscViewer_Matlab *)viewer->data;

102:   PetscFunctionBegin;
103:   vmatlab->btype = type;
104:   PetscFunctionReturn(PETSC_SUCCESS);
105: }

107: /*
108:         Actually opens the file
109: */
110: static PetscErrorCode PetscViewerFileSetName_Matlab(PetscViewer viewer, const char name[])
111: {
112:   PetscViewer_Matlab *vmatlab = (PetscViewer_Matlab *)viewer->data;
113:   PetscFileMode       type    = vmatlab->btype;

115:   PetscFunctionBegin;
116:   PetscCheck(type != (PetscFileMode)-1, PETSC_COMM_SELF, PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
117:   if (vmatlab->ep) matClose(vmatlab->ep);

119:   /* only first processor opens file */
120:   if (!vmatlab->rank) {
121:     if (type == FILE_MODE_READ) vmatlab->ep = matOpen(name, "r");
122:     else if (type == FILE_MODE_WRITE) vmatlab->ep = matOpen(name, "w");
123:     else {
124:       PetscCheck(type != FILE_MODE_UNDEFINED, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ORDER, "Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
125:       SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unsupported file mode %s", PetscFileModes[type]);
126:     }
127:   }
128:   PetscFunctionReturn(PETSC_SUCCESS);
129: }

131: static PetscErrorCode PetscViewerDestroy_Matlab(PetscViewer v)
132: {
133:   PetscViewer_Matlab *vf = (PetscViewer_Matlab *)v->data;

135:   PetscFunctionBegin;
136:   if (vf->ep) matClose(vf->ep);
137:   PetscCall(PetscFree(vf));
138:   PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetName_C", NULL));
139:   PetscCall(PetscObjectComposeFunction((PetscObject)v, "PetscViewerFileSetMode_C", NULL));
140:   PetscFunctionReturn(PETSC_SUCCESS);
141: }

143: /*MC
144:    PETSCVIEWERMATLAB - A viewer that saves the variables into a MATLAB .mat file that may be read into MATLAB
145:        with load('filename').

147:    Level: intermediate

149:        Notes:
150:              Currently can only save PETSc vectors to .mat files, not matrices (use the `PETSCVIEWERBINARY` and
151:              ${PETSC_DIR}/share/petsc/matlab/PetscBinaryRead.m to read matrices into MATLAB).

153:              For parallel vectors obtained with `DMCreateGlobalVector()` or `DMGetGlobalVector()` the vectors are saved to
154:              the .mat file in natural ordering. You can use DMView() to save the `DMDA` information to the .mat file
155:              the fields in the MATLAB loaded da variable give the array dimensions so you can reshape the MATLAB
156:              vector to the same multidimensional shape as it had in PETSc for plotting etc. For example,

158:              In your PETSc C/C++ code (assuming a two dimensional `DMDA` with one degree of freedom per node)
159: .vb
160:                 PetscObjectSetName((PetscObject)x,"x");
161:                 VecView(x,PETSC_VIEWER_MATLAB_WORLD);
162:                 PetscObjectSetName((PetscObject)da,"da");
163:                 DMView(x,PETSC_VIEWER_MATLAB_WORLD);
164: .ve
165:              Then from MATLAB
166: .vb
167:                 load('matlaboutput.mat')   % matlaboutput.mat is the default filename
168:                 xnew = zeros(da.n,da.m);
169:                 xnew(:) = x;    % reshape one dimensional vector back to two dimensions
170: .ve

172:               If you wish to put the same variable into the .mat file several times you need to give it a new
173:               name before each call to view.

175:               Use `PetscViewerMatlabPutArray()` to just put an array of doubles into the .mat file

177: .seealso: `PETSC_VIEWER_MATLAB_()`, `PETSC_VIEWER_MATLAB_SELF`, `PETSC_VIEWER_MATLAB_WORLD`, `PetscViewerCreate()`,
178:           `PetscViewerMatlabOpen()`, `VecView()`, `DMView()`, `PetscViewerMatlabPutArray()`, `PETSCVIEWERBINARY`, `PETSCVIEWERASCII`, `PETSCVIEWERDRAW`,
179:           `PETSC_VIEWER_STDOUT_()`, `PetscViewerFileSetName()`, `PetscViewerFileSetMode()`, `PetscViewerFormat`, `PetscMatlabEngine`
180: M*/
181: PETSC_EXTERN PetscErrorCode PetscViewerCreate_Matlab(PetscViewer viewer)
182: {
183:   PetscViewer_Matlab *e;

185:   PetscFunctionBegin;
186:   PetscCall(PetscNew(&e));
187:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &e->rank));
188:   e->btype     = FILE_MODE_UNDEFINED;
189:   viewer->data = (void *)e;

191:   PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetName_C", PetscViewerFileSetName_Matlab));
192:   PetscCall(PetscObjectComposeFunction((PetscObject)viewer, "PetscViewerFileSetMode_C", PetscViewerFileSetMode_Matlab));

194:   viewer->ops->destroy = PetscViewerDestroy_Matlab;
195:   PetscFunctionReturn(PETSC_SUCCESS);
196: }

198: /*@
199:   PetscViewerMatlabOpen - Opens a MATLAB .mat file for output

201:   Collective

203:   Input Parameters:
204: + comm - MPI communicator
205: . name - name of file
206: - type - type of file
207: .vb
208:     FILE_MODE_WRITE - create new file for MATLAB output
209:     FILE_MODE_READ - open existing file for MATLAB input
210:     FILE_MODE_WRITE - open existing file for MATLAB output
211: .ve

213:   Output Parameter:
214: . binv - PetscViewer for MATLAB output to use with the specified file

216:   Level: beginner

218:   Notes:
219:   This `PetscViewer` should be destroyed with `PetscViewerDestroy()`.

221:   For writing files it only opens the file on processor 0 in the communicator.

223:   This only saves `Vec`s it cannot be used to save `Mat`s. We recommend using the `PETSCVIEWERBINARY` to save objects to be loaded into MATLAB
224:   instead of this routine.

226:   PETSc must be configured with the option `--with-matlab` for this functionality

228: .seealso: `PETSCVIEWERMATLAB`, `PetscViewerASCIIOpen()`, `PetscViewerPushFormat()`, `PetscViewerDestroy()`, `PETSCVIEWERBINARY`, `PetscViewerBinaryOpen()`
229:           `VecView()`, `MatView()`, `VecLoad()`, `MatLoad()`
230: @*/
231: PetscErrorCode PetscViewerMatlabOpen(MPI_Comm comm, const char name[], PetscFileMode type, PetscViewer *binv)
232: {
233:   PetscFunctionBegin;
234:   PetscCall(PetscViewerCreate(comm, binv));
235:   PetscCall(PetscViewerSetType(*binv, PETSCVIEWERMATLAB));
236:   PetscCall(PetscViewerFileSetMode(*binv, type));
237:   PetscCall(PetscViewerFileSetName(*binv, name));
238:   PetscFunctionReturn(PETSC_SUCCESS);
239: }

241: static PetscMPIInt Petsc_Viewer_Matlab_keyval = MPI_KEYVAL_INVALID;

243: /*@C
244:    PETSC_VIEWER_MATLAB_ - Creates a `PETSCVIEWERMATLAB` `PetscViewer` shared by all processors
245:    in a communicator.

247:    Collective

249:    Input Parameter:
250: .  comm - the MPI communicator to share the MATLAB `PetscViewer`

252:    Options Database Key:
253: .    -viewer_matlab_filename <name> - name of the MATLAB file

255:    Environmental variable:
256: .  `PETSC_VIEWER_MATLAB_FILENAME` - name of the MATLAB file

258:    Level: intermediate

260:    Notes:
261:    This object is destroyed in `PetscFinalize()`, `PetscViewerDestroy()` should never be called on it

263:    Unlike almost all other PETSc routines, `PETSC_VIEWER_MATLAB_()` does not return
264:    an error code.  The MATLAB `PetscViewer` is usually used in the form `XXXView(XXX object, PETSC_VIEWER_MATLAB_(comm))`

266:    Use `PETSC_VIEWER_SOCKET_()` or `PetscViewerSocketOpen()` to communicator with an interactive MATLAB session.

268: .seealso: `PETSC_VIEWER_MATLAB_WORLD`, `PETSC_VIEWER_MATLAB_SELF`, `PetscViewerMatlabOpen()`, `PetscViewerCreate()`,
269:           `PetscViewerDestroy()`
270: @*/
271: PetscViewer PETSC_VIEWER_MATLAB_(MPI_Comm comm)
272: {
273:   PetscBool   flg;
274:   PetscViewer viewer;
275:   char        fname[PETSC_MAX_PATH_LEN];
276:   MPI_Comm    ncomm;

278:   PetscFunctionBegin;
279:   PetscCallNull(PetscCommDuplicate(comm, &ncomm, NULL));
280:   if (Petsc_Viewer_Matlab_keyval == MPI_KEYVAL_INVALID) { PetscCallMPINull(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, MPI_COMM_NULL_DELETE_FN, &Petsc_Viewer_Matlab_keyval, 0)); }
281:   PetscCallMPINull(MPI_Comm_get_attr(ncomm, Petsc_Viewer_Matlab_keyval, (void **)&viewer, (int *)&flg));
282:   if (!flg) { /* PetscViewer not yet created */
283:     PetscCallNull(PetscOptionsGetenv(ncomm, "PETSC_VIEWER_MATLAB_FILENAME", fname, PETSC_MAX_PATH_LEN, &flg));
284:     if (!flg) { PetscCallNull(PetscStrncpy(fname, "matlaboutput.mat", sizeof(fname))); }
285:     PetscCallNull(PetscViewerMatlabOpen(ncomm, fname, FILE_MODE_WRITE, &viewer));
286:     PetscCallNull(PetscObjectRegisterDestroy((PetscObject)viewer));
287:     PetscCallMPINull(MPI_Comm_set_attr(ncomm, Petsc_Viewer_Matlab_keyval, (void *)viewer));
288:   }
289:   PetscCallNull(PetscCommDestroy(&ncomm));
290:   PetscFunctionReturn(viewer);
291: }