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