Actual source code: vmatlab.c
petsc-3.9.4 2018-09-11
2: #include <petsc/private/viewerimpl.h>
3: #include <mat.h>
6: typedef struct {
7: MATFile *ep;
8: PetscMPIInt rank;
9: PetscFileMode btype;
10: } PetscViewer_Matlab;
12: /*@C
13: PetscViewerMatlabPutArray - Puts an array into the MATLAB viewer.
15: Not collective: only processor zero saves the array
17: Input Parameters:
18: + mfile - the viewer
19: . m,n - the dimensions of the array
20: . array - the array (represented in one dimension)
21: - name - the name of the array
23: Level: advanced
25: Notes: Only writes array values on processor 0.
27: @*/
28: PetscErrorCode PetscViewerMatlabPutArray(PetscViewer mfile,int m,int n,const PetscScalar *array,const char *name)
29: {
30: PetscErrorCode ierr;
31: PetscViewer_Matlab *ml = (PetscViewer_Matlab*)mfile->data;
32: mxArray *mat;
35: if (!ml->rank) {
36: PetscInfo1(mfile,"Putting MATLAB array %s\n",name);
37: #if !defined(PETSC_USE_COMPLEX)
38: mat = mxCreateDoubleMatrix(m,n,mxREAL);
39: #else
40: mat = mxCreateDoubleMatrix(m,n,mxCOMPLEX);
41: #endif
42: PetscMemcpy(mxGetPr(mat),array,m*n*sizeof(PetscScalar));
43: matPutVariable(ml->ep,name,mat);
45: PetscInfo1(mfile,"Put MATLAB array %s\n",name);
46: }
47: return(0);
48: }
50: PetscErrorCode PetscViewerMatlabPutVariable(PetscViewer viewer,const char *name,void *mat)
51: {
52: PetscViewer_Matlab *ml = (PetscViewer_Matlab*)viewer->data;
55: matPutVariable(ml->ep,name,(mxArray*)mat);
56: return(0);
57: }
59: /*@C
60: PetscViewerMatlabGetArray - Gets a variable from a MATLAB viewer into an array
62: Not Collective; only processor zero reads in the array
64: Input Parameters:
65: + mfile - the MATLAB file viewer
66: . m,n - the dimensions of the array
67: . array - the array (represented in one dimension)
68: - name - the name of the array
70: Level: advanced
72: Notes: Only reads in array values on processor 0.
74: @*/
75: PetscErrorCode PetscViewerMatlabGetArray(PetscViewer mfile,int m,int n,PetscScalar *array,const char *name)
76: {
77: PetscErrorCode ierr;
78: PetscViewer_Matlab *ml = (PetscViewer_Matlab*)mfile->data;
79: mxArray *mat;
82: if (!ml->rank) {
83: PetscInfo1(mfile,"Getting MATLAB array %s\n",name);
84: mat = matGetVariable(ml->ep,name);
85: if (!mat) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Unable to get array %s from matlab",name);
86: PetscMemcpy(array,mxGetPr(mat),m*n*sizeof(PetscScalar));
87: PetscInfo1(mfile,"Got MATLAB array %s\n",name);
88: }
89: return(0);
90: }
92: PetscErrorCode PetscViewerFileSetMode_Matlab(PetscViewer viewer,PetscFileMode type)
93: {
94: PetscViewer_Matlab *vmatlab = (PetscViewer_Matlab*)viewer->data;
97: vmatlab->btype = type;
98: return(0);
99: }
101: /*
102: Actually opens the file
103: */
104: PetscErrorCode PetscViewerFileSetName_Matlab(PetscViewer viewer,const char name[])
105: {
106: PetscViewer_Matlab *vmatlab = (PetscViewer_Matlab*)viewer->data;
107: PetscFileMode type = vmatlab->btype;
110: if (type == (PetscFileMode) -1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Must call PetscViewerFileSetMode() before PetscViewerFileSetName()");
111: if (vmatlab->ep) matClose(vmatlab->ep);
113: /* only first processor opens file */
114: if (!vmatlab->rank) {
115: if (type == FILE_MODE_READ) vmatlab->ep = matOpen(name,"r");
116: else if (type == FILE_MODE_WRITE || type == FILE_MODE_WRITE) vmatlab->ep = matOpen(name,"w");
117: else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Unknown file type");
118: }
119: return(0);
120: }
122: PetscErrorCode PetscViewerDestroy_Matlab(PetscViewer v)
123: {
124: PetscErrorCode ierr;
125: PetscViewer_Matlab *vf = (PetscViewer_Matlab*)v->data;
128: if (vf->ep) matClose(vf->ep);
129: PetscFree(vf);
130: return(0);
131: }
133: /*MC
134: PETSCVIEWERMATLAB - A viewer that saves the variables into a MATLAB .mat file that may be read into MATLAB
135: with load('filename').
137: Level: intermediate
139: Note: Currently can only save PETSc vectors to .mat files, not matrices (use the PETSCVIEWERBINARY and
140: ${PETSC_DIR}/share/petsc/matlab/PetscBinaryRead.m to read matrices into MATLAB).
142: For parallel vectors obtained with DMCreateGlobalVector() or DMGetGlobalVector() the vectors are saved to
143: the .mat file in natural ordering. You can use DMView() to save the DMDA information to the .mat file
144: the fields in the MATLAB loaded da variable give the array dimensions so you can reshape the MATLAB
145: vector to the same multidimensional shape as it had in PETSc for plotting etc. For example,
147: $ In your PETSc C/C++ code (assuming a two dimensional DMDA with one degree of freedom per node)
148: $ PetscObjectSetName((PetscObject)x,"x");
149: $ VecView(x,PETSC_VIEWER_MATLAB_WORLD);
150: $ PetscObjectSetName((PetscObject)da,"da");
151: $ DMView(x,PETSC_VIEWER_MATLAB_WORLD);
152: $ Then from MATLAB
153: $ load('matlaboutput.mat') % matlaboutput.mat is the default filename
154: $ xnew = zeros(da.n,da.m);
155: $ xnew(:) = x; % reshape one dimensional vector back to two dimensions
157: If you wish to put the same variable into the .mat file several times you need to give it a new
158: name before each call to view.
160: Use PetscViewerMatlabPutArray() to just put an array of doubles into the .mat file
162: .seealso: PETSC_VIEWER_MATLAB_(),PETSC_VIEWER_MATLAB_SELF, PETSC_VIEWER_MATLAB_WORLD,PetscViewerCreate(),
163: PetscViewerMatlabOpen(), VecView(), DMView(), PetscViewerMatlabPutArray(), PETSCVIEWERBINARY, PETSCVIEWERASCII, PETSCVIEWERDRAW,
164: PETSC_VIEWER_STDOUT_(), PetscViewerFileSetName(), PetscViewerFileSetMode(), PetscViewerFormat
166: M*/
167: PETSC_EXTERN PetscErrorCode PetscViewerCreate_Matlab(PetscViewer viewer)
168: {
169: PetscErrorCode ierr;
170: PetscViewer_Matlab *e;
173: PetscNewLog(viewer,&e);
174: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&e->rank);
175: e->btype = (PetscFileMode)-1;
176: viewer->data = (void*) e;
178: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",PetscViewerFileSetName_Matlab);
179: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetMode_C",PetscViewerFileSetMode_Matlab);
181: viewer->ops->destroy = PetscViewerDestroy_Matlab;
182: return(0);
183: }
185: /*@C
186: PetscViewerMatlabOpen - Opens a Matlab .mat file for output
188: Collective on MPI_Comm
190: Input Parameters:
191: + comm - MPI communicator
192: . name - name of file
193: - type - type of file
194: $ FILE_MODE_WRITE - create new file for MATLAB output
195: $ FILE_MODE_READ - open existing file for MATLAB input
196: $ FILE_MODE_WRITE - open existing file for MATLAB output
198: Output Parameter:
199: . binv - PetscViewer for MATLAB output to use with the specified file
201: Level: beginner
203: Note: This PetscViewer should be destroyed with PetscViewerDestroy().
205: For writing files it only opens the file on processor 0 in the communicator.
207: This only saves Vecs it cannot be used to save Mats. We recommend using the PETSCVIEWERBINARY to save objects to be loaded into MATLAB
208: instead of this routine.
210: Concepts: MATLAB .mat files
211: Concepts: PetscViewerMatlab^creating
213: .seealso: PetscViewerASCIIOpen(), PetscViewerPushFormat(), PetscViewerDestroy(), PETSCVIEWERBINARY, PetscViewerBinaryOpen()
214: VecView(), MatView(), VecLoad(), MatLoad()
215: @*/
216: PetscErrorCode PetscViewerMatlabOpen(MPI_Comm comm,const char name[],PetscFileMode type,PetscViewer *binv)
217: {
221: PetscViewerCreate(comm,binv);
222: PetscViewerSetType(*binv,PETSCVIEWERMATLAB);
223: PetscViewerFileSetMode(*binv,type);
224: PetscViewerFileSetName(*binv,name);
225: return(0);
226: }
228: static PetscMPIInt Petsc_Viewer_Matlab_keyval = MPI_KEYVAL_INVALID;
230: /*@C
231: PETSC_VIEWER_MATLAB_ - Creates a Matlab PetscViewer shared by all processors
232: in a communicator.
234: Collective on MPI_Comm
236: Input Parameter:
237: . comm - the MPI communicator to share the Matlab PetscViewer
239: Level: intermediate
241: Options Database Keys:
242: . -viewer_matlab_filename <name>
244: Environmental variables:
245: . PETSC_VIEWER_MATLAB_FILENAME
247: Notes:
248: Unlike almost all other PETSc routines, PETSC_VIEWER_MATLAB_ does not return
249: an error code. The matlab PetscViewer is usually used in the form
250: $ XXXView(XXX object,PETSC_VIEWER_MATLAB_(comm));
252: Use PETSC_VIEWER_SOCKET_() or PetscViewerSocketOpen() to communicator with an interactive MATLAB session.
254: .seealso: PETSC_VIEWER_MATLAB_WORLD, PETSC_VIEWER_MATLAB_SELF, PetscViewerMatlabOpen(), PetscViewerCreate(),
255: PetscViewerDestroy()
256: @*/
257: PetscViewer PETSC_VIEWER_MATLAB_(MPI_Comm comm)
258: {
260: PetscBool flg;
261: PetscViewer viewer;
262: char fname[PETSC_MAX_PATH_LEN];
263: MPI_Comm ncomm;
266: PetscCommDuplicate(comm,&ncomm,NULL);if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
267: if (Petsc_Viewer_Matlab_keyval == MPI_KEYVAL_INVALID) {
268: MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN,MPI_COMM_NULL_DELETE_FN,&Petsc_Viewer_Matlab_keyval,0);
269: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
270: }
271: MPI_Comm_get_attr(ncomm,Petsc_Viewer_Matlab_keyval,(void**)&viewer,(int*)&flg);
272: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
273: if (!flg) { /* PetscViewer not yet created */
274: PetscOptionsGetenv(ncomm,"PETSC_VIEWER_MATLAB_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
275: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
276: if (!flg) {
277: PetscStrcpy(fname,"matlaboutput.mat");
278: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
279: }
280: PetscViewerMatlabOpen(ncomm,fname,FILE_MODE_WRITE,&viewer);
281: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
282: PetscObjectRegisterDestroy((PetscObject)viewer);
283: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
284: MPI_Comm_set_attr(ncomm,Petsc_Viewer_Matlab_keyval,(void*)viewer);
285: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
286: }
287: PetscCommDestroy(&ncomm);
288: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_MATLAB_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
289: PetscFunctionReturn(viewer);
290: }