xref: /petsc/src/ksp/pc/impls/factor/cholesky/cholesky.c (revision a4fd02ac6712e7c4a27c9ffdf1171591c41a72f0)
1dba47a55SKris Buschelman #define PETSCKSP_DLL
29b54502bSHong Zhang 
39b54502bSHong Zhang /*
49b54502bSHong Zhang    Defines a direct factorization preconditioner for any Mat implementation
59b54502bSHong Zhang    Note: this need not be consided a preconditioner since it supplies
69b54502bSHong Zhang          a direct solver.
79b54502bSHong Zhang */
86356e834SBarry Smith #include "private/pcimpl.h"                /*I "petscpc.h" I*/
99b54502bSHong Zhang 
109b54502bSHong Zhang typedef struct {
119b54502bSHong Zhang   Mat             fact;             /* factored matrix */
129b54502bSHong Zhang   PetscReal       actualfill;       /* actual fill in factor */
139b54502bSHong Zhang   PetscTruth      inplace;          /* flag indicating in-place factorization */
149b54502bSHong Zhang   IS              row,col;          /* index sets used for reordering */
159b54502bSHong Zhang   MatOrderingType ordering;         /* matrix ordering */
169b54502bSHong Zhang   PetscTruth      reuseordering;    /* reuses previous reordering computed */
179b54502bSHong Zhang   PetscTruth      reusefill;        /* reuse fill from previous Cholesky */
189b54502bSHong Zhang   MatFactorInfo   info;
199b54502bSHong Zhang } PC_Cholesky;
209b54502bSHong Zhang 
219b54502bSHong Zhang EXTERN_C_BEGIN
229b54502bSHong Zhang #undef __FUNCT__
23afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetZeroPivot_Cholesky"
24dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_Cholesky(PC pc,PetscReal z)
25afaefe49SHong Zhang {
26afaefe49SHong Zhang   PC_Cholesky *ch;
27afaefe49SHong Zhang 
28afaefe49SHong Zhang   PetscFunctionBegin;
29afaefe49SHong Zhang   ch                 = (PC_Cholesky*)pc->data;
30afaefe49SHong Zhang   ch->info.zeropivot = z;
31afaefe49SHong Zhang   PetscFunctionReturn(0);
32afaefe49SHong Zhang }
33afaefe49SHong Zhang EXTERN_C_END
34afaefe49SHong Zhang 
35afaefe49SHong Zhang EXTERN_C_BEGIN
36afaefe49SHong Zhang #undef __FUNCT__
37afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetShiftNonzero_Cholesky"
38dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_Cholesky(PC pc,PetscReal shift)
39afaefe49SHong Zhang {
40afaefe49SHong Zhang   PC_Cholesky *dir;
41afaefe49SHong Zhang 
42afaefe49SHong Zhang   PetscFunctionBegin;
43afaefe49SHong Zhang   dir = (PC_Cholesky*)pc->data;
44afaefe49SHong Zhang   if (shift == (PetscReal) PETSC_DECIDE) {
45afaefe49SHong Zhang     dir->info.shiftnz = 1.e-12;
46afaefe49SHong Zhang   } else {
47afaefe49SHong Zhang     dir->info.shiftnz = shift;
48afaefe49SHong Zhang   }
49afaefe49SHong Zhang   PetscFunctionReturn(0);
50afaefe49SHong Zhang }
51afaefe49SHong Zhang EXTERN_C_END
52afaefe49SHong Zhang 
53afaefe49SHong Zhang EXTERN_C_BEGIN
54afaefe49SHong Zhang #undef __FUNCT__
55afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetShiftPd_Cholesky"
56dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_Cholesky(PC pc,PetscTruth shift)
57afaefe49SHong Zhang {
58afaefe49SHong Zhang   PC_Cholesky *dir;
59afaefe49SHong Zhang 
60afaefe49SHong Zhang   PetscFunctionBegin;
61afaefe49SHong Zhang   dir = (PC_Cholesky*)pc->data;
62fbf22428SSatish Balay   if (shift) {
63fbf22428SSatish Balay     dir->info.shift_fraction = 0.0;
64fbf22428SSatish Balay     dir->info.shiftpd = 1.0;
65fbf22428SSatish Balay   } else {
66fbf22428SSatish Balay     dir->info.shiftpd = 0.0;
67fbf22428SSatish Balay   }
68afaefe49SHong Zhang   PetscFunctionReturn(0);
69afaefe49SHong Zhang }
70afaefe49SHong Zhang EXTERN_C_END
71afaefe49SHong Zhang 
72afaefe49SHong Zhang EXTERN_C_BEGIN
73afaefe49SHong Zhang #undef __FUNCT__
742401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering_Cholesky"
752401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_Cholesky(PC pc,PetscTruth flag)
769b54502bSHong Zhang {
779b54502bSHong Zhang   PC_Cholesky *lu;
789b54502bSHong Zhang 
799b54502bSHong Zhang   PetscFunctionBegin;
809b54502bSHong Zhang   lu               = (PC_Cholesky*)pc->data;
819b54502bSHong Zhang   lu->reuseordering = flag;
829b54502bSHong Zhang   PetscFunctionReturn(0);
839b54502bSHong Zhang }
849b54502bSHong Zhang EXTERN_C_END
859b54502bSHong Zhang 
869b54502bSHong Zhang EXTERN_C_BEGIN
879b54502bSHong Zhang #undef __FUNCT__
882401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill_Cholesky"
892401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_Cholesky(PC pc,PetscTruth flag)
909b54502bSHong Zhang {
919b54502bSHong Zhang   PC_Cholesky *lu;
929b54502bSHong Zhang 
939b54502bSHong Zhang   PetscFunctionBegin;
949b54502bSHong Zhang   lu = (PC_Cholesky*)pc->data;
959b54502bSHong Zhang   lu->reusefill = flag;
969b54502bSHong Zhang   PetscFunctionReturn(0);
979b54502bSHong Zhang }
989b54502bSHong Zhang EXTERN_C_END
999b54502bSHong Zhang 
1009b54502bSHong Zhang #undef __FUNCT__
1019b54502bSHong Zhang #define __FUNCT__ "PCSetFromOptions_Cholesky"
1029b54502bSHong Zhang static PetscErrorCode PCSetFromOptions_Cholesky(PC pc)
1039b54502bSHong Zhang {
1049b54502bSHong Zhang   PC_Cholesky    *lu = (PC_Cholesky*)pc->data;
1059b54502bSHong Zhang   PetscErrorCode ierr;
1069b54502bSHong Zhang   PetscTruth     flg;
1079b54502bSHong Zhang   char           tname[256];
1089b54502bSHong Zhang   PetscFList     ordlist;
1099b54502bSHong Zhang 
1109b54502bSHong Zhang   PetscFunctionBegin;
1119b54502bSHong Zhang   ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr);
1129b54502bSHong Zhang   ierr = PetscOptionsHead("Cholesky options");CHKERRQ(ierr);
1132401956bSBarry Smith   ierr = PetscOptionsName("-pc_factor_in_place","Form Cholesky in the same memory as the matrix","PCFactorSetUseInPlace",&flg);CHKERRQ(ierr);
1149b54502bSHong Zhang   if (flg) {
1152401956bSBarry Smith     ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr);
1169b54502bSHong Zhang   }
11755ba2a51SBarry Smith   ierr = PetscOptionsReal("-pc_factor_fill","Expected non-zeros in Cholesky/non-zeros in matrix","PCFactorSetFill",lu->info.fill,&lu->info.fill,0);CHKERRQ(ierr);
1189b54502bSHong Zhang 
1192401956bSBarry Smith   ierr = PetscOptionsName("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",&flg);CHKERRQ(ierr);
1209b54502bSHong Zhang   if (flg) {
1212401956bSBarry Smith     ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr);
1229b54502bSHong Zhang   }
1232401956bSBarry Smith   ierr = PetscOptionsName("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",&flg);CHKERRQ(ierr);
1249b54502bSHong Zhang   if (flg) {
1252401956bSBarry Smith     ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr);
1269b54502bSHong Zhang   }
1279b54502bSHong Zhang 
1289b54502bSHong Zhang   ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr);
129e5a9bf91SBarry Smith   ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in Cholesky","PCFactorSetMatOrderingType",ordlist,lu->ordering,tname,256,&flg);CHKERRQ(ierr);
1309b54502bSHong Zhang   if (flg) {
131e5a9bf91SBarry Smith     ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr);
1329b54502bSHong Zhang   }
1339f95998fSHong Zhang   ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr);
1349b54502bSHong Zhang   if (flg) {
135afaefe49SHong Zhang     ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr);
1369b54502bSHong Zhang   }
1379f95998fSHong Zhang   ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",lu->info.shiftnz,&lu->info.shiftnz,0);CHKERRQ(ierr);
1389f95998fSHong Zhang   ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",&flg);CHKERRQ(ierr);
1399b54502bSHong Zhang   if (flg) {
140afaefe49SHong Zhang     ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr);
1419b54502bSHong Zhang   }
142ee45ca4aSHong Zhang   ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",lu->info.zeropivot,&lu->info.zeropivot,0);CHKERRQ(ierr);
1439b54502bSHong Zhang 
1449b54502bSHong Zhang   ierr = PetscOptionsTail();CHKERRQ(ierr);
1459b54502bSHong Zhang   PetscFunctionReturn(0);
1469b54502bSHong Zhang }
1479b54502bSHong Zhang 
1489b54502bSHong Zhang #undef __FUNCT__
1499b54502bSHong Zhang #define __FUNCT__ "PCView_Cholesky"
1509b54502bSHong Zhang static PetscErrorCode PCView_Cholesky(PC pc,PetscViewer viewer)
1519b54502bSHong Zhang {
1529b54502bSHong Zhang   PC_Cholesky    *lu = (PC_Cholesky*)pc->data;
1539b54502bSHong Zhang   PetscErrorCode ierr;
1549b54502bSHong Zhang   PetscTruth     iascii,isstring;
1559b54502bSHong Zhang 
1569b54502bSHong Zhang   PetscFunctionBegin;
1579b54502bSHong Zhang   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr);
1589b54502bSHong Zhang   ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr);
1599b54502bSHong Zhang   if (iascii) {
1609b54502bSHong Zhang 
1619b54502bSHong Zhang     if (lu->inplace) {ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: in-place factorization\n");CHKERRQ(ierr);}
1629b54502bSHong Zhang     else             {ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: out-of-place factorization\n");CHKERRQ(ierr);}
1639b54502bSHong Zhang     ierr = PetscViewerASCIIPrintf(viewer,"    matrix ordering: %s\n",lu->ordering);CHKERRQ(ierr);
1649b54502bSHong Zhang     if (lu->fact) {
165f3a39becSBarry Smith       ierr = PetscViewerASCIIPrintf(viewer,"  Cholesky: factor fill ratio needed %G\n",lu->actualfill);CHKERRQ(ierr);
166f3a39becSBarry Smith       ierr = PetscViewerASCIIPrintf(viewer,"       Factored matrix follows\n");CHKERRQ(ierr);
167f3a39becSBarry Smith       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
168f3a39becSBarry Smith       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
169f3a39becSBarry Smith       ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
170f3a39becSBarry Smith       ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr);
1719b54502bSHong Zhang       ierr = MatView(lu->fact,viewer);CHKERRQ(ierr);
1729b54502bSHong Zhang       ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr);
173f3a39becSBarry Smith       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
174f3a39becSBarry Smith       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
175f3a39becSBarry Smith       ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
1769b54502bSHong Zhang     }
1779b54502bSHong Zhang     if (lu->reusefill)    {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing fill from past factorization\n");CHKERRQ(ierr);}
1789b54502bSHong Zhang     if (lu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer,"       Reusing reordering from past factorization\n");CHKERRQ(ierr);}
1799b54502bSHong Zhang   } else if (isstring) {
1809b54502bSHong Zhang     ierr = PetscViewerStringSPrintf(viewer," order=%s",lu->ordering);CHKERRQ(ierr);CHKERRQ(ierr);
1819b54502bSHong Zhang   } else {
1822401956bSBarry Smith     SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCCHOLESKY",((PetscObject)viewer)->type_name);
1839b54502bSHong Zhang   }
1849b54502bSHong Zhang   PetscFunctionReturn(0);
1859b54502bSHong Zhang }
1869b54502bSHong Zhang 
1879b54502bSHong Zhang #undef __FUNCT__
188*a4fd02acSBarry Smith #define __FUNCT__ "PCFactorGetMatrix_Cholesky"
189*a4fd02acSBarry Smith static PetscErrorCode PCFactorGetMatrix_Cholesky(PC pc,Mat *mat)
1909b54502bSHong Zhang {
1919b54502bSHong Zhang   PC_Cholesky *dir = (PC_Cholesky*)pc->data;
1929b54502bSHong Zhang 
1939b54502bSHong Zhang   PetscFunctionBegin;
1949b54502bSHong Zhang   if (!dir->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()");
1959b54502bSHong Zhang   *mat = dir->fact;
1969b54502bSHong Zhang   PetscFunctionReturn(0);
1979b54502bSHong Zhang }
1989b54502bSHong Zhang 
1999b54502bSHong Zhang #undef __FUNCT__
2009b54502bSHong Zhang #define __FUNCT__ "PCSetUp_Cholesky"
2019b54502bSHong Zhang static PetscErrorCode PCSetUp_Cholesky(PC pc)
2029b54502bSHong Zhang {
2039b54502bSHong Zhang   PetscErrorCode ierr;
2049b54502bSHong Zhang   PetscTruth     flg;
2059b54502bSHong Zhang   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
2069b54502bSHong Zhang 
2079b54502bSHong Zhang   PetscFunctionBegin;
2089b54502bSHong Zhang   if (dir->reusefill && pc->setupcalled) dir->info.fill = dir->actualfill;
2099b54502bSHong Zhang 
2109b54502bSHong Zhang   if (dir->inplace) {
2119b54502bSHong Zhang     if (dir->row && dir->col && (dir->row != dir->col)) {
2129b54502bSHong Zhang       ierr = ISDestroy(dir->row);CHKERRQ(ierr);
2139b54502bSHong Zhang       dir->row = 0;
2149b54502bSHong Zhang     }
2159b54502bSHong Zhang     if (dir->col) {
2169b54502bSHong Zhang       ierr = ISDestroy(dir->col);CHKERRQ(ierr);
2179b54502bSHong Zhang       dir->col = 0;
2189b54502bSHong Zhang     }
2199b54502bSHong Zhang     ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
2209b54502bSHong Zhang     if (dir->col && (dir->row != dir->col)) {  /* only use row ordering for SBAIJ */
2219b54502bSHong Zhang       ierr = ISDestroy(dir->col);CHKERRQ(ierr);
2229b54502bSHong Zhang       dir->col=0;
2239b54502bSHong Zhang     }
224efee365bSSatish Balay     if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
2259b54502bSHong Zhang     ierr = MatCholeskyFactor(pc->pmat,dir->row,&dir->info);CHKERRQ(ierr);
2269b54502bSHong Zhang     dir->fact = pc->pmat;
2279b54502bSHong Zhang   } else {
2289b54502bSHong Zhang     MatInfo info;
2299b54502bSHong Zhang     if (!pc->setupcalled) {
2309b54502bSHong Zhang       ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
2319bfd6278SHong Zhang       /* check if dir->row == dir->col */
2329bfd6278SHong Zhang       ierr = ISEqual(dir->row,dir->col,&flg);CHKERRQ(ierr);
2339bfd6278SHong Zhang       if (!flg) SETERRQ(PETSC_ERR_ARG_INCOMP,"row and column permutations must equal");
2349bfd6278SHong Zhang       ierr = ISDestroy(dir->col);CHKERRQ(ierr); /* only pass one ordering into CholeskyFactor */
2359b54502bSHong Zhang       dir->col=0;
2369bfd6278SHong Zhang 
2372401956bSBarry Smith       ierr = PetscOptionsHasName(pc->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr);
2389b54502bSHong Zhang       if (flg) {
2399b54502bSHong Zhang         PetscReal tol = 1.e-10;
2402401956bSBarry Smith         ierr = PetscOptionsGetReal(pc->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr);
2419b54502bSHong Zhang         ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr);
2429b54502bSHong Zhang       }
243efee365bSSatish Balay       if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
2449b54502bSHong Zhang       ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr);
2459b54502bSHong Zhang       ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
2469b54502bSHong Zhang       dir->actualfill = info.fill_ratio_needed;
24752e6d16bSBarry Smith       ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr);
2489b54502bSHong Zhang     } else if (pc->flag != SAME_NONZERO_PATTERN) {
2499b54502bSHong Zhang       if (!dir->reuseordering) {
2509b54502bSHong Zhang         if (dir->row && dir->col && (dir->row != dir->col)) {
2519b54502bSHong Zhang           ierr = ISDestroy(dir->row);CHKERRQ(ierr);
2529b54502bSHong Zhang           dir->row = 0;
2539b54502bSHong Zhang         }
2549b54502bSHong Zhang         if (dir->col) {
2559b54502bSHong Zhang           ierr = ISDestroy(dir->col);CHKERRQ(ierr);
2569b54502bSHong Zhang           dir->col =0;
2579b54502bSHong Zhang         }
2589b54502bSHong Zhang         ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr);
2599b54502bSHong Zhang         if (dir->col && (dir->row != dir->col)) {  /* only use row ordering for SBAIJ */
2609b54502bSHong Zhang           ierr = ISDestroy(dir->col);CHKERRQ(ierr);
2619b54502bSHong Zhang           dir->col=0;
2629b54502bSHong Zhang         }
2632401956bSBarry Smith         ierr = PetscOptionsHasName(pc->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr);
2649b54502bSHong Zhang         if (flg) {
2659b54502bSHong Zhang           PetscReal tol = 1.e-10;
2662401956bSBarry Smith           ierr = PetscOptionsGetReal(pc->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr);
2679b54502bSHong Zhang           ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr);
2689b54502bSHong Zhang         }
269efee365bSSatish Balay         if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);}
2709b54502bSHong Zhang       }
2719b54502bSHong Zhang       ierr = MatDestroy(dir->fact);CHKERRQ(ierr);
2729b54502bSHong Zhang       ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr);
2739b54502bSHong Zhang       ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr);
2749b54502bSHong Zhang       dir->actualfill = info.fill_ratio_needed;
27552e6d16bSBarry Smith       ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr);
2769b54502bSHong Zhang     }
2779b54502bSHong Zhang     ierr = MatCholeskyFactorNumeric(pc->pmat,&dir->info,&dir->fact);CHKERRQ(ierr);
2789b54502bSHong Zhang   }
2799b54502bSHong Zhang   PetscFunctionReturn(0);
2809b54502bSHong Zhang }
2819b54502bSHong Zhang 
2829b54502bSHong Zhang #undef __FUNCT__
2839b54502bSHong Zhang #define __FUNCT__ "PCDestroy_Cholesky"
2849b54502bSHong Zhang static PetscErrorCode PCDestroy_Cholesky(PC pc)
2859b54502bSHong Zhang {
2869b54502bSHong Zhang   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
2879b54502bSHong Zhang   PetscErrorCode ierr;
2889b54502bSHong Zhang 
2899b54502bSHong Zhang   PetscFunctionBegin;
2909b54502bSHong Zhang   if (!dir->inplace && dir->fact) {ierr = MatDestroy(dir->fact);CHKERRQ(ierr);}
2919b54502bSHong Zhang   if (dir->row) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);}
2929b54502bSHong Zhang   if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);}
2939b54502bSHong Zhang   ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr);
2949b54502bSHong Zhang   ierr = PetscFree(dir);CHKERRQ(ierr);
2959b54502bSHong Zhang   PetscFunctionReturn(0);
2969b54502bSHong Zhang }
2979b54502bSHong Zhang 
2989b54502bSHong Zhang #undef __FUNCT__
2999b54502bSHong Zhang #define __FUNCT__ "PCApply_Cholesky"
3009b54502bSHong Zhang static PetscErrorCode PCApply_Cholesky(PC pc,Vec x,Vec y)
3019b54502bSHong Zhang {
3029b54502bSHong Zhang   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
3039b54502bSHong Zhang   PetscErrorCode ierr;
3049b54502bSHong Zhang 
3059b54502bSHong Zhang   PetscFunctionBegin;
3069b54502bSHong Zhang   if (dir->inplace) {ierr = MatSolve(pc->pmat,x,y);CHKERRQ(ierr);}
3079b54502bSHong Zhang   else              {ierr = MatSolve(dir->fact,x,y);CHKERRQ(ierr);}
3089b54502bSHong Zhang   PetscFunctionReturn(0);
3099b54502bSHong Zhang }
3109b54502bSHong Zhang 
3119b54502bSHong Zhang #undef __FUNCT__
3129b54502bSHong Zhang #define __FUNCT__ "PCApplyTranspose_Cholesky"
3139b54502bSHong Zhang static PetscErrorCode PCApplyTranspose_Cholesky(PC pc,Vec x,Vec y)
3149b54502bSHong Zhang {
3159b54502bSHong Zhang   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
3169b54502bSHong Zhang   PetscErrorCode ierr;
3179b54502bSHong Zhang 
3189b54502bSHong Zhang   PetscFunctionBegin;
3199b54502bSHong Zhang   if (dir->inplace) {ierr = MatSolveTranspose(pc->pmat,x,y);CHKERRQ(ierr);}
3209b54502bSHong Zhang   else              {ierr = MatSolveTranspose(dir->fact,x,y);CHKERRQ(ierr);}
3219b54502bSHong Zhang   PetscFunctionReturn(0);
3229b54502bSHong Zhang }
3239b54502bSHong Zhang 
3249b54502bSHong Zhang /* -----------------------------------------------------------------------------------*/
3259b54502bSHong Zhang 
3269b54502bSHong Zhang EXTERN_C_BEGIN
3279b54502bSHong Zhang #undef __FUNCT__
32855ba2a51SBarry Smith #define __FUNCT__ "PCFactorSetFill_Cholesky"
32955ba2a51SBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Cholesky(PC pc,PetscReal fill)
3309b54502bSHong Zhang {
3319b54502bSHong Zhang   PC_Cholesky *dir;
3329b54502bSHong Zhang 
3339b54502bSHong Zhang   PetscFunctionBegin;
3349b54502bSHong Zhang   dir = (PC_Cholesky*)pc->data;
3359b54502bSHong Zhang   dir->info.fill = fill;
3369b54502bSHong Zhang   PetscFunctionReturn(0);
3379b54502bSHong Zhang }
3389b54502bSHong Zhang EXTERN_C_END
3399b54502bSHong Zhang 
3409b54502bSHong Zhang EXTERN_C_BEGIN
3419b54502bSHong Zhang #undef __FUNCT__
3422401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseInPlace_Cholesky"
3432401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_Cholesky(PC pc)
3449b54502bSHong Zhang {
3459b54502bSHong Zhang   PC_Cholesky *dir;
3469b54502bSHong Zhang 
3479b54502bSHong Zhang   PetscFunctionBegin;
3489b54502bSHong Zhang   dir = (PC_Cholesky*)pc->data;
3499b54502bSHong Zhang   dir->inplace = PETSC_TRUE;
3509b54502bSHong Zhang   PetscFunctionReturn(0);
3519b54502bSHong Zhang }
3529b54502bSHong Zhang EXTERN_C_END
3539b54502bSHong Zhang 
3549b54502bSHong Zhang EXTERN_C_BEGIN
3559b54502bSHong Zhang #undef __FUNCT__
356e5a9bf91SBarry Smith #define __FUNCT__ "PCFactorSetMatOrderingType_Cholesky"
357e5a9bf91SBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Cholesky(PC pc,MatOrderingType ordering)
3589b54502bSHong Zhang {
3599b54502bSHong Zhang   PC_Cholesky    *dir = (PC_Cholesky*)pc->data;
3609b54502bSHong Zhang   PetscErrorCode ierr;
3619b54502bSHong Zhang 
3629b54502bSHong Zhang   PetscFunctionBegin;
3639b54502bSHong Zhang   ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr);
3649b54502bSHong Zhang   ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr);
3659b54502bSHong Zhang   PetscFunctionReturn(0);
3669b54502bSHong Zhang }
3679b54502bSHong Zhang EXTERN_C_END
3689b54502bSHong Zhang 
3699b54502bSHong Zhang /* -----------------------------------------------------------------------------------*/
3709b54502bSHong Zhang 
3719b54502bSHong Zhang #undef __FUNCT__
3722401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering"
3739b54502bSHong Zhang /*@
3742401956bSBarry Smith    PCFactorSetReuseOrdering - When similar matrices are factored, this
3759b54502bSHong Zhang    causes the ordering computed in the first factor to be used for all
3769b54502bSHong Zhang    following factors.
3779b54502bSHong Zhang 
3789b54502bSHong Zhang    Collective on PC
3799b54502bSHong Zhang 
3809b54502bSHong Zhang    Input Parameters:
3819b54502bSHong Zhang +  pc - the preconditioner context
3829b54502bSHong Zhang -  flag - PETSC_TRUE to reuse else PETSC_FALSE
3839b54502bSHong Zhang 
3849b54502bSHong Zhang    Options Database Key:
3852401956bSBarry Smith .  -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering()
3869b54502bSHong Zhang 
3879b54502bSHong Zhang    Level: intermediate
3889b54502bSHong Zhang 
3899b54502bSHong Zhang .keywords: PC, levels, reordering, factorization, incomplete, LU
3909b54502bSHong Zhang 
3912401956bSBarry Smith .seealso: PCFactorSetReuseFill()
3929b54502bSHong Zhang @*/
3932401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering(PC pc,PetscTruth flag)
3949b54502bSHong Zhang {
3959b54502bSHong Zhang   PetscErrorCode ierr,(*f)(PC,PetscTruth);
3969b54502bSHong Zhang 
3979b54502bSHong Zhang   PetscFunctionBegin;
3989b54502bSHong Zhang   PetscValidHeaderSpecific(pc,PC_COOKIE,1);
3992401956bSBarry Smith   ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",(void (**)(void))&f);CHKERRQ(ierr);
4009b54502bSHong Zhang   if (f) {
4019b54502bSHong Zhang     ierr = (*f)(pc,flag);CHKERRQ(ierr);
4029b54502bSHong Zhang   }
4039b54502bSHong Zhang   PetscFunctionReturn(0);
4049b54502bSHong Zhang }
4059b54502bSHong Zhang 
4069b54502bSHong Zhang #undef __FUNCT__
4072401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill"
4089b54502bSHong Zhang /*@
4092401956bSBarry Smith    PCFactorSetReuseFill - When matrices with same different nonzero structure are factored,
4102401956bSBarry Smith    this causes later ones to use the fill ratio computed in the initial factorization.
4119b54502bSHong Zhang 
4129b54502bSHong Zhang    Collective on PC
4139b54502bSHong Zhang 
4149b54502bSHong Zhang    Input Parameters:
4159b54502bSHong Zhang +  pc - the preconditioner context
4169b54502bSHong Zhang -  flag - PETSC_TRUE to reuse else PETSC_FALSE
4179b54502bSHong Zhang 
4189b54502bSHong Zhang    Options Database Key:
4192401956bSBarry Smith .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
4209b54502bSHong Zhang 
4219b54502bSHong Zhang    Level: intermediate
4229b54502bSHong Zhang 
4239b54502bSHong Zhang .keywords: PC, levels, reordering, factorization, incomplete, Cholesky
4249b54502bSHong Zhang 
4252401956bSBarry Smith .seealso: PCFactorSetReuseOrdering()
4269b54502bSHong Zhang @*/
4272401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill(PC pc,PetscTruth flag)
4289b54502bSHong Zhang {
4299b54502bSHong Zhang   PetscErrorCode ierr,(*f)(PC,PetscTruth);
4309b54502bSHong Zhang 
4319b54502bSHong Zhang   PetscFunctionBegin;
4329b54502bSHong Zhang   PetscValidHeaderSpecific(pc,PC_COOKIE,2);
4332401956bSBarry Smith   ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseFill_C",(void (**)(void))&f);CHKERRQ(ierr);
4349b54502bSHong Zhang   if (f) {
4359b54502bSHong Zhang     ierr = (*f)(pc,flag);CHKERRQ(ierr);
4369b54502bSHong Zhang   }
4379b54502bSHong Zhang   PetscFunctionReturn(0);
4389b54502bSHong Zhang }
4399b54502bSHong Zhang 
4409b54502bSHong Zhang /*MC
44196fc60bcSBarry Smith    PCCHOLESKY - Uses a direct solver, based on Cholesky factorization, as a preconditioner
4429b54502bSHong Zhang 
4439b54502bSHong Zhang    Options Database Keys:
4442401956bSBarry Smith +  -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering()
4452401956bSBarry Smith .  -pc_factor_reuse_fill - Activates PCFactorSetReuseFill()
44655ba2a51SBarry Smith .  -pc_factor_fill <fill> - Sets fill amount
4472401956bSBarry Smith .  -pc_factor_in_place - Activates in-place factorization
4482401956bSBarry Smith .  -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine
449f251bdbdSHong Zhang .  -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default
450f251bdbdSHong Zhang -  -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value
451f251bdbdSHong Zhang    is optional with PETSC_TRUE being the default
4529b54502bSHong Zhang 
4539b54502bSHong Zhang    Notes: Not all options work for all matrix formats
4549b54502bSHong Zhang 
4559b54502bSHong Zhang    Level: beginner
4569b54502bSHong Zhang 
4579b54502bSHong Zhang    Concepts: Cholesky factorization, direct solver
4589b54502bSHong Zhang 
4599b54502bSHong Zhang    Notes: Usually this will compute an "exact" solution in one iteration and does
4609b54502bSHong Zhang           not need a Krylov method (i.e. you can use -ksp_type preonly, or
4619b54502bSHong Zhang           KSPSetType(ksp,KSPPREONLY) for the Krylov method
4629b54502bSHong Zhang 
4639b54502bSHong Zhang .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
464*a4fd02acSBarry Smith            PCILU, PCLU, PCICC, PCFactorSetReuseOrdering(), PCFactorSetReuseFill(), PCFactorGetMatrix(),
46555ba2a51SBarry Smith            PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(),
466e5a9bf91SBarry Smith 	   PCFactorSetUseInPlace(), PCFactorSetMatOrderingType(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd()
4679b54502bSHong Zhang 
4689b54502bSHong Zhang M*/
4699b54502bSHong Zhang 
4709b54502bSHong Zhang EXTERN_C_BEGIN
4719b54502bSHong Zhang #undef __FUNCT__
4729b54502bSHong Zhang #define __FUNCT__ "PCCreate_Cholesky"
473dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Cholesky(PC pc)
4749b54502bSHong Zhang {
4759b54502bSHong Zhang   PetscErrorCode ierr;
4769b54502bSHong Zhang   PC_Cholesky    *dir;
4779b54502bSHong Zhang 
4789b54502bSHong Zhang   PetscFunctionBegin;
47938f2d2fdSLisandro Dalcin   ierr = PetscNewLog(pc,PC_Cholesky,&dir);CHKERRQ(ierr);
4809b54502bSHong Zhang 
4819b54502bSHong Zhang   dir->fact                   = 0;
4829b54502bSHong Zhang   dir->inplace                = PETSC_FALSE;
4839b54502bSHong Zhang   ierr = MatFactorInfoInitialize(&dir->info);CHKERRQ(ierr);
4849b54502bSHong Zhang   dir->info.fill              = 5.0;
4850a29876aSHong Zhang   dir->info.shiftnz           = 0.0;
486fbf22428SSatish Balay   dir->info.shiftpd           = 0.0; /* false */
4879b54502bSHong Zhang   dir->info.shift_fraction    = 0.0;
4889b54502bSHong Zhang   dir->info.pivotinblocks     = 1.0;
4899b54502bSHong Zhang   dir->col                    = 0;
4909b54502bSHong Zhang   dir->row                    = 0;
4919b54502bSHong Zhang   ierr = PetscStrallocpy(MATORDERING_NATURAL,&dir->ordering);CHKERRQ(ierr);
4929b54502bSHong Zhang   dir->reusefill        = PETSC_FALSE;
4939b54502bSHong Zhang   dir->reuseordering    = PETSC_FALSE;
4949b54502bSHong Zhang   pc->data              = (void*)dir;
4959b54502bSHong Zhang 
4969b54502bSHong Zhang   pc->ops->destroy           = PCDestroy_Cholesky;
4979b54502bSHong Zhang   pc->ops->apply             = PCApply_Cholesky;
4989b54502bSHong Zhang   pc->ops->applytranspose    = PCApplyTranspose_Cholesky;
4999b54502bSHong Zhang   pc->ops->setup             = PCSetUp_Cholesky;
5009b54502bSHong Zhang   pc->ops->setfromoptions    = PCSetFromOptions_Cholesky;
5019b54502bSHong Zhang   pc->ops->view              = PCView_Cholesky;
5029b54502bSHong Zhang   pc->ops->applyrichardson   = 0;
503*a4fd02acSBarry Smith   pc->ops->getfactoredmatrix = PCFactorGetMatrix_Cholesky;
5049b54502bSHong Zhang 
505afaefe49SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Cholesky",
506afaefe49SHong Zhang                     PCFactorSetZeroPivot_Cholesky);CHKERRQ(ierr);
507afaefe49SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Cholesky",
508afaefe49SHong Zhang                     PCFactorSetShiftNonzero_Cholesky);CHKERRQ(ierr);
509afaefe49SHong Zhang   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Cholesky",
510afaefe49SHong Zhang                     PCFactorSetShiftPd_Cholesky);CHKERRQ(ierr);
511afaefe49SHong Zhang 
51255ba2a51SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Cholesky",
51355ba2a51SBarry Smith                     PCFactorSetFill_Cholesky);CHKERRQ(ierr);
5142401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_Cholesky",
5152401956bSBarry Smith                     PCFactorSetUseInPlace_Cholesky);CHKERRQ(ierr);
516e5a9bf91SBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Cholesky",
517e5a9bf91SBarry Smith                     PCFactorSetMatOrderingType_Cholesky);CHKERRQ(ierr);
5182401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_Cholesky",
5192401956bSBarry Smith                     PCFactorSetReuseOrdering_Cholesky);CHKERRQ(ierr);
5202401956bSBarry Smith   ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_Cholesky",
5212401956bSBarry Smith                     PCFactorSetReuseFill_Cholesky);CHKERRQ(ierr);
5229b54502bSHong Zhang   PetscFunctionReturn(0);
5239b54502bSHong Zhang }
5249b54502bSHong Zhang EXTERN_C_END
525