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