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