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