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