xref: /petsc/src/ksp/pc/impls/factor/ilu/ilu.c (revision ac793be5d4bbb0eafaea99146873fee4f01ef996)
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 EXTERN_C_BEGIN
9afaefe49SHong Zhang #undef __FUNCT__
102401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill_ILU"
117087cfbeSBarry Smith PetscErrorCode  PCFactorSetReuseFill_ILU(PC pc,PetscBool  flag)
122401956bSBarry Smith {
13075768bcSBarry Smith   PC_ILU *lu = (PC_ILU*)pc->data;
142401956bSBarry Smith 
152401956bSBarry Smith   PetscFunctionBegin;
162401956bSBarry Smith   lu->reusefill = flag;
172401956bSBarry Smith   PetscFunctionReturn(0);
182401956bSBarry Smith }
192401956bSBarry Smith EXTERN_C_END
202401956bSBarry Smith 
212401956bSBarry Smith EXTERN_C_BEGIN
222401956bSBarry Smith #undef __FUNCT__
232401956bSBarry Smith #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal_ILU"
247087cfbeSBarry Smith PetscErrorCode  PCFactorReorderForNonzeroDiagonal_ILU(PC pc,PetscReal z)
259b54502bSHong Zhang {
269b54502bSHong Zhang   PC_ILU *ilu = (PC_ILU*)pc->data;
279b54502bSHong Zhang 
289b54502bSHong Zhang   PetscFunctionBegin;
299b54502bSHong Zhang   ilu->nonzerosalongdiagonal = PETSC_TRUE;
309b54502bSHong Zhang   if (z == PETSC_DECIDE) {
319b54502bSHong Zhang     ilu->nonzerosalongdiagonaltol = 1.e-10;
329b54502bSHong Zhang   } else {
339b54502bSHong Zhang     ilu->nonzerosalongdiagonaltol = z;
349b54502bSHong Zhang   }
359b54502bSHong Zhang   PetscFunctionReturn(0);
369b54502bSHong Zhang }
379b54502bSHong Zhang EXTERN_C_END
389b54502bSHong Zhang 
399b54502bSHong Zhang #undef __FUNCT__
40574deadeSBarry Smith #define __FUNCT__ "PCReset_ILU"
41574deadeSBarry Smith PetscErrorCode PCReset_ILU(PC pc)
429b54502bSHong Zhang {
439b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
449b54502bSHong Zhang   PetscErrorCode ierr;
459b54502bSHong Zhang 
469b54502bSHong Zhang   PetscFunctionBegin;
47298cc208SBarry Smith   if (!ilu->inplace) {ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);}
486bf464f9SBarry Smith   if (ilu->row && ilu->col && ilu->row != ilu->col) {ierr = ISDestroy(&ilu->row);CHKERRQ(ierr);}
49fcfd50ebSBarry Smith   ierr = ISDestroy(&ilu->col);CHKERRQ(ierr);
509b54502bSHong Zhang   PetscFunctionReturn(0);
519b54502bSHong Zhang }
529b54502bSHong Zhang 
539b54502bSHong Zhang EXTERN_C_BEGIN
549b54502bSHong Zhang #undef __FUNCT__
55b7c853c4SBarry Smith #define __FUNCT__ "PCFactorSetDropTolerance_ILU"
567087cfbeSBarry Smith PetscErrorCode  PCFactorSetDropTolerance_ILU(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount)
579b54502bSHong Zhang {
58075768bcSBarry Smith   PC_ILU         *ilu = (PC_ILU*)pc->data;
599b54502bSHong Zhang 
609b54502bSHong Zhang   PetscFunctionBegin;
614c9036c7SBarry Smith   if (pc->setupcalled && (((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) {
62e7e72b3dSBarry Smith     SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_SUP,"Cannot change drop tolerance after using PC");
639b54502bSHong Zhang   }
64075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dt      = dt;
65075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcol   = dtcol;
66075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcount = dtcount;
674c9036c7SBarry Smith   ((PC_Factor*)ilu)->info.usedt   = 1.0;
689b54502bSHong Zhang   PetscFunctionReturn(0);
699b54502bSHong Zhang }
709b54502bSHong Zhang EXTERN_C_END
719b54502bSHong Zhang 
729b54502bSHong Zhang EXTERN_C_BEGIN
739b54502bSHong Zhang #undef __FUNCT__
742401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering_ILU"
757087cfbeSBarry Smith PetscErrorCode  PCFactorSetReuseOrdering_ILU(PC pc,PetscBool  flag)
769b54502bSHong Zhang {
77075768bcSBarry Smith   PC_ILU *ilu = (PC_ILU*)pc->data;
789b54502bSHong Zhang 
799b54502bSHong Zhang   PetscFunctionBegin;
809b54502bSHong Zhang   ilu->reuseordering = flag;
819b54502bSHong Zhang   PetscFunctionReturn(0);
829b54502bSHong Zhang }
839b54502bSHong Zhang EXTERN_C_END
849b54502bSHong Zhang 
859b54502bSHong Zhang EXTERN_C_BEGIN
869b54502bSHong Zhang #undef __FUNCT__
872401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseInPlace_ILU"
887087cfbeSBarry Smith PetscErrorCode  PCFactorSetUseInPlace_ILU(PC pc)
899b54502bSHong Zhang {
90075768bcSBarry Smith   PC_ILU *dir = (PC_ILU*)pc->data;
919b54502bSHong Zhang 
929b54502bSHong Zhang   PetscFunctionBegin;
939b54502bSHong Zhang   dir->inplace = PETSC_TRUE;
949b54502bSHong Zhang   PetscFunctionReturn(0);
959b54502bSHong Zhang }
969b54502bSHong Zhang EXTERN_C_END
979b54502bSHong Zhang 
989b54502bSHong Zhang #undef __FUNCT__
999b54502bSHong Zhang #define __FUNCT__ "PCSetFromOptions_ILU"
1009b54502bSHong Zhang static PetscErrorCode PCSetFromOptions_ILU(PC pc)
1019b54502bSHong Zhang {
1029b54502bSHong Zhang   PetscErrorCode ierr;
10378fc6b22SHong Zhang   PetscInt       itmp;
104ace3abfcSBarry Smith   PetscBool      flg;
1055a586d82SBarry Smith   /* PetscReal      dt[3]; */
1069b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
1079b54502bSHong Zhang   PetscReal      tol;
1089b54502bSHong Zhang 
1099b54502bSHong Zhang   PetscFunctionBegin;
1109b54502bSHong Zhang   ierr = PetscOptionsHead("ILU Options");CHKERRQ(ierr);
1118ff23777SHong Zhang     ierr = PCSetFromOptions_Factor(pc);CHKERRQ(ierr);
1128ff23777SHong Zhang 
113075768bcSBarry Smith     ierr = PetscOptionsInt("-pc_factor_levels","levels of fill","PCFactorSetLevels",(PetscInt)((PC_Factor*)ilu)->info.levels,&itmp,&flg);CHKERRQ(ierr);
114075768bcSBarry Smith     if (flg) ((PC_Factor*)ilu)->info.levels = itmp;
11590d69ab7SBarry Smith     flg  = PETSC_FALSE;
116acfcf0e5SJed Brown     ierr = PetscOptionsBool("-pc_factor_diagonal_fill","Allow fill into empty diagonal entry","PCFactorSetAllowDiagonalFill",flg,&flg,PETSC_NULL);CHKERRQ(ierr);
117075768bcSBarry Smith     ((PC_Factor*)ilu)->info.diagonal_fill = (double) flg;
1185a586d82SBarry Smith     /*
119075768bcSBarry Smith     dt[0] = ((PC_Factor*)ilu)->info.dt;
120075768bcSBarry Smith     dt[1] = ((PC_Factor*)ilu)->info.dtcol;
121075768bcSBarry Smith     dt[2] = ((PC_Factor*)ilu)->info.dtcount;
1225a586d82SBarry Smith 
12378fc6b22SHong Zhang     PetscInt       dtmax = 3;
12478fc6b22SHong Zhang     ierr = PetscOptionsRealArray("-pc_factor_drop_tolerance,","<dt,dtcol,maxrowcount>","PCFactorSetDropTolerance",dt,&dtmax,&flg);CHKERRQ(ierr);
1259b54502bSHong Zhang     if (flg) {
126b7c853c4SBarry Smith       ierr = PCFactorSetDropTolerance(pc,dt[0],dt[1],(PetscInt)dt[2]);CHKERRQ(ierr);
1279b54502bSHong Zhang     }
12878fc6b22SHong Zhang     */
1292401956bSBarry Smith     ierr = PetscOptionsName("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",&flg);CHKERRQ(ierr);
1309b54502bSHong Zhang     if (flg) {
1319b54502bSHong Zhang       tol = PETSC_DECIDE;
1322401956bSBarry Smith       ierr = PetscOptionsReal("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",ilu->nonzerosalongdiagonaltol,&tol,0);CHKERRQ(ierr);
1332401956bSBarry Smith       ierr = PCFactorReorderForNonzeroDiagonal(pc,tol);CHKERRQ(ierr);
1349b54502bSHong Zhang     }
1359b54502bSHong Zhang 
1369b54502bSHong Zhang   ierr = PetscOptionsTail();CHKERRQ(ierr);
1379b54502bSHong Zhang   PetscFunctionReturn(0);
1389b54502bSHong Zhang }
1399b54502bSHong Zhang 
1409b54502bSHong Zhang #undef __FUNCT__
1419b54502bSHong Zhang #define __FUNCT__ "PCView_ILU"
1429b54502bSHong Zhang static PetscErrorCode PCView_ILU(PC pc,PetscViewer viewer)
1439b54502bSHong Zhang {
1449b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
1459b54502bSHong Zhang   PetscErrorCode ierr;
146ace3abfcSBarry Smith   PetscBool      iascii;
1479b54502bSHong Zhang 
1489b54502bSHong Zhang   PetscFunctionBegin;
1492692d6eeSBarry Smith   ierr = PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
1509b54502bSHong Zhang   if (iascii) {
151914a5d51SHong Zhang     if (ilu->inplace) {
152914a5d51SHong Zhang       ierr = PetscViewerASCIIPrintf(viewer,"  ILU: in-place factorization\n");CHKERRQ(ierr);
1539b54502bSHong Zhang     } else {
154914a5d51SHong Zhang       ierr = PetscViewerASCIIPrintf(viewer,"  ILU: out-of-place factorization\n");CHKERRQ(ierr);
155145b9266SHong Zhang     }
156145b9266SHong Zhang 
157914a5d51SHong Zhang     if (ilu->reusefill)     {ierr = PetscViewerASCIIPrintf(viewer,"  ILU: Reusing fill from past factorization\n");CHKERRQ(ierr);}
158914a5d51SHong Zhang     if (ilu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer,"  ILU: Reusing reordering from past factorization\n");CHKERRQ(ierr);}
1599b54502bSHong Zhang   }
160914a5d51SHong Zhang   ierr = PCView_Factor(pc,viewer);CHKERRQ(ierr);
1619b54502bSHong Zhang   PetscFunctionReturn(0);
1629b54502bSHong Zhang }
1639b54502bSHong Zhang 
1649b54502bSHong Zhang #undef __FUNCT__
1659b54502bSHong Zhang #define __FUNCT__ "PCSetUp_ILU"
1669b54502bSHong Zhang static PetscErrorCode PCSetUp_ILU(PC pc)
1679b54502bSHong Zhang {
1689b54502bSHong Zhang   PetscErrorCode ierr;
1699b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
170f3a39becSBarry Smith   MatInfo        info;
171ace3abfcSBarry Smith   PetscBool      flg;
1729b54502bSHong Zhang 
1739b54502bSHong Zhang   PetscFunctionBegin;
17492927226SBarry Smith   /* ugly hack to change default, since it is not support by some matrix types */
17592927226SBarry Smith   if (((PC_Factor*)ilu)->info.shifttype == (PetscReal)MAT_SHIFT_NONZERO) {
17692927226SBarry Smith     ierr = PetscTypeCompare((PetscObject)pc->pmat,MATSEQAIJ,&flg);
17792927226SBarry Smith     if (!flg) {
17892927226SBarry Smith       ierr = PetscTypeCompare((PetscObject)pc->pmat,MATMPIAIJ,&flg);
17992927226SBarry Smith       if (!flg) {
18092927226SBarry Smith         ((PC_Factor*)ilu)->info.shifttype = (PetscReal)MAT_SHIFT_INBLOCKS;
18192927226SBarry Smith         PetscInfo(pc,"Changing shift type from NONZERO to INBLOCKS because block matrices do not support NONZERO");
18292927226SBarry Smith       }
18392927226SBarry Smith     }
18492927226SBarry Smith   }
18592927226SBarry Smith 
1869b54502bSHong Zhang   if (ilu->inplace) {
187075768bcSBarry Smith     CHKMEMQ;
1889b54502bSHong Zhang     if (!pc->setupcalled) {
1899b54502bSHong Zhang 
1909b54502bSHong Zhang       /* In-place factorization only makes sense with the natural ordering,
1919b54502bSHong Zhang          so we only need to get the ordering once, even if nonzero structure changes */
192075768bcSBarry Smith       ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
193efee365bSSatish Balay       if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
194efee365bSSatish Balay       if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
1959b54502bSHong Zhang     }
1969b54502bSHong Zhang 
1979b54502bSHong Zhang     /* In place ILU only makes sense with fill factor of 1.0 because
1989b54502bSHong Zhang        cannot have levels of fill */
199075768bcSBarry Smith     ((PC_Factor*)ilu)->info.fill          = 1.0;
20075567043SBarry Smith     ((PC_Factor*)ilu)->info.diagonal_fill = 0.0;
201fe97e370SBarry Smith     ierr = MatILUFactor(pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);CHKMEMQ;
202075768bcSBarry Smith     ((PC_Factor*)ilu)->fact = pc->pmat;
2039b54502bSHong Zhang   } else {
2049b54502bSHong Zhang     if (!pc->setupcalled) {
2059b54502bSHong Zhang       /* first time in so compute reordering and symbolic factorization */
206075768bcSBarry Smith       ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
207efee365bSSatish Balay       if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
208efee365bSSatish Balay       if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
2099b54502bSHong Zhang       /*  Remove zeros along diagonal?     */
2109b54502bSHong Zhang       if (ilu->nonzerosalongdiagonal) {
2119b54502bSHong Zhang         ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr);
2129b54502bSHong Zhang       }
21386daa42cSBarry Smith       if (!((PC_Factor*)ilu)->fact){
21411bfe483SHong Zhang         ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
215a1f19f5aSHong Zhang       }
216075768bcSBarry Smith       ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
217075768bcSBarry Smith       ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
218f3a39becSBarry Smith       ilu->actualfill = info.fill_ratio_needed;
219075768bcSBarry Smith       ierr = PetscLogObjectParent(pc,((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
2209b54502bSHong Zhang     } else if (pc->flag != SAME_NONZERO_PATTERN) {
2219b54502bSHong Zhang       if (!ilu->reuseordering) {
2229b54502bSHong Zhang         /* compute a new ordering for the ILU */
223fcfd50ebSBarry Smith         ierr = ISDestroy(&ilu->row);CHKERRQ(ierr);
224fcfd50ebSBarry Smith         ierr = ISDestroy(&ilu->col);CHKERRQ(ierr);
225075768bcSBarry Smith         ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
226efee365bSSatish Balay         if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
227efee365bSSatish Balay         if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
2289b54502bSHong Zhang         /*  Remove zeros along diagonal?     */
2299b54502bSHong Zhang         if (ilu->nonzerosalongdiagonal) {
2309b54502bSHong Zhang           ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr);
2319b54502bSHong Zhang         }
2329b54502bSHong Zhang       }
2336bf464f9SBarry Smith       ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
2342692d6eeSBarry Smith       ierr = MatGetFactor(pc->pmat,MATSOLVERPETSC,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
235075768bcSBarry Smith       ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
236075768bcSBarry Smith       ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
237f3a39becSBarry Smith       ilu->actualfill = info.fill_ratio_needed;
238075768bcSBarry Smith       ierr = PetscLogObjectParent(pc,((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
2399b54502bSHong Zhang     }
240075768bcSBarry Smith     CHKMEMQ;
241075768bcSBarry Smith     ierr = MatLUFactorNumeric(((PC_Factor*)ilu)->fact,pc->pmat,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
242075768bcSBarry Smith     CHKMEMQ;
2439b54502bSHong Zhang   }
2449b54502bSHong Zhang   PetscFunctionReturn(0);
2459b54502bSHong Zhang }
2469b54502bSHong Zhang 
2479b54502bSHong Zhang #undef __FUNCT__
2489b54502bSHong Zhang #define __FUNCT__ "PCDestroy_ILU"
2499b54502bSHong Zhang static PetscErrorCode PCDestroy_ILU(PC pc)
2509b54502bSHong Zhang {
2519b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2529b54502bSHong Zhang   PetscErrorCode ierr;
2539b54502bSHong Zhang 
2549b54502bSHong Zhang   PetscFunctionBegin;
255574deadeSBarry Smith   ierr = PCReset_ILU(pc);CHKERRQ(ierr);
256503cfb0cSBarry Smith   ierr = PetscFree(((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr);
257503cfb0cSBarry Smith   ierr = PetscFree(((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);
258c31cb41cSBarry Smith   ierr = PetscFree(pc->data);CHKERRQ(ierr);
2599b54502bSHong Zhang   PetscFunctionReturn(0);
2609b54502bSHong Zhang }
2619b54502bSHong Zhang 
2629b54502bSHong Zhang #undef __FUNCT__
2639b54502bSHong Zhang #define __FUNCT__ "PCApply_ILU"
2649b54502bSHong Zhang static PetscErrorCode PCApply_ILU(PC pc,Vec x,Vec y)
2659b54502bSHong Zhang {
2669b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2679b54502bSHong Zhang   PetscErrorCode ierr;
2689b54502bSHong Zhang 
2699b54502bSHong Zhang   PetscFunctionBegin;
270075768bcSBarry Smith   ierr = MatSolve(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr);
2719b54502bSHong Zhang   PetscFunctionReturn(0);
2729b54502bSHong Zhang }
2739b54502bSHong Zhang 
2749b54502bSHong Zhang #undef __FUNCT__
2759b54502bSHong Zhang #define __FUNCT__ "PCApplyTranspose_ILU"
2769b54502bSHong Zhang static PetscErrorCode PCApplyTranspose_ILU(PC pc,Vec x,Vec y)
2779b54502bSHong Zhang {
2789b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2799b54502bSHong Zhang   PetscErrorCode ierr;
2809b54502bSHong Zhang 
2819b54502bSHong Zhang   PetscFunctionBegin;
282075768bcSBarry Smith   ierr = MatSolveTranspose(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr);
2839b54502bSHong Zhang   PetscFunctionReturn(0);
2849b54502bSHong Zhang }
2859b54502bSHong Zhang 
286f0b9ad6cSBarry Smith #undef __FUNCT__
287f0b9ad6cSBarry Smith #define __FUNCT__ "PCApplySymmetricLeft_ILU"
288f0b9ad6cSBarry Smith static PetscErrorCode PCApplySymmetricLeft_ILU(PC pc,Vec x,Vec y)
289f0b9ad6cSBarry Smith {
290f0b9ad6cSBarry Smith   PetscErrorCode ierr;
291f0b9ad6cSBarry Smith   PC_ILU         *icc = (PC_ILU*)pc->data;
292f0b9ad6cSBarry Smith 
293f0b9ad6cSBarry Smith   PetscFunctionBegin;
294f0b9ad6cSBarry Smith   ierr = MatForwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr);
295f0b9ad6cSBarry Smith   PetscFunctionReturn(0);
296f0b9ad6cSBarry Smith }
297f0b9ad6cSBarry Smith 
298f0b9ad6cSBarry Smith #undef __FUNCT__
299f0b9ad6cSBarry Smith #define __FUNCT__ "PCApplySymmetricRight_ILU"
300f0b9ad6cSBarry Smith static PetscErrorCode PCApplySymmetricRight_ILU(PC pc,Vec x,Vec y)
301f0b9ad6cSBarry Smith {
302f0b9ad6cSBarry Smith   PetscErrorCode ierr;
303f0b9ad6cSBarry Smith   PC_ILU         *icc = (PC_ILU*)pc->data;
304f0b9ad6cSBarry Smith 
305f0b9ad6cSBarry Smith   PetscFunctionBegin;
306f0b9ad6cSBarry Smith   ierr = MatBackwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr);
307f0b9ad6cSBarry Smith   PetscFunctionReturn(0);
308f0b9ad6cSBarry Smith }
309f0b9ad6cSBarry Smith 
3109b54502bSHong Zhang /*MC
3119b54502bSHong Zhang      PCILU - Incomplete factorization preconditioners.
3129b54502bSHong Zhang 
3139b54502bSHong Zhang    Options Database Keys:
3142401956bSBarry Smith +  -pc_factor_levels <k> - number of levels of fill for ILU(k)
3152401956bSBarry Smith .  -pc_factor_in_place - only for ILU(0) with natural ordering, reuses the space of the matrix for
3169b54502bSHong Zhang                       its factorization (overwrites original matrix)
3172401956bSBarry Smith .  -pc_factor_diagonal_fill - fill in a zero diagonal even if levels of fill indicate it wouldn't be fill
3182401956bSBarry Smith .  -pc_factor_reuse_ordering - reuse ordering of factorized matrix from previous factorization
31955ba2a51SBarry Smith .  -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1
3202401956bSBarry Smith .  -pc_factor_nonzeros_along_diagonal - reorder the matrix before factorization to remove zeros from the diagonal,
3219b54502bSHong Zhang                                    this decreases the chance of getting a zero pivot
3222401956bSBarry Smith .  -pc_factor_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix
3232401956bSBarry Smith .  -pc_factor_pivot_in_blocks - for block ILU(k) factorization, i.e. with BAIJ matrices with block size larger
3249b54502bSHong Zhang                              than 1 the diagonal blocks are factored with partial pivoting (this increases the
3259b54502bSHong Zhang                              stability of the ILU factorization
3269b54502bSHong Zhang 
3279b54502bSHong Zhang    Level: beginner
3289b54502bSHong Zhang 
3299b54502bSHong Zhang   Concepts: incomplete factorization
3309b54502bSHong Zhang 
331d9ba8547SBarry Smith    Notes: Only implemented for some matrix formats. (for parallel see PCHYPRE for hypre's ILU)
3329b54502bSHong Zhang 
3339b54502bSHong Zhang           For BAIJ matrices this implements a point block ILU
3349b54502bSHong Zhang 
335f0b9ad6cSBarry Smith           The "symmetric" application of this preconditioner is not actually symmetric since L is not transpose(U)
336f0b9ad6cSBarry Smith           even when the matrix is not symmetric since the U stores the diagonals of the factorization.
337f0b9ad6cSBarry Smith 
338*ac793be5SBarry Smith           If you are using MATSEQAIJCUSP matrices (or MATMPIAIJCUSP matrices with block Jacobi) you must have ./configured
339*ac793be5SBarry Smith           PETSc with also --download-txpetscgpu to have the triangular solves performed on the GPU (factorization is never
340*ac793be5SBarry Smith           done on the GPU).
341*ac793be5SBarry Smith 
342c582cd25SBarry Smith    References:
343c582cd25SBarry Smith    T. Dupont, R. Kendall, and H. Rachford. An approximate factorization procedure for solving
344c582cd25SBarry Smith    self-adjoint elliptic difference equations. SIAM J. Numer. Anal., 5:559--573, 1968.
345c582cd25SBarry Smith 
346c582cd25SBarry Smith    T.A. Oliphant. An implicit numerical method for solving two-dimensional time-dependent dif-
347c582cd25SBarry Smith    fusion problems. Quart. Appl. Math., 19:221--229, 1961.
348c582cd25SBarry Smith 
349c582cd25SBarry Smith    Review article: APPROXIMATE AND INCOMPLETE FACTORIZATIONS, TONY F. CHAN AND HENK A. VAN DER VORST
350c582cd25SBarry Smith       http://igitur-archive.library.uu.nl/math/2001-0621-115821/proc.pdf chapter in Parallel Numerical
351c582cd25SBarry Smith       Algorithms, edited by D. Keyes, A. Semah, V. Venkatakrishnan, ICASE/LaRC Interdisciplinary Series in
352c582cd25SBarry Smith       Science and Engineering, Kluwer, pp. 167--202.
353c582cd25SBarry Smith 
354c582cd25SBarry Smith 
3559b54502bSHong Zhang .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType,
356145b9266SHong Zhang            PCFactorSetZeroPivot(), PCFactorSetShiftSetType(), PCFactorSetAmount(),
357b7c853c4SBarry Smith            PCFactorSetDropTolerance(),PCFactorSetFill(), PCFactorSetMatOrderingType(), PCFactorSetReuseOrdering(),
3588ff23777SHong Zhang            PCFactorSetLevels(), PCFactorSetUseInPlace(), PCFactorSetAllowDiagonalFill(), PCFactorSetPivotInBlocks()
3599b54502bSHong Zhang 
3609b54502bSHong Zhang M*/
3619b54502bSHong Zhang 
3629b54502bSHong Zhang EXTERN_C_BEGIN
3639b54502bSHong Zhang #undef __FUNCT__
3649b54502bSHong Zhang #define __FUNCT__ "PCCreate_ILU"
3657087cfbeSBarry Smith PetscErrorCode  PCCreate_ILU(PC pc)
3669b54502bSHong Zhang {
3679b54502bSHong Zhang   PetscErrorCode ierr;
3689b54502bSHong Zhang   PC_ILU         *ilu;
3699b54502bSHong Zhang 
3709b54502bSHong Zhang   PetscFunctionBegin;
37138f2d2fdSLisandro Dalcin   ierr = PetscNewLog(pc,PC_ILU,&ilu);CHKERRQ(ierr);
3729b54502bSHong Zhang 
373075768bcSBarry Smith   ((PC_Factor*)ilu)->fact                    = 0;
374075768bcSBarry Smith   ierr = MatFactorInfoInitialize(&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
375879e8a4dSBarry Smith   ((PC_Factor*)ilu)->factortype              = MAT_FACTOR_ILU;
37675567043SBarry Smith   ((PC_Factor*)ilu)->info.levels             = 0.;
377075768bcSBarry Smith   ((PC_Factor*)ilu)->info.fill               = 1.0;
3789b54502bSHong Zhang   ilu->col                                   = 0;
3799b54502bSHong Zhang   ilu->row                                   = 0;
3809b54502bSHong Zhang   ilu->inplace                               = PETSC_FALSE;
3812692d6eeSBarry Smith   ierr = PetscStrallocpy(MATSOLVERPETSC,&((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr);
3822692d6eeSBarry Smith   ierr = PetscStrallocpy(MATORDERINGNATURAL,&((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);
3839b54502bSHong Zhang   ilu->reuseordering                         = PETSC_FALSE;
384075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dt                 = PETSC_DEFAULT;
385075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcount            = PETSC_DEFAULT;
386075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcol              = PETSC_DEFAULT;
387357fed3dSBarry Smith   ((PC_Factor*)ilu)->info.shifttype          = (PetscReal)MAT_SHIFT_NONZERO;
3880ed735ceSHong Zhang   ((PC_Factor*)ilu)->info.shiftamount        = 100.0*PETSC_MACHINE_EPSILON;
3890ed735ceSHong Zhang   ((PC_Factor*)ilu)->info.zeropivot          = 100.0*PETSC_MACHINE_EPSILON;
390075768bcSBarry Smith   ((PC_Factor*)ilu)->info.pivotinblocks      = 1.0;
3919b54502bSHong Zhang   ilu->reusefill                             = PETSC_FALSE;
3920ed735ceSHong Zhang   ((PC_Factor*)ilu)->info.diagonal_fill      = 0.0;
3939b54502bSHong Zhang   pc->data                     = (void*)ilu;
3949b54502bSHong Zhang 
395574deadeSBarry Smith   pc->ops->reset               = PCReset_ILU;
3969b54502bSHong Zhang   pc->ops->destroy             = PCDestroy_ILU;
3979b54502bSHong Zhang   pc->ops->apply               = PCApply_ILU;
3989b54502bSHong Zhang   pc->ops->applytranspose      = PCApplyTranspose_ILU;
3999b54502bSHong Zhang   pc->ops->setup               = PCSetUp_ILU;
4009b54502bSHong Zhang   pc->ops->setfromoptions      = PCSetFromOptions_ILU;
40185317021SBarry Smith   pc->ops->getfactoredmatrix   = PCFactorGetMatrix_Factor;
4029b54502bSHong Zhang   pc->ops->view                = PCView_ILU;
403f0b9ad6cSBarry Smith   pc->ops->applysymmetricleft  = PCApplySymmetricLeft_ILU;
404f0b9ad6cSBarry Smith   pc->ops->applysymmetricright = PCApplySymmetricRight_ILU;
4059b54502bSHong Zhang   pc->ops->applyrichardson     = 0;
4069b54502bSHong Zhang 
40785317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Factor",
40885317021SBarry Smith                     PCFactorSetZeroPivot_Factor);CHKERRQ(ierr);
409d90ac83dSHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftType_C","PCFactorSetShiftType_Factor",
410d90ac83dSHong Zhang                     PCFactorSetShiftType_Factor);CHKERRQ(ierr);
411d90ac83dSHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftAmount_C","PCFactorSetShiftAmount_Factor",
412d90ac83dSHong Zhang                     PCFactorSetShiftAmount_Factor);CHKERRQ(ierr);
4137112b564SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorGetMatSolverPackage_C","PCFactorGetMatSolverPackage_Factor",
4147112b564SBarry Smith                     PCFactorGetMatSolverPackage_Factor);CHKERRQ(ierr);
41511bfe483SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatSolverPackage_C","PCFactorSetMatSolverPackage_Factor",
41611bfe483SHong Zhang                     PCFactorSetMatSolverPackage_Factor);CHKERRQ(ierr);
417f8260c8fSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUpMatSolverPackage_C","PCFactorSetUpMatSolverPackage_Factor",
418f8260c8fSBarry Smith                     PCFactorSetUpMatSolverPackage_Factor);CHKERRQ(ierr);
419b7c853c4SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetDropTolerance_C","PCFactorSetDropTolerance_ILU",
420b7c853c4SBarry Smith                     PCFactorSetDropTolerance_ILU);CHKERRQ(ierr);
42185317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Factor",
42285317021SBarry Smith                     PCFactorSetFill_Factor);CHKERRQ(ierr);
42385317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Factor",
42485317021SBarry Smith                     PCFactorSetMatOrderingType_Factor);CHKERRQ(ierr);
4252401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_ILU",
4262401956bSBarry Smith                     PCFactorSetReuseOrdering_ILU);CHKERRQ(ierr);
4272401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_ILU",
4282401956bSBarry Smith                     PCFactorSetReuseFill_ILU);CHKERRQ(ierr);
42985317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetLevels_C","PCFactorSetLevels_Factor",
43085317021SBarry Smith                     PCFactorSetLevels_Factor);CHKERRQ(ierr);
4312401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_ILU",
4322401956bSBarry Smith                     PCFactorSetUseInPlace_ILU);CHKERRQ(ierr);
43385317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C","PCFactorSetAllowDiagonalFill_Factor",
43485317021SBarry Smith                     PCFactorSetAllowDiagonalFill_Factor);CHKERRQ(ierr);
43585317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetPivotInBlocks_C","PCFactorSetPivotInBlocks_Factor",
43685317021SBarry Smith                     PCFactorSetPivotInBlocks_Factor);CHKERRQ(ierr);
4372401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorReorderForNonzeroDiagonal_C","PCFactorReorderForNonzeroDiagonal_ILU",
4382401956bSBarry Smith                     PCFactorReorderForNonzeroDiagonal_ILU);CHKERRQ(ierr);
4399b54502bSHong Zhang   PetscFunctionReturn(0);
4409b54502bSHong Zhang }
4419b54502bSHong Zhang EXTERN_C_END
442