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