xref: /petsc/src/ksp/pc/impls/factor/ilu/ilu.c (revision fe97e370694949207a72f2d229e0cb264fc7a055)
1dba47a55SKris Buschelman #define PETSCKSP_DLL
2dba47a55SKris Buschelman 
39b54502bSHong Zhang /*
49b54502bSHong Zhang    Defines a ILU factorization preconditioner for any Mat implementation
59b54502bSHong Zhang */
67c4f633dSBarry Smith #include "../src/ksp/pc/impls/factor/ilu/ilu.h"     /*I "petscpc.h"  I*/
79b54502bSHong Zhang 
89b54502bSHong Zhang /* ------------------------------------------------------------------------------------------*/
9afaefe49SHong Zhang EXTERN_C_BEGIN
10afaefe49SHong Zhang #undef __FUNCT__
112401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill_ILU"
122401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_ILU(PC pc,PetscTruth flag)
132401956bSBarry Smith {
14075768bcSBarry Smith   PC_ILU *lu = (PC_ILU*)pc->data;
152401956bSBarry Smith 
162401956bSBarry Smith   PetscFunctionBegin;
172401956bSBarry Smith   lu->reusefill = flag;
182401956bSBarry Smith   PetscFunctionReturn(0);
192401956bSBarry Smith }
202401956bSBarry Smith EXTERN_C_END
212401956bSBarry Smith 
222401956bSBarry Smith EXTERN_C_BEGIN
232401956bSBarry Smith #undef __FUNCT__
242401956bSBarry Smith #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal_ILU"
252401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorReorderForNonzeroDiagonal_ILU(PC pc,PetscReal z)
269b54502bSHong Zhang {
279b54502bSHong Zhang   PC_ILU *ilu = (PC_ILU*)pc->data;
289b54502bSHong Zhang 
299b54502bSHong Zhang   PetscFunctionBegin;
309b54502bSHong Zhang   ilu->nonzerosalongdiagonal = PETSC_TRUE;
319b54502bSHong Zhang   if (z == PETSC_DECIDE) {
329b54502bSHong Zhang     ilu->nonzerosalongdiagonaltol = 1.e-10;
339b54502bSHong Zhang   } else {
349b54502bSHong Zhang     ilu->nonzerosalongdiagonaltol = z;
359b54502bSHong Zhang   }
369b54502bSHong Zhang   PetscFunctionReturn(0);
379b54502bSHong Zhang }
389b54502bSHong Zhang EXTERN_C_END
399b54502bSHong Zhang 
409b54502bSHong Zhang #undef __FUNCT__
419b54502bSHong Zhang #define __FUNCT__ "PCDestroy_ILU_Internal"
429b54502bSHong Zhang PetscErrorCode PCDestroy_ILU_Internal(PC pc)
439b54502bSHong Zhang {
449b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
459b54502bSHong Zhang   PetscErrorCode ierr;
469b54502bSHong Zhang 
479b54502bSHong Zhang   PetscFunctionBegin;
48075768bcSBarry Smith   if (!ilu->inplace && ((PC_Factor*)ilu)->fact) {ierr = MatDestroy(((PC_Factor*)ilu)->fact);CHKERRQ(ierr);}
499b54502bSHong Zhang   if (ilu->row && ilu->col && ilu->row != ilu->col) {ierr = ISDestroy(ilu->row);CHKERRQ(ierr);}
509b54502bSHong Zhang   if (ilu->col) {ierr = ISDestroy(ilu->col);CHKERRQ(ierr);}
519b54502bSHong Zhang   PetscFunctionReturn(0);
529b54502bSHong Zhang }
539b54502bSHong Zhang 
549b54502bSHong Zhang EXTERN_C_BEGIN
559b54502bSHong Zhang #undef __FUNCT__
562401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseDropTolerance_ILU"
572401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseDropTolerance_ILU(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount)
589b54502bSHong Zhang {
59075768bcSBarry Smith   PC_ILU         *ilu = (PC_ILU*)pc->data;
609b54502bSHong Zhang   PetscErrorCode ierr;
619b54502bSHong Zhang 
629b54502bSHong Zhang   PetscFunctionBegin;
63075768bcSBarry Smith   if (pc->setupcalled && (!ilu->usedt || ((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) {
649b54502bSHong Zhang     pc->setupcalled   = 0;
659b54502bSHong Zhang     ierr = PCDestroy_ILU_Internal(pc);CHKERRQ(ierr);
669b54502bSHong Zhang   }
679b54502bSHong Zhang   ilu->usedt                      = PETSC_TRUE;
68075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dt      = dt;
69075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcol   = dtcol;
70075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcount = dtcount;
71075768bcSBarry Smith   ((PC_Factor*)ilu)->info.fill    = PETSC_DEFAULT;
729b54502bSHong Zhang   PetscFunctionReturn(0);
739b54502bSHong Zhang }
749b54502bSHong Zhang EXTERN_C_END
759b54502bSHong Zhang 
769b54502bSHong Zhang EXTERN_C_BEGIN
779b54502bSHong Zhang #undef __FUNCT__
782401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering_ILU"
792401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_ILU(PC pc,PetscTruth flag)
809b54502bSHong Zhang {
81075768bcSBarry Smith   PC_ILU *ilu = (PC_ILU*)pc->data;
829b54502bSHong Zhang 
839b54502bSHong Zhang   PetscFunctionBegin;
849b54502bSHong Zhang   ilu->reuseordering = flag;
859b54502bSHong Zhang   PetscFunctionReturn(0);
869b54502bSHong Zhang }
879b54502bSHong Zhang EXTERN_C_END
889b54502bSHong Zhang 
899b54502bSHong Zhang EXTERN_C_BEGIN
909b54502bSHong Zhang #undef __FUNCT__
912401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseInPlace_ILU"
922401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_ILU(PC pc)
939b54502bSHong Zhang {
94075768bcSBarry Smith   PC_ILU *dir = (PC_ILU*)pc->data;
959b54502bSHong Zhang 
969b54502bSHong Zhang   PetscFunctionBegin;
979b54502bSHong Zhang   dir->inplace = PETSC_TRUE;
989b54502bSHong Zhang   PetscFunctionReturn(0);
999b54502bSHong Zhang }
1009b54502bSHong Zhang EXTERN_C_END
1019b54502bSHong Zhang 
1029b54502bSHong Zhang #undef __FUNCT__
1039b54502bSHong Zhang #define __FUNCT__ "PCSetFromOptions_ILU"
1049b54502bSHong Zhang static PetscErrorCode PCSetFromOptions_ILU(PC pc)
1059b54502bSHong Zhang {
1069b54502bSHong Zhang   PetscErrorCode ierr;
1079b54502bSHong Zhang   PetscInt       dtmax = 3,itmp;
1088ff23777SHong Zhang   PetscTruth     flg;
1099b54502bSHong Zhang   PetscReal      dt[3];
1109b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
1119b54502bSHong Zhang   PetscReal      tol;
1129b54502bSHong Zhang 
1139b54502bSHong Zhang   PetscFunctionBegin;
1149b54502bSHong Zhang   ierr = PetscOptionsHead("ILU Options");CHKERRQ(ierr);
1158ff23777SHong Zhang     ierr = PCSetFromOptions_Factor(pc);CHKERRQ(ierr);
1168ff23777SHong Zhang 
117075768bcSBarry Smith     ierr = PetscOptionsInt("-pc_factor_levels","levels of fill","PCFactorSetLevels",(PetscInt)((PC_Factor*)ilu)->info.levels,&itmp,&flg);CHKERRQ(ierr);
118075768bcSBarry Smith     if (flg) ((PC_Factor*)ilu)->info.levels = itmp;
11990d69ab7SBarry Smith     flg  = PETSC_FALSE;
12090d69ab7SBarry Smith     ierr = PetscOptionsTruth("-pc_factor_diagonal_fill","Allow fill into empty diagonal entry","PCFactorSetAllowDiagonalFill",flg,&flg,PETSC_NULL);CHKERRQ(ierr);
121075768bcSBarry Smith     ((PC_Factor*)ilu)->info.diagonal_fill = (double) flg;
1229b54502bSHong Zhang 
123075768bcSBarry Smith     dt[0] = ((PC_Factor*)ilu)->info.dt;
124075768bcSBarry Smith     dt[1] = ((PC_Factor*)ilu)->info.dtcol;
125075768bcSBarry Smith     dt[2] = ((PC_Factor*)ilu)->info.dtcount;
1262401956bSBarry Smith     ierr = PetscOptionsRealArray("-pc_factor_use_drop_tolerance","<dt,dtcol,maxrowcount>","PCFactorSetUseDropTolerance",dt,&dtmax,&flg);CHKERRQ(ierr);
1279b54502bSHong Zhang     if (flg) {
1282401956bSBarry Smith       ierr = PCFactorSetUseDropTolerance(pc,dt[0],dt[1],(PetscInt)dt[2]);CHKERRQ(ierr);
1299b54502bSHong Zhang     }
1308ff23777SHong Zhang 
1312401956bSBarry Smith     ierr = PetscOptionsName("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",&flg);CHKERRQ(ierr);
1329b54502bSHong Zhang     if (flg) {
1339b54502bSHong Zhang       tol = PETSC_DECIDE;
1342401956bSBarry Smith       ierr = PetscOptionsReal("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",ilu->nonzerosalongdiagonaltol,&tol,0);CHKERRQ(ierr);
1352401956bSBarry Smith       ierr = PCFactorReorderForNonzeroDiagonal(pc,tol);CHKERRQ(ierr);
1369b54502bSHong Zhang     }
1379b54502bSHong Zhang 
1389b54502bSHong Zhang   ierr = PetscOptionsTail();CHKERRQ(ierr);
1399b54502bSHong Zhang   PetscFunctionReturn(0);
1409b54502bSHong Zhang }
1419b54502bSHong Zhang 
1429b54502bSHong Zhang #undef __FUNCT__
1439b54502bSHong Zhang #define __FUNCT__ "PCView_ILU"
1449b54502bSHong Zhang static PetscErrorCode PCView_ILU(PC pc,PetscViewer viewer)
1459b54502bSHong Zhang {
1469b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
1479b54502bSHong Zhang   PetscErrorCode ierr;
1489b54502bSHong Zhang   PetscTruth     isstring,iascii;
1499b54502bSHong Zhang 
1509b54502bSHong Zhang   PetscFunctionBegin;
1519b54502bSHong Zhang   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr);
1529b54502bSHong Zhang   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr);
1539b54502bSHong Zhang   if (iascii) {
1549b54502bSHong Zhang     if (ilu->usedt) {
155075768bcSBarry Smith         ierr = PetscViewerASCIIPrintf(viewer,"  ILU: drop tolerance %G\n",((PC_Factor*)ilu)->info.dt);CHKERRQ(ierr);
156075768bcSBarry Smith         ierr = PetscViewerASCIIPrintf(viewer,"  ILU: max nonzeros per row %D\n",(PetscInt)((PC_Factor*)ilu)->info.dtcount);CHKERRQ(ierr);
157075768bcSBarry Smith         ierr = PetscViewerASCIIPrintf(viewer,"  ILU: column permutation tolerance %G\n",((PC_Factor*)ilu)->info.dtcol);CHKERRQ(ierr);
158075768bcSBarry Smith     } else if (((PC_Factor*)ilu)->info.levels == 1) {
159075768bcSBarry Smith         ierr = PetscViewerASCIIPrintf(viewer,"  ILU: %D level of fill\n",(PetscInt)((PC_Factor*)ilu)->info.levels);CHKERRQ(ierr);
1609b54502bSHong Zhang     } else {
161075768bcSBarry Smith         ierr = PetscViewerASCIIPrintf(viewer,"  ILU: %D levels of fill\n",(PetscInt)((PC_Factor*)ilu)->info.levels);CHKERRQ(ierr);
1629b54502bSHong Zhang     }
163075768bcSBarry Smith     ierr = PetscViewerASCIIPrintf(viewer,"  ILU: factor fill ratio allocated %G\n",((PC_Factor*)ilu)->info.fill);CHKERRQ(ierr);
164075768bcSBarry Smith     ierr = PetscViewerASCIIPrintf(viewer,"  ILU: tolerance for zero pivot %G\n",((PC_Factor*)ilu)->info.zeropivot);CHKERRQ(ierr);
165075768bcSBarry Smith     if (((PC_Factor*)ilu)->info.shiftpd) {ierr = PetscViewerASCIIPrintf(viewer,"  ILU: using Manteuffel shift\n");CHKERRQ(ierr);}
166075768bcSBarry Smith     if (((PC_Factor*)ilu)->info.shiftnz) {ierr = PetscViewerASCIIPrintf(viewer,"  ILU: using diagonal shift to prevent zero pivot\n");CHKERRQ(ierr);}
167075768bcSBarry Smith     if (((PC_Factor*)ilu)->info.shiftinblocks) {ierr = PetscViewerASCIIPrintf(viewer,"  ILU: using diagonal shift on blocks to prevent zero pivot\n");CHKERRQ(ierr);}
1689b54502bSHong Zhang     if (ilu->inplace) {ierr = PetscViewerASCIIPrintf(viewer,"       in-place factorization\n");CHKERRQ(ierr);}
1699b54502bSHong Zhang     else              {ierr = PetscViewerASCIIPrintf(viewer,"       out-of-place factorization\n");CHKERRQ(ierr);}
170075768bcSBarry Smith     ierr = PetscViewerASCIIPrintf(viewer,"       matrix ordering: %s\n",((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);
1719b54502bSHong Zhang     if (ilu->reusefill)     {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing fill from past factorization\n");CHKERRQ(ierr);}
1729b54502bSHong Zhang     if (ilu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing reordering from past factorization\n");CHKERRQ(ierr);}
173075768bcSBarry Smith     if (((PC_Factor*)ilu)->fact) {
174f3a39becSBarry Smith       ierr = PetscViewerASCIIPrintf(viewer,"  ILU: factor fill ratio needed %G\n",ilu->actualfill);CHKERRQ(ierr);
1759b54502bSHong Zhang       ierr = PetscViewerASCIIPrintf(viewer,"       Factored matrix follows\n");CHKERRQ(ierr);
1769b54502bSHong Zhang       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
177f3a39becSBarry Smith       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
178f3a39becSBarry Smith       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
1799b54502bSHong Zhang       ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr);
180075768bcSBarry Smith       ierr = MatView(((PC_Factor*)ilu)->fact,viewer);CHKERRQ(ierr);
1819b54502bSHong Zhang       ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr);
1829b54502bSHong Zhang       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
183f3a39becSBarry Smith       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
184f3a39becSBarry Smith       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
1859b54502bSHong Zhang     }
1869b54502bSHong Zhang   } else if (isstring) {
187075768bcSBarry Smith     ierr = PetscViewerStringSPrintf(viewer," lvls=%D,order=%s",(PetscInt)((PC_Factor*)ilu)->info.levels,((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);CHKERRQ(ierr);
1889b54502bSHong Zhang   } else {
1899b54502bSHong Zhang     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCILU",((PetscObject)viewer)->type_name);
1909b54502bSHong Zhang   }
1919b54502bSHong Zhang   PetscFunctionReturn(0);
1929b54502bSHong Zhang }
1939b54502bSHong Zhang 
1949b54502bSHong Zhang #undef __FUNCT__
1959b54502bSHong Zhang #define __FUNCT__ "PCSetUp_ILU"
1969b54502bSHong Zhang static PetscErrorCode PCSetUp_ILU(PC pc)
1979b54502bSHong Zhang {
1989b54502bSHong Zhang   PetscErrorCode ierr;
1999b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
200f3a39becSBarry Smith   MatInfo        info;
2019b54502bSHong Zhang 
2029b54502bSHong Zhang   PetscFunctionBegin;
2039b54502bSHong Zhang   if (ilu->inplace) {
204075768bcSBarry Smith     CHKMEMQ;
2059b54502bSHong Zhang     if (!pc->setupcalled) {
2069b54502bSHong Zhang 
2079b54502bSHong Zhang       /* In-place factorization only makes sense with the natural ordering,
2089b54502bSHong Zhang          so we only need to get the ordering once, even if nonzero structure changes */
209075768bcSBarry Smith       ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
210efee365bSSatish Balay       if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
211efee365bSSatish Balay       if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
2129b54502bSHong Zhang     }
2139b54502bSHong Zhang 
2149b54502bSHong Zhang     /* In place ILU only makes sense with fill factor of 1.0 because
2159b54502bSHong Zhang        cannot have levels of fill */
216075768bcSBarry Smith     ((PC_Factor*)ilu)->info.fill          = 1.0;
21775567043SBarry Smith     ((PC_Factor*)ilu)->info.diagonal_fill = 0.0;
218*fe97e370SBarry Smith     ierr = MatILUFactor(pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);CHKMEMQ;
219075768bcSBarry Smith     ((PC_Factor*)ilu)->fact = pc->pmat;
2209b54502bSHong Zhang   } else {
2219b54502bSHong Zhang     if (!pc->setupcalled) {
2229b54502bSHong Zhang       /* first time in so compute reordering and symbolic factorization */
223075768bcSBarry Smith       ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
224efee365bSSatish Balay       if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
225efee365bSSatish Balay       if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
2269b54502bSHong Zhang       /*  Remove zeros along diagonal?     */
2279b54502bSHong Zhang       if (ilu->nonzerosalongdiagonal) {
2289b54502bSHong Zhang         ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr);
2299b54502bSHong Zhang       }
230075768bcSBarry Smith     CHKMEMQ;
23111bfe483SHong Zhang       ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
232075768bcSBarry Smith     CHKMEMQ;
233075768bcSBarry Smith       ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
234075768bcSBarry Smith       ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
235f3a39becSBarry Smith       ilu->actualfill = info.fill_ratio_needed;
236075768bcSBarry Smith       ierr = PetscLogObjectParent(pc,((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
2379b54502bSHong Zhang     } else if (pc->flag != SAME_NONZERO_PATTERN) {
2389b54502bSHong Zhang       if (!ilu->reuseordering) {
2399b54502bSHong Zhang         /* compute a new ordering for the ILU */
2409b54502bSHong Zhang         ierr = ISDestroy(ilu->row);CHKERRQ(ierr);
2419b54502bSHong Zhang         ierr = ISDestroy(ilu->col);CHKERRQ(ierr);
242075768bcSBarry Smith         ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr);
243efee365bSSatish Balay         if (ilu->row) {ierr = PetscLogObjectParent(pc,ilu->row);CHKERRQ(ierr);}
244efee365bSSatish Balay         if (ilu->col) {ierr = PetscLogObjectParent(pc,ilu->col);CHKERRQ(ierr);}
2459b54502bSHong Zhang         /*  Remove zeros along diagonal?     */
2469b54502bSHong Zhang         if (ilu->nonzerosalongdiagonal) {
2479b54502bSHong Zhang           ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr);
2489b54502bSHong Zhang         }
2499b54502bSHong Zhang       }
250075768bcSBarry Smith       ierr = MatDestroy(((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
251075768bcSBarry Smith       ierr = MatGetFactor(pc->pmat,MAT_SOLVER_PETSC,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
252075768bcSBarry Smith       ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
253075768bcSBarry Smith       ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
254f3a39becSBarry Smith       ilu->actualfill = info.fill_ratio_needed;
255075768bcSBarry Smith       ierr = PetscLogObjectParent(pc,((PC_Factor*)ilu)->fact);CHKERRQ(ierr);
2569b54502bSHong Zhang     }
257075768bcSBarry Smith     CHKMEMQ;
258075768bcSBarry Smith     ierr = MatLUFactorNumeric(((PC_Factor*)ilu)->fact,pc->pmat,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
259075768bcSBarry Smith     CHKMEMQ;
2609b54502bSHong Zhang   }
2619b54502bSHong Zhang   PetscFunctionReturn(0);
2629b54502bSHong Zhang }
2639b54502bSHong Zhang 
2649b54502bSHong Zhang #undef __FUNCT__
2659b54502bSHong Zhang #define __FUNCT__ "PCDestroy_ILU"
2669b54502bSHong Zhang static PetscErrorCode PCDestroy_ILU(PC pc)
2679b54502bSHong Zhang {
2689b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2699b54502bSHong Zhang   PetscErrorCode ierr;
2709b54502bSHong Zhang 
2719b54502bSHong Zhang   PetscFunctionBegin;
2729b54502bSHong Zhang   ierr = PCDestroy_ILU_Internal(pc);CHKERRQ(ierr);
27362ba1fd3SLisandro Dalcin   ierr = PetscStrfree(((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr);
274075768bcSBarry Smith   ierr = PetscStrfree(((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);
2759b54502bSHong Zhang   ierr = PetscFree(ilu);CHKERRQ(ierr);
2769b54502bSHong Zhang   PetscFunctionReturn(0);
2779b54502bSHong Zhang }
2789b54502bSHong Zhang 
2799b54502bSHong Zhang #undef __FUNCT__
2809b54502bSHong Zhang #define __FUNCT__ "PCApply_ILU"
2819b54502bSHong Zhang static PetscErrorCode PCApply_ILU(PC pc,Vec x,Vec y)
2829b54502bSHong Zhang {
2839b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2849b54502bSHong Zhang   PetscErrorCode ierr;
2859b54502bSHong Zhang 
2869b54502bSHong Zhang   PetscFunctionBegin;
287075768bcSBarry Smith   ierr = MatSolve(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr);
2889b54502bSHong Zhang   PetscFunctionReturn(0);
2899b54502bSHong Zhang }
2909b54502bSHong Zhang 
2919b54502bSHong Zhang #undef __FUNCT__
2929b54502bSHong Zhang #define __FUNCT__ "PCApplyTranspose_ILU"
2939b54502bSHong Zhang static PetscErrorCode PCApplyTranspose_ILU(PC pc,Vec x,Vec y)
2949b54502bSHong Zhang {
2959b54502bSHong Zhang   PC_ILU         *ilu = (PC_ILU*)pc->data;
2969b54502bSHong Zhang   PetscErrorCode ierr;
2979b54502bSHong Zhang 
2989b54502bSHong Zhang   PetscFunctionBegin;
299075768bcSBarry Smith   ierr = MatSolveTranspose(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr);
3009b54502bSHong Zhang   PetscFunctionReturn(0);
3019b54502bSHong Zhang }
3029b54502bSHong Zhang 
3039b54502bSHong Zhang /*MC
3049b54502bSHong Zhang      PCILU - Incomplete factorization preconditioners.
3059b54502bSHong Zhang 
3069b54502bSHong Zhang    Options Database Keys:
3072401956bSBarry Smith +  -pc_factor_levels <k> - number of levels of fill for ILU(k)
3082401956bSBarry Smith .  -pc_factor_in_place - only for ILU(0) with natural ordering, reuses the space of the matrix for
3099b54502bSHong Zhang                       its factorization (overwrites original matrix)
3102401956bSBarry Smith .  -pc_factor_diagonal_fill - fill in a zero diagonal even if levels of fill indicate it wouldn't be fill
3112401956bSBarry Smith .  -pc_factor_reuse_ordering - reuse ordering of factorized matrix from previous factorization
3122401956bSBarry Smith .  -pc_factor_use_drop_tolerance <dt,dtcol,maxrowcount> - use Saad's drop tolerance ILUdt
31355ba2a51SBarry Smith .  -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1
3142401956bSBarry Smith .  -pc_factor_nonzeros_along_diagonal - reorder the matrix before factorization to remove zeros from the diagonal,
3159b54502bSHong Zhang                                    this decreases the chance of getting a zero pivot
3162401956bSBarry Smith .  -pc_factor_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix
3172401956bSBarry Smith .  -pc_factor_pivot_in_blocks - for block ILU(k) factorization, i.e. with BAIJ matrices with block size larger
3189b54502bSHong Zhang                              than 1 the diagonal blocks are factored with partial pivoting (this increases the
3199b54502bSHong Zhang                              stability of the ILU factorization
32062bba022SBarry Smith .  -pc_factor_shift_in_blocks - adds a small diagonal to any block if it is singular during ILU factorization
321f251bdbdSHong Zhang .  -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default
32253ae6812SBarry Smith -  -pc_factor_shift_positive_definite true or false - Activate/Deactivate PCFactorSetShiftPd(); the value
32353ae6812SBarry Smith    is optional with true being the default
3249b54502bSHong Zhang 
3259b54502bSHong Zhang    Level: beginner
3269b54502bSHong Zhang 
3279b54502bSHong Zhang   Concepts: incomplete factorization
3289b54502bSHong Zhang 
329d9ba8547SBarry Smith    Notes: Only implemented for some matrix formats. (for parallel see PCHYPRE for hypre's ILU)
3309b54502bSHong Zhang 
3319b54502bSHong Zhang           For BAIJ matrices this implements a point block ILU
3329b54502bSHong Zhang 
333c582cd25SBarry Smith    References:
334c582cd25SBarry Smith    T. Dupont, R. Kendall, and H. Rachford. An approximate factorization procedure for solving
335c582cd25SBarry Smith    self-adjoint elliptic difference equations. SIAM J. Numer. Anal., 5:559--573, 1968.
336c582cd25SBarry Smith 
337c582cd25SBarry Smith    T.A. Oliphant. An implicit numerical method for solving two-dimensional time-dependent dif-
338c582cd25SBarry Smith    fusion problems. Quart. Appl. Math., 19:221--229, 1961.
339c582cd25SBarry Smith 
340c582cd25SBarry Smith    Review article: APPROXIMATE AND INCOMPLETE FACTORIZATIONS, TONY F. CHAN AND HENK A. VAN DER VORST
341c582cd25SBarry Smith       http://igitur-archive.library.uu.nl/math/2001-0621-115821/proc.pdf chapter in Parallel Numerical
342c582cd25SBarry Smith       Algorithms, edited by D. Keyes, A. Semah, V. Venkatakrishnan, ICASE/LaRC Interdisciplinary Series in
343c582cd25SBarry Smith       Science and Engineering, Kluwer, pp. 167--202.
344c582cd25SBarry Smith 
345c582cd25SBarry Smith 
3469b54502bSHong Zhang .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType,
3478ff23777SHong Zhang            PCFactorSetZeroPivot(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(), PCFactorSetShiftInBlocks(),
3488ff23777SHong Zhang            PCFactorSetUseDropTolerance(),PCFactorSetFill(), PCFactorSetMatOrderingType(), PCFactorSetReuseOrdering(),
3498ff23777SHong Zhang            PCFactorSetLevels(), PCFactorSetUseInPlace(), PCFactorSetAllowDiagonalFill(), PCFactorSetPivotInBlocks()
3509b54502bSHong Zhang 
3519b54502bSHong Zhang M*/
3529b54502bSHong Zhang 
3539b54502bSHong Zhang EXTERN_C_BEGIN
3549b54502bSHong Zhang #undef __FUNCT__
3559b54502bSHong Zhang #define __FUNCT__ "PCCreate_ILU"
356dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_ILU(PC pc)
3579b54502bSHong Zhang {
3589b54502bSHong Zhang   PetscErrorCode ierr;
3599b54502bSHong Zhang   PC_ILU         *ilu;
3609b54502bSHong Zhang 
3619b54502bSHong Zhang   PetscFunctionBegin;
36238f2d2fdSLisandro Dalcin   ierr = PetscNewLog(pc,PC_ILU,&ilu);CHKERRQ(ierr);
3639b54502bSHong Zhang 
364075768bcSBarry Smith   ((PC_Factor*)ilu)->fact                    = 0;
365075768bcSBarry Smith   ierr = MatFactorInfoInitialize(&((PC_Factor*)ilu)->info);CHKERRQ(ierr);
36675567043SBarry Smith   ((PC_Factor*)ilu)->info.levels             = 0.;
367075768bcSBarry Smith   ((PC_Factor*)ilu)->info.fill               = 1.0;
3689b54502bSHong Zhang   ilu->col                     = 0;
3699b54502bSHong Zhang   ilu->row                     = 0;
3709b54502bSHong Zhang   ilu->inplace                 = PETSC_FALSE;
37111bfe483SHong Zhang   ierr = PetscStrallocpy(MAT_SOLVER_PETSC,&((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr);
372075768bcSBarry Smith   ierr = PetscStrallocpy(MATORDERING_NATURAL,&((PC_Factor*)ilu)->ordering);CHKERRQ(ierr);
3739b54502bSHong Zhang   ilu->reuseordering           = PETSC_FALSE;
3749b54502bSHong Zhang   ilu->usedt                   = PETSC_FALSE;
375075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dt                 = PETSC_DEFAULT;
376075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcount            = PETSC_DEFAULT;
377075768bcSBarry Smith   ((PC_Factor*)ilu)->info.dtcol              = PETSC_DEFAULT;
378075768bcSBarry Smith   ((PC_Factor*)ilu)->info.shiftnz            = 1.e-12;
379075768bcSBarry Smith   ((PC_Factor*)ilu)->info.shiftpd            = 0.0; /* false */
380075768bcSBarry Smith   ((PC_Factor*)ilu)->info.zeropivot          = 1.e-12;
381075768bcSBarry Smith   ((PC_Factor*)ilu)->info.pivotinblocks      = 1.0;
382075768bcSBarry Smith   ((PC_Factor*)ilu)->info.shiftinblocks      = 1.e-12;
3839b54502bSHong Zhang   ilu->reusefill               = PETSC_FALSE;
38475567043SBarry Smith   ((PC_Factor*)ilu)->info.diagonal_fill      = 0.;
3859b54502bSHong Zhang   pc->data                     = (void*)ilu;
3869b54502bSHong Zhang 
3879b54502bSHong Zhang   pc->ops->destroy             = PCDestroy_ILU;
3889b54502bSHong Zhang   pc->ops->apply               = PCApply_ILU;
3899b54502bSHong Zhang   pc->ops->applytranspose      = PCApplyTranspose_ILU;
3909b54502bSHong Zhang   pc->ops->setup               = PCSetUp_ILU;
3919b54502bSHong Zhang   pc->ops->setfromoptions      = PCSetFromOptions_ILU;
39285317021SBarry Smith   pc->ops->getfactoredmatrix   = PCFactorGetMatrix_Factor;
3939b54502bSHong Zhang   pc->ops->view                = PCView_ILU;
3949b54502bSHong Zhang   pc->ops->applyrichardson     = 0;
3959b54502bSHong Zhang 
39685317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Factor",
39785317021SBarry Smith                     PCFactorSetZeroPivot_Factor);CHKERRQ(ierr);
39885317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Factor",
39985317021SBarry Smith                     PCFactorSetShiftNonzero_Factor);CHKERRQ(ierr);
40085317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Factor",
40185317021SBarry Smith                     PCFactorSetShiftPd_Factor);CHKERRQ(ierr);
4028ff23777SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftInBlocks_C","PCFactorSetShiftInBlocks_Factor",
4038ff23777SHong Zhang                     PCFactorSetShiftInBlocks_Factor);CHKERRQ(ierr);
404afaefe49SHong Zhang 
4057112b564SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorGetMatSolverPackage_C","PCFactorGetMatSolverPackage_Factor",
4067112b564SBarry Smith                     PCFactorGetMatSolverPackage_Factor);CHKERRQ(ierr);
40711bfe483SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatSolverPackage_C","PCFactorSetMatSolverPackage_Factor",
40811bfe483SHong Zhang                     PCFactorSetMatSolverPackage_Factor);CHKERRQ(ierr);
4092401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseDropTolerance_C","PCFactorSetUseDropTolerance_ILU",
4102401956bSBarry Smith                     PCFactorSetUseDropTolerance_ILU);CHKERRQ(ierr);
41185317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Factor",
41285317021SBarry Smith                     PCFactorSetFill_Factor);CHKERRQ(ierr);
41385317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Factor",
41485317021SBarry Smith                     PCFactorSetMatOrderingType_Factor);CHKERRQ(ierr);
4152401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_ILU",
4162401956bSBarry Smith                     PCFactorSetReuseOrdering_ILU);CHKERRQ(ierr);
4172401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_ILU",
4182401956bSBarry Smith                     PCFactorSetReuseFill_ILU);CHKERRQ(ierr);
41985317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetLevels_C","PCFactorSetLevels_Factor",
42085317021SBarry Smith                     PCFactorSetLevels_Factor);CHKERRQ(ierr);
4212401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_ILU",
4222401956bSBarry Smith                     PCFactorSetUseInPlace_ILU);CHKERRQ(ierr);
42385317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C","PCFactorSetAllowDiagonalFill_Factor",
42485317021SBarry Smith                     PCFactorSetAllowDiagonalFill_Factor);CHKERRQ(ierr);
42585317021SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetPivotInBlocks_C","PCFactorSetPivotInBlocks_Factor",
42685317021SBarry Smith                     PCFactorSetPivotInBlocks_Factor);CHKERRQ(ierr);
4272401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorReorderForNonzeroDiagonal_C","PCFactorReorderForNonzeroDiagonal_ILU",
4282401956bSBarry Smith                     PCFactorReorderForNonzeroDiagonal_ILU);CHKERRQ(ierr);
4299b54502bSHong Zhang   PetscFunctionReturn(0);
4309b54502bSHong Zhang }
4319b54502bSHong Zhang EXTERN_C_END
432