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