xref: /petsc/src/mat/interface/matrix.c (revision 2715565ab77787a22c210452f37a58b3f059f08a)
1 #ifndef lint
2 static char vcid[] = "$Id: matrix.c,v 1.235 1997/03/26 01:35:37 bsmith Exp balay $";
3 #endif
4 
5 /*
6    This is where the abstract matrix operations are defined
7 */
8 
9 #include "src/mat/matimpl.h"        /*I "mat.h" I*/
10 #include "src/vec/vecimpl.h"
11 #include "pinclude/pviewer.h"
12 
13 
14 
15 #undef __FUNC__
16 #define __FUNC__ "MatGetRow"
17 /*@C
18    MatGetRow - Gets a row of a matrix.  You MUST call MatRestoreRow()
19    for each row that you get to ensure that your application does
20    not bleed memory.
21 
22    Input Parameters:
23 .  mat - the matrix
24 .  row - the row to get
25 
26    Output Parameters:
27 .  ncols -  the number of nonzeros in the row
28 .  cols - if nonzero, the column numbers
29 .  vals - if nonzero, the values
30 
31    Notes:
32    This routine is provided for people who need to have direct access
33    to the structure of a matrix.  We hope that we provide enough
34    high-level matrix routines that few users will need it.
35 
36    For better efficiency, set cols and/or vals to PETSC_NULL if you do
37    not wish to extract these quantities.
38 
39    The user can only examine the values extracted with MatGetRow();
40    the values cannot be altered.  To change the matrix entries, one
41    must use MatSetValues().
42 
43    Caution:
44    Do not try to change the contents of the output arrays (cols and vals).
45    In some cases, this may corrupt the matrix.
46 
47 .keywords: matrix, row, get, extract
48 
49 .seealso: MatRestoreRow(), MatSetValues()
50 @*/
51 int MatGetRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals)
52 {
53   int   ierr;
54   PetscValidHeaderSpecific(mat,MAT_COOKIE);
55   PetscValidIntPointer(ncols);
56   if (!mat->assembled) SETERRQ(1,0,"Not for unassembled matrix");
57   if (mat->factor) SETERRQ(1,0,"Not for factored matrix");
58   if (!mat->ops.getrow) SETERRQ(PETSC_ERR_SUP,0,"");
59   PLogEventBegin(MAT_GetRow,mat,0,0,0);
60   ierr = (*mat->ops.getrow)(mat,row,ncols,cols,vals); CHKERRQ(ierr);
61   PLogEventEnd(MAT_GetRow,mat,0,0,0);
62   return 0;
63 }
64 
65 #undef __FUNC__
66 #define __FUNC__ "MatRestoreRow" /* ADIC Ignore */
67 /*@C
68    MatRestoreRow - Frees any temporary space allocated by MatGetRow().
69 
70    Input Parameters:
71 .  mat - the matrix
72 .  row - the row to get
73 .  ncols, cols - the number of nonzeros and their columns
74 .  vals - if nonzero the column values
75 
76 .keywords: matrix, row, restore
77 
78 .seealso:  MatGetRow()
79 @*/
80 int MatRestoreRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals)
81 {
82   PetscValidHeaderSpecific(mat,MAT_COOKIE);
83   PetscValidIntPointer(ncols);
84   if (!mat->assembled) SETERRQ(1,0,"Not for unassembled matrix");
85   if (!mat->ops.restorerow) return 0;
86   return (*mat->ops.restorerow)(mat,row,ncols,cols,vals);
87 }
88 
89 #undef __FUNC__
90 #define __FUNC__ "MatView" /* ADIC Ignore */
91 /*@C
92    MatView - Visualizes a matrix object.
93 
94    Input Parameters:
95 .  mat - the matrix
96 .  ptr - visualization context
97 
98    Notes:
99    The available visualization contexts include
100 $     VIEWER_STDOUT_SELF - standard output (default)
101 $     VIEWER_STDOUT_WORLD - synchronized standard
102 $       output where only the first processor opens
103 $       the file.  All other processors send their
104 $       data to the first processor to print.
105 $     VIEWER_DRAWX_WORLD - graphical display of nonzero structure
106 
107    The user can open alternative vistualization contexts with
108 $    ViewerFileOpenASCII() - output to a specified file
109 $    ViewerFileOpenBinary() - output in binary to a
110 $         specified file; corresponding input uses MatLoad()
111 $    ViewerDrawOpenX() - output nonzero matrix structure to
112 $         an X window display
113 $    ViewerMatlabOpen() - output matrix to Matlab viewer.
114 $         Currently only the sequential dense and AIJ
115 $         matrix types support the Matlab viewer.
116 
117    The user can call ViewerSetFormat() to specify the output
118    format of ASCII printed objects (when using VIEWER_STDOUT_SELF,
119    VIEWER_STDOUT_WORLD and ViewerFileOpenASCII).  Available formats include
120 $    VIEWER_FORMAT_ASCII_DEFAULT - default, prints matrix contents
121 $    VIEWER_FORMAT_ASCII_MATLAB - Matlab format
122 $    VIEWER_FORMAT_ASCII_IMPL - implementation-specific format
123 $      (which is in many cases the same as the default)
124 $    VIEWER_FORMAT_ASCII_INFO - basic information about the matrix
125 $      size and structure (not the matrix entries)
126 $    VIEWER_FORMAT_ASCII_INFO_LONG - more detailed information about the
127 $      matrix structure
128 
129 .keywords: matrix, view, visualize, output, print, write, draw
130 
131 .seealso: ViewerSetFormat(), ViewerFileOpenASCII(), ViewerDrawOpenX(),
132           ViewerMatlabOpen(), ViewerFileOpenBinary(), MatLoad()
133 @*/
134 int MatView(Mat mat,Viewer viewer)
135 {
136   int          format, ierr, rows, cols;
137   FILE         *fd;
138   char         *cstr;
139   ViewerType   vtype;
140   MPI_Comm     comm = mat->comm;
141 
142   PetscValidHeaderSpecific(mat,MAT_COOKIE);
143   if (viewer) PetscValidHeaderSpecific(viewer,VIEWER_COOKIE);
144   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
145 
146   if (!viewer) {
147     viewer = VIEWER_STDOUT_SELF;
148   }
149 
150   ierr = ViewerGetType(viewer,&vtype);
151   if (vtype == ASCII_FILE_VIEWER || vtype == ASCII_FILES_VIEWER) {
152     ierr = ViewerGetFormat(viewer,&format); CHKERRQ(ierr);
153     ierr = ViewerASCIIGetPointer(viewer,&fd); CHKERRQ(ierr);
154     if (format == VIEWER_FORMAT_ASCII_INFO || format == VIEWER_FORMAT_ASCII_INFO_LONG) {
155       PetscFPrintf(comm,fd,"Matrix Object:\n");
156       ierr = MatGetType(mat,PETSC_NULL,&cstr); CHKERRQ(ierr);
157       ierr = MatGetSize(mat,&rows,&cols); CHKERRQ(ierr);
158       PetscFPrintf(comm,fd,"  type=%s, rows=%d, cols=%d\n",cstr,rows,cols);
159       if (mat->ops.getinfo) {
160         MatInfo info;
161         ierr = MatGetInfo(mat,MAT_GLOBAL_SUM,&info); CHKERRQ(ierr);
162         PetscFPrintf(comm,fd,"  total: nonzeros=%d, allocated nonzeros=%d\n",
163                      (int)info.nz_used,(int)info.nz_allocated);
164       }
165     }
166   }
167   if (mat->view) {ierr = (*mat->view)((PetscObject)mat,viewer); CHKERRQ(ierr);}
168   return 0;
169 }
170 
171 #undef __FUNC__
172 #define __FUNC__ "MatDestroy" /* ADIC Ignore */
173 /*@C
174    MatDestroy - Frees space taken by a matrix.
175 
176    Input Parameter:
177 .  mat - the matrix
178 
179 .keywords: matrix, destroy
180 @*/
181 int MatDestroy(Mat mat)
182 {
183   int ierr;
184   PetscValidHeaderSpecific(mat,MAT_COOKIE);
185   ierr = (*mat->destroy)((PetscObject)mat); CHKERRQ(ierr);
186   return 0;
187 }
188 
189 #undef __FUNC__
190 #define __FUNC__ "MatValid" /* ADIC Ignore */
191 /*@
192    MatValid - Checks whether a matrix object is valid.
193 
194    Input Parameter:
195 .  m - the matrix to check
196 
197    Output Parameter:
198    flg - flag indicating matrix status, either
199 $     PETSC_TRUE if matrix is valid;
200 $     PETSC_FALSE otherwise.
201 
202 .keywords: matrix, valid
203 @*/
204 int MatValid(Mat m,PetscTruth *flg)
205 {
206   PetscValidIntPointer(flg);
207   if (!m)                           *flg = PETSC_FALSE;
208   else if (m->cookie != MAT_COOKIE) *flg = PETSC_FALSE;
209   else                              *flg = PETSC_TRUE;
210   return 0;
211 }
212 
213 #undef __FUNC__
214 #define __FUNC__ "MatSetValues"
215 /*@
216    MatSetValues - Inserts or adds a block of values into a matrix.
217    These values may be cached, so MatAssemblyBegin() and MatAssemblyEnd()
218    MUST be called after all calls to MatSetValues() have been completed.
219 
220    Input Parameters:
221 .  mat - the matrix
222 .  v - a logically two-dimensional array of values
223 .  m, idxm - the number of rows and their global indices
224 .  n, idxn - the number of columns and their global indices
225 .  addv - either ADD_VALUES or INSERT_VALUES, where
226 $     ADD_VALUES - adds values to any existing entries
227 $     INSERT_VALUES - replaces existing entries with new values
228 
229    Notes:
230    By default the values, v, are row-oriented and unsorted.
231    See MatSetOptions() for other options.
232 
233    Calls to MatSetValues() with the INSERT_VALUES and ADD_VALUES
234    options cannot be mixed without intervening calls to the assembly
235    routines.
236 
237    MatSetValues() uses 0-based row and column numbers in Fortran
238    as well as in C.
239 
240 .keywords: matrix, insert, add, set, values
241 
242 .seealso: MatSetOptions(), MatAssemblyBegin(), MatAssemblyEnd()
243 @*/
244 int MatSetValues(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v,InsertMode addv)
245 {
246   int ierr;
247   PetscValidHeaderSpecific(mat,MAT_COOKIE);
248   if (!m || !n) return 0; /* no values to insert */
249   PetscValidIntPointer(idxm);
250   PetscValidIntPointer(idxn);
251   PetscValidScalarPointer(v);
252   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
253   if (mat->insertmode == NOT_SET_VALUES) {
254     mat->insertmode = addv;
255   } else if (mat->insertmode != addv) {
256     SETERRQ(1,1,"Cannot mix add values and insert values");
257   }
258 
259   if (mat->assembled) {
260     mat->was_assembled = PETSC_TRUE;
261     mat->assembled     = PETSC_FALSE;
262   }
263   PLogEventBegin(MAT_SetValues,mat,0,0,0);
264   ierr = (*mat->ops.setvalues)(mat,m,idxm,n,idxn,v,addv);CHKERRQ(ierr);
265   PLogEventEnd(MAT_SetValues,mat,0,0,0);
266   return 0;
267 }
268 
269 #undef __FUNC__
270 #define __FUNC__ "MatSetValuesBlocked"
271 /*@
272    MatSetValuesBlocked - Inserts or adds a block of values into a matrix.
273 
274    Input Parameters:
275 .  mat - the matrix
276 .  v - a logically two-dimensional array of values
277 .  m, idxm - the number of block rows and their global block indices
278 .  n, idxn - the number of block columns and their global block indices
279 .  addv - either ADD_VALUES or INSERT_VALUES, where
280 $     ADD_VALUES - adds values to any existing entries
281 $     INSERT_VALUES - replaces existing entries with new values
282 
283    Notes:
284    By default the values, v, are row-oriented and unsorted.
285    See MatSetOptions() for other options.
286 
287    Calls to MatSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
288    options cannot be mixed without intervening calls to the assembly
289    routines.
290 
291    MatSetValuesBlocked() uses 0-based row and column numbers in Fortran
292    as well as in C.
293 
294    Restrictions:
295    MatSetValuesBlocked() is currently supported only for the block AIJ
296    matrix format (MATSEQBAIJ and MATMPIBAIJ).
297 
298 .keywords: matrix, insert, add, set, values
299 
300 .seealso: MatSetOptions(), MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues()
301 @*/
302 int MatSetValuesBlocked(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v,InsertMode addv)
303 {
304   int ierr;
305   PetscValidHeaderSpecific(mat,MAT_COOKIE);
306   if (!m || !n) return 0; /* no values to insert */
307   PetscValidIntPointer(idxm);
308   PetscValidIntPointer(idxn);
309   PetscValidScalarPointer(v);
310   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
311   if (mat->insertmode == NOT_SET_VALUES) {
312     mat->insertmode = addv;
313   } else if (mat->insertmode != addv) {
314     SETERRQ(1,1,"Cannot mix add values and insert values");
315   }
316 
317   if (mat->assembled) {
318     mat->was_assembled = PETSC_TRUE;
319     mat->assembled     = PETSC_FALSE;
320   }
321   PLogEventBegin(MAT_SetValues,mat,0,0,0);
322   ierr = (*mat->ops.setvaluesblocked)(mat,m,idxm,n,idxn,v,addv);CHKERRQ(ierr);
323   PLogEventEnd(MAT_SetValues,mat,0,0,0);
324   return 0;
325 }
326 
327 /*MC
328    MatSetValue - Set a single entry into a matrix.
329 
330    Input Parameters:
331 .  m - the matrix
332 .  row - the row location of the entry
333 .  col - the column location of the entry
334 .  value - the value to insert
335 .  mode - either INSERT_VALUES or ADD_VALUES
336 
337    Synopsis:
338    void MatSetValue(Mat m,int row,int col,Scalar value,InsertMode mode);
339 
340    Notes: For efficiency one should use MatSetValues() and set
341 several or many values simultaneously.
342 
343 .seealso: MatSetValues()
344 M*/
345 
346 #undef __FUNC__
347 #define __FUNC__ "MatGetValues"
348 /*@
349    MatGetValues - Gets a block of values from a matrix.
350 
351    Input Parameters:
352 .  mat - the matrix
353 .  v - a logically two-dimensional array for storing the values
354 .  m, idxm - the number of rows and their global indices
355 .  n, idxn - the number of columns and their global indices
356 
357    Notes:
358    The user must allocate space (m*n Scalars) for the values, v.
359    The values, v, are then returned in a row-oriented format,
360    analogous to that used by default in MatSetValues().
361 
362    MatGetValues() uses 0-based row and column numbers in
363    Fortran as well as in C.
364 
365    MatGetValues() requires that the matrix has been assembled
366    with MatAssemblyBegin()/MatAssemblyEnd().  Thus, calls to
367    MatSetValues() and MatGetValues() CANNOT be made in succession
368    without intermediate matrix assembly.
369 
370 .keywords: matrix, get, values
371 
372 .seealso: MatGetRow(), MatGetSubmatrices(), MatSetValues()
373 @*/
374 int MatGetValues(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v)
375 {
376   int ierr;
377 
378   PetscValidHeaderSpecific(mat,MAT_COOKIE);
379   PetscValidIntPointer(idxm);
380   PetscValidIntPointer(idxn);
381   PetscValidScalarPointer(v);
382   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
383   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
384   if (!mat->ops.getvalues) SETERRQ(PETSC_ERR_SUP,0,"");
385 
386   PLogEventBegin(MAT_GetValues,mat,0,0,0);
387   ierr = (*mat->ops.getvalues)(mat,m,idxm,n,idxn,v); CHKERRQ(ierr);
388   PLogEventEnd(MAT_GetValues,mat,0,0,0);
389   return 0;
390 }
391 
392 #undef __FUNC__
393 #define __FUNC__ "MatSetLocalToGlobalMapping" /* ADIC Ignore */
394 /*@
395    MatSetLocalToGlobalMapping - Sets a local numbering to global numbering used
396    by the routine MatSetValuesLocal() to allow users to insert matrices entries
397    using a local (per-processor) numbering.
398 
399    Input Parameters:
400 .  x - the matrix
401 .  n - number of local indices
402 .  indices - global index for each local index
403 
404 .keywords: matrix, set, values, local ordering
405 
406 .seealso:  MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetValuesLocal()
407 @*/
408 int MatSetLocalToGlobalMapping(Mat x, int n,int *indices)
409 {
410   int ierr;
411   PetscValidHeaderSpecific(x,MAT_COOKIE);
412   PetscValidIntPointer(indices);
413 
414   if (x->mapping) {
415     SETERRQ(1,0,"Mapping already set for matrix");
416   }
417 
418   ierr = ISLocalToGlobalMappingCreate(n,indices,&x->mapping);CHKERRQ(ierr);
419   return 0;
420 }
421 
422 #undef __FUNC__
423 #define __FUNC__ "MatSetValuesLocal"
424 /*@
425    MatSetValuesLocal - Inserts or adds values into certain locations of a matrix,
426         using a local ordering of the nodes.
427 
428    Input Parameters:
429 .  x - matrix to insert in
430 .  nrow - number of row elements to add
431 .  irow - row indices where to add
432 .  ncol - number of column elements to add
433 .  icol - column indices where to add
434 .  y - array of values
435 .  iora - either INSERT_VALUES or ADD_VALUES
436 
437    Notes:
438    Calls to MatSetValuesLocal() with the INSERT_VALUES and ADD_VALUES
439    options cannot be mixed without intervening calls to the assembly
440    routines.
441    These values may be cached, so MatAssemblyBegin() and MatAssemblyEnd()
442    MUST be called after all calls to MatSetValuesLocal() have been completed.
443 
444 .keywords: matrix, set, values, local ordering
445 
446 .seealso:  MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetLocalToGlobalMapping()
447 @*/
448 int MatSetValuesLocal(Mat mat,int nrow,int *irow,int ncol, int *icol,Scalar *y,InsertMode addv)
449 {
450   int ierr,irowm[128],icolm[128];
451 
452   PetscValidHeaderSpecific(mat,MAT_COOKIE);
453   PetscValidIntPointer(irow);
454   PetscValidIntPointer(icol);
455   PetscValidScalarPointer(y);
456   if (!mat->mapping) {
457     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Local to global never set with MatSetLocalToGlobalMapping");
458   }
459   if (nrow > 128 || ncol > 128) {
460     SETERRQ(PETSC_ERR_SUP,0,"Number indices must be <= 128");
461   }
462   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
463   if (mat->insertmode == NOT_SET_VALUES) {
464     mat->insertmode = addv;
465   } else if (mat->insertmode != addv) {
466     SETERRQ(1,1,"Cannot mix add values and insert values");
467   }
468 
469   if (mat->assembled) {
470     mat->was_assembled = PETSC_TRUE;
471     mat->assembled     = PETSC_FALSE;
472   }
473   PLogEventBegin(MAT_SetValues,mat,0,0,0);
474   ISLocalToGlobalMappingApply(mat->mapping,nrow,irow,irowm);
475   ISLocalToGlobalMappingApply(mat->mapping,ncol,icol,icolm);
476   ierr = (*mat->ops.setvalues)(mat,nrow,irowm,ncol,icolm,y,addv);CHKERRQ(ierr);
477   PLogEventEnd(MAT_SetValues,mat,0,0,0);
478   return 0;
479 }
480 
481 /* --------------------------------------------------------*/
482 #undef __FUNC__
483 #define __FUNC__ "MatMult"
484 /*@
485    MatMult - Computes the matrix-vector product, y = Ax.
486 
487    Input Parameters:
488 .  mat - the matrix
489 .  x   - the vector to be multilplied
490 
491    Output Parameters:
492 .  y - the result
493 
494    Notes:
495    The vectors x and y cannot be the same.  I.e., one cannot
496    call MatMult(A,y,y).
497 
498 .keywords: matrix, multiply, matrix-vector product
499 
500 .seealso: MatMultTrans(), MatMultAdd(), MatMultTransAdd()
501 @*/
502 int MatMult(Mat mat,Vec x,Vec y)
503 {
504   int ierr;
505   PetscValidHeaderSpecific(mat,MAT_COOKIE);
506   PetscValidHeaderSpecific(x,VEC_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE);
507   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
508   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
509   if (x == y) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"x and y must be different vectors");
510   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
511   if (mat->M != y->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim");
512   if (mat->m != y->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: local dim");
513 
514   PLogEventBegin(MAT_Mult,mat,x,y,0);
515   ierr = (*mat->ops.mult)(mat,x,y); CHKERRQ(ierr);
516   PLogEventEnd(MAT_Mult,mat,x,y,0);
517 
518   return 0;
519 }
520 
521 #undef __FUNC__
522 #define __FUNC__ "MatMultTrans"
523 /*@
524    MatMultTrans - Computes matrix transpose times a vector.
525 
526    Input Parameters:
527 .  mat - the matrix
528 .  x   - the vector to be multilplied
529 
530    Output Parameters:
531 .  y - the result
532 
533    Notes:
534    The vectors x and y cannot be the same.  I.e., one cannot
535    call MatMultTrans(A,y,y).
536 
537 .keywords: matrix, multiply, matrix-vector product, transpose
538 
539 .seealso: MatMult(), MatMultAdd(), MatMultTransAdd()
540 @*/
541 int MatMultTrans(Mat mat,Vec x,Vec y)
542 {
543   int ierr;
544   PetscValidHeaderSpecific(mat,MAT_COOKIE);
545   PetscValidHeaderSpecific(x,VEC_COOKIE); PetscValidHeaderSpecific(y,VEC_COOKIE);
546   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
547   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
548   if (x == y) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"x and y must be different vectors");
549   if (mat->M != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
550   if (mat->N != y->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim");
551   PLogEventBegin(MAT_MultTrans,mat,x,y,0);
552   ierr = (*mat->ops.multtrans)(mat,x,y); CHKERRQ(ierr);
553   PLogEventEnd(MAT_MultTrans,mat,x,y,0);
554   return 0;
555 }
556 
557 #undef __FUNC__
558 #define __FUNC__ "MatMultAdd"
559 /*@
560     MatMultAdd -  Computes v3 = v2 + A * v1.
561 
562     Input Parameters:
563 .   mat - the matrix
564 .   v1, v2 - the vectors
565 
566     Output Parameters:
567 .   v3 - the result
568 
569    Notes:
570    The vectors v1 and v3 cannot be the same.  I.e., one cannot
571    call MatMultAdd(A,v1,v2,v1).
572 
573 .keywords: matrix, multiply, matrix-vector product, add
574 
575 .seealso: MatMultTrans(), MatMult(), MatMultTransAdd()
576 @*/
577 int MatMultAdd(Mat mat,Vec v1,Vec v2,Vec v3)
578 {
579   int ierr;
580   PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v1,VEC_COOKIE);
581   PetscValidHeaderSpecific(v2,VEC_COOKIE); PetscValidHeaderSpecific(v3,VEC_COOKIE);
582   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
583   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
584   if (mat->N != v1->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim");
585   if (mat->M != v2->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim");
586   if (mat->M != v3->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim");
587   if (mat->m != v3->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: local dim");
588   if (mat->m != v2->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: local dim");
589   if (v1 == v3) SETERRQ(PETSC_ERR_ARG_IDN,0,"v1 and v3 must be different vectors");
590 
591   PLogEventBegin(MAT_MultAdd,mat,v1,v2,v3);
592   ierr = (*mat->ops.multadd)(mat,v1,v2,v3); CHKERRQ(ierr);
593   PLogEventEnd(MAT_MultAdd,mat,v1,v2,v3);
594   return 0;
595 }
596 
597 #undef __FUNC__
598 #define __FUNC__ "MatMultTransAdd"
599 /*@
600    MatMultTransAdd - Computes v3 = v2 + A' * v1.
601 
602    Input Parameters:
603 .  mat - the matrix
604 .  v1, v2 - the vectors
605 
606    Output Parameters:
607 .  v3 - the result
608 
609    Notes:
610    The vectors v1 and v3 cannot be the same.  I.e., one cannot
611    call MatMultTransAdd(A,v1,v2,v1).
612 
613 .keywords: matrix, multiply, matrix-vector product, transpose, add
614 
615 .seealso: MatMultTrans(), MatMultAdd(), MatMult()
616 @*/
617 int MatMultTransAdd(Mat mat,Vec v1,Vec v2,Vec v3)
618 {
619   int ierr;
620   PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v1,VEC_COOKIE);
621   PetscValidHeaderSpecific(v2,VEC_COOKIE);PetscValidHeaderSpecific(v3,VEC_COOKIE);
622   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
623   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
624   if (!mat->ops.multtransadd) SETERRQ(PETSC_ERR_SUP,0,"");
625   if (v1 == v3) SETERRQ(PETSC_ERR_ARG_IDN,0,"v1 and v3 must be different vectors");
626   if (mat->M != v1->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim");
627   if (mat->N != v2->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim");
628   if (mat->N != v3->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim");
629 
630   PLogEventBegin(MAT_MultTransAdd,mat,v1,v2,v3);
631   ierr = (*mat->ops.multtransadd)(mat,v1,v2,v3); CHKERRQ(ierr);
632   PLogEventEnd(MAT_MultTransAdd,mat,v1,v2,v3);
633   return 0;
634 }
635 /* ------------------------------------------------------------*/
636 #undef __FUNC__
637 #define __FUNC__ "MatGetInfo"  /* ADIC Ignore */
638 /*@C
639    MatGetInfo - Returns information about matrix storage (number of
640    nonzeros, memory, etc.).
641 
642    Input Parameters:
643 .  mat - the matrix
644 
645    Output Parameters:
646 .  flag - flag indicating the type of parameters to be returned
647 $    flag = MAT_LOCAL: local matrix
648 $    flag = MAT_GLOBAL_MAX: maximum over all processors
649 $    flag = MAT_GLOBAL_SUM: sum over all processors
650 .  info - matrix information context
651 
652    Notes:
653    The MatInfo context contains a variety of matrix data, including
654    number of nonzeros allocated and used, number of mallocs during
655    matrix assembly, etc.  Additional information for factored matrices
656    is provided (such as the fill ratio, number of mallocs during
657    factorization, etc.).  Much of this info is printed to STDOUT
658    when using the runtime options
659 $       -log_info -mat_view_info
660 
661    Example for C/C++ Users:
662    See the file $(PETSC_DIR)/include/mat.h for a complete list of
663    data within the MatInfo context.  For example,
664 $
665 $      MatInfo *info;
666 $      Mat     A;
667 $      double  mal, nz_a, nz_u;
668 $
669 $      MatGetInfo(A,MAT_LOCAL,&info);
670 $      mal  = info->mallocs;
671 $      nz_a = info->nz_allocated;
672 $
673 
674    Example for Fortran Users:
675    Fortran users should declare info as a double precision
676    array of dimension MAT_INFO_SIZE, and then extract the parameters
677    of interest.  See the file $(PETSC_DIR)/include/FINCLUDE/mat.h
678    a complete list of parameter names.
679 $
680 $      double  precision info(MAT_INFO_SIZE)
681 $      double  precision mal, nz_a
682 $      Mat     A
683 $      integer ierr
684 $
685 $      call MatGetInfo(A,MAT_LOCAL,info,ierr)
686 $      mal = info(MAT_INFO_MALLOCS)
687 $      nz_a = info(MAT_INFO_NZ_ALLOCATED)
688 $
689 
690 .keywords: matrix, get, info, storage, nonzeros, memory, fill
691 @*/
692 int MatGetInfo(Mat mat,MatInfoType flag,MatInfo *info)
693 {
694   PetscValidHeaderSpecific(mat,MAT_COOKIE);
695   PetscValidPointer(info);
696   if (!mat->ops.getinfo) SETERRQ(PETSC_ERR_SUP,0,"");
697   return  (*mat->ops.getinfo)(mat,flag,info);
698 }
699 
700 /* ----------------------------------------------------------*/
701 #undef __FUNC__
702 #define __FUNC__ "MatILUDTFactor"
703 /*@
704    MatILUDTFactor - Performs a drop tolerance ILU factorization.
705 
706    Input Parameters:
707 .  mat - the matrix
708 .  dt  - the drop tolerance
709 .  maxnz - the maximum number of nonzeros per row allowed?
710 .  row - row permutation
711 .  col - column permutation
712 
713    Output Parameters:
714 .  fact - the factored matrix
715 
716 .keywords: matrix, factor, LU, in-place
717 
718 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor()
719 @*/
720 int MatILUDTFactor(Mat mat,double dt,int maxnz,IS row,IS col,Mat *fact)
721 {
722   int ierr;
723   PetscValidHeaderSpecific(mat,MAT_COOKIE);
724   PetscValidPointer(fact);
725   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
726   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
727   if (!mat->ops.iludtfactor) SETERRQ(PETSC_ERR_SUP,0,"");
728 
729   PLogEventBegin(MAT_ILUFactor,mat,row,col,0);
730   ierr = (*mat->ops.iludtfactor)(mat,dt,maxnz,row,col,fact); CHKERRQ(ierr);
731   PLogEventEnd(MAT_ILUFactor,mat,row,col,0);
732 
733   return 0;
734 }
735 
736 #undef __FUNC__
737 #define __FUNC__ "MatLUFactor"
738 /*@
739    MatLUFactor - Performs in-place LU factorization of matrix.
740 
741    Input Parameters:
742 .  mat - the matrix
743 .  row - row permutation
744 .  col - column permutation
745 .  f - expected fill as ratio of original fill.
746 
747 .keywords: matrix, factor, LU, in-place
748 
749 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor()
750 @*/
751 int MatLUFactor(Mat mat,IS row,IS col,double f)
752 {
753   int ierr;
754   PetscValidHeaderSpecific(mat,MAT_COOKIE);
755   if (mat->M != mat->N) SETERRQ(1,0,"matrix must be square");
756   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
757   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
758   if (!mat->ops.lufactor) SETERRQ(PETSC_ERR_SUP,0,"");
759 
760   PLogEventBegin(MAT_LUFactor,mat,row,col,0);
761   ierr = (*mat->ops.lufactor)(mat,row,col,f); CHKERRQ(ierr);
762   PLogEventEnd(MAT_LUFactor,mat,row,col,0);
763   return 0;
764 }
765 
766 #undef __FUNC__
767 #define __FUNC__ "MatLUFactorSymbolic"
768 /*@
769    MatILUFactor - Performs in-place ILU factorization of matrix.
770 
771    Input Parameters:
772 .  mat - the matrix
773 .  row - row permutation
774 .  col - column permutation
775 .  f - expected fill as ratio of original fill.
776 .  level - number of levels of fill.
777 
778    Note: probably really only in-place when level is zero.
779 .keywords: matrix, factor, ILU, in-place
780 
781 .seealso: MatILUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor()
782 @*/
783 int MatILUFactor(Mat mat,IS row,IS col,double f,int level)
784 {
785   int ierr;
786   PetscValidHeaderSpecific(mat,MAT_COOKIE);
787   if (mat->M != mat->N) SETERRQ(1,0,"matrix must be square");
788   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
789   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
790   if (!mat->ops.ilufactor) SETERRQ(PETSC_ERR_SUP,0,"");
791 
792   PLogEventBegin(MAT_ILUFactor,mat,row,col,0);
793   ierr = (*mat->ops.ilufactor)(mat,row,col,f,level); CHKERRQ(ierr);
794   PLogEventEnd(MAT_ILUFactor,mat,row,col,0);
795   return 0;
796 }
797 
798 #undef __FUNC__
799 #define __FUNC__ "MatLUFactorSymbolic"
800 /*@
801    MatLUFactorSymbolic - Performs symbolic LU factorization of matrix.
802    Call this routine before calling MatLUFactorNumeric().
803 
804    Input Parameters:
805 .  mat - the matrix
806 .  row, col - row and column permutations
807 .  f - expected fill as ratio of the original number of nonzeros,
808        for example 3.0; choosing this parameter well can result in
809        more efficient use of time and space.
810 
811    Output Parameter:
812 .  fact - new matrix that has been symbolically factored
813 
814    Options Database Key:
815 $     -mat_lu_fill <f>, where f is the fill ratio
816 
817    Notes:
818    See the file $(PETSC_DIR)/Performace for additional information about
819    choosing the fill factor for better efficiency.
820 
821 .keywords: matrix, factor, LU, symbolic, fill
822 
823 .seealso: MatLUFactor(), MatLUFactorNumeric(), MatCholeskyFactor()
824 @*/
825 int MatLUFactorSymbolic(Mat mat,IS row,IS col,double f,Mat *fact)
826 {
827   int ierr,flg;
828   PetscValidHeaderSpecific(mat,MAT_COOKIE);
829   if (mat->M != mat->N) SETERRQ(1,0,"matrix must be square");
830   PetscValidPointer(fact);
831   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
832   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
833   if (!mat->ops.lufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"");
834 
835   ierr = OptionsGetDouble(PETSC_NULL,"-mat_lu_fill",&f,&flg); CHKERRQ(ierr);
836   PLogEventBegin(MAT_LUFactorSymbolic,mat,row,col,0);
837   ierr = (*mat->ops.lufactorsymbolic)(mat,row,col,f,fact); CHKERRQ(ierr);
838   PLogEventEnd(MAT_LUFactorSymbolic,mat,row,col,0);
839   return 0;
840 }
841 
842 #undef __FUNC__
843 #define __FUNC__ "MatLUFactorNumeric"
844 /*@
845    MatLUFactorNumeric - Performs numeric LU factorization of a matrix.
846    Call this routine after first calling MatLUFactorSymbolic().
847 
848    Input Parameters:
849 .  mat - the matrix
850 .  row, col - row and  column permutations
851 
852    Output Parameters:
853 .  fact - symbolically factored matrix that must have been generated
854           by MatLUFactorSymbolic()
855 
856    Notes:
857    See MatLUFactor() for in-place factorization.  See
858    MatCholeskyFactorNumeric() for the symmetric, positive definite case.
859 
860 .keywords: matrix, factor, LU, numeric
861 
862 .seealso: MatLUFactorSymbolic(), MatLUFactor(), MatCholeskyFactor()
863 @*/
864 int MatLUFactorNumeric(Mat mat,Mat *fact)
865 {
866   int ierr,flg;
867 
868   PetscValidHeaderSpecific(mat,MAT_COOKIE);
869   PetscValidPointer(fact);
870   PetscValidHeaderSpecific(*fact,MAT_COOKIE);
871   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
872   if (mat->M != (*fact)->M || mat->N != (*fact)->N)
873     SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dim");
874   if (!mat->ops.lufactornumeric) SETERRQ(PETSC_ERR_SUP,0,"");
875 
876   PLogEventBegin(MAT_LUFactorNumeric,mat,*fact,0,0);
877   ierr = (*mat->ops.lufactornumeric)(mat,fact); CHKERRQ(ierr);
878   PLogEventEnd(MAT_LUFactorNumeric,mat,*fact,0,0);
879   ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr);
880   if (flg) {
881     ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr);
882     if (flg) {
883       ViewerPushFormat(VIEWER_DRAWX_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr);
884     }
885     ierr = MatView(*fact,VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr);
886     ierr = ViewerFlush(VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr);
887     if (flg) {
888       ViewerPopFormat(VIEWER_DRAWX_(mat->comm));CHKERRQ(ierr);
889     }
890   }
891   return 0;
892 }
893 
894 #undef __FUNC__
895 #define __FUNC__ "MatCholeskyFactor"
896 /*@
897    MatCholeskyFactor - Performs in-place Cholesky factorization of a
898    symmetric matrix.
899 
900    Input Parameters:
901 .  mat - the matrix
902 .  perm - row and column permutations
903 .  f - expected fill as ratio of original fill
904 
905    Notes:
906    See MatLUFactor() for the nonsymmetric case.  See also
907    MatCholeskyFactorSymbolic(), and MatCholeskyFactorNumeric().
908 
909 .keywords: matrix, factor, in-place, Cholesky
910 
911 .seealso: MatLUFactor(), MatCholeskyFactorSymbolic(), MatCholeskyFactorNumeric()
912 @*/
913 int MatCholeskyFactor(Mat mat,IS perm,double f)
914 {
915   int ierr;
916   PetscValidHeaderSpecific(mat,MAT_COOKIE);
917   if (mat->M != mat->N) SETERRQ(1,0,"matrix must be square");
918   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
919   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
920   if (!mat->ops.choleskyfactor) SETERRQ(PETSC_ERR_SUP,0,"");
921 
922   PLogEventBegin(MAT_CholeskyFactor,mat,perm,0,0);
923   ierr = (*mat->ops.choleskyfactor)(mat,perm,f); CHKERRQ(ierr);
924   PLogEventEnd(MAT_CholeskyFactor,mat,perm,0,0);
925   return 0;
926 }
927 
928 #undef __FUNC__
929 #define __FUNC__ "MatCholeskyFactorSymbolic"
930 /*@
931    MatCholeskyFactorSymbolic - Performs symbolic Cholesky factorization
932    of a symmetric matrix.
933 
934    Input Parameters:
935 .  mat - the matrix
936 .  perm - row and column permutations
937 .  f - expected fill as ratio of original
938 
939    Output Parameter:
940 .  fact - the factored matrix
941 
942    Notes:
943    See MatLUFactorSymbolic() for the nonsymmetric case.  See also
944    MatCholeskyFactor() and MatCholeskyFactorNumeric().
945 
946 .keywords: matrix, factor, factorization, symbolic, Cholesky
947 
948 .seealso: MatLUFactorSymbolic(), MatCholeskyFactor(), MatCholeskyFactorNumeric()
949 @*/
950 int MatCholeskyFactorSymbolic(Mat mat,IS perm,double f,Mat *fact)
951 {
952   int ierr;
953   PetscValidHeaderSpecific(mat,MAT_COOKIE);
954   PetscValidPointer(fact);
955   if (mat->M != mat->N) SETERRQ(1,0,"matrix must be square");
956   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
957   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
958   if (!mat->ops.choleskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"");
959 
960   PLogEventBegin(MAT_CholeskyFactorSymbolic,mat,perm,0,0);
961   ierr = (*mat->ops.choleskyfactorsymbolic)(mat,perm,f,fact); CHKERRQ(ierr);
962   PLogEventEnd(MAT_CholeskyFactorSymbolic,mat,perm,0,0);
963   return 0;
964 }
965 
966 #undef __FUNC__
967 #define __FUNC__ "MatCholeskyFactorNumeric"
968 /*@
969    MatCholeskyFactorNumeric - Performs numeric Cholesky factorization
970    of a symmetric matrix. Call this routine after first calling
971    MatCholeskyFactorSymbolic().
972 
973    Input Parameter:
974 .  mat - the initial matrix
975 
976    Output Parameter:
977 .  fact - the factored matrix
978 
979 .keywords: matrix, factor, numeric, Cholesky
980 
981 .seealso: MatCholeskyFactorSymbolic(), MatCholeskyFactor(), MatLUFactorNumeric()
982 @*/
983 int MatCholeskyFactorNumeric(Mat mat,Mat *fact)
984 {
985   int ierr;
986   PetscValidHeaderSpecific(mat,MAT_COOKIE);
987   PetscValidPointer(fact);
988   if (!mat->ops.choleskyfactornumeric) SETERRQ(PETSC_ERR_SUP,0,"");
989   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
990   if (mat->M != (*fact)->M || mat->N != (*fact)->N)
991     SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dim");
992 
993   PLogEventBegin(MAT_CholeskyFactorNumeric,mat,*fact,0,0);
994   ierr = (*mat->ops.choleskyfactornumeric)(mat,fact); CHKERRQ(ierr);
995   PLogEventEnd(MAT_CholeskyFactorNumeric,mat,*fact,0,0);
996   return 0;
997 }
998 
999 /* ----------------------------------------------------------------*/
1000 #undef __FUNC__
1001 #define __FUNC__ "MatSolve"
1002 /*@
1003    MatSolve - Solves A x = b, given a factored matrix.
1004 
1005    Input Parameters:
1006 .  mat - the factored matrix
1007 .  b - the right-hand-side vector
1008 
1009    Output Parameter:
1010 .  x - the result vector
1011 
1012    Notes:
1013    The vectors b and x cannot be the same.  I.e., one cannot
1014    call MatSolve(A,x,x).
1015 
1016 .keywords: matrix, linear system, solve, LU, Cholesky, triangular solve
1017 
1018 .seealso: MatSolveAdd(), MatSolveTrans(), MatSolveTransAdd()
1019 @*/
1020 int MatSolve(Mat mat,Vec b,Vec x)
1021 {
1022   int ierr;
1023   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1024   PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE);
1025   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1026   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1027   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1028   if (mat->M != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1029   if (mat->m != b->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim");
1030 
1031   if (!mat->ops.solve) SETERRQ(PETSC_ERR_SUP,0,"");
1032   PLogEventBegin(MAT_Solve,mat,b,x,0);
1033   ierr = (*mat->ops.solve)(mat,b,x); CHKERRQ(ierr);
1034   PLogEventEnd(MAT_Solve,mat,b,x,0);
1035   return 0;
1036 }
1037 
1038 #undef __FUNC__
1039 #define __FUNC__ "MatForwardSolve"
1040 /* @
1041    MatForwardSolve - Solves L x = b, given a factored matrix, A = LU.
1042 
1043    Input Parameters:
1044 .  mat - the factored matrix
1045 .  b - the right-hand-side vector
1046 
1047    Output Parameter:
1048 .  x - the result vector
1049 
1050    Notes:
1051    MatSolve() should be used for most applications, as it performs
1052    a forward solve followed by a backward solve.
1053 
1054    The vectors b and x cannot be the same.  I.e., one cannot
1055    call MatForwardSolve(A,x,x).
1056 
1057 .keywords: matrix, forward, LU, Cholesky, triangular solve
1058 
1059 .seealso: MatSolve(), MatBackwardSolve()
1060 @ */
1061 int MatForwardSolve(Mat mat,Vec b,Vec x)
1062 {
1063   int ierr;
1064   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1065   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1066   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1067   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1068   if (!mat->ops.forwardsolve) SETERRQ(PETSC_ERR_SUP,0,"");
1069   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1070   if (mat->M != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1071   if (mat->m != b->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim");
1072 
1073   PLogEventBegin(MAT_ForwardSolve,mat,b,x,0);
1074   ierr = (*mat->ops.forwardsolve)(mat,b,x); CHKERRQ(ierr);
1075   PLogEventEnd(MAT_ForwardSolve,mat,b,x,0);
1076   return 0;
1077 }
1078 
1079 #undef __FUNC__
1080 #define __FUNC__ "MatBackwardSolve"
1081 /* @
1082    MatBackwardSolve - Solves U x = b, given a factored matrix, A = LU.
1083 
1084    Input Parameters:
1085 .  mat - the factored matrix
1086 .  b - the right-hand-side vector
1087 
1088    Output Parameter:
1089 .  x - the result vector
1090 
1091    Notes:
1092    MatSolve() should be used for most applications, as it performs
1093    a forward solve followed by a backward solve.
1094 
1095    The vectors b and x cannot be the same.  I.e., one cannot
1096    call MatBackwardSolve(A,x,x).
1097 
1098 .keywords: matrix, backward, LU, Cholesky, triangular solve
1099 
1100 .seealso: MatSolve(), MatForwardSolve()
1101 @ */
1102 int MatBackwardSolve(Mat mat,Vec b,Vec x)
1103 {
1104   int ierr;
1105   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1106   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1107   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1108   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1109   if (!mat->ops.backwardsolve) SETERRQ(PETSC_ERR_SUP,0,"");
1110   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1111   if (mat->M != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1112   if (mat->m != b->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim");
1113 
1114   PLogEventBegin(MAT_BackwardSolve,mat,b,x,0);
1115   ierr = (*mat->ops.backwardsolve)(mat,b,x); CHKERRQ(ierr);
1116   PLogEventEnd(MAT_BackwardSolve,mat,b,x,0);
1117   return 0;
1118 }
1119 
1120 #undef __FUNC__
1121 #define __FUNC__ "MatSolveAdd"
1122 /*@
1123    MatSolveAdd - Computes x = y + inv(A)*b, given a factored matrix.
1124 
1125    Input Parameters:
1126 .  mat - the factored matrix
1127 .  b - the right-hand-side vector
1128 .  y - the vector to be added to
1129 
1130    Output Parameter:
1131 .  x - the result vector
1132 
1133    Notes:
1134    The vectors b and x cannot be the same.  I.e., one cannot
1135    call MatSolveAdd(A,x,y,x).
1136 
1137 .keywords: matrix, linear system, solve, LU, Cholesky, add
1138 
1139 .seealso: MatSolve(), MatSolveTrans(), MatSolveTransAdd()
1140 @*/
1141 int MatSolveAdd(Mat mat,Vec b,Vec y,Vec x)
1142 {
1143   Scalar one = 1.0;
1144   Vec    tmp;
1145   int    ierr;
1146   PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE);
1147   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1148   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1149   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1150   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1151   if (mat->M != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1152   if (mat->M != y->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim");
1153   if (mat->m != b->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim");
1154   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim");
1155 
1156   PLogEventBegin(MAT_SolveAdd,mat,b,x,y);
1157   if (mat->ops.solveadd)  {
1158     ierr = (*mat->ops.solveadd)(mat,b,y,x); CHKERRQ(ierr);
1159   }
1160   else {
1161     /* do the solve then the add manually */
1162     if (x != y) {
1163       ierr = MatSolve(mat,b,x); CHKERRQ(ierr);
1164       ierr = VecAXPY(&one,y,x); CHKERRQ(ierr);
1165     }
1166     else {
1167       ierr = VecDuplicate(x,&tmp); CHKERRQ(ierr);
1168       PLogObjectParent(mat,tmp);
1169       ierr = VecCopy(x,tmp); CHKERRQ(ierr);
1170       ierr = MatSolve(mat,b,x); CHKERRQ(ierr);
1171       ierr = VecAXPY(&one,tmp,x); CHKERRQ(ierr);
1172       ierr = VecDestroy(tmp); CHKERRQ(ierr);
1173     }
1174   }
1175   PLogEventEnd(MAT_SolveAdd,mat,b,x,y);
1176   return 0;
1177 }
1178 
1179 #undef __FUNC__
1180 #define __FUNC__ "MatSolveTrans"
1181 /*@
1182    MatSolveTrans - Solves A' x = b, given a factored matrix.
1183 
1184    Input Parameters:
1185 .  mat - the factored matrix
1186 .  b - the right-hand-side vector
1187 
1188    Output Parameter:
1189 .  x - the result vector
1190 
1191    Notes:
1192    The vectors b and x cannot be the same.  I.e., one cannot
1193    call MatSolveTrans(A,x,x).
1194 
1195 .keywords: matrix, linear system, solve, LU, Cholesky, transpose
1196 
1197 .seealso: MatSolve(), MatSolveAdd(), MatSolveTransAdd()
1198 @*/
1199 int MatSolveTrans(Mat mat,Vec b,Vec x)
1200 {
1201   int ierr;
1202   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1203   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1204   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1205   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1206   if (!mat->ops.solvetrans) SETERRQ(PETSC_ERR_SUP,0,"");
1207   if (mat->M != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1208   if (mat->N != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1209 
1210   PLogEventBegin(MAT_SolveTrans,mat,b,x,0);
1211   ierr = (*mat->ops.solvetrans)(mat,b,x); CHKERRQ(ierr);
1212   PLogEventEnd(MAT_SolveTrans,mat,b,x,0);
1213   return 0;
1214 }
1215 
1216 #undef __FUNC__
1217 #define __FUNC__ "MatSolveTransAdd"
1218 /*@
1219    MatSolveTransAdd - Computes x = y + inv(trans(A)) b, given a
1220                       factored matrix.
1221 
1222    Input Parameters:
1223 .  mat - the factored matrix
1224 .  b - the right-hand-side vector
1225 .  y - the vector to be added to
1226 
1227    Output Parameter:
1228 .  x - the result vector
1229 
1230    Notes:
1231    The vectors b and x cannot be the same.  I.e., one cannot
1232    call MatSolveTransAdd(A,x,y,x).
1233 
1234 .keywords: matrix, linear system, solve, LU, Cholesky, transpose, add
1235 
1236 .seealso: MatSolve(), MatSolveAdd(), MatSolveTrans()
1237 @*/
1238 int MatSolveTransAdd(Mat mat,Vec b,Vec y,Vec x)
1239 {
1240   Scalar one = 1.0;
1241   int    ierr;
1242   Vec    tmp;
1243   PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE);
1244   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1245   if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors");
1246   if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix");
1247   if (mat->M != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1248   if (mat->N != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1249   if (mat->N != y->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim");
1250   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim");
1251 
1252   PLogEventBegin(MAT_SolveTransAdd,mat,b,x,y);
1253   if (mat->ops.solvetransadd) {
1254     ierr = (*mat->ops.solvetransadd)(mat,b,y,x); CHKERRQ(ierr);
1255   }
1256   else {
1257     /* do the solve then the add manually */
1258     if (x != y) {
1259       ierr = MatSolveTrans(mat,b,x); CHKERRQ(ierr);
1260       ierr = VecAXPY(&one,y,x); CHKERRQ(ierr);
1261     }
1262     else {
1263       ierr = VecDuplicate(x,&tmp); CHKERRQ(ierr);
1264       PLogObjectParent(mat,tmp);
1265       ierr = VecCopy(x,tmp); CHKERRQ(ierr);
1266       ierr = MatSolveTrans(mat,b,x); CHKERRQ(ierr);
1267       ierr = VecAXPY(&one,tmp,x); CHKERRQ(ierr);
1268       ierr = VecDestroy(tmp); CHKERRQ(ierr);
1269     }
1270   }
1271   PLogEventEnd(MAT_SolveTransAdd,mat,b,x,y);
1272   return 0;
1273 }
1274 /* ----------------------------------------------------------------*/
1275 
1276 #undef __FUNC__
1277 #define __FUNC__ "MatRelax"
1278 /*@
1279    MatRelax - Computes one relaxation sweep.
1280 
1281    Input Parameters:
1282 .  mat - the matrix
1283 .  b - the right hand side
1284 .  omega - the relaxation factor
1285 .  flag - flag indicating the type of SOR, one of
1286 $     SOR_FORWARD_SWEEP
1287 $     SOR_BACKWARD_SWEEP
1288 $     SOR_SYMMETRIC_SWEEP (SSOR method)
1289 $     SOR_LOCAL_FORWARD_SWEEP
1290 $     SOR_LOCAL_BACKWARD_SWEEP
1291 $     SOR_LOCAL_SYMMETRIC_SWEEP (local SSOR)
1292 $     SOR_APPLY_UPPER, SOR_APPLY_LOWER - applies
1293 $       upper/lower triangular part of matrix to
1294 $       vector (with omega)
1295 $     SOR_ZERO_INITIAL_GUESS - zero initial guess
1296 .  shift -  diagonal shift
1297 .  its - the number of iterations
1298 
1299    Output Parameters:
1300 .  x - the solution (can contain an initial guess)
1301 
1302    Notes:
1303    SOR_LOCAL_FORWARD_SWEEP, SOR_LOCAL_BACKWARD_SWEEP, and
1304    SOR_LOCAL_SYMMETRIC_SWEEP perform seperate independent smoothings
1305    on each processor.
1306 
1307    Application programmers will not generally use MatRelax() directly,
1308    but instead will employ the SLES/PC interface.
1309 
1310    Notes for Advanced Users:
1311    The flags are implemented as bitwise inclusive or operations.
1312    For example, use (SOR_ZERO_INITIAL_GUESS | SOR_SYMMETRIC_SWEEP)
1313    to specify a zero initial guess for SSOR.
1314 
1315 .keywords: matrix, relax, relaxation, sweep
1316 @*/
1317 int MatRelax(Mat mat,Vec b,double omega,MatSORType flag,double shift,
1318              int its,Vec x)
1319 {
1320   int ierr;
1321   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1322   PetscValidHeaderSpecific(b,VEC_COOKIE);  PetscValidHeaderSpecific(x,VEC_COOKIE);
1323   if (!mat->ops.relax) SETERRQ(PETSC_ERR_SUP,0,"");
1324   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1325   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1326   if (mat->N != x->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim");
1327   if (mat->M != b->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim");
1328   if (mat->m != b->n) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim");
1329 
1330   PLogEventBegin(MAT_Relax,mat,b,x,0);
1331   ierr =(*mat->ops.relax)(mat,b,omega,flag,shift,its,x); CHKERRQ(ierr);
1332   PLogEventEnd(MAT_Relax,mat,b,x,0);
1333   return 0;
1334 }
1335 
1336 #undef __FUNC__
1337 #define __FUNC__ "MatCopy_Basic"
1338 /*
1339       Default matrix copy routine.
1340 */
1341 int MatCopy_Basic(Mat A,Mat B)
1342 {
1343   int    ierr,i,rstart,rend,nz,*cwork;
1344   Scalar *vwork;
1345 
1346   ierr = MatZeroEntries(B); CHKERRQ(ierr);
1347   ierr = MatGetOwnershipRange(A,&rstart,&rend); CHKERRQ(ierr);
1348   for (i=rstart; i<rend; i++) {
1349     ierr = MatGetRow(A,i,&nz,&cwork,&vwork); CHKERRQ(ierr);
1350     ierr = MatSetValues(B,1,&i,nz,cwork,vwork,INSERT_VALUES); CHKERRQ(ierr);
1351     ierr = MatRestoreRow(A,i,&nz,&cwork,&vwork); CHKERRQ(ierr);
1352   }
1353   ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
1354   ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
1355   return 0;
1356 }
1357 
1358 #undef __FUNC__
1359 #define __FUNC__ "MatCopy"
1360 /*@C
1361    MatCopy - Copys a matrix to another matrix.
1362 
1363    Input Parameters:
1364 .  A - the matrix
1365 
1366    Output Parameter:
1367 .  B - where the copy is put
1368 
1369    Notes:
1370    MatCopy() copies the matrix entries of a matrix to another existing
1371    matrix (after first zeroing the second matrix).  A related routine is
1372    MatConvert(), which first creates a new matrix and then copies the data.
1373 
1374 .keywords: matrix, copy, convert
1375 
1376 .seealso: MatConvert()
1377 @*/
1378 int MatCopy(Mat A,Mat B)
1379 {
1380   int ierr;
1381   PetscValidHeaderSpecific(A,MAT_COOKIE); PetscValidHeaderSpecific(B,MAT_COOKIE);
1382   if (!A->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1383   if (A->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1384   if (A->M != B->M || A->N != B->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim");
1385 
1386   PLogEventBegin(MAT_Copy,A,B,0,0);
1387   if (A->ops.copy) {
1388     ierr = (*A->ops.copy)(A,B); CHKERRQ(ierr);
1389   }
1390   else { /* generic conversion */
1391     ierr = MatCopy_Basic(A,B); CHKERRQ(ierr);
1392   }
1393   PLogEventEnd(MAT_Copy,A,B,0,0);
1394   return 0;
1395 }
1396 
1397 static int MatConvertersSet = 0;
1398 static int (*MatConverters[MAX_MATRIX_TYPES][MAX_MATRIX_TYPES])(Mat,MatType,Mat*) =
1399            {{0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0},
1400             {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0},
1401             {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0},
1402             {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0},
1403             {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0},
1404             {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}};
1405 
1406 #undef __FUNC__
1407 #define __FUNC__ "MatConvertRegister" /* ADIC Ignore */
1408 /*@C
1409     MatConvertRegister - Allows one to register a routine that converts between
1410         two matrix types.
1411 
1412   Input Parameters:
1413 .   intype - the type of matrix (defined in include/mat.h), for example, MATSEQAIJ.
1414 .   outtype - new matrix type, or MATSAME
1415 
1416 .seealso: MatConvertRegisterAll()
1417 
1418 @*/
1419 int MatConvertRegister(MatType intype,MatType outtype,int (*converter)(Mat,MatType,Mat*))
1420 {
1421   MatConverters[intype][outtype] = converter;
1422   MatConvertersSet               = 1;
1423   return 0;
1424 }
1425 
1426 #undef __FUNC__
1427 #define __FUNC__ "MatConvert"
1428 /*@C
1429    MatConvert - Converts a matrix to another matrix, either of the same
1430    or different type.
1431 
1432    Input Parameters:
1433 .  mat - the matrix
1434 .  newtype - new matrix type.  Use MATSAME to create a new matrix of the
1435    same type as the original matrix.
1436 
1437    Output Parameter:
1438 .  M - pointer to place new matrix
1439 
1440    Notes:
1441    MatConvert() first creates a new matrix and then copies the data from
1442    the first matrix.  A related routine is MatCopy(), which copies the matrix
1443    entries of one matrix to another already existing matrix context.
1444 
1445 .keywords: matrix, copy, convert
1446 
1447 .seealso: MatCopy()
1448 @*/
1449 int MatConvert(Mat mat,MatType newtype,Mat *M)
1450 {
1451   int ierr;
1452   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1453   PetscValidPointer(M);
1454   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1455   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1456 
1457   if (newtype > MAX_MATRIX_TYPES || newtype < -1) {
1458     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,1,"Not a valid matrix type");
1459   }
1460   *M  = 0;
1461 
1462   if (!MatConvertersSet) {
1463     ierr = MatLoadRegisterAll(); CHKERRQ(ierr);
1464   }
1465 
1466   PLogEventBegin(MAT_Convert,mat,0,0,0);
1467   if ((newtype == mat->type || newtype == MATSAME) && mat->ops.convertsametype) {
1468     ierr = (*mat->ops.convertsametype)(mat,M,COPY_VALUES); CHKERRQ(ierr);
1469   } else {
1470     if (!MatConvertersSet) {
1471       ierr = MatConvertRegisterAll(); CHKERRQ(ierr);
1472     }
1473     if (!MatConverters[mat->type][newtype]) {
1474       SETERRQ(PETSC_ERR_ARG_WRONG,1,"Invalid matrix type, or matrix converter not registered");
1475     }
1476     ierr = (*MatConverters[mat->type][newtype])(mat,newtype,M); CHKERRQ(ierr);
1477   }
1478   PLogEventEnd(MAT_Convert,mat,0,0,0);
1479   return 0;
1480 }
1481 
1482 #undef __FUNC__
1483 #define __FUNC__ "MatGetDiagonal" /* ADIC Ignore */
1484 /*@
1485    MatGetDiagonal - Gets the diagonal of a matrix.
1486 
1487    Input Parameters:
1488 .  mat - the matrix
1489 .  v - the vector for storing the diagonal
1490 
1491    Output Parameter:
1492 .  v - the diagonal of the matrix
1493 
1494    Notes: For the SeqAIJ matrix format, this routine may also be called
1495     on a LU factored matrix; in that case it routines the reciprocal of
1496     the diagonal entries in U. It returns the entries permuted by the
1497     row and column permutation used during the symbolic factorization.
1498 
1499 .keywords: matrix, get, diagonal
1500 @*/
1501 int MatGetDiagonal(Mat mat,Vec v)
1502 {
1503   PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v,VEC_COOKIE);
1504   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1505   /*
1506      The error checking for a factored matrix is handled inside
1507     each matrix type, since MatGetDiagonal() is supported by
1508     factored AIJ matrices
1509   */
1510   /* if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");  */
1511   if (!mat->ops.getdiagonal) SETERRQ(PETSC_ERR_SUP,0,"");
1512   return (*mat->ops.getdiagonal)(mat,v);
1513 }
1514 
1515 #undef __FUNC__
1516 #define __FUNC__ "MatTranspose"
1517 /*@C
1518    MatTranspose - Computes an in-place or out-of-place transpose of a matrix.
1519 
1520    Input Parameter:
1521 .  mat - the matrix to transpose
1522 
1523    Output Parameters:
1524 .  B - the transpose (or pass in PETSC_NULL for an in-place transpose)
1525 
1526 .keywords: matrix, transpose
1527 
1528 .seealso: MatMultTrans(), MatMultTransAdd()
1529 @*/
1530 int MatTranspose(Mat mat,Mat *B)
1531 {
1532   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1533   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1534   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1535   if (!mat->ops.transpose) SETERRQ(PETSC_ERR_SUP,0,"");
1536   return (*mat->ops.transpose)(mat,B);
1537 }
1538 
1539 #undef __FUNC__
1540 #define __FUNC__ "MatPermute"
1541 /*@C
1542    MatPermute - Creates a new matrix with rows and columns permuted from the
1543        original.
1544 
1545    Input Parameter:
1546 .  mat - the matrix to permute
1547 .  row - row permutation
1548 .  col - column permutation
1549 
1550    Output Parameters:
1551 .  B - the permuted matrix
1552 
1553 .keywords: matrix, transpose
1554 
1555 .seealso: MatGetReordering()
1556 @*/
1557 int MatPermute(Mat mat,IS row,IS col,Mat *B)
1558 {
1559   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1560   PetscValidHeaderSpecific(row,IS_COOKIE);
1561   PetscValidHeaderSpecific(col,IS_COOKIE);
1562   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1563   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1564   if (!mat->ops.permute) SETERRQ(PETSC_ERR_SUP,0,"");
1565   return (*mat->ops.permute)(mat,row,col,B);
1566 }
1567 
1568 #undef __FUNC__
1569 #define __FUNC__ "MatEqual"
1570 /*@
1571    MatEqual - Compares two matrices.
1572 
1573    Input Parameters:
1574 .  A - the first matrix
1575 .  B - the second matrix
1576 
1577    Output Parameter:
1578 .  flg : PETSC_TRUE if the matrices are equal;
1579          PETSC_FALSE otherwise.
1580 
1581 .keywords: matrix, equal, equivalent
1582 @*/
1583 int MatEqual(Mat A,Mat B,PetscTruth *flg)
1584 {
1585   PetscValidHeaderSpecific(A,MAT_COOKIE); PetscValidHeaderSpecific(B,MAT_COOKIE);
1586   PetscValidIntPointer(flg);
1587   if (!A->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1588   if (!B->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1589   if (A->M != B->M || A->N != B->N) SETERRQ(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim");
1590   if (!A->ops.equal) SETERRQ(PETSC_ERR_SUP,0,"");
1591   return (*A->ops.equal)(A,B,flg);
1592 }
1593 
1594 #undef __FUNC__
1595 #define __FUNC__ "MatDiagonalScale"
1596 /*@
1597    MatDiagonalScale - Scales a matrix on the left and right by diagonal
1598    matrices that are stored as vectors.  Either of the two scaling
1599    matrices can be PETSC_NULL.
1600 
1601    Input Parameters:
1602 .  mat - the matrix to be scaled
1603 .  l - the left scaling vector (or PETSC_NULL)
1604 .  r - the right scaling vector (or PETSC_NULL)
1605 
1606    Notes:
1607    MatDiagonalScale() computes A <- LAR, where
1608 $      L = a diagonal matrix
1609 $      R = a diagonal matrix
1610 
1611 .keywords: matrix, diagonal, scale
1612 
1613 .seealso: MatDiagonalScale()
1614 @*/
1615 int MatDiagonalScale(Mat mat,Vec l,Vec r)
1616 {
1617   int ierr;
1618   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1619   if (!mat->ops.diagonalscale) SETERRQ(PETSC_ERR_SUP,0,"");
1620   if (l) PetscValidHeaderSpecific(l,VEC_COOKIE);
1621   if (r) PetscValidHeaderSpecific(r,VEC_COOKIE);
1622   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1623   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1624 
1625   PLogEventBegin(MAT_Scale,mat,0,0,0);
1626   ierr = (*mat->ops.diagonalscale)(mat,l,r); CHKERRQ(ierr);
1627   PLogEventEnd(MAT_Scale,mat,0,0,0);
1628   return 0;
1629 }
1630 
1631 #undef __FUNC__
1632 #define __FUNC__ "MatScale"
1633 /*@
1634     MatScale - Scales all elements of a matrix by a given number.
1635 
1636     Input Parameters:
1637 .   mat - the matrix to be scaled
1638 .   a  - the scaling value
1639 
1640     Output Parameter:
1641 .   mat - the scaled matrix
1642 
1643 .keywords: matrix, scale
1644 
1645 .seealso: MatDiagonalScale()
1646 @*/
1647 int MatScale(Scalar *a,Mat mat)
1648 {
1649   int ierr;
1650   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1651   PetscValidScalarPointer(a);
1652   if (!mat->ops.scale) SETERRQ(PETSC_ERR_SUP,0,"");
1653   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1654   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1655 
1656   PLogEventBegin(MAT_Scale,mat,0,0,0);
1657   ierr = (*mat->ops.scale)(a,mat); CHKERRQ(ierr);
1658   PLogEventEnd(MAT_Scale,mat,0,0,0);
1659   return 0;
1660 }
1661 
1662 #undef __FUNC__
1663 #define __FUNC__ "MatNorm"
1664 /*@
1665    MatNorm - Calculates various norms of a matrix.
1666 
1667    Input Parameters:
1668 .  mat - the matrix
1669 .  type - the type of norm, NORM_1, NORM_2, NORM_FROBENIUS, NORM_INFINITY
1670 
1671    Output Parameters:
1672 .  norm - the resulting norm
1673 
1674 .keywords: matrix, norm, Frobenius
1675 @*/
1676 int MatNorm(Mat mat,NormType type,double *norm)
1677 {
1678   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1679   PetscValidScalarPointer(norm);
1680 
1681   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1682   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1683   if (!mat->ops.norm) SETERRQ(PETSC_ERR_SUP,0,"Not for this matrix type");
1684   return (*mat->ops.norm)(mat,type,norm);
1685 }
1686 
1687 /*
1688      This variable is used to prevent counting of MatAssemblyBegin() that
1689    are called from within a MatAssemblyEnd().
1690 */
1691 static int MatAssemblyEnd_InUse = 0;
1692 #undef __FUNC__
1693 #define __FUNC__ "MatAssemblyBegin"
1694 /*@
1695    MatAssemblyBegin - Begins assembling the matrix.  This routine should
1696    be called after completing all calls to MatSetValues().
1697 
1698    Input Parameters:
1699 .  mat - the matrix
1700 .  type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY
1701 
1702    Notes:
1703    MatSetValues() generally caches the values.  The matrix is ready to
1704    use only after MatAssemblyBegin() and MatAssemblyEnd() have been called.
1705    Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES
1706    in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before
1707    using the matrix.
1708 
1709 .keywords: matrix, assembly, assemble, begin
1710 
1711 .seealso: MatAssemblyEnd(), MatSetValues()
1712 @*/
1713 int MatAssemblyBegin(Mat mat,MatAssemblyType type)
1714 {
1715   int ierr;
1716   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1717   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1718   if (mat->assembled) {
1719     mat->was_assembled = PETSC_TRUE;
1720     mat->assembled     = PETSC_FALSE;
1721   }
1722   if (!MatAssemblyEnd_InUse) {
1723     PLogEventBegin(MAT_AssemblyBegin,mat,0,0,0);
1724     if (mat->ops.assemblybegin){ierr = (*mat->ops.assemblybegin)(mat,type);CHKERRQ(ierr);}
1725     PLogEventEnd(MAT_AssemblyBegin,mat,0,0,0);
1726   } else {
1727     if (mat->ops.assemblybegin){ierr = (*mat->ops.assemblybegin)(mat,type);CHKERRQ(ierr);}
1728   }
1729   return 0;
1730 }
1731 
1732 
1733 #undef __FUNC__
1734 #define __FUNC__ "MatView_Private"
1735 /*
1736     Processes command line options to determine if/how a matrix
1737   is to be viewed.
1738 */
1739 int MatView_Private(Mat mat)
1740 {
1741   int ierr,flg;
1742 
1743   ierr = OptionsHasName(PETSC_NULL,"-mat_view_info",&flg); CHKERRQ(ierr);
1744   if (flg) {
1745     Viewer viewer;
1746     ierr = ViewerFileOpenASCII(mat->comm,"stdout",&viewer);CHKERRQ(ierr);
1747     ierr = ViewerSetFormat(viewer,VIEWER_FORMAT_ASCII_INFO,0);CHKERRQ(ierr);
1748     ierr = MatView(mat,viewer); CHKERRQ(ierr);
1749     ierr = ViewerDestroy(viewer); CHKERRQ(ierr);
1750   }
1751   ierr = OptionsHasName(PETSC_NULL,"-mat_view_info_detailed",&flg);CHKERRQ(ierr);
1752   if (flg) {
1753     Viewer viewer;
1754     ierr = ViewerFileOpenASCII(mat->comm,"stdout",&viewer);CHKERRQ(ierr);
1755     ierr = ViewerSetFormat(viewer,VIEWER_FORMAT_ASCII_INFO_LONG,0);CHKERRQ(ierr);
1756     ierr = MatView(mat,viewer); CHKERRQ(ierr);
1757     ierr = ViewerDestroy(viewer); CHKERRQ(ierr);
1758   }
1759   ierr = OptionsHasName(PETSC_NULL,"-mat_view",&flg); CHKERRQ(ierr);
1760   if (flg) {
1761     Viewer viewer;
1762     ierr = ViewerFileOpenASCII(mat->comm,"stdout",&viewer);CHKERRQ(ierr);
1763     ierr = MatView(mat,viewer); CHKERRQ(ierr);
1764     ierr = ViewerDestroy(viewer); CHKERRQ(ierr);
1765   }
1766   ierr = OptionsHasName(PETSC_NULL,"-mat_view_matlab",&flg); CHKERRQ(ierr);
1767   if (flg) {
1768     Viewer viewer;
1769     ierr = ViewerFileOpenASCII(mat->comm,"stdout",&viewer);CHKERRQ(ierr);
1770     ierr = ViewerSetFormat(viewer,VIEWER_FORMAT_ASCII_MATLAB,"M");CHKERRQ(ierr);
1771     ierr = MatView(mat,viewer); CHKERRQ(ierr);
1772     ierr = ViewerDestroy(viewer); CHKERRQ(ierr);
1773   }
1774   ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr);
1775   if (flg) {
1776     ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr);
1777     if (flg) {
1778       ViewerPushFormat(VIEWER_DRAWX_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr);
1779     }
1780     ierr = MatView(mat,VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr);
1781     ierr = ViewerFlush(VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr);
1782     if (flg) {
1783       ViewerPopFormat(VIEWER_DRAWX_(mat->comm));CHKERRQ(ierr);
1784     }
1785   }
1786   return 0;
1787 }
1788 
1789 #undef __FUNC__
1790 #define __FUNC__ "MatAssemblyEnd"
1791 /*@
1792    MatAssemblyEnd - Completes assembling the matrix.  This routine should
1793    be called after MatAssemblyBegin().
1794 
1795    Input Parameters:
1796 .  mat - the matrix
1797 .  type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY
1798 
1799    Options Database Keys:
1800 $  -mat_view_info : Prints info on matrix at
1801 $      conclusion of MatEndAssembly()
1802 $  -mat_view_info_detailed: Prints more detailed info.
1803 $  -mat_view : Prints matrix in ASCII format.
1804 $  -mat_view_matlab : Prints matrix in Matlab format.
1805 $  -mat_view_draw : Draws nonzero structure of matrix,
1806 $      using MatView() and DrawOpenX().
1807 $  -display <name> : Set display name (default is host)
1808 $  -draw_pause <sec> : Set number of seconds to pause after display
1809 
1810    Notes:
1811    MatSetValues() generally caches the values.  The matrix is ready to
1812    use only after MatAssemblyBegin() and MatAssemblyEnd() have been called.
1813    Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES
1814    in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before
1815    using the matrix.
1816 
1817 .keywords: matrix, assembly, assemble, end
1818 
1819 .seealso: MatAssemblyBegin(), MatSetValues(), DrawOpenX(), MatView()
1820 @*/
1821 int MatAssemblyEnd(Mat mat,MatAssemblyType type)
1822 {
1823   int        ierr;
1824   static int inassm = 0;
1825 
1826   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1827 
1828   inassm++;
1829   MatAssemblyEnd_InUse++;
1830   PLogEventBegin(MAT_AssemblyEnd,mat,0,0,0);
1831   if (mat->ops.assemblyend) {
1832     ierr = (*mat->ops.assemblyend)(mat,type); CHKERRQ(ierr);
1833   }
1834 
1835   /* Flush assembly is not a true assembly */
1836   if (type != MAT_FLUSH_ASSEMBLY) {
1837     mat->assembled  = PETSC_TRUE; mat->num_ass++;
1838   }
1839   mat->insertmode = NOT_SET_VALUES;
1840   PLogEventEnd(MAT_AssemblyEnd,mat,0,0,0);
1841   MatAssemblyEnd_InUse--;
1842 
1843   if (inassm == 1 && type != MAT_FLUSH_ASSEMBLY) {
1844     ierr = MatView_Private(mat); CHKERRQ(ierr);
1845   }
1846   inassm--;
1847   return 0;
1848 }
1849 
1850 
1851 #undef __FUNC__
1852 #define __FUNC__ "MatCompress"
1853 /*@
1854    MatCompress - Tries to store the matrix in as little space as
1855    possible.  May fail if memory is already fully used, since it
1856    tries to allocate new space.
1857 
1858    Input Parameters:
1859 .  mat - the matrix
1860 
1861 .keywords: matrix, compress
1862 @*/
1863 int MatCompress(Mat mat)
1864 {
1865   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1866   if (mat->ops.compress) return (*mat->ops.compress)(mat);
1867   return 0;
1868 }
1869 
1870 #undef __FUNC__
1871 #define __FUNC__ "MatSetOption" /* ADIC Ignore */
1872 /*@
1873    MatSetOption - Sets a parameter option for a matrix. Some options
1874    may be specific to certain storage formats.  Some options
1875    determine how values will be inserted (or added). Sorted,
1876    row-oriented input will generally assemble the fastest. The default
1877    is row-oriented, nonsorted input.
1878 
1879    Input Parameters:
1880 .  mat - the matrix
1881 .  option - the option, one of the following:
1882 $    MAT_ROW_ORIENTED
1883 $    MAT_COLUMN_ORIENTED,
1884 $    MAT_ROWS_SORTED,
1885 $    MAT_ROWS_UNSORTED,
1886 $    MAT_COLUMNS_SORTED,
1887 $    MAT_COLUMNS_UNSORTED,
1888 $    MAT_NO_NEW_NONZERO_LOCATIONS,
1889 $    MAT_YES_NEW_NONZERO_LOCATIONS,
1890 $    MAT_SYMMETRIC,
1891 $    MAT_STRUCTURALLY_SYMMETRIC,
1892 $    MAT_NO_NEW_DIAGONALS,
1893 $    MAT_YES_NEW_DIAGONALS,
1894 $    MAT_IGNORE_OFF_PROCESSOR_ENTRIES
1895 $    MAT_NEW_NONZERO_LOCATION_ERROR
1896 $    and possibly others.
1897 
1898    Notes:
1899    Some options are relevant only for particular matrix types and
1900    are thus ignored by others.  Other options are not supported by
1901    certain matrix types and will generate an error message if set.
1902 
1903    If using a Fortran 77 module to compute a matrix, one may need to
1904    use the column-oriented option (or convert to the row-oriented
1905    format).
1906 
1907    MAT_NO_NEW_NONZERO_LOCATIONS indicates that any add or insertion
1908    that will generate a new entry in the nonzero structure is ignored.
1909    Thus, if memory has not alredy been allocated for this particular
1910    data, then the insertion is ignored. For dense matrices, in which
1911    the entire array is allocated, no entries are ever ignored.
1912 
1913    MAT_NEW_NONZERO_LOCATION_ERROR indicates any add or insertion
1914    that will generate a new entry in the nonzero structure generates
1915    an error. Supported for AIJ and BAIJ formats.
1916 
1917    MAT_IGNORE_OFF_PROCESSOR_ENTRIES indicates entries destined for
1918    other processors are dropped, rather than stashed.
1919 
1920 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero
1921 @*/
1922 int MatSetOption(Mat mat,MatOption op)
1923 {
1924   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1925   if (mat->ops.setoption) return (*mat->ops.setoption)(mat,op);
1926   return 0;
1927 }
1928 
1929 #undef __FUNC__
1930 #define __FUNC__ "MatZeroEntries"
1931 /*@
1932    MatZeroEntries - Zeros all entries of a matrix.  For sparse matrices
1933    this routine retains the old nonzero structure.
1934 
1935    Input Parameters:
1936 .  mat - the matrix
1937 
1938 .keywords: matrix, zero, entries
1939 
1940 .seealso: MatZeroRows()
1941 @*/
1942 int MatZeroEntries(Mat mat)
1943 {
1944   int ierr;
1945   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1946   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1947   if (!mat->ops.zeroentries) SETERRQ(PETSC_ERR_SUP,0,"");
1948 
1949   PLogEventBegin(MAT_ZeroEntries,mat,0,0,0);
1950   ierr = (*mat->ops.zeroentries)(mat); CHKERRQ(ierr);
1951   PLogEventEnd(MAT_ZeroEntries,mat,0,0,0);
1952   return 0;
1953 }
1954 
1955 #undef __FUNC__
1956 #define __FUNC__ "MatZeroRows"
1957 /*@
1958    MatZeroRows - Zeros all entries (except possibly the main diagonal)
1959    of a set of rows of a matrix.
1960 
1961    Input Parameters:
1962 .  mat - the matrix
1963 .  is - index set of rows to remove
1964 .  diag - pointer to value put in all diagonals of eliminated rows.
1965           Note that diag is not a pointer to an array, but merely a
1966           pointer to a single value.
1967 
1968    Notes:
1969    For the AIJ matrix formats this removes the old nonzero structure,
1970    but does not release memory.  For the dense and block diagonal
1971    formats this does not alter the nonzero structure.
1972 
1973    The user can set a value in the diagonal entry (or for the AIJ and
1974    row formats can optionally remove the main diagonal entry from the
1975    nonzero structure as well, by passing a null pointer as the final
1976    argument).
1977 
1978    For the parallel case, all processes that share the matrix (i.e.,
1979    those in the communicator used for matrix creation) MUST call this
1980    routine, regardless of whether any rows being zeroed are owned by
1981    them.
1982 
1983 .keywords: matrix, zero, rows, boundary conditions
1984 
1985 .seealso: MatZeroEntries(),
1986 @*/
1987 int MatZeroRows(Mat mat,IS is, Scalar *diag)
1988 {
1989   int ierr;
1990 
1991   PetscValidHeaderSpecific(mat,MAT_COOKIE);
1992   PetscValidHeaderSpecific(is,IS_COOKIE);
1993   if (diag) PetscValidScalarPointer(diag);
1994   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
1995   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
1996   if (!mat->ops.zerorows) SETERRQ(PETSC_ERR_SUP,0,"");
1997 
1998   ierr = (*mat->ops.zerorows)(mat,is,diag); CHKERRQ(ierr);
1999   ierr = MatView_Private(mat); CHKERRQ(ierr);
2000   return 0;
2001 }
2002 
2003 #undef __FUNC__
2004 #define __FUNC__ "MatZeroRowsLocal"
2005 /*@
2006    MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal)
2007    of a set of rows of a matrix; using local numbering of rows.
2008 
2009    Input Parameters:
2010 .  mat - the matrix
2011 .  is - index set of rows to remove
2012 .  diag - pointer to value put in all diagonals of eliminated rows.
2013           Note that diag is not a pointer to an array, but merely a
2014           pointer to a single value.
2015 
2016    Notes:
2017    For the AIJ matrix formats this removes the old nonzero structure,
2018    but does not release memory.  For the dense and block diagonal
2019    formats this does not alter the nonzero structure.
2020 
2021    The user can set a value in the diagonal entry (or for the AIJ and
2022    row formats can optionally remove the main diagonal entry from the
2023    nonzero structure as well, by passing a null pointer as the final
2024    argument).
2025 
2026 .keywords: matrix, zero, rows, boundary conditions
2027 
2028 .seealso: MatZeroEntries(),
2029 @*/
2030 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag)
2031 {
2032   int ierr;
2033   IS  newis;
2034 
2035   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2036   PetscValidHeaderSpecific(is,IS_COOKIE);
2037   if (diag) PetscValidScalarPointer(diag);
2038   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
2039   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
2040   if (!mat->ops.zerorows) SETERRQ(PETSC_ERR_SUP,0,"");
2041 
2042   ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis); CHKERRQ(ierr);
2043   ierr =  (*mat->ops.zerorows)(mat,newis,diag); CHKERRQ(ierr);
2044   ierr = ISDestroy(newis);
2045   return 0;
2046 }
2047 
2048 #undef __FUNC__
2049 #define __FUNC__ "MatGetSize" /* ADIC Ignore */
2050 /*@
2051    MatGetSize - Returns the numbers of rows and columns in a matrix.
2052 
2053    Input Parameter:
2054 .  mat - the matrix
2055 
2056    Output Parameters:
2057 .  m - the number of global rows
2058 .  n - the number of global columns
2059 
2060 .keywords: matrix, dimension, size, rows, columns, global, get
2061 
2062 .seealso: MatGetLocalSize()
2063 @*/
2064 int MatGetSize(Mat mat,int *m,int* n)
2065 {
2066   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2067   PetscValidIntPointer(m);
2068   PetscValidIntPointer(n);
2069   return (*mat->ops.getsize)(mat,m,n);
2070 }
2071 
2072 #undef __FUNC__
2073 #define __FUNC__ "MatGetLocalSize" /* ADIC Ignore */
2074 /*@
2075    MatGetLocalSize - Returns the number of rows and columns in a matrix
2076    stored locally.  This information may be implementation dependent, so
2077    use with care.
2078 
2079    Input Parameters:
2080 .  mat - the matrix
2081 
2082    Output Parameters:
2083 .  m - the number of local rows
2084 .  n - the number of local columns
2085 
2086 .keywords: matrix, dimension, size, local, rows, columns, get
2087 
2088 .seealso: MatGetSize()
2089 @*/
2090 int MatGetLocalSize(Mat mat,int *m,int* n)
2091 {
2092   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2093   PetscValidIntPointer(m);
2094   PetscValidIntPointer(n);
2095   return (*mat->ops.getlocalsize)(mat,m,n);
2096 }
2097 
2098 #undef __FUNC__
2099 #define __FUNC__ "MatGetOwnershipRange" /* ADIC Ignore */
2100 /*@
2101    MatGetOwnershipRange - Returns the range of matrix rows owned by
2102    this processor, assuming that the matrix is laid out with the first
2103    n1 rows on the first processor, the next n2 rows on the second, etc.
2104    For certain parallel layouts this range may not be well defined.
2105 
2106    Input Parameters:
2107 .  mat - the matrix
2108 
2109    Output Parameters:
2110 .  m - the first local row
2111 .  n - one more then the last local row
2112 
2113 .keywords: matrix, get, range, ownership
2114 @*/
2115 int MatGetOwnershipRange(Mat mat,int *m,int* n)
2116 {
2117   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2118   PetscValidIntPointer(m);
2119   PetscValidIntPointer(n);
2120   if (!mat->ops.getownershiprange) SETERRQ(PETSC_ERR_SUP,0,"");
2121   return (*mat->ops.getownershiprange)(mat,m,n);
2122 }
2123 
2124 #undef __FUNC__
2125 #define __FUNC__ "MatILUFactorSymbolic"
2126 /*@
2127    MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix.
2128    Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric()
2129    to complete the factorization.
2130 
2131    Input Parameters:
2132 .  mat - the matrix
2133 .  row - row permutation
2134 .  column - column permutation
2135 .  fill - number of levels of fill
2136 .  f - expected fill as ratio of the original number of nonzeros,
2137        for example 3.0; choosing this parameter well can result in
2138        more efficient use of time and space.
2139 
2140    Output Parameters:
2141 .  fact - new matrix that has been symbolically factored
2142 
2143    Options Database Key:
2144 $   -mat_ilu_fill <f>, where f is the fill ratio
2145 
2146    Notes:
2147    See the file $(PETSC_DIR)/Performace for additional information about
2148    choosing the fill factor for better efficiency.
2149 
2150 .keywords: matrix, factor, incomplete, ILU, symbolic, fill
2151 
2152 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric()
2153 @*/
2154 int MatILUFactorSymbolic(Mat mat,IS row,IS col,double f,int fill,Mat *fact)
2155 {
2156   int ierr,flg;
2157 
2158   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2159   PetscValidPointer(fact);
2160   if (fill < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative");
2161   if (!mat->ops.ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"");
2162   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
2163   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
2164 
2165   ierr = OptionsGetDouble(PETSC_NULL,"-mat_ilu_fill",&f,&flg); CHKERRQ(ierr);
2166   PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0);
2167   ierr = (*mat->ops.ilufactorsymbolic)(mat,row,col,f,fill,fact); CHKERRQ(ierr);
2168   PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0);
2169   return 0;
2170 }
2171 
2172 #undef __FUNC__
2173 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic"
2174 /*@
2175    MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete
2176    Cholesky factorization for a symmetric matrix.  Use
2177    MatCholeskyFactorNumeric() to complete the factorization.
2178 
2179    Input Parameters:
2180 .  mat - the matrix
2181 .  perm - row and column permutation
2182 .  fill - levels of fill
2183 .  f - expected fill as ratio of original fill
2184 
2185    Output Parameter:
2186 .  fact - the factored matrix
2187 
2188    Note:  Currently only no-fill factorization is supported.
2189 
2190 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill
2191 
2192 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor()
2193 @*/
2194 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill,
2195                                         Mat *fact)
2196 {
2197   int ierr;
2198   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2199   PetscValidPointer(fact);
2200   if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
2201   if (fill < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative");
2202   if (!mat->ops.incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"");
2203   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
2204 
2205   PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0);
2206   ierr = (*mat->ops.incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr);
2207   PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0);
2208   return 0;
2209 }
2210 
2211 #undef __FUNC__
2212 #define __FUNC__ "MatGetArray" /* ADIC Ignore */
2213 /*@C
2214    MatGetArray - Returns a pointer to the element values in the matrix.
2215    This routine  is implementation dependent, and may not even work for
2216    certain matrix types. You MUST call MatRestoreArray() when you no
2217    longer need to access the array.
2218 
2219    Input Parameter:
2220 .  mat - the matrix
2221 
2222    Output Parameter:
2223 .  v - the location of the values
2224 
2225    Fortran Note:
2226    The Fortran interface is slightly different from that given below.
2227    See the Fortran chapter of the users manual and
2228    petsc/src/mat/examples for details.
2229 
2230 .keywords: matrix, array, elements, values
2231 
2232 .seealso: MatRestoreArray()
2233 @*/
2234 int MatGetArray(Mat mat,Scalar **v)
2235 {
2236   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2237   PetscValidPointer(v);
2238   if (!mat->ops.getarray) SETERRQ(PETSC_ERR_SUP,0,"");
2239   return (*mat->ops.getarray)(mat,v);
2240 }
2241 
2242 #undef __FUNC__
2243 #define __FUNC__ "MatRestoreArray" /* ADIC Ignore */
2244 /*@C
2245    MatRestoreArray - Restores the matrix after MatGetArray has been called.
2246 
2247    Input Parameter:
2248 .  mat - the matrix
2249 .  v - the location of the values
2250 
2251    Fortran Note:
2252    The Fortran interface is slightly different from that given below.
2253    See the users manual and petsc/src/mat/examples for details.
2254 
2255 .keywords: matrix, array, elements, values, restore
2256 
2257 .seealso: MatGetArray()
2258 @*/
2259 int MatRestoreArray(Mat mat,Scalar **v)
2260 {
2261   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2262   PetscValidPointer(v);
2263   if (!mat->ops.restorearray) SETERRQ(PETSC_ERR_SUP,0,"");
2264   return (*mat->ops.restorearray)(mat,v);
2265 }
2266 
2267 #undef __FUNC__
2268 #define __FUNC__ "MatGetSubMatrices"
2269 /*@C
2270    MatGetSubMatrices - Extracts several submatrices from a matrix. If submat
2271    points to an array of valid matrices, it may be reused.
2272 
2273    Input Parameters:
2274 .  mat - the matrix
2275 .  n   - the number of submatrixes to be extracted
2276 .  irow, icol - index sets of rows and columns to extract
2277 .  scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX
2278 
2279    Output Parameter:
2280 .  submat - the array of submatrices
2281 
2282    Limitations:
2283    Currently, MatGetSubMatrices() can extract only sequential submatrices
2284    (from both sequential and parallel matrices).
2285 
2286    Notes:
2287    When extracting submatrices from a parallel matrix, each processor can
2288    form a different submatrix by setting the rows and columns of its
2289    individual index sets according to the local submatrix desired.
2290 
2291    When finished using the submatrices, the user should destroy
2292    them with MatDestroySubMatrices().
2293 
2294 .keywords: matrix, get, submatrix, submatrices
2295 
2296 .seealso: MatDestroyMatrices()
2297 @*/
2298 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatGetSubMatrixCall scall,
2299                       Mat **submat)
2300 {
2301   int ierr;
2302   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2303   if (!mat->ops.getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,"");
2304   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
2305 
2306   PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0);
2307   ierr = (*mat->ops.getsubmatrices)(mat,n,irow,icol,scall,submat); CHKERRQ(ierr);
2308   PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0);
2309 
2310   return 0;
2311 }
2312 
2313 #undef __FUNC__
2314 #define __FUNC__ "MatDestroyMatrices" /* ADIC Ignore */
2315 /*@C
2316    MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices().
2317 
2318    Input Parameters:
2319 .  n - the number of local matrices
2320 .  mat - the matrices
2321 
2322 .keywords: matrix, destroy, submatrix, submatrices
2323 
2324 .seealso: MatGetSubMatrices()
2325 @*/
2326 int MatDestroyMatrices(int n,Mat **mat)
2327 {
2328   int ierr,i;
2329 
2330   if (n < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices");
2331   PetscValidPointer(mat);
2332   for ( i=0; i<n; i++ ) {
2333     ierr = MatDestroy((*mat)[i]); CHKERRQ(ierr);
2334   }
2335   if (n) PetscFree(*mat);
2336   return 0;
2337 }
2338 
2339 #undef __FUNC__
2340 #define __FUNC__ "MatIncreaseOverlap"
2341 /*@
2342    MatIncreaseOverlap - Given a set of submatrices indicated by index sets,
2343    replaces the index sets by larger ones that represent submatrices with
2344    additional overlap.
2345 
2346    Input Parameters:
2347 .  mat - the matrix
2348 .  n   - the number of index sets
2349 .  is  - the array of pointers to index sets
2350 .  ov  - the additional overlap requested
2351 
2352 .keywords: matrix, overlap, Schwarz
2353 
2354 .seealso: MatGetSubMatrices()
2355 @*/
2356 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov)
2357 {
2358   int ierr;
2359   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2360   if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix");
2361   if (mat->factor)     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix");
2362 
2363   if (ov == 0) return 0;
2364   if (!mat->ops.increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,"");
2365   PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0);
2366   ierr = (*mat->ops.increaseoverlap)(mat,n,is,ov); CHKERRQ(ierr);
2367   PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0);
2368   return 0;
2369 }
2370 
2371 #undef __FUNC__
2372 #define __FUNC__ "MatPrintHelp" /* ADIC Ignore */
2373 /*@
2374    MatPrintHelp - Prints all the options for the matrix.
2375 
2376    Input Parameter:
2377 .  mat - the matrix
2378 
2379    Options Database Keys:
2380 $  -help, -h
2381 
2382 .keywords: mat, help
2383 
2384 .seealso: MatCreate(), MatCreateXXX()
2385 @*/
2386 int MatPrintHelp(Mat mat)
2387 {
2388   static int called = 0;
2389   MPI_Comm   comm;
2390   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2391 
2392   comm = mat->comm;
2393   if (!called) {
2394     PetscPrintf(comm,"General matrix options:\n");
2395     PetscPrintf(comm,"  -mat_view_info: view basic matrix info during MatAssemblyEnd()\n");
2396     PetscPrintf(comm,"  -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n");
2397     PetscPrintf(comm,"  -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n");
2398     PetscPrintf(comm,"      -draw_pause <sec>: set seconds of display pause\n");
2399     PetscPrintf(comm,"      -display <name>: set alternate display\n");
2400     called = 1;
2401   }
2402   if (mat->ops.printhelp) (*mat->ops.printhelp)(mat);
2403   return 0;
2404 }
2405 
2406 #undef __FUNC__
2407 #define __FUNC__ "MatGetBlockSize" /* ADIC Ignore */
2408 /*@
2409    MatGetBlockSize - Returns the matrix block size; useful especially for the
2410    block row and block diagonal formats.
2411 
2412    Input Parameter:
2413 .  mat - the matrix
2414 
2415    Output Parameter:
2416 .  bs - block size
2417 
2418    Notes:
2419 $  block diagonal formats: MATSEQBDIAG, MATMPIBDIAG
2420 $  block row formats: MATSEQBAIJ, MATMPIBAIJ
2421 
2422 .keywords: matrix, get, block, size
2423 
2424 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag()
2425 @*/
2426 int MatGetBlockSize(Mat mat,int *bs)
2427 {
2428   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2429   PetscValidIntPointer(bs);
2430   if (!mat->ops.getblocksize) SETERRQ(PETSC_ERR_SUP,0,"");
2431   return (*mat->ops.getblocksize)(mat,bs);
2432 }
2433 
2434 #undef __FUNC__
2435 #define __FUNC__ "MatGetRowIJ"
2436 /*@C
2437       MatGetRowIJ - Returns the compress row storage i and j indices for sequential matrices.
2438                  EXPERTS ONLY.
2439 
2440   Input Parameters:
2441 .   mat - the matrix
2442 .   shift - 1 or zero indicating we want the indices starting at 0 or 1
2443 .   symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be
2444                 symmetrized
2445 
2446   Output Parameters:
2447 .   n - number of rows and columns in the (possibly compressed) matrix
2448 .   ia - the row indices
2449 .   ja - the column indices
2450 .   done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned
2451 @*/
2452 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done)
2453 {
2454   int ierr;
2455   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2456   if (ia) PetscValidIntPointer(ia);
2457   if (ja) PetscValidIntPointer(ja);
2458   PetscValidIntPointer(done);
2459   if (!mat->ops.getrowij) *done = PETSC_FALSE;
2460   else {
2461     *done = PETSC_TRUE;
2462     ierr  = (*mat->ops.getrowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr);
2463   }
2464   return 0;
2465 }
2466 
2467 #undef __FUNC__
2468 #define __FUNC__ "MatGetColumnIJ" /* ADIC Ignore */
2469 /*@C
2470       MatGetColumnIJ - Returns the compress Column storage i and j indices for sequential matrices.
2471                  EXPERTS ONLY.
2472 
2473   Input Parameters:
2474 .   mat - the matrix
2475 .   shift - 1 or zero indicating we want the indices starting at 0 or 1
2476 .   symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be
2477                 symmetrized
2478 
2479   Output Parameters:
2480 .   n - number of Columns and columns in the (possibly compressed) matrix
2481 .   ia - the Column indices
2482 .   ja - the column indices
2483 .   done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned
2484 @*/
2485 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done)
2486 {
2487   int ierr;
2488   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2489   if (ia) PetscValidIntPointer(ia);
2490   if (ja) PetscValidIntPointer(ja);
2491   PetscValidIntPointer(done);
2492 
2493   if (!mat->ops.getcolumnij) *done = PETSC_FALSE;
2494   else {
2495     *done = PETSC_TRUE;
2496     ierr  = (*mat->ops.getcolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr);
2497   }
2498   return 0;
2499 }
2500 
2501 #undef __FUNC__
2502 #define __FUNC__ "MatRestoreRowIJ" /* ADIC Ignore */
2503 /*@C
2504       MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with
2505                      MatGetRowIJ(). EXPERTS ONLY.
2506 
2507   Input Parameters:
2508 .   mat - the matrix
2509 .   shift - 1 or zero indicating we want the indices starting at 0 or 1
2510 .   symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be
2511                 symmetrized
2512 
2513   Output Parameters:
2514 .   n - size of (possibly compressed) matrix
2515 .   ia - the row indices
2516 .   ja - the column indices
2517 .   done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned
2518 @*/
2519 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done)
2520 {
2521   int ierr;
2522   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2523   if (ia) PetscValidIntPointer(ia);
2524   if (ja) PetscValidIntPointer(ja);
2525   PetscValidIntPointer(done);
2526 
2527   if (!mat->ops.restorerowij) *done = PETSC_FALSE;
2528   else {
2529     *done = PETSC_TRUE;
2530     ierr  = (*mat->ops.restorerowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr);
2531   }
2532   return 0;
2533 }
2534 
2535 #undef __FUNC__
2536 #define __FUNC__ "MatRestoreColumnIJ" /* ADIC Ignore */
2537 /*@C
2538       MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with
2539                      MatGetColumnIJ(). EXPERTS ONLY.
2540 
2541   Input Parameters:
2542 .   mat - the matrix
2543 .   shift - 1 or zero indicating we want the indices starting at 0 or 1
2544 .   symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be
2545                 symmetrized
2546 
2547   Output Parameters:
2548 .   n - size of (possibly compressed) matrix
2549 .   ia - the Column indices
2550 .   ja - the column indices
2551 .   done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned
2552 @*/
2553 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done)
2554 {
2555   int ierr;
2556   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2557   if (ia) PetscValidIntPointer(ia);
2558   if (ja) PetscValidIntPointer(ja);
2559   PetscValidIntPointer(done);
2560 
2561   if (!mat->ops.restorecolumnij) *done = PETSC_FALSE;
2562   else {
2563     *done = PETSC_TRUE;
2564     ierr  = (*mat->ops.restorecolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr);
2565   }
2566   return 0;
2567 }
2568 
2569 #undef __FUNC__
2570 #define __FUNC__ "MatColoringPatch"
2571 /*@C
2572       MatColoringPatch - EXPERTS ONLY, used inside matrix coloring routines that
2573           use matGetRowIJ() and/or MatGetColumnIJ().
2574 
2575   Input Parameters:
2576 .   mat - the matrix
2577 .   n   - number of colors
2578 .   colorarray - array indicating color for each column
2579 
2580   Output Parameters:
2581 .   iscoloring - coloring generated using colorarray information
2582 
2583 @*/
2584 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring)
2585 {
2586   int ierr;
2587   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2588   PetscValidIntPointer(colorarray);
2589 
2590   if (!mat->ops.coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");}
2591   else {
2592     ierr  = (*mat->ops.coloringpatch)(mat,n,colorarray,iscoloring); CHKERRQ(ierr);
2593   }
2594   return 0;
2595 }
2596 
2597 
2598 /*@
2599    MatSetUnfactored - Resets a factored matrix to be treated as unfactored.
2600 
2601    Input Paramter:
2602 .  mat - the factored matrix to be reset
2603 
2604    Notes:
2605    This routine should be used only with factored matrices formed by in-place
2606    factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE
2607    format).  This option can save memory, for example, when solving nonlinear
2608    systems with a matrix-free Newton-Krylov method and a matrix-based, in-place
2609    ILU(0) preconditioner.
2610 
2611   Note that one can specify in-place ILU(0) factorization by calling
2612 $     PCType(pc,PCILU);
2613 $     PCILUSeUseInPlace(pc);
2614   or by using the options -pc_type ilu -pc_ilu_in_place
2615 
2616   In-place factorization ILU(0) can also be used as a local
2617   solver for the blocks within the block Jacobi or additive Schwarz
2618   methods (runtime option: -sub_pc_ilu_in_place).  See the discussion
2619   of these preconditioners in the users manual for details on setting
2620   local solver options.
2621 
2622 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace()
2623 
2624 .keywords: matrix-free, in-place ILU, in-place LU
2625 @*/
2626 int MatSetUnfactored(Mat mat)
2627 {
2628   int ierr;
2629 
2630   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2631   mat->factor = 0;
2632   if (!mat->ops.setunfactored) return 0;
2633   ierr = (*mat->ops.setunfactored)(mat); CHKERRQ(ierr);
2634   return 0;
2635 }
2636 
2637 #undef __FUNC__
2638 #define __FUNC__ "MatGetType" /* ADIC Ignore */
2639 /*@C
2640    MatGetType - Gets the matrix type and name (as a string) from the matrix.
2641 
2642    Input Parameter:
2643 .  mat - the matrix
2644 
2645    Output Parameter:
2646 .  type - the matrix type (or use PETSC_NULL)
2647 .  name - name of matrix type (or use PETSC_NULL)
2648 
2649 .keywords: matrix, get, type, name
2650 @*/
2651 int MatGetType(Mat mat,MatType *type,char **name)
2652 {
2653   int  itype = (int)mat->type;
2654   char *matname[10];
2655 
2656   PetscValidHeaderSpecific(mat,MAT_COOKIE);
2657 
2658   if (type) *type = (MatType) mat->type;
2659   if (name) {
2660     /* Note:  Be sure that this list corresponds to the enum in mat.h */
2661     matname[0] = "MATSEQDENSE";
2662     matname[1] = "MATSEQAIJ";
2663     matname[2] = "MATMPIAIJ";
2664     matname[3] = "MATSHELL";
2665     matname[4] = "MATMPIROWBS";
2666     matname[5] = "MATSEQBDIAG";
2667     matname[6] = "MATMPIBDIAG";
2668     matname[7] = "MATMPIDENSE";
2669     matname[8] = "MATSEQBAIJ";
2670     matname[9] = "MATMPIBAIJ";
2671 
2672     if (itype < 0 || itype > 9) *name = "Unknown matrix type";
2673     else                        *name = matname[itype];
2674   }
2675   return 0;
2676 }
2677