1dba47a55SKris Buschelman 29b54502bSHong Zhang /* 39b54502bSHong Zhang Defines a ILU factorization preconditioner for any Mat implementation 49b54502bSHong Zhang */ 5c6db04a5SJed Brown #include <../src/ksp/pc/impls/factor/ilu/ilu.h> /*I "petscpc.h" I*/ 69b54502bSHong Zhang 79b54502bSHong Zhang /* ------------------------------------------------------------------------------------------*/ 8afaefe49SHong Zhang #undef __FUNCT__ 92401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill_ILU" 107087cfbeSBarry Smith PetscErrorCode PCFactorSetReuseFill_ILU(PC pc,PetscBool flag) 112401956bSBarry Smith { 12075768bcSBarry Smith PC_ILU *lu = (PC_ILU*)pc->data; 132401956bSBarry Smith 142401956bSBarry Smith PetscFunctionBegin; 152401956bSBarry Smith lu->reusefill = flag; 162401956bSBarry Smith PetscFunctionReturn(0); 172401956bSBarry Smith } 182401956bSBarry Smith 192401956bSBarry Smith #undef __FUNCT__ 202401956bSBarry Smith #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal_ILU" 217087cfbeSBarry Smith PetscErrorCode PCFactorReorderForNonzeroDiagonal_ILU(PC pc,PetscReal z) 229b54502bSHong Zhang { 239b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 249b54502bSHong Zhang 259b54502bSHong Zhang PetscFunctionBegin; 269b54502bSHong Zhang ilu->nonzerosalongdiagonal = PETSC_TRUE; 272fa5cd67SKarl Rupp if (z == PETSC_DECIDE) ilu->nonzerosalongdiagonaltol = 1.e-10; 282fa5cd67SKarl Rupp else ilu->nonzerosalongdiagonaltol = z; 299b54502bSHong Zhang PetscFunctionReturn(0); 309b54502bSHong Zhang } 319b54502bSHong Zhang 329b54502bSHong Zhang #undef __FUNCT__ 33574deadeSBarry Smith #define __FUNCT__ "PCReset_ILU" 34574deadeSBarry Smith PetscErrorCode PCReset_ILU(PC pc) 359b54502bSHong Zhang { 369b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 379b54502bSHong Zhang PetscErrorCode ierr; 389b54502bSHong Zhang 399b54502bSHong Zhang PetscFunctionBegin; 40298cc208SBarry Smith if (!ilu->inplace) {ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);} 416bf464f9SBarry Smith if (ilu->row && ilu->col && ilu->row != ilu->col) {ierr = ISDestroy(&ilu->row);CHKERRQ(ierr);} 42fcfd50ebSBarry Smith ierr = ISDestroy(&ilu->col);CHKERRQ(ierr); 439b54502bSHong Zhang PetscFunctionReturn(0); 449b54502bSHong Zhang } 459b54502bSHong Zhang 469b54502bSHong Zhang #undef __FUNCT__ 47b7c853c4SBarry Smith #define __FUNCT__ "PCFactorSetDropTolerance_ILU" 487087cfbeSBarry Smith PetscErrorCode PCFactorSetDropTolerance_ILU(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount) 499b54502bSHong Zhang { 50075768bcSBarry Smith PC_ILU *ilu = (PC_ILU*)pc->data; 519b54502bSHong Zhang 529b54502bSHong Zhang PetscFunctionBegin; 534c9036c7SBarry Smith if (pc->setupcalled && (((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) { 54ce94432eSBarry Smith SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot change drop tolerance after using PC"); 559b54502bSHong Zhang } 56075768bcSBarry Smith ((PC_Factor*)ilu)->info.dt = dt; 57075768bcSBarry Smith ((PC_Factor*)ilu)->info.dtcol = dtcol; 58075768bcSBarry Smith ((PC_Factor*)ilu)->info.dtcount = dtcount; 594c9036c7SBarry Smith ((PC_Factor*)ilu)->info.usedt = 1.0; 609b54502bSHong Zhang PetscFunctionReturn(0); 619b54502bSHong Zhang } 629b54502bSHong Zhang 639b54502bSHong Zhang #undef __FUNCT__ 642401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering_ILU" 657087cfbeSBarry Smith PetscErrorCode PCFactorSetReuseOrdering_ILU(PC pc,PetscBool flag) 669b54502bSHong Zhang { 67075768bcSBarry Smith PC_ILU *ilu = (PC_ILU*)pc->data; 689b54502bSHong Zhang 699b54502bSHong Zhang PetscFunctionBegin; 709b54502bSHong Zhang ilu->reuseordering = flag; 719b54502bSHong Zhang PetscFunctionReturn(0); 729b54502bSHong Zhang } 739b54502bSHong Zhang 749b54502bSHong Zhang #undef __FUNCT__ 752401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseInPlace_ILU" 768e37d05fSBarry Smith PetscErrorCode PCFactorSetUseInPlace_ILU(PC pc,PetscBool flg) 779b54502bSHong Zhang { 78075768bcSBarry Smith PC_ILU *dir = (PC_ILU*)pc->data; 799b54502bSHong Zhang 809b54502bSHong Zhang PetscFunctionBegin; 818e37d05fSBarry Smith dir->inplace = flg; 828e37d05fSBarry Smith PetscFunctionReturn(0); 838e37d05fSBarry Smith } 848e37d05fSBarry Smith 858e37d05fSBarry Smith #undef __FUNCT__ 868e37d05fSBarry Smith #define __FUNCT__ "PCFactorGetUseInPlace_ILU" 878e37d05fSBarry Smith PetscErrorCode PCFactorGetUseInPlace_ILU(PC pc,PetscBool *flg) 888e37d05fSBarry Smith { 898e37d05fSBarry Smith PC_ILU *dir = (PC_ILU*)pc->data; 908e37d05fSBarry Smith 918e37d05fSBarry Smith PetscFunctionBegin; 928e37d05fSBarry Smith *flg = dir->inplace; 939b54502bSHong Zhang PetscFunctionReturn(0); 949b54502bSHong Zhang } 959b54502bSHong Zhang 969b54502bSHong Zhang #undef __FUNCT__ 979b54502bSHong Zhang #define __FUNCT__ "PCSetFromOptions_ILU" 984416b707SBarry Smith static PetscErrorCode PCSetFromOptions_ILU(PetscOptionItems *PetscOptionsObject,PC pc) 999b54502bSHong Zhang { 1009b54502bSHong Zhang PetscErrorCode ierr; 10178fc6b22SHong Zhang PetscInt itmp; 1028afaa268SBarry Smith PetscBool flg,set; 1039b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 1049b54502bSHong Zhang PetscReal tol; 1059b54502bSHong Zhang 1069b54502bSHong Zhang PetscFunctionBegin; 107e55864a3SBarry Smith ierr = PetscOptionsHead(PetscOptionsObject,"ILU Options");CHKERRQ(ierr); 108e55864a3SBarry Smith ierr = PCSetFromOptions_Factor(PetscOptionsObject,pc);CHKERRQ(ierr); 1098ff23777SHong Zhang 110075768bcSBarry Smith ierr = PetscOptionsInt("-pc_factor_levels","levels of fill","PCFactorSetLevels",(PetscInt)((PC_Factor*)ilu)->info.levels,&itmp,&flg);CHKERRQ(ierr); 111075768bcSBarry Smith if (flg) ((PC_Factor*)ilu)->info.levels = itmp; 1122fa5cd67SKarl Rupp 1138afaa268SBarry Smith ierr = PetscOptionsBool("-pc_factor_diagonal_fill","Allow fill into empty diagonal entry","PCFactorSetAllowDiagonalFill",((PC_Factor*)ilu)->info.diagonal_fill ? PETSC_TRUE : PETSC_FALSE,&flg,&set);CHKERRQ(ierr); 1148afaa268SBarry Smith if (set) ((PC_Factor*)ilu)->info.diagonal_fill = (PetscReal) flg; 1152401956bSBarry Smith ierr = PetscOptionsName("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",&flg);CHKERRQ(ierr); 1169b54502bSHong Zhang if (flg) { 1179b54502bSHong Zhang tol = PETSC_DECIDE; 1182401956bSBarry Smith ierr = PetscOptionsReal("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",ilu->nonzerosalongdiagonaltol,&tol,0);CHKERRQ(ierr); 1192401956bSBarry Smith ierr = PCFactorReorderForNonzeroDiagonal(pc,tol);CHKERRQ(ierr); 1209b54502bSHong Zhang } 1219b54502bSHong Zhang 1229b54502bSHong Zhang ierr = PetscOptionsTail();CHKERRQ(ierr); 1239b54502bSHong Zhang PetscFunctionReturn(0); 1249b54502bSHong Zhang } 1259b54502bSHong Zhang 1269b54502bSHong Zhang #undef __FUNCT__ 1279b54502bSHong Zhang #define __FUNCT__ "PCView_ILU" 1289b54502bSHong Zhang static PetscErrorCode PCView_ILU(PC pc,PetscViewer viewer) 1299b54502bSHong Zhang { 1309b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 1319b54502bSHong Zhang PetscErrorCode ierr; 132ace3abfcSBarry Smith PetscBool iascii; 1339b54502bSHong Zhang 1349b54502bSHong Zhang PetscFunctionBegin; 135251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 1369b54502bSHong Zhang if (iascii) { 137914a5d51SHong Zhang if (ilu->inplace) { 138914a5d51SHong Zhang ierr = PetscViewerASCIIPrintf(viewer," ILU: in-place factorization\n");CHKERRQ(ierr); 1399b54502bSHong Zhang } else { 140914a5d51SHong Zhang ierr = PetscViewerASCIIPrintf(viewer," ILU: out-of-place factorization\n");CHKERRQ(ierr); 141145b9266SHong Zhang } 142145b9266SHong Zhang 143914a5d51SHong Zhang if (ilu->reusefill) {ierr = PetscViewerASCIIPrintf(viewer," ILU: Reusing fill from past factorization\n");CHKERRQ(ierr);} 144914a5d51SHong Zhang if (ilu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer," ILU: Reusing reordering from past factorization\n");CHKERRQ(ierr);} 1459b54502bSHong Zhang } 146914a5d51SHong Zhang ierr = PCView_Factor(pc,viewer);CHKERRQ(ierr); 1479b54502bSHong Zhang PetscFunctionReturn(0); 1489b54502bSHong Zhang } 1499b54502bSHong Zhang 1509b54502bSHong Zhang #undef __FUNCT__ 1519b54502bSHong Zhang #define __FUNCT__ "PCSetUp_ILU" 1529b54502bSHong Zhang static PetscErrorCode PCSetUp_ILU(PC pc) 1539b54502bSHong Zhang { 1549b54502bSHong Zhang PetscErrorCode ierr; 1559b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 156f3a39becSBarry Smith MatInfo info; 157ace3abfcSBarry Smith PetscBool flg; 15800c67f3bSHong Zhang const MatSolverPackage stype; 1599b54502bSHong Zhang 1609b54502bSHong Zhang PetscFunctionBegin; 161*c6e4fdc6SHong Zhang pc->failedreason = PC_NOERROR; 16292927226SBarry Smith /* ugly hack to change default, since it is not support by some matrix types */ 16392927226SBarry Smith if (((PC_Factor*)ilu)->info.shifttype == (PetscReal)MAT_SHIFT_NONZERO) { 16422d28d08SBarry Smith ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATSEQAIJ,&flg);CHKERRQ(ierr); 16592927226SBarry Smith if (!flg) { 16622d28d08SBarry Smith ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATMPIAIJ,&flg);CHKERRQ(ierr); 16792927226SBarry Smith if (!flg) { 16892927226SBarry Smith ((PC_Factor*)ilu)->info.shifttype = (PetscReal)MAT_SHIFT_INBLOCKS; 169955c1f14SBarry Smith PetscInfo(pc,"Changing shift type from NONZERO to INBLOCKS because block matrices do not support NONZERO\n");CHKERRQ(ierr); 17092927226SBarry Smith } 17192927226SBarry Smith } 17292927226SBarry Smith } 17392927226SBarry Smith 17484d44b13SHong Zhang ierr = MatSetErrorIfFailure(pc->pmat,pc->erroriffailure);CHKERRQ(ierr); 1759b54502bSHong Zhang if (ilu->inplace) { 1769b54502bSHong Zhang if (!pc->setupcalled) { 1779b54502bSHong Zhang 1789b54502bSHong Zhang /* In-place factorization only makes sense with the natural ordering, 1799b54502bSHong Zhang so we only need to get the ordering once, even if nonzero structure changes */ 180075768bcSBarry Smith ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 1813bb1ff40SBarry Smith if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 1823bb1ff40SBarry Smith if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 1839b54502bSHong Zhang } 1849b54502bSHong Zhang 1859b54502bSHong Zhang /* In place ILU only makes sense with fill factor of 1.0 because 1869b54502bSHong Zhang cannot have levels of fill */ 187075768bcSBarry Smith ((PC_Factor*)ilu)->info.fill = 1.0; 18875567043SBarry Smith ((PC_Factor*)ilu)->info.diagonal_fill = 0.0; 1892fa5cd67SKarl Rupp 190365a8a9eSBarry Smith ierr = MatILUFactor(pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);CHKERRQ(ierr); 1914cccfbddSHong Zhang if (pc->pmat->errortype) { /* Factor() fails */ 1924cccfbddSHong Zhang pc->failedreason = (PCFailedReason)pc->pmat->errortype; 1936baea169SHong Zhang PetscFunctionReturn(0); 1946baea169SHong Zhang } 1956baea169SHong Zhang 196075768bcSBarry Smith ((PC_Factor*)ilu)->fact = pc->pmat; 197b33856dcSBarry Smith /* must update the pc record of the matrix state or the PC will attempt to run PCSetUp() yet again */ 198b33856dcSBarry Smith ierr = PetscObjectStateGet((PetscObject)pc->pmat,&pc->matstate);CHKERRQ(ierr); 1999b54502bSHong Zhang } else { 2004cccfbddSHong Zhang Mat F; 2019b54502bSHong Zhang if (!pc->setupcalled) { 2029b54502bSHong Zhang /* first time in so compute reordering and symbolic factorization */ 203075768bcSBarry Smith ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 2043bb1ff40SBarry Smith if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 2053bb1ff40SBarry Smith if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 2069b54502bSHong Zhang /* Remove zeros along diagonal? */ 2079b54502bSHong Zhang if (ilu->nonzerosalongdiagonal) { 2089b54502bSHong Zhang ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr); 2099b54502bSHong Zhang } 21086daa42cSBarry Smith if (!((PC_Factor*)ilu)->fact) { 21111bfe483SHong Zhang ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 212a1f19f5aSHong Zhang } 213075768bcSBarry Smith ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 214075768bcSBarry Smith ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 2152fa5cd67SKarl Rupp 216f3a39becSBarry Smith ilu->actualfill = info.fill_ratio_needed; 2172fa5cd67SKarl Rupp 2183bb1ff40SBarry Smith ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 2199b54502bSHong Zhang } else if (pc->flag != SAME_NONZERO_PATTERN) { 2209b54502bSHong Zhang if (!ilu->reuseordering) { 2219b54502bSHong Zhang /* compute a new ordering for the ILU */ 222fcfd50ebSBarry Smith ierr = ISDestroy(&ilu->row);CHKERRQ(ierr); 223fcfd50ebSBarry Smith ierr = ISDestroy(&ilu->col);CHKERRQ(ierr); 224075768bcSBarry Smith ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 2253bb1ff40SBarry Smith if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 2263bb1ff40SBarry Smith if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 2279b54502bSHong Zhang /* Remove zeros along diagonal? */ 2289b54502bSHong Zhang if (ilu->nonzerosalongdiagonal) { 2299b54502bSHong Zhang ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr); 2309b54502bSHong Zhang } 2319b54502bSHong Zhang } 2326bf464f9SBarry Smith ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 23314798fb4SJed Brown ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 234075768bcSBarry Smith ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 235075768bcSBarry Smith ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 2362fa5cd67SKarl Rupp 237f3a39becSBarry Smith ilu->actualfill = info.fill_ratio_needed; 2382fa5cd67SKarl Rupp 2393bb1ff40SBarry Smith ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 2409b54502bSHong Zhang } 2414cccfbddSHong Zhang F = ((PC_Factor*)ilu)->fact; 2424cccfbddSHong Zhang if (F->errortype) { /* FactorSymbolic() fails */ 2434cccfbddSHong Zhang pc->failedreason = (PCFailedReason)F->errortype; 2446baea169SHong Zhang PetscFunctionReturn(0); 2456baea169SHong Zhang } 2466baea169SHong Zhang 247075768bcSBarry Smith ierr = MatLUFactorNumeric(((PC_Factor*)ilu)->fact,pc->pmat,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 2484cccfbddSHong Zhang if (F->errortype) { /* FactorNumeric() fails */ 2494cccfbddSHong Zhang pc->failedreason = (PCFailedReason)F->errortype; 2506baea169SHong Zhang } 2519b54502bSHong Zhang } 25200c67f3bSHong Zhang 25300c67f3bSHong Zhang ierr = PCFactorGetMatSolverPackage(pc,&stype);CHKERRQ(ierr); 25400c67f3bSHong Zhang if (!stype) { 25500c67f3bSHong Zhang ierr = PCFactorSetMatSolverPackage(pc,((PC_Factor*)ilu)->fact->solvertype);CHKERRQ(ierr); 25600c67f3bSHong Zhang } 2579b54502bSHong Zhang PetscFunctionReturn(0); 2589b54502bSHong Zhang } 2599b54502bSHong Zhang 2609b54502bSHong Zhang #undef __FUNCT__ 2619b54502bSHong Zhang #define __FUNCT__ "PCDestroy_ILU" 2629b54502bSHong Zhang static PetscErrorCode PCDestroy_ILU(PC pc) 2639b54502bSHong Zhang { 2649b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 2659b54502bSHong Zhang PetscErrorCode ierr; 2669b54502bSHong Zhang 2679b54502bSHong Zhang PetscFunctionBegin; 268574deadeSBarry Smith ierr = PCReset_ILU(pc);CHKERRQ(ierr); 269503cfb0cSBarry Smith ierr = PetscFree(((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr); 270503cfb0cSBarry Smith ierr = PetscFree(((PC_Factor*)ilu)->ordering);CHKERRQ(ierr); 271c31cb41cSBarry Smith ierr = PetscFree(pc->data);CHKERRQ(ierr); 2729b54502bSHong Zhang PetscFunctionReturn(0); 2739b54502bSHong Zhang } 2749b54502bSHong Zhang 2759b54502bSHong Zhang #undef __FUNCT__ 2769b54502bSHong Zhang #define __FUNCT__ "PCApply_ILU" 2779b54502bSHong Zhang static PetscErrorCode PCApply_ILU(PC pc,Vec x,Vec y) 2789b54502bSHong Zhang { 2799b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 2809b54502bSHong Zhang PetscErrorCode ierr; 2819b54502bSHong Zhang 2829b54502bSHong Zhang PetscFunctionBegin; 283075768bcSBarry Smith ierr = MatSolve(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr); 2849b54502bSHong Zhang PetscFunctionReturn(0); 2859b54502bSHong Zhang } 2869b54502bSHong Zhang 2879b54502bSHong Zhang #undef __FUNCT__ 2889b54502bSHong Zhang #define __FUNCT__ "PCApplyTranspose_ILU" 2899b54502bSHong Zhang static PetscErrorCode PCApplyTranspose_ILU(PC pc,Vec x,Vec y) 2909b54502bSHong Zhang { 2919b54502bSHong Zhang PC_ILU *ilu = (PC_ILU*)pc->data; 2929b54502bSHong Zhang PetscErrorCode ierr; 2939b54502bSHong Zhang 2949b54502bSHong Zhang PetscFunctionBegin; 295075768bcSBarry Smith ierr = MatSolveTranspose(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr); 2969b54502bSHong Zhang PetscFunctionReturn(0); 2979b54502bSHong Zhang } 2989b54502bSHong Zhang 299f0b9ad6cSBarry Smith #undef __FUNCT__ 300f0b9ad6cSBarry Smith #define __FUNCT__ "PCApplySymmetricLeft_ILU" 301f0b9ad6cSBarry Smith static PetscErrorCode PCApplySymmetricLeft_ILU(PC pc,Vec x,Vec y) 302f0b9ad6cSBarry Smith { 303f0b9ad6cSBarry Smith PetscErrorCode ierr; 304f0b9ad6cSBarry Smith PC_ILU *icc = (PC_ILU*)pc->data; 305f0b9ad6cSBarry Smith 306f0b9ad6cSBarry Smith PetscFunctionBegin; 307f0b9ad6cSBarry Smith ierr = MatForwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr); 308f0b9ad6cSBarry Smith PetscFunctionReturn(0); 309f0b9ad6cSBarry Smith } 310f0b9ad6cSBarry Smith 311f0b9ad6cSBarry Smith #undef __FUNCT__ 312f0b9ad6cSBarry Smith #define __FUNCT__ "PCApplySymmetricRight_ILU" 313f0b9ad6cSBarry Smith static PetscErrorCode PCApplySymmetricRight_ILU(PC pc,Vec x,Vec y) 314f0b9ad6cSBarry Smith { 315f0b9ad6cSBarry Smith PetscErrorCode ierr; 316f0b9ad6cSBarry Smith PC_ILU *icc = (PC_ILU*)pc->data; 317f0b9ad6cSBarry Smith 318f0b9ad6cSBarry Smith PetscFunctionBegin; 319f0b9ad6cSBarry Smith ierr = MatBackwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr); 320f0b9ad6cSBarry Smith PetscFunctionReturn(0); 321f0b9ad6cSBarry Smith } 322f0b9ad6cSBarry Smith 3239b54502bSHong Zhang /*MC 3249b54502bSHong Zhang PCILU - Incomplete factorization preconditioners. 3259b54502bSHong Zhang 3269b54502bSHong Zhang Options Database Keys: 3272401956bSBarry Smith + -pc_factor_levels <k> - number of levels of fill for ILU(k) 3282401956bSBarry Smith . -pc_factor_in_place - only for ILU(0) with natural ordering, reuses the space of the matrix for 3299b54502bSHong Zhang its factorization (overwrites original matrix) 3302401956bSBarry Smith . -pc_factor_diagonal_fill - fill in a zero diagonal even if levels of fill indicate it wouldn't be fill 3312401956bSBarry Smith . -pc_factor_reuse_ordering - reuse ordering of factorized matrix from previous factorization 33255ba2a51SBarry Smith . -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1 3332401956bSBarry Smith . -pc_factor_nonzeros_along_diagonal - reorder the matrix before factorization to remove zeros from the diagonal, 3349b54502bSHong Zhang this decreases the chance of getting a zero pivot 3352401956bSBarry Smith . -pc_factor_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix 336967c93d3SBarry Smith - -pc_factor_pivot_in_blocks - for block ILU(k) factorization, i.e. with BAIJ matrices with block size larger 3379b54502bSHong Zhang than 1 the diagonal blocks are factored with partial pivoting (this increases the 3389b54502bSHong Zhang stability of the ILU factorization 3399b54502bSHong Zhang 3409b54502bSHong Zhang Level: beginner 3419b54502bSHong Zhang 3429b54502bSHong Zhang Concepts: incomplete factorization 3439b54502bSHong Zhang 344d9ba8547SBarry Smith Notes: Only implemented for some matrix formats. (for parallel see PCHYPRE for hypre's ILU) 3459b54502bSHong Zhang 3469b54502bSHong Zhang For BAIJ matrices this implements a point block ILU 3479b54502bSHong Zhang 348f0b9ad6cSBarry Smith The "symmetric" application of this preconditioner is not actually symmetric since L is not transpose(U) 349f0b9ad6cSBarry Smith even when the matrix is not symmetric since the U stores the diagonals of the factorization. 350f0b9ad6cSBarry Smith 351cd0a26f6SPaul Mullowney If you are using MATSEQAIJCUSPARSE matrices (or MATMPIAIJCUSPARESE matrices with block Jacobi), factorization 352cd0a26f6SPaul Mullowney is never done on the GPU). 353ac793be5SBarry Smith 354c582cd25SBarry Smith References: 35596a0c994SBarry Smith + 1. - T. Dupont, R. Kendall, and H. Rachford. An approximate factorization procedure for solving 35696a0c994SBarry Smith self adjoint elliptic difference equations. SIAM J. Numer. Anal., 5, 1968. 35796a0c994SBarry Smith . 2. - T.A. Oliphant. An implicit numerical method for solving two dimensional timedependent diffusion problems. Quart. Appl. Math., 19, 1961. 35896a0c994SBarry Smith - 3. - TONY F. CHAN AND HENK A. VAN DER VORST, APPROXIMATE AND INCOMPLETE FACTORIZATIONS, 35996a0c994SBarry Smith Chapter in Parallel Numerical 360c582cd25SBarry Smith Algorithms, edited by D. Keyes, A. Semah, V. Venkatakrishnan, ICASE/LaRC Interdisciplinary Series in 36196a0c994SBarry Smith Science and Engineering, Kluwer. 362c582cd25SBarry Smith 363c582cd25SBarry Smith 3649b54502bSHong Zhang .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType, 365145b9266SHong Zhang PCFactorSetZeroPivot(), PCFactorSetShiftSetType(), PCFactorSetAmount(), 366b7c853c4SBarry Smith PCFactorSetDropTolerance(),PCFactorSetFill(), PCFactorSetMatOrderingType(), PCFactorSetReuseOrdering(), 36792e9c092SBarry Smith PCFactorSetLevels(), PCFactorSetUseInPlace(), PCFactorSetAllowDiagonalFill(), PCFactorSetPivotInBlocks(), 36892e9c092SBarry Smith PCFactorGetAllowDiagonalFill(), PCFactorGetUseInPlace() 3699b54502bSHong Zhang 3709b54502bSHong Zhang M*/ 3719b54502bSHong Zhang 3729b54502bSHong Zhang #undef __FUNCT__ 3739b54502bSHong Zhang #define __FUNCT__ "PCCreate_ILU" 3748cc058d9SJed Brown PETSC_EXTERN PetscErrorCode PCCreate_ILU(PC pc) 3759b54502bSHong Zhang { 3769b54502bSHong Zhang PetscErrorCode ierr; 3779b54502bSHong Zhang PC_ILU *ilu; 3789b54502bSHong Zhang 3799b54502bSHong Zhang PetscFunctionBegin; 380b00a9115SJed Brown ierr = PetscNewLog(pc,&ilu);CHKERRQ(ierr); 3819b54502bSHong Zhang 382075768bcSBarry Smith ((PC_Factor*)ilu)->fact = 0; 383075768bcSBarry Smith ierr = MatFactorInfoInitialize(&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 384879e8a4dSBarry Smith ((PC_Factor*)ilu)->factortype = MAT_FACTOR_ILU; 38575567043SBarry Smith ((PC_Factor*)ilu)->info.levels = 0.; 386075768bcSBarry Smith ((PC_Factor*)ilu)->info.fill = 1.0; 3879b54502bSHong Zhang ilu->col = 0; 3889b54502bSHong Zhang ilu->row = 0; 3899b54502bSHong Zhang ilu->inplace = PETSC_FALSE; 39019fd82e9SBarry Smith ierr = PetscStrallocpy(MATORDERINGNATURAL,(char**)&((PC_Factor*)ilu)->ordering);CHKERRQ(ierr); 3919b54502bSHong Zhang ilu->reuseordering = PETSC_FALSE; 392075768bcSBarry Smith ((PC_Factor*)ilu)->info.dt = PETSC_DEFAULT; 393075768bcSBarry Smith ((PC_Factor*)ilu)->info.dtcount = PETSC_DEFAULT; 394075768bcSBarry Smith ((PC_Factor*)ilu)->info.dtcol = PETSC_DEFAULT; 395c378d667SBarry Smith ((PC_Factor*)ilu)->info.shifttype = (PetscReal)MAT_SHIFT_NONE; 3960ed735ceSHong Zhang ((PC_Factor*)ilu)->info.shiftamount = 100.0*PETSC_MACHINE_EPSILON; 3970ed735ceSHong Zhang ((PC_Factor*)ilu)->info.zeropivot = 100.0*PETSC_MACHINE_EPSILON; 398075768bcSBarry Smith ((PC_Factor*)ilu)->info.pivotinblocks = 1.0; 3999b54502bSHong Zhang ilu->reusefill = PETSC_FALSE; 4000ed735ceSHong Zhang ((PC_Factor*)ilu)->info.diagonal_fill = 0.0; 4019b54502bSHong Zhang pc->data = (void*)ilu; 4029b54502bSHong Zhang 403574deadeSBarry Smith pc->ops->reset = PCReset_ILU; 4049b54502bSHong Zhang pc->ops->destroy = PCDestroy_ILU; 4059b54502bSHong Zhang pc->ops->apply = PCApply_ILU; 4069b54502bSHong Zhang pc->ops->applytranspose = PCApplyTranspose_ILU; 4079b54502bSHong Zhang pc->ops->setup = PCSetUp_ILU; 4089b54502bSHong Zhang pc->ops->setfromoptions = PCSetFromOptions_ILU; 40985317021SBarry Smith pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor; 4109b54502bSHong Zhang pc->ops->view = PCView_ILU; 411f0b9ad6cSBarry Smith pc->ops->applysymmetricleft = PCApplySymmetricLeft_ILU; 412f0b9ad6cSBarry Smith pc->ops->applysymmetricright = PCApplySymmetricRight_ILU; 4139b54502bSHong Zhang pc->ops->applyrichardson = 0; 4149b54502bSHong Zhang 415bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);CHKERRQ(ierr); 416c7f610a1SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetZeroPivot_C",PCFactorGetZeroPivot_Factor);CHKERRQ(ierr); 417bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);CHKERRQ(ierr); 418c7f610a1SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftType_C",PCFactorGetShiftType_Factor);CHKERRQ(ierr); 419bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);CHKERRQ(ierr); 420c7f610a1SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetShiftAmount_C",PCFactorGetShiftAmount_Factor);CHKERRQ(ierr); 421bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverPackage_C",PCFactorGetMatSolverPackage_Factor);CHKERRQ(ierr); 422bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverPackage_C",PCFactorSetMatSolverPackage_Factor);CHKERRQ(ierr); 423bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverPackage_C",PCFactorSetUpMatSolverPackage_Factor);CHKERRQ(ierr); 424bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetDropTolerance_C",PCFactorSetDropTolerance_ILU);CHKERRQ(ierr); 425bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);CHKERRQ(ierr); 426bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);CHKERRQ(ierr); 427bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_ILU);CHKERRQ(ierr); 428bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_ILU);CHKERRQ(ierr); 429bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);CHKERRQ(ierr); 4302591b318SToby Isaac ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);CHKERRQ(ierr); 431bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_ILU);CHKERRQ(ierr); 4328e37d05fSBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_ILU);CHKERRQ(ierr); 433bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);CHKERRQ(ierr); 43492e9c092SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);CHKERRQ(ierr); 435bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);CHKERRQ(ierr); 436bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorReorderForNonzeroDiagonal_C",PCFactorReorderForNonzeroDiagonal_ILU);CHKERRQ(ierr); 4379b54502bSHong Zhang PetscFunctionReturn(0); 4389b54502bSHong Zhang } 439