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 (*F)->ops->lufactornumeric = MatLUFactorNumeric_UMFPACK; 253 (*F)->ops->solve = MatSolve_UMFPACK; 254 PetscFunctionReturn(0); 255 } 256 257 #undef __FUNCT__ 258 #define __FUNCT__ "MatFactorInfo_UMFPACK" 259 PetscErrorCode MatFactorInfo_UMFPACK(Mat A,PetscViewer viewer) 260 { 261 Mat_UMFPACK *lu= (Mat_UMFPACK*)A->spptr; 262 PetscErrorCode ierr; 263 264 PetscFunctionBegin; 265 /* check if matrix is UMFPACK type */ 266 if (A->ops->solve != MatSolve_UMFPACK) PetscFunctionReturn(0); 267 268 ierr = PetscViewerASCIIPrintf(viewer,"UMFPACK run parameters:\n");CHKERRQ(ierr); 269 /* Control parameters used by reporting routiones */ 270 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_PRL]: %g\n",lu->Control[UMFPACK_PRL]);CHKERRQ(ierr); 271 272 /* Control parameters used by symbolic factorization */ 273 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_STRATEGY]: %g\n",lu->Control[UMFPACK_STRATEGY]);CHKERRQ(ierr); 274 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_DENSE_COL]: %g\n",lu->Control[UMFPACK_DENSE_COL]);CHKERRQ(ierr); 275 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_DENSE_ROW]: %g\n",lu->Control[UMFPACK_DENSE_ROW]);CHKERRQ(ierr); 276 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_AMD_DENSE]: %g\n",lu->Control[UMFPACK_AMD_DENSE]);CHKERRQ(ierr); 277 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_BLOCK_SIZE]: %g\n",lu->Control[UMFPACK_BLOCK_SIZE]);CHKERRQ(ierr); 278 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_2BY2_TOLERANCE]: %g\n",lu->Control[UMFPACK_2BY2_TOLERANCE]);CHKERRQ(ierr); 279 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_FIXQ]: %g\n",lu->Control[UMFPACK_FIXQ]);CHKERRQ(ierr); 280 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_AGGRESSIVE]: %g\n",lu->Control[UMFPACK_AGGRESSIVE]);CHKERRQ(ierr); 281 282 /* Control parameters used by numeric factorization */ 283 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_PIVOT_TOLERANCE]);CHKERRQ(ierr); 284 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_SYM_PIVOT_TOLERANCE]: %g\n",lu->Control[UMFPACK_SYM_PIVOT_TOLERANCE]);CHKERRQ(ierr); 285 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_SCALE]: %g\n",lu->Control[UMFPACK_SCALE]);CHKERRQ(ierr); 286 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_ALLOC_INIT]: %g\n",lu->Control[UMFPACK_ALLOC_INIT]);CHKERRQ(ierr); 287 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_DROPTOL]: %g\n",lu->Control[UMFPACK_DROPTOL]);CHKERRQ(ierr); 288 289 /* Control parameters used by solve */ 290 ierr = PetscViewerASCIIPrintf(viewer," Control[UMFPACK_IRSTEP]: %g\n",lu->Control[UMFPACK_IRSTEP]);CHKERRQ(ierr); 291 292 /* mat ordering */ 293 if(!lu->PetscMatOdering) ierr = PetscViewerASCIIPrintf(viewer," UMFPACK default matrix ordering is used (not the PETSc matrix ordering) \n");CHKERRQ(ierr); 294 295 PetscFunctionReturn(0); 296 } 297 298 #undef __FUNCT__ 299 #define __FUNCT__ "MatView_UMFPACK" 300 PetscErrorCode MatView_UMFPACK(Mat A,PetscViewer viewer) 301 { 302 PetscErrorCode ierr; 303 PetscTruth iascii; 304 PetscViewerFormat format; 305 306 PetscFunctionBegin; 307 ierr = MatView_SeqAIJ(A,viewer);CHKERRQ(ierr); 308 309 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 310 if (iascii) { 311 ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr); 312 if (format == PETSC_VIEWER_ASCII_INFO) { 313 ierr = MatFactorInfo_UMFPACK(A,viewer);CHKERRQ(ierr); 314 } 315 } 316 PetscFunctionReturn(0); 317 } 318 319 /*MC 320 MAT_SOLVER_UMFPACK = "umfpack" - A matrix type providing direct solvers (LU) for sequential matrices 321 via the external package UMFPACK. 322 323 config/configure.py --download-umfpack to install PETSc to use UMFPACK 324 325 Consult UMFPACK documentation for more information about the Control parameters 326 which correspond to the options database keys below. 327 328 Options Database Keys: 329 + -mat_umfpack_prl - UMFPACK print level: Control[UMFPACK_PRL] 330 . -mat_umfpack_strategy <AUTO> (choose one of) AUTO UNSYMMETRIC SYMMETRIC 2BY2 331 . -mat_umfpack_dense_col <alpha_c> - UMFPACK dense column threshold: Control[UMFPACK_DENSE_COL] 332 . -mat_umfpack_dense_row <0.2>: Control[UMFPACK_DENSE_ROW] 333 . -mat_umfpack_amd_dense <10>: Control[UMFPACK_AMD_DENSE] 334 . -mat_umfpack_block_size <bs> - UMFPACK block size for BLAS-Level 3 calls: Control[UMFPACK_BLOCK_SIZE] 335 . -mat_umfpack_2by2_tolerance <0.01>: Control[UMFPACK_2BY2_TOLERANCE] 336 . -mat_umfpack_fixq <0>: Control[UMFPACK_FIXQ] 337 . -mat_umfpack_aggressive <1>: Control[UMFPACK_AGGRESSIVE] 338 . -mat_umfpack_pivot_tolerance <delta> - UMFPACK partial pivot tolerance: Control[UMFPACK_PIVOT_TOLERANCE] 339 . -mat_umfpack_sym_pivot_tolerance <0.001>: Control[UMFPACK_SYM_PIVOT_TOLERANCE] 340 . -mat_umfpack_scale <NONE> (choose one of) NONE SUM MAX 341 . -mat_umfpack_alloc_init <delta> - UMFPACK factorized matrix allocation modifier: Control[UMFPACK_ALLOC_INIT] 342 . -mat_umfpack_droptol <0>: Control[UMFPACK_DROPTOL] 343 - -mat_umfpack_irstep <maxit> - UMFPACK maximum number of iterative refinement steps: Control[UMFPACK_IRSTEP] 344 345 Level: beginner 346 347 .seealso: PCLU, MAT_SOLVER_SUPERLU, MAT_SOLVER_MUMPS, MAT_SOLVER_SPOOLES 348 M*/ 349 EXTERN_C_BEGIN 350 #undef __FUNCT__ 351 #define __FUNCT__ "MatGetFactor_seqaij_umfpack" 352 PetscErrorCode MatGetFactor_seqaij_umfpack(Mat A,MatFactorType ftype,Mat *F) 353 { 354 Mat B; 355 Mat_UMFPACK *lu; 356 PetscErrorCode ierr; 357 int m=A->rmap->n,n=A->cmap->n,idx; 358 359 const char *strategy[]={"AUTO","UNSYMMETRIC","SYMMETRIC","2BY2"}, 360 *scale[]={"NONE","SUM","MAX"}; 361 PetscTruth flg; 362 363 PetscFunctionBegin; 364 /* Create the factorization matrix F */ 365 ierr = MatCreate(((PetscObject)A)->comm,&B);CHKERRQ(ierr); 366 ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr); 367 ierr = MatSetType(B,((PetscObject)A)->type_name);CHKERRQ(ierr); 368 ierr = MatSeqAIJSetPreallocation(B,0,PETSC_NULL);CHKERRQ(ierr); 369 ierr = PetscNewLog(B,Mat_UMFPACK,&lu);CHKERRQ(ierr); 370 B->spptr = lu; 371 B->ops->lufactorsymbolic = MatLUFactorSymbolic_UMFPACK; 372 B->ops->destroy = MatDestroy_UMFPACK; 373 B->ops->view = MatView_UMFPACK; 374 B->factor = MAT_FACTOR_LU; 375 B->assembled = PETSC_TRUE; /* required by -ksp_view */ 376 B->preallocated = PETSC_TRUE; 377 378 /* initializations */ 379 /* ------------------------------------------------*/ 380 /* get the default control parameters */ 381 #if defined(PETSC_USE_COMPLEX) 382 umfpack_zl_defaults(lu->Control); 383 #else 384 umfpack_dl_defaults(lu->Control); 385 #endif 386 lu->perm_c = PETSC_NULL; /* use defaul UMFPACK col permutation */ 387 lu->Control[UMFPACK_IRSTEP] = 0; /* max num of iterative refinement steps to attempt */ 388 389 ierr = PetscOptionsBegin(((PetscObject)A)->comm,((PetscObject)A)->prefix,"UMFPACK Options","Mat");CHKERRQ(ierr); 390 /* Control parameters used by reporting routiones */ 391 ierr = PetscOptionsReal("-mat_umfpack_prl","Control[UMFPACK_PRL]","None",lu->Control[UMFPACK_PRL],&lu->Control[UMFPACK_PRL],PETSC_NULL);CHKERRQ(ierr); 392 393 /* Control parameters for symbolic factorization */ 394 ierr = PetscOptionsEList("-mat_umfpack_strategy","ordering and pivoting strategy","None",strategy,4,strategy[0],&idx,&flg);CHKERRQ(ierr); 395 if (flg) { 396 switch (idx){ 397 case 0: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_AUTO; break; 398 case 1: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_UNSYMMETRIC; break; 399 case 2: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_SYMMETRIC; break; 400 case 3: lu->Control[UMFPACK_STRATEGY] = UMFPACK_STRATEGY_2BY2; break; 401 } 402 } 403 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); 404 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); 405 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); 406 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); 407 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); 408 ierr = PetscOptionsReal("-mat_umfpack_fixq","Control[UMFPACK_FIXQ]","None",lu->Control[UMFPACK_FIXQ],&lu->Control[UMFPACK_FIXQ],PETSC_NULL);CHKERRQ(ierr); 409 ierr = PetscOptionsReal("-mat_umfpack_aggressive","Control[UMFPACK_AGGRESSIVE]","None",lu->Control[UMFPACK_AGGRESSIVE],&lu->Control[UMFPACK_AGGRESSIVE],PETSC_NULL);CHKERRQ(ierr); 410 411 /* Control parameters used by numeric factorization */ 412 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); 413 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); 414 ierr = PetscOptionsEList("-mat_umfpack_scale","Control[UMFPACK_SCALE]","None",scale,3,scale[0],&idx,&flg);CHKERRQ(ierr); 415 if (flg) { 416 switch (idx){ 417 case 0: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_NONE; break; 418 case 1: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_SUM; break; 419 case 2: lu->Control[UMFPACK_SCALE] = UMFPACK_SCALE_MAX; break; 420 } 421 } 422 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); 423 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); 424 ierr = PetscOptionsReal("-mat_umfpack_droptol","Control[UMFPACK_DROPTOL]","None",lu->Control[UMFPACK_DROPTOL],&lu->Control[UMFPACK_DROPTOL],PETSC_NULL);CHKERRQ(ierr); 425 426 /* Control parameters used by solve */ 427 ierr = PetscOptionsReal("-mat_umfpack_irstep","Control[UMFPACK_IRSTEP]","None",lu->Control[UMFPACK_IRSTEP],&lu->Control[UMFPACK_IRSTEP],PETSC_NULL);CHKERRQ(ierr); 428 429 /* use Petsc mat ordering (note: size is for the transpose, and PETSc r = Umfpack perm_c) */ 430 ierr = PetscOptionsHasName(PETSC_NULL,"-pc_factor_mat_ordering_type",&lu->PetscMatOdering);CHKERRQ(ierr); 431 PetscOptionsEnd(); 432 *F = B; 433 PetscFunctionReturn(0); 434 } 435 EXTERN_C_END 436 437