xref: /petsc/src/mat/impls/aij/seq/umfpack/umfpack.c (revision d32f9abdbc052d6e1fd06679b17a55415c3aae30)
1 #define PETSCMAT_DLL
2 
3 /*
4    Provides an interface to the UMFPACKv5.1 sparse solver
5 
6    This interface uses the "UF_long version" of the UMFPACK API
7    (*_dl_* and *_zl_* routines, instead of *_di_* and *_zi_* routines)
8    so that UMFPACK can address more than 2Gb of memory on 64 bit
9    machines.
10 
11    If sizeof(UF_long) == 32 the interface only allocates the memory
12    necessary for UMFPACK's working arrays (W, Wi and possibly
13    perm_c). If sizeof(UF_long) == 64, in addition to allocating the
14    working arrays, the interface also has to re-allocate the matrix
15    index arrays (ai and aj, which must be stored as UF_long).
16 
17    The interface is implemented for both real and complex
18    arithmetic. Complex numbers are assumed to be in packed format,
19    which requires UMFPACK >= 4.4.
20 
21    We thank Christophe Geuzaine <geuzaine@acm.caltech.edu> for upgrading this interface to the UMFPACKv5.1
22 */
23 
24 #include "src/mat/impls/aij/seq/aij.h"
25 
26 EXTERN_C_BEGIN
27 #include "umfpack.h"
28 EXTERN_C_END
29 
30 typedef struct {
31   void         *Symbolic, *Numeric;
32   double       Info[UMFPACK_INFO], Control[UMFPACK_CONTROL],*W;
33   UF_long      *Wi,*ai,*aj,*perm_c;
34   PetscScalar  *av;
35   MatStructure flg;
36   PetscTruth   PetscMatOdering;
37 
38   /* Flag to clean up UMFPACK objects during Destroy */
39   PetscTruth CleanUpUMFPACK;
40 } Mat_UMFPACK;
41 
42 #undef __FUNCT__
43 #define __FUNCT__ "MatDestroy_UMFPACK"
44 PetscErrorCode MatDestroy_UMFPACK(Mat A)
45 {
46   PetscErrorCode ierr;
47   Mat_UMFPACK    *lu=(Mat_UMFPACK*)A->spptr;
48 
49   PetscFunctionBegin;
50   if (lu->CleanUpUMFPACK) {
51 #if defined(PETSC_USE_COMPLEX)
52     umfpack_zl_free_symbolic(&lu->Symbolic);
53     umfpack_zl_free_numeric(&lu->Numeric);
54 #else
55     umfpack_dl_free_symbolic(&lu->Symbolic);
56     umfpack_dl_free_numeric(&lu->Numeric);
57 #endif
58     ierr = PetscFree(lu->Wi);CHKERRQ(ierr);
59     ierr = PetscFree(lu->W);CHKERRQ(ierr);
60     if(sizeof(UF_long) != sizeof(int)){
61       ierr = PetscFree(lu->ai);CHKERRQ(ierr);
62       ierr = PetscFree(lu->aj);CHKERRQ(ierr);
63     }
64     if (lu->PetscMatOdering) {
65       ierr = PetscFree(lu->perm_c);CHKERRQ(ierr);
66     }
67   }
68   ierr = MatDestroy_SeqAIJ(A);CHKERRQ(ierr);
69   PetscFunctionReturn(0);
70 }
71 
72 #undef __FUNCT__
73 #define __FUNCT__ "MatSolve_UMFPACK"
74 PetscErrorCode MatSolve_UMFPACK(Mat A,Vec b,Vec x)
75 {
76   Mat_UMFPACK *lu = (Mat_UMFPACK*)A->spptr;
77   PetscScalar *av=lu->av,*ba,*xa;
78   PetscErrorCode ierr;
79   UF_long     *ai=lu->ai,*aj=lu->aj,status;
80 
81   PetscFunctionBegin;
82   /* solve Ax = b by umfpack_*_wsolve */
83   /* ----------------------------------*/
84 
85 #if defined(PETSC_USE_COMPLEX)
86   ierr = VecConjugate(b);
87 #endif
88 
89   ierr = VecGetArray(b,&ba);
90   ierr = VecGetArray(x,&xa);
91 
92 #if defined(PETSC_USE_COMPLEX)
93   status = umfpack_zl_wsolve(UMFPACK_At,ai,aj,(double*)av,NULL,(double*)xa,NULL,(double*)ba,NULL,
94 			     lu->Numeric,lu->Control,lu->Info,lu->Wi,lu->W);
95   umfpack_zl_report_info(lu->Control, lu->Info);
96   if (status < 0){
97     umfpack_zl_report_status(lu->Control, status);
98     SETERRQ(PETSC_ERR_LIB,"umfpack_zl_wsolve failed");
99   }
100 #else
101   status = umfpack_dl_wsolve(UMFPACK_At,ai,aj,av,xa,ba,
102 			     lu->Numeric,lu->Control,lu->Info,lu->Wi,lu->W);
103   umfpack_dl_report_info(lu->Control, lu->Info);
104   if (status < 0){
105     umfpack_dl_report_status(lu->Control, status);
106     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_wsolve failed");
107   }
108 #endif
109 
110   ierr = VecRestoreArray(b,&ba);
111   ierr = VecRestoreArray(x,&xa);
112 
113 #if defined(PETSC_USE_COMPLEX)
114   ierr = VecConjugate(b);
115   ierr = VecConjugate(x);
116 #endif
117 
118   PetscFunctionReturn(0);
119 }
120 
121 #undef __FUNCT__
122 #define __FUNCT__ "MatLUFactorNumeric_UMFPACK"
123 PetscErrorCode MatLUFactorNumeric_UMFPACK(Mat A,MatFactorInfo *info,Mat *F)
124 {
125   Mat_UMFPACK *lu=(Mat_UMFPACK*)(*F)->spptr;
126   PetscErrorCode ierr;
127   UF_long     *ai=lu->ai,*aj=lu->aj,m=A->rmap.n,status;
128   PetscScalar *av=lu->av;
129 
130   PetscFunctionBegin;
131   /* numeric factorization of A' */
132   /* ----------------------------*/
133 
134 #if defined(PETSC_USE_COMPLEX)
135   if (lu->flg == SAME_NONZERO_PATTERN && lu->Numeric){
136     umfpack_zl_free_numeric(&lu->Numeric);
137   }
138   status = umfpack_zl_numeric(ai,aj,(double*)av,NULL,lu->Symbolic,&lu->Numeric,lu->Control,lu->Info);
139   if (status < 0) {
140     umfpack_zl_report_status(lu->Control, status);
141     SETERRQ(PETSC_ERR_LIB,"umfpack_zl_numeric failed");
142   }
143   /* report numeric factorization of A' when Control[PRL] > 3 */
144   (void) umfpack_zl_report_numeric(lu->Numeric, lu->Control);
145 #else
146   if (lu->flg == SAME_NONZERO_PATTERN && lu->Numeric){
147     umfpack_dl_free_numeric(&lu->Numeric);
148   }
149   status = umfpack_dl_numeric(ai,aj,av,lu->Symbolic,&lu->Numeric,lu->Control,lu->Info);
150   if (status < 0) {
151     umfpack_zl_report_status(lu->Control, status);
152     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_numeric failed");
153   }
154   /* report numeric factorization of A' when Control[PRL] > 3 */
155   (void) umfpack_dl_report_numeric(lu->Numeric, lu->Control);
156 #endif
157 
158   if (lu->flg == DIFFERENT_NONZERO_PATTERN){  /* first numeric factorization */
159     /* allocate working space to be used by Solve */
160     ierr = PetscMalloc(m * sizeof(UF_long), &lu->Wi);CHKERRQ(ierr);
161 #if defined(PETSC_USE_COMPLEX)
162     ierr = PetscMalloc(2*5*m * sizeof(double), &lu->W);CHKERRQ(ierr);
163 #else
164     ierr = PetscMalloc(5*m * sizeof(double), &lu->W);CHKERRQ(ierr);
165 #endif
166   }
167 
168   lu->flg = SAME_NONZERO_PATTERN;
169   lu->CleanUpUMFPACK = PETSC_TRUE;
170   PetscFunctionReturn(0);
171 }
172 
173 /*
174    Note the r permutation is ignored
175 */
176 #undef __FUNCT__
177 #define __FUNCT__ "MatLUFactorSymbolic_UMFPACK"
178 PetscErrorCode MatLUFactorSymbolic_UMFPACK(Mat A,IS r,IS c,MatFactorInfo *info,Mat *F)
179 {
180   Mat_SeqAIJ     *mat=(Mat_SeqAIJ*)A->data;
181   Mat_UMFPACK    *lu = (Mat_UMFPACK*)((*F)->spptr);
182   PetscErrorCode ierr;
183   int            i,m=A->rmap.n,n=A->cmap.n,*ra;
184   UF_long        status;
185   PetscScalar    *av=mat->a;
186 
187   PetscFunctionBegin;
188   if (lu->PetscMatOdering) {
189     ierr = ISGetIndices(r,&ra);CHKERRQ(ierr);
190     ierr = PetscMalloc(m*sizeof(UF_long),&lu->perm_c);CHKERRQ(ierr);
191     /* we cannot simply memcpy on 64 bit archs */
192     for(i = 0; i < m; i++) lu->perm_c[i] = ra[i];
193     ierr = ISRestoreIndices(r,&ra);CHKERRQ(ierr);
194   }
195 
196   if(sizeof(UF_long) != sizeof(int)){
197     /* we cannot directly use mat->i and mat->j on 64 bit archs */
198     ierr = PetscMalloc((m+1)*sizeof(UF_long),&lu->ai);CHKERRQ(ierr);
199     ierr = PetscMalloc(mat->nz*sizeof(UF_long),&lu->aj);CHKERRQ(ierr);
200     for(i = 0; i < m + 1; i++) lu->ai[i] = mat->i[i];
201     for(i = 0; i < mat->nz; i++) lu->aj[i] = mat->j[i];
202   }
203   else{
204     lu->ai = (UF_long*)mat->i;
205     lu->aj = (UF_long*)mat->j;
206   }
207 
208   /* print the control parameters */
209 #if defined(PETSC_USE_COMPLEX)
210   if(lu->Control[UMFPACK_PRL] > 1) umfpack_zl_report_control(lu->Control);
211 #else
212   if(lu->Control[UMFPACK_PRL] > 1) umfpack_dl_report_control(lu->Control);
213 #endif
214 
215   /* symbolic factorization of A' */
216   /* ---------------------------------------------------------------------- */
217 #if defined(PETSC_USE_COMPLEX)
218   if (lu->PetscMatOdering) { /* use Petsc row ordering */
219     status = umfpack_zl_qsymbolic(n,m,lu->ai,lu->aj,(double*)av,NULL,
220 				  lu->perm_c,&lu->Symbolic,lu->Control,lu->Info);
221   } else { /* use Umfpack col ordering */
222     status = umfpack_zl_symbolic(n,m,lu->ai,lu->aj,(double*)av,NULL,
223 				 &lu->Symbolic,lu->Control,lu->Info);
224   }
225   if (status < 0){
226     umfpack_zl_report_info(lu->Control, lu->Info);
227     umfpack_zl_report_status(lu->Control, status);
228     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_symbolic failed");
229   }
230   /* report sumbolic factorization of A' when Control[PRL] > 3 */
231   (void) umfpack_zl_report_symbolic(lu->Symbolic, lu->Control);
232 #else
233   if (lu->PetscMatOdering) { /* use Petsc row ordering */
234     status = umfpack_dl_qsymbolic(n,m,lu->ai,lu->aj,av,
235 				  lu->perm_c,&lu->Symbolic,lu->Control,lu->Info);
236   } else { /* use Umfpack col ordering */
237     status = umfpack_dl_symbolic(n,m,lu->ai,lu->aj,av,
238 				 &lu->Symbolic,lu->Control,lu->Info);
239   }
240   if (status < 0){
241     umfpack_dl_report_info(lu->Control, lu->Info);
242     umfpack_dl_report_status(lu->Control, status);
243     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_symbolic failed");
244   }
245   /* report sumbolic factorization of A' when Control[PRL] > 3 */
246   (void) umfpack_dl_report_symbolic(lu->Symbolic, lu->Control);
247 #endif
248 
249   lu->flg = DIFFERENT_NONZERO_PATTERN;
250   lu->av  = av;
251   lu->CleanUpUMFPACK = PETSC_TRUE;
252   PetscFunctionReturn(0);
253 }
254 
255 #undef __FUNCT__
256 #define __FUNCT__ "MatGetFactor_seqaij_umfpack"
257 PetscErrorCode MatGetFactor_seqaij_umfpack(Mat A,MatFactorType ftype,Mat *F)
258 {
259   Mat            B;
260   Mat_UMFPACK    *lu;
261   PetscErrorCode ierr;
262   int            m=A->rmap.n,n=A->cmap.n,idx;
263 
264   const char     *strategy[]={"AUTO","UNSYMMETRIC","SYMMETRIC","2BY2"},
265                  *scale[]={"NONE","SUM","MAX"};
266   PetscTruth     flg;
267 
268   PetscFunctionBegin;
269   /* Create the factorization matrix F */
270   ierr = MatCreate(((PetscObject)A)->comm,&B);CHKERRQ(ierr);
271   ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr);
272   ierr = MatSetType(B,((PetscObject)A)->type_name);CHKERRQ(ierr);
273   ierr = MatSeqAIJSetPreallocation(B,0,PETSC_NULL);CHKERRQ(ierr);
274   ierr = PetscNewLog(B,Mat_UMFPACK,&lu);CHKERRQ(ierr);
275   B->spptr                 = lu;
276   B->ops->lufactornumeric  = MatLUFactorNumeric_UMFPACK;
277   B->ops->lufactorsymbolic = MatLUFactorSymbolic_UMFPACK;
278   B->ops->solve            = MatSolve_UMFPACK;
279   B->ops->destroy          = MatDestroy_UMFPACK;
280   B->factor                = MAT_FACTOR_LU;
281   B->assembled             = PETSC_TRUE;  /* required by -ksp_view */
282   B->preallocated          = PETSC_TRUE;
283 
284   /* initializations */
285   /* ------------------------------------------------*/
286   /* get the default control parameters */
287 #if defined(PETSC_USE_COMPLEX)
288   umfpack_zl_defaults(lu->Control);
289 #else
290   umfpack_dl_defaults(lu->Control);
291 #endif
292   lu->perm_c = PETSC_NULL;  /* use defaul UMFPACK col permutation */
293   lu->Control[UMFPACK_IRSTEP] = 0; /* max num of iterative refinement steps to attempt */
294 
295   ierr = PetscOptionsBegin(((PetscObject)A)->comm,((PetscObject)A)->prefix,"UMFPACK Options","Mat");CHKERRQ(ierr);
296   /* Control parameters used by reporting routiones */
297   ierr = PetscOptionsReal("-mat_umfpack_prl","Control[UMFPACK_PRL]","None",lu->Control[UMFPACK_PRL],&lu->Control[UMFPACK_PRL],PETSC_NULL);CHKERRQ(ierr);
298 
299   /* Control parameters for symbolic factorization */
300   ierr = PetscOptionsEList("-mat_umfpack_strategy","ordering and pivoting strategy","None",strategy,4,strategy[0],&idx,&flg);CHKERRQ(ierr);
301   if (flg) {
302     switch (idx){
303     case 0: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_AUTO; break;
304     case 1: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_UNSYMMETRIC; break;
305     case 2: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_SYMMETRIC; break;
306     case 3: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_2BY2; break;
307     }
308   }
309   ierr = PetscOptionsReal("-mat_umfpack_dense_col","Control[UMFPACK_DENSE_COL]","None",lu->Control[UMFPACK_DENSE_COL],&lu->Control[UMFPACK_DENSE_COL],PETSC_NULL);CHKERRQ(ierr);
310   ierr = PetscOptionsReal("-mat_umfpack_dense_row","Control[UMFPACK_DENSE_ROW]","None",lu->Control[UMFPACK_DENSE_ROW],&lu->Control[UMFPACK_DENSE_ROW],PETSC_NULL);CHKERRQ(ierr);
311   ierr = PetscOptionsReal("-mat_umfpack_amd_dense","Control[UMFPACK_AMD_DENSE]","None",lu->Control[UMFPACK_AMD_DENSE],&lu->Control[UMFPACK_AMD_DENSE],PETSC_NULL);CHKERRQ(ierr);
312   ierr = PetscOptionsReal("-mat_umfpack_block_size","Control[UMFPACK_BLOCK_SIZE]","None",lu->Control[UMFPACK_BLOCK_SIZE],&lu->Control[UMFPACK_BLOCK_SIZE],PETSC_NULL);CHKERRQ(ierr);
313   ierr = PetscOptionsReal("-mat_umfpack_2by2_tolerance","Control[UMFPACK_2BY2_TOLERANCE]","None",lu->Control[UMFPACK_2BY2_TOLERANCE],&lu->Control[UMFPACK_2BY2_TOLERANCE],PETSC_NULL);CHKERRQ(ierr);
314   ierr = PetscOptionsReal("-mat_umfpack_fixq","Control[UMFPACK_FIXQ]","None",lu->Control[UMFPACK_FIXQ],&lu->Control[UMFPACK_FIXQ],PETSC_NULL);CHKERRQ(ierr);
315   ierr = PetscOptionsReal("-mat_umfpack_aggressive","Control[UMFPACK_AGGRESSIVE]","None",lu->Control[UMFPACK_AGGRESSIVE],&lu->Control[UMFPACK_AGGRESSIVE],PETSC_NULL);CHKERRQ(ierr);
316 
317   /* Control parameters used by numeric factorization */
318   ierr = PetscOptionsReal("-mat_umfpack_pivot_tolerance","Control[UMFPACK_PIVOT_TOLERANCE]","None",lu->Control[UMFPACK_PIVOT_TOLERANCE],&lu->Control[UMFPACK_PIVOT_TOLERANCE],PETSC_NULL);CHKERRQ(ierr);
319   ierr = PetscOptionsReal("-mat_umfpack_sym_pivot_tolerance","Control[UMFPACK_SYM_PIVOT_TOLERANCE]","None",lu->Control[UMFPACK_SYM_PIVOT_TOLERANCE],&lu->Control[UMFPACK_SYM_PIVOT_TOLERANCE],PETSC_NULL);CHKERRQ(ierr);
320   ierr = PetscOptionsEList("-mat_umfpack_scale","Control[UMFPACK_SCALE]","None",scale,3,scale[0],&idx,&flg);CHKERRQ(ierr);
321   if (flg) {
322     switch (idx){
323     case 0: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_NONE; break;
324     case 1: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_SUM; break;
325     case 2: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_MAX; break;
326     }
327   }
328   ierr = PetscOptionsReal("-mat_umfpack_alloc_init","Control[UMFPACK_ALLOC_INIT]","None",lu->Control[UMFPACK_ALLOC_INIT],&lu->Control[UMFPACK_ALLOC_INIT],PETSC_NULL);CHKERRQ(ierr);
329   ierr = PetscOptionsReal("-mat_umfpack_front_alloc_init","Control[UMFPACK_FRONT_ALLOC_INIT]","None",lu->Control[UMFPACK_FRONT_ALLOC_INIT],&lu->Control[UMFPACK_ALLOC_INIT],PETSC_NULL);CHKERRQ(ierr);
330   ierr = PetscOptionsReal("-mat_umfpack_droptol","Control[UMFPACK_DROPTOL]","None",lu->Control[UMFPACK_DROPTOL],&lu->Control[UMFPACK_DROPTOL],PETSC_NULL);CHKERRQ(ierr);
331 
332   /* Control parameters used by solve */
333   ierr = PetscOptionsReal("-mat_umfpack_irstep","Control[UMFPACK_IRSTEP]","None",lu->Control[UMFPACK_IRSTEP],&lu->Control[UMFPACK_IRSTEP],PETSC_NULL);CHKERRQ(ierr);
334 
335   /* use Petsc mat ordering (note: size is for the transpose, and PETSc r = Umfpack perm_c) */
336   ierr = PetscOptionsHasName(PETSC_NULL,"-pc_factor_mat_ordering_type",&lu->PetscMatOdering);CHKERRQ(ierr);
337   PetscOptionsEnd();
338   *F = B;
339   PetscFunctionReturn(0);
340 }
341 
342 #undef __FUNCT__
343 #define __FUNCT__ "MatFactorInfo_UMFPACK"
344 PetscErrorCode MatFactorInfo_UMFPACK(Mat A,PetscViewer viewer)
345 {
346   Mat_UMFPACK    *lu= (Mat_UMFPACK*)A->spptr;
347   PetscErrorCode ierr;
348 
349   PetscFunctionBegin;
350   /* check if matrix is UMFPACK type */
351   if (A->ops->solve != MatSolve_UMFPACK) PetscFunctionReturn(0);
352 
353   ierr = PetscViewerASCIIPrintf(viewer,"UMFPACK run parameters:\n");CHKERRQ(ierr);
354   /* Control parameters used by reporting routiones */
355   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_PRL]: %g\n",lu->Control[UMFPACK_PRL]);CHKERRQ(ierr);
356 
357   /* Control parameters used by symbolic factorization */
358   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_STRATEGY]: %g\n",lu->Control[UMFPACK_STRATEGY]);CHKERRQ(ierr);
359   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DENSE_COL]: %g\n",lu->Control[UMFPACK_DENSE_COL]);CHKERRQ(ierr);
360   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DENSE_ROW]: %g\n",lu->Control[UMFPACK_DENSE_ROW]);CHKERRQ(ierr);
361   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_AMD_DENSE]: %g\n",lu->Control[UMFPACK_AMD_DENSE]);CHKERRQ(ierr);
362   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_BLOCK_SIZE]: %g\n",lu->Control[UMFPACK_BLOCK_SIZE]);CHKERRQ(ierr);
363   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_2BY2_TOLERANCE]: %g\n",lu->Control[UMFPACK_2BY2_TOLERANCE]);CHKERRQ(ierr);
364   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_FIXQ]: %g\n",lu->Control[UMFPACK_FIXQ]);CHKERRQ(ierr);
365   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_AGGRESSIVE]: %g\n",lu->Control[UMFPACK_AGGRESSIVE]);CHKERRQ(ierr);
366 
367   /* Control parameters used by numeric factorization */
368   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_PIVOT_TOLERANCE]);CHKERRQ(ierr);
369   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_SYM_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_SYM_PIVOT_TOLERANCE]);CHKERRQ(ierr);
370   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_SCALE]: %g\n",lu->Control[UMFPACK_SCALE]);CHKERRQ(ierr);
371   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_ALLOC_INIT]: %g\n",lu->Control[UMFPACK_ALLOC_INIT]);CHKERRQ(ierr);
372   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DROPTOL]: %g\n",lu->Control[UMFPACK_DROPTOL]);CHKERRQ(ierr);
373 
374   /* Control parameters used by solve */
375   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_IRSTEP]: %g\n",lu->Control[UMFPACK_IRSTEP]);CHKERRQ(ierr);
376 
377   /* mat ordering */
378   if(!lu->PetscMatOdering) ierr = PetscViewerASCIIPrintf(viewer,"  UMFPACK default matrix ordering is used (not the PETSc matrix ordering) \n");CHKERRQ(ierr);
379 
380   PetscFunctionReturn(0);
381 }
382 
383 #undef __FUNCT__
384 #define __FUNCT__ "MatView_UMFPACK"
385 PetscErrorCode MatView_UMFPACK(Mat A,PetscViewer viewer)
386 {
387   PetscErrorCode    ierr;
388   PetscTruth        iascii;
389   PetscViewerFormat format;
390 
391   PetscFunctionBegin;
392   ierr = MatView_SeqAIJ(A,viewer);CHKERRQ(ierr);
393 
394   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr);
395   if (iascii) {
396     ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr);
397     if (format == PETSC_VIEWER_ASCII_INFO) {
398       ierr = MatFactorInfo_UMFPACK(A,viewer);CHKERRQ(ierr);
399     }
400   }
401   PetscFunctionReturn(0);
402 }
403 
404 /*MC
405   MAT_SOLVER_UMFPACK = "umfpack" - A matrix type providing direct solvers (LU) for sequential matrices
406   via the external package UMFPACK.
407 
408   config/configure.py --download-umfpack to install PETSc to use UMFPACK
409 
410   Consult UMFPACK documentation for more information about the Control parameters
411   which correspond to the options database keys below.
412 
413   Options Database Keys:
414 + -mat_umfpack_prl - UMFPACK print level: Control[UMFPACK_PRL]
415 . -mat_umfpack_dense_col <alpha_c> - UMFPACK dense column threshold: Control[UMFPACK_DENSE_COL]
416 . -mat_umfpack_block_size <bs> - UMFPACK block size for BLAS-Level 3 calls: Control[UMFPACK_BLOCK_SIZE]
417 . -mat_umfpack_pivot_tolerance <delta> - UMFPACK partial pivot tolerance: Control[UMFPACK_PIVOT_TOLERANCE]
418 . -mat_umfpack_alloc_init <delta> - UMFPACK factorized matrix allocation modifier: Control[UMFPACK_ALLOC_INIT]
419 - -mat_umfpack_irstep <maxit> - UMFPACK maximum number of iterative refinement steps: Control[UMFPACK_IRSTEP]
420 
421    Level: beginner
422 
423 .seealso: PCLU
424 M*/
425 
426 
427 
428