xref: /petsc/src/mat/impls/aij/seq/umfpack/umfpack.c (revision 5dffd610ba0f7d1b98cf004897cd9ae601ebd0a5)
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            B;
181   Mat_SeqAIJ     *mat=(Mat_SeqAIJ*)A->data;
182   Mat_UMFPACK    *lu = (Mat_UMFPACK*)((*F)->spptr);
183   PetscErrorCode ierr;
184   int            i,m=A->rmap.n,n=A->cmap.n,*ra,idx;
185   UF_long        status;
186   PetscScalar    *av=mat->a;
187 
188   PetscFunctionBegin;
189   if (lu->PetscMatOdering) {
190     ierr = ISGetIndices(r,&ra);CHKERRQ(ierr);
191     ierr = PetscMalloc(m*sizeof(UF_long),&lu->perm_c);CHKERRQ(ierr);
192     /* we cannot simply memcpy on 64 bit archs */
193     for(i = 0; i < m; i++) lu->perm_c[i] = ra[i];
194     ierr = ISRestoreIndices(r,&ra);CHKERRQ(ierr);
195   }
196 
197   if(sizeof(UF_long) != sizeof(int)){
198     /* we cannot directly use mat->i and mat->j on 64 bit archs */
199     ierr = PetscMalloc((m+1)*sizeof(UF_long),&lu->ai);CHKERRQ(ierr);
200     ierr = PetscMalloc(mat->nz*sizeof(UF_long),&lu->aj);CHKERRQ(ierr);
201     for(i = 0; i < m + 1; i++) lu->ai[i] = mat->i[i];
202     for(i = 0; i < mat->nz; i++) lu->aj[i] = mat->j[i];
203   }
204   else{
205     lu->ai = (UF_long*)mat->i;
206     lu->aj = (UF_long*)mat->j;
207   }
208 
209   /* print the control parameters */
210 #if defined(PETSC_USE_COMPLEX)
211   if(lu->Control[UMFPACK_PRL] > 1) umfpack_zl_report_control(lu->Control);
212 #else
213   if(lu->Control[UMFPACK_PRL] > 1) umfpack_dl_report_control(lu->Control);
214 #endif
215 
216   /* symbolic factorization of A' */
217   /* ---------------------------------------------------------------------- */
218 #if defined(PETSC_USE_COMPLEX)
219   if (lu->PetscMatOdering) { /* use Petsc row ordering */
220     status = umfpack_zl_qsymbolic(n,m,lu->ai,lu->aj,(double*)av,NULL,
221 				  lu->perm_c,&lu->Symbolic,lu->Control,lu->Info);
222   } else { /* use Umfpack col ordering */
223     status = umfpack_zl_symbolic(n,m,lu->ai,lu->aj,(double*)av,NULL,
224 				 &lu->Symbolic,lu->Control,lu->Info);
225   }
226   if (status < 0){
227     umfpack_zl_report_info(lu->Control, lu->Info);
228     umfpack_zl_report_status(lu->Control, status);
229     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_symbolic failed");
230   }
231   /* report sumbolic factorization of A' when Control[PRL] > 3 */
232   (void) umfpack_zl_report_symbolic(lu->Symbolic, lu->Control);
233 #else
234   if (lu->PetscMatOdering) { /* use Petsc row ordering */
235     status = umfpack_dl_qsymbolic(n,m,lu->ai,lu->aj,av,
236 				  lu->perm_c,&lu->Symbolic,lu->Control,lu->Info);
237   } else { /* use Umfpack col ordering */
238     status = umfpack_dl_symbolic(n,m,lu->ai,lu->aj,av,
239 				 &lu->Symbolic,lu->Control,lu->Info);
240   }
241   if (status < 0){
242     umfpack_dl_report_info(lu->Control, lu->Info);
243     umfpack_dl_report_status(lu->Control, status);
244     SETERRQ(PETSC_ERR_LIB,"umfpack_dl_symbolic failed");
245   }
246   /* report sumbolic factorization of A' when Control[PRL] > 3 */
247   (void) umfpack_dl_report_symbolic(lu->Symbolic, lu->Control);
248 #endif
249 
250   lu->flg = DIFFERENT_NONZERO_PATTERN;
251   lu->av  = av;
252   lu->CleanUpUMFPACK = PETSC_TRUE;
253   PetscFunctionReturn(0);
254 }
255 
256 #undef __FUNCT__
257 #define __FUNCT__ "MatGetFactor_seqaij_umfpack"
258 PetscErrorCode MatGetFactor_seqaij_umfpack(Mat A,MatFactorType ftype,Mat *F)
259 {
260   Mat            B;
261   Mat_SeqAIJ     *mat=(Mat_SeqAIJ*)A->data;
262   Mat_UMFPACK    *lu;
263   PetscErrorCode ierr;
264   int            i,m=A->rmap.n,n=A->cmap.n,*ra,idx;
265   UF_long        status;
266 
267   PetscScalar    *av=mat->a;
268   const char     *strategy[]={"AUTO","UNSYMMETRIC","SYMMETRIC","2BY2"},
269                  *scale[]={"NONE","SUM","MAX"};
270   PetscTruth     flg;
271 
272   PetscFunctionBegin;
273   /* Create the factorization matrix F */
274   ierr = MatCreate(((PetscObject)A)->comm,&B);CHKERRQ(ierr);
275   ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr);
276   ierr = MatSetType(B,((PetscObject)A)->type_name);CHKERRQ(ierr);
277   ierr = MatSeqAIJSetPreallocation(B,0,PETSC_NULL);CHKERRQ(ierr);
278   ierr = PetscNewLog(B,Mat_UMFPACK,&lu);CHKERRQ(ierr);
279   B->spptr                 = lu;
280   B->ops->lufactornumeric  = MatLUFactorNumeric_UMFPACK;
281   B->ops->lufactorsymbolic = MatLUFactorSymbolic_UMFPACK;
282   B->ops->solve            = MatSolve_UMFPACK;
283   B->ops->destroy          = MatDestroy_UMFPACK;
284   B->factor                = MAT_FACTOR_LU;
285   B->assembled             = PETSC_TRUE;  /* required by -ksp_view */
286   B->preallocated          = PETSC_TRUE;
287 
288   /* initializations */
289   /* ------------------------------------------------*/
290   /* get the default control parameters */
291 #if defined(PETSC_USE_COMPLEX)
292   umfpack_zl_defaults(lu->Control);
293 #else
294   umfpack_dl_defaults(lu->Control);
295 #endif
296   lu->perm_c = PETSC_NULL;  /* use defaul UMFPACK col permutation */
297   lu->Control[UMFPACK_IRSTEP] = 0; /* max num of iterative refinement steps to attempt */
298 
299   ierr = PetscOptionsBegin(((PetscObject)A)->comm,((PetscObject)A)->prefix,"UMFPACK Options","Mat");CHKERRQ(ierr);
300   /* Control parameters used by reporting routiones */
301   ierr = PetscOptionsReal("-mat_umfpack_prl","Control[UMFPACK_PRL]","None",lu->Control[UMFPACK_PRL],&lu->Control[UMFPACK_PRL],PETSC_NULL);CHKERRQ(ierr);
302 
303   /* Control parameters for symbolic factorization */
304   ierr = PetscOptionsEList("-mat_umfpack_strategy","ordering and pivoting strategy","None",strategy,4,strategy[0],&idx,&flg);CHKERRQ(ierr);
305   if (flg) {
306     switch (idx){
307     case 0: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_AUTO; break;
308     case 1: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_UNSYMMETRIC; break;
309     case 2: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_SYMMETRIC; break;
310     case 3: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_2BY2; break;
311     }
312   }
313   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);
314   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);
315   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);
316   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);
317   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);
318   ierr = PetscOptionsReal("-mat_umfpack_fixq","Control[UMFPACK_FIXQ]","None",lu->Control[UMFPACK_FIXQ],&lu->Control[UMFPACK_FIXQ],PETSC_NULL);CHKERRQ(ierr);
319   ierr = PetscOptionsReal("-mat_umfpack_aggressive","Control[UMFPACK_AGGRESSIVE]","None",lu->Control[UMFPACK_AGGRESSIVE],&lu->Control[UMFPACK_AGGRESSIVE],PETSC_NULL);CHKERRQ(ierr);
320 
321   /* Control parameters used by numeric factorization */
322   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);
323   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);
324   ierr = PetscOptionsEList("-mat_umfpack_scale","Control[UMFPACK_SCALE]","None",scale,3,scale[0],&idx,&flg);CHKERRQ(ierr);
325   if (flg) {
326     switch (idx){
327     case 0: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_NONE; break;
328     case 1: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_SUM; break;
329     case 2: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_MAX; break;
330     }
331   }
332   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);
333   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);
334   ierr = PetscOptionsReal("-mat_umfpack_droptol","Control[UMFPACK_DROPTOL]","None",lu->Control[UMFPACK_DROPTOL],&lu->Control[UMFPACK_DROPTOL],PETSC_NULL);CHKERRQ(ierr);
335 
336   /* Control parameters used by solve */
337   ierr = PetscOptionsReal("-mat_umfpack_irstep","Control[UMFPACK_IRSTEP]","None",lu->Control[UMFPACK_IRSTEP],&lu->Control[UMFPACK_IRSTEP],PETSC_NULL);CHKERRQ(ierr);
338 
339   /* use Petsc mat ordering (note: size is for the transpose, and PETSc r = Umfpack perm_c) */
340   ierr = PetscOptionsHasName(PETSC_NULL,"-pc_factor_mat_ordering_type",&lu->PetscMatOdering);CHKERRQ(ierr);
341   PetscOptionsEnd();
342   *F = B;
343   PetscFunctionReturn(0);
344 }
345 
346 #undef __FUNCT__
347 #define __FUNCT__ "MatFactorInfo_UMFPACK"
348 PetscErrorCode MatFactorInfo_UMFPACK(Mat A,PetscViewer viewer)
349 {
350   Mat_UMFPACK    *lu= (Mat_UMFPACK*)A->spptr;
351   PetscErrorCode ierr;
352 
353   PetscFunctionBegin;
354   /* check if matrix is UMFPACK type */
355   if (A->ops->solve != MatSolve_UMFPACK) PetscFunctionReturn(0);
356 
357   ierr = PetscViewerASCIIPrintf(viewer,"UMFPACK run parameters:\n");CHKERRQ(ierr);
358   /* Control parameters used by reporting routiones */
359   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_PRL]: %g\n",lu->Control[UMFPACK_PRL]);CHKERRQ(ierr);
360 
361   /* Control parameters used by symbolic factorization */
362   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_STRATEGY]: %g\n",lu->Control[UMFPACK_STRATEGY]);CHKERRQ(ierr);
363   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DENSE_COL]: %g\n",lu->Control[UMFPACK_DENSE_COL]);CHKERRQ(ierr);
364   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DENSE_ROW]: %g\n",lu->Control[UMFPACK_DENSE_ROW]);CHKERRQ(ierr);
365   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_AMD_DENSE]: %g\n",lu->Control[UMFPACK_AMD_DENSE]);CHKERRQ(ierr);
366   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_BLOCK_SIZE]: %g\n",lu->Control[UMFPACK_BLOCK_SIZE]);CHKERRQ(ierr);
367   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_2BY2_TOLERANCE]: %g\n",lu->Control[UMFPACK_2BY2_TOLERANCE]);CHKERRQ(ierr);
368   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_FIXQ]: %g\n",lu->Control[UMFPACK_FIXQ]);CHKERRQ(ierr);
369   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_AGGRESSIVE]: %g\n",lu->Control[UMFPACK_AGGRESSIVE]);CHKERRQ(ierr);
370 
371   /* Control parameters used by numeric factorization */
372   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_PIVOT_TOLERANCE]);CHKERRQ(ierr);
373   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_SYM_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_SYM_PIVOT_TOLERANCE]);CHKERRQ(ierr);
374   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_SCALE]: %g\n",lu->Control[UMFPACK_SCALE]);CHKERRQ(ierr);
375   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_ALLOC_INIT]: %g\n",lu->Control[UMFPACK_ALLOC_INIT]);CHKERRQ(ierr);
376   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_DROPTOL]: %g\n",lu->Control[UMFPACK_DROPTOL]);CHKERRQ(ierr);
377 
378   /* Control parameters used by solve */
379   ierr = PetscViewerASCIIPrintf(viewer,"  Control[UMFPACK_IRSTEP]: %g\n",lu->Control[UMFPACK_IRSTEP]);CHKERRQ(ierr);
380 
381   /* mat ordering */
382   if(!lu->PetscMatOdering) ierr = PetscViewerASCIIPrintf(viewer,"  UMFPACK default matrix ordering is used (not the PETSc matrix ordering) \n");CHKERRQ(ierr);
383 
384   PetscFunctionReturn(0);
385 }
386 
387 #undef __FUNCT__
388 #define __FUNCT__ "MatView_UMFPACK"
389 PetscErrorCode MatView_UMFPACK(Mat A,PetscViewer viewer)
390 {
391   PetscErrorCode    ierr;
392   PetscTruth        iascii;
393   PetscViewerFormat format;
394   Mat_UMFPACK       *lu=(Mat_UMFPACK*)(A->spptr);
395 
396   PetscFunctionBegin;
397   ierr = MatView_SeqAIJ(A,viewer);CHKERRQ(ierr);
398 
399   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr);
400   if (iascii) {
401     ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr);
402     if (format == PETSC_VIEWER_ASCII_INFO) {
403       ierr = MatFactorInfo_UMFPACK(A,viewer);CHKERRQ(ierr);
404     }
405   }
406   PetscFunctionReturn(0);
407 }
408 
409 /*MC
410   MAT_SOLVER_UMFPACK = "umfpack" - A matrix type providing direct solvers (LU) for sequential matrices
411   via the external package UMFPACK.
412 
413   config/configure.py --download-umfpack to install PETSc to use UMFPACK
414 
415   Consult UMFPACK documentation for more information about the Control parameters
416   which correspond to the options database keys below.
417 
418   Options Database Keys:
419 + -mat_umfpack_prl - UMFPACK print level: Control[UMFPACK_PRL]
420 . -mat_umfpack_dense_col <alpha_c> - UMFPACK dense column threshold: Control[UMFPACK_DENSE_COL]
421 . -mat_umfpack_block_size <bs> - UMFPACK block size for BLAS-Level 3 calls: Control[UMFPACK_BLOCK_SIZE]
422 . -mat_umfpack_pivot_tolerance <delta> - UMFPACK partial pivot tolerance: Control[UMFPACK_PIVOT_TOLERANCE]
423 . -mat_umfpack_alloc_init <delta> - UMFPACK factorized matrix allocation modifier: Control[UMFPACK_ALLOC_INIT]
424 - -mat_umfpack_irstep <maxit> - UMFPACK maximum number of iterative refinement steps: Control[UMFPACK_IRSTEP]
425 
426    Level: beginner
427 
428 .seealso: PCLU
429 M*/
430 
431 
432 
433