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; 19*c7393fdbSBarry Smith MatSolverPackage solvertype; 209b54502bSHong Zhang } PC_Cholesky; 219b54502bSHong Zhang 229b54502bSHong Zhang EXTERN_C_BEGIN 239b54502bSHong Zhang #undef __FUNCT__ 24*c7393fdbSBarry Smith #define __FUNCT__ "PCFactorSetMatSolverPackage_LU" 25*c7393fdbSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatSolverPackage_Cholesky(PC pc,const MatSolverPackage stype) 26c62fe420SBarry Smith { 27c62fe420SBarry Smith PetscErrorCode ierr; 28c62fe420SBarry Smith PC_Cholesky *choleksy = (PC_Cholesky*)pc->data; 29c62fe420SBarry Smith 30c62fe420SBarry Smith PetscFunctionBegin; 31c62fe420SBarry Smith ierr = PetscStrallocpy(stype,&choleksy->solvertype);CHKERRQ(ierr); 32c62fe420SBarry Smith PetscFunctionReturn(0); 33c62fe420SBarry Smith } 34c62fe420SBarry Smith EXTERN_C_END 35c62fe420SBarry Smith 36c62fe420SBarry Smith EXTERN_C_BEGIN 37c62fe420SBarry Smith #undef __FUNCT__ 38afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetZeroPivot_Cholesky" 39dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_Cholesky(PC pc,PetscReal z) 40afaefe49SHong Zhang { 41afaefe49SHong Zhang PC_Cholesky *ch; 42afaefe49SHong Zhang 43afaefe49SHong Zhang PetscFunctionBegin; 44afaefe49SHong Zhang ch = (PC_Cholesky*)pc->data; 45afaefe49SHong Zhang ch->info.zeropivot = z; 46afaefe49SHong Zhang PetscFunctionReturn(0); 47afaefe49SHong Zhang } 48afaefe49SHong Zhang EXTERN_C_END 49afaefe49SHong Zhang 50afaefe49SHong Zhang EXTERN_C_BEGIN 51afaefe49SHong Zhang #undef __FUNCT__ 52afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetShiftNonzero_Cholesky" 53dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_Cholesky(PC pc,PetscReal shift) 54afaefe49SHong Zhang { 55afaefe49SHong Zhang PC_Cholesky *dir; 56afaefe49SHong Zhang 57afaefe49SHong Zhang PetscFunctionBegin; 58afaefe49SHong Zhang dir = (PC_Cholesky*)pc->data; 59afaefe49SHong Zhang if (shift == (PetscReal) PETSC_DECIDE) { 60afaefe49SHong Zhang dir->info.shiftnz = 1.e-12; 61afaefe49SHong Zhang } else { 62afaefe49SHong Zhang dir->info.shiftnz = shift; 63afaefe49SHong Zhang } 64afaefe49SHong Zhang PetscFunctionReturn(0); 65afaefe49SHong Zhang } 66afaefe49SHong Zhang EXTERN_C_END 67afaefe49SHong Zhang 68afaefe49SHong Zhang EXTERN_C_BEGIN 69afaefe49SHong Zhang #undef __FUNCT__ 70afaefe49SHong Zhang #define __FUNCT__ "PCFactorSetShiftPd_Cholesky" 71dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_Cholesky(PC pc,PetscTruth shift) 72afaefe49SHong Zhang { 73afaefe49SHong Zhang PC_Cholesky *dir; 74afaefe49SHong Zhang 75afaefe49SHong Zhang PetscFunctionBegin; 76afaefe49SHong Zhang dir = (PC_Cholesky*)pc->data; 77fbf22428SSatish Balay if (shift) { 78fbf22428SSatish Balay dir->info.shift_fraction = 0.0; 79fbf22428SSatish Balay dir->info.shiftpd = 1.0; 80fbf22428SSatish Balay } else { 81fbf22428SSatish Balay dir->info.shiftpd = 0.0; 82fbf22428SSatish Balay } 83afaefe49SHong Zhang PetscFunctionReturn(0); 84afaefe49SHong Zhang } 85afaefe49SHong Zhang EXTERN_C_END 86afaefe49SHong Zhang 87afaefe49SHong Zhang EXTERN_C_BEGIN 88afaefe49SHong Zhang #undef __FUNCT__ 892401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering_Cholesky" 902401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_Cholesky(PC pc,PetscTruth flag) 919b54502bSHong Zhang { 929b54502bSHong Zhang PC_Cholesky *lu; 939b54502bSHong Zhang 949b54502bSHong Zhang PetscFunctionBegin; 959b54502bSHong Zhang lu = (PC_Cholesky*)pc->data; 969b54502bSHong Zhang lu->reuseordering = flag; 979b54502bSHong Zhang PetscFunctionReturn(0); 989b54502bSHong Zhang } 999b54502bSHong Zhang EXTERN_C_END 1009b54502bSHong Zhang 1019b54502bSHong Zhang EXTERN_C_BEGIN 1029b54502bSHong Zhang #undef __FUNCT__ 1032401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill_Cholesky" 1042401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_Cholesky(PC pc,PetscTruth flag) 1059b54502bSHong Zhang { 1069b54502bSHong Zhang PC_Cholesky *lu; 1079b54502bSHong Zhang 1089b54502bSHong Zhang PetscFunctionBegin; 1099b54502bSHong Zhang lu = (PC_Cholesky*)pc->data; 1109b54502bSHong Zhang lu->reusefill = flag; 1119b54502bSHong Zhang PetscFunctionReturn(0); 1129b54502bSHong Zhang } 1139b54502bSHong Zhang EXTERN_C_END 1149b54502bSHong Zhang 1159b54502bSHong Zhang #undef __FUNCT__ 1169b54502bSHong Zhang #define __FUNCT__ "PCSetFromOptions_Cholesky" 1179b54502bSHong Zhang static PetscErrorCode PCSetFromOptions_Cholesky(PC pc) 1189b54502bSHong Zhang { 1199b54502bSHong Zhang PC_Cholesky *lu = (PC_Cholesky*)pc->data; 1209b54502bSHong Zhang PetscErrorCode ierr; 1219b54502bSHong Zhang PetscTruth flg; 1225c9eb25fSBarry Smith char tname[256], solvertype[64]; 1239b54502bSHong Zhang PetscFList ordlist; 1249b54502bSHong Zhang 1259b54502bSHong Zhang PetscFunctionBegin; 1269b54502bSHong Zhang ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr); 1279b54502bSHong Zhang ierr = PetscOptionsHead("Cholesky options");CHKERRQ(ierr); 1282401956bSBarry Smith ierr = PetscOptionsName("-pc_factor_in_place","Form Cholesky in the same memory as the matrix","PCFactorSetUseInPlace",&flg);CHKERRQ(ierr); 1299b54502bSHong Zhang if (flg) { 1302401956bSBarry Smith ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr); 1319b54502bSHong Zhang } 13255ba2a51SBarry 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); 1339b54502bSHong Zhang 1342401956bSBarry Smith ierr = PetscOptionsName("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",&flg);CHKERRQ(ierr); 1359b54502bSHong Zhang if (flg) { 1362401956bSBarry Smith ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr); 1379b54502bSHong Zhang } 1382401956bSBarry Smith ierr = PetscOptionsName("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",&flg);CHKERRQ(ierr); 1399b54502bSHong Zhang if (flg) { 1402401956bSBarry Smith ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr); 1419b54502bSHong Zhang } 1429b54502bSHong Zhang 1439b54502bSHong Zhang ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 144e5a9bf91SBarry Smith ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in Cholesky","PCFactorSetMatOrderingType",ordlist,lu->ordering,tname,256,&flg);CHKERRQ(ierr); 1459b54502bSHong Zhang if (flg) { 146e5a9bf91SBarry Smith ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr); 1479b54502bSHong Zhang } 1485c9eb25fSBarry Smith 1495c9eb25fSBarry Smith /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */ 150*c7393fdbSBarry Smith ierr = PetscOptionsString("-pc_factor_mat_solver_package","Specific Cholesky solver to use","MatGetFactor",lu->solvertype,solvertype,64,&flg);CHKERRQ(ierr); 1515c9eb25fSBarry Smith if (flg) { 152*c7393fdbSBarry Smith ierr = PCFactorSetMatSolverPackage(pc,solvertype);CHKERRQ(ierr); 1535c9eb25fSBarry Smith } 1545c9eb25fSBarry Smith 1559f95998fSHong Zhang ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 1569b54502bSHong Zhang if (flg) { 157afaefe49SHong Zhang ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr); 1589b54502bSHong Zhang } 1599f95998fSHong Zhang ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",lu->info.shiftnz,&lu->info.shiftnz,0);CHKERRQ(ierr); 1609f95998fSHong Zhang ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",&flg);CHKERRQ(ierr); 1619b54502bSHong Zhang if (flg) { 162afaefe49SHong Zhang ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 1639b54502bSHong Zhang } 164ee45ca4aSHong Zhang ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",lu->info.zeropivot,&lu->info.zeropivot,0);CHKERRQ(ierr); 1659b54502bSHong Zhang 1669b54502bSHong Zhang ierr = PetscOptionsTail();CHKERRQ(ierr); 1679b54502bSHong Zhang PetscFunctionReturn(0); 1689b54502bSHong Zhang } 1699b54502bSHong Zhang 1709b54502bSHong Zhang #undef __FUNCT__ 1719b54502bSHong Zhang #define __FUNCT__ "PCView_Cholesky" 1729b54502bSHong Zhang static PetscErrorCode PCView_Cholesky(PC pc,PetscViewer viewer) 1739b54502bSHong Zhang { 1749b54502bSHong Zhang PC_Cholesky *lu = (PC_Cholesky*)pc->data; 1759b54502bSHong Zhang PetscErrorCode ierr; 1769b54502bSHong Zhang PetscTruth iascii,isstring; 1779b54502bSHong Zhang 1789b54502bSHong Zhang PetscFunctionBegin; 1799b54502bSHong Zhang ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 1809b54502bSHong Zhang ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 1819b54502bSHong Zhang if (iascii) { 1829b54502bSHong Zhang 1839b54502bSHong Zhang if (lu->inplace) {ierr = PetscViewerASCIIPrintf(viewer," Cholesky: in-place factorization\n");CHKERRQ(ierr);} 1849b54502bSHong Zhang else {ierr = PetscViewerASCIIPrintf(viewer," Cholesky: out-of-place factorization\n");CHKERRQ(ierr);} 1859b54502bSHong Zhang ierr = PetscViewerASCIIPrintf(viewer," matrix ordering: %s\n",lu->ordering);CHKERRQ(ierr); 1869b54502bSHong Zhang if (lu->fact) { 187f3a39becSBarry Smith ierr = PetscViewerASCIIPrintf(viewer," Cholesky: factor fill ratio needed %G\n",lu->actualfill);CHKERRQ(ierr); 188f3a39becSBarry Smith ierr = PetscViewerASCIIPrintf(viewer," Factored matrix follows\n");CHKERRQ(ierr); 189f3a39becSBarry Smith ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 190f3a39becSBarry Smith ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 191f3a39becSBarry Smith ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 192f3a39becSBarry Smith ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr); 1939b54502bSHong Zhang ierr = MatView(lu->fact,viewer);CHKERRQ(ierr); 1949b54502bSHong Zhang ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 195f3a39becSBarry Smith ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 196f3a39becSBarry Smith ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 197f3a39becSBarry Smith ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1989b54502bSHong Zhang } 1999b54502bSHong Zhang if (lu->reusefill) {ierr = PetscViewerASCIIPrintf(viewer," Reusing fill from past factorization\n");CHKERRQ(ierr);} 2009b54502bSHong Zhang if (lu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer," Reusing reordering from past factorization\n");CHKERRQ(ierr);} 2019b54502bSHong Zhang } else if (isstring) { 2029b54502bSHong Zhang ierr = PetscViewerStringSPrintf(viewer," order=%s",lu->ordering);CHKERRQ(ierr);CHKERRQ(ierr); 2039b54502bSHong Zhang } else { 2042401956bSBarry Smith SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCCHOLESKY",((PetscObject)viewer)->type_name); 2059b54502bSHong Zhang } 2069b54502bSHong Zhang PetscFunctionReturn(0); 2079b54502bSHong Zhang } 2089b54502bSHong Zhang 2099b54502bSHong Zhang #undef __FUNCT__ 210a4fd02acSBarry Smith #define __FUNCT__ "PCFactorGetMatrix_Cholesky" 211a4fd02acSBarry Smith static PetscErrorCode PCFactorGetMatrix_Cholesky(PC pc,Mat *mat) 2129b54502bSHong Zhang { 2139b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 2149b54502bSHong Zhang 2159b54502bSHong Zhang PetscFunctionBegin; 2169b54502bSHong Zhang if (!dir->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()"); 2179b54502bSHong Zhang *mat = dir->fact; 2189b54502bSHong Zhang PetscFunctionReturn(0); 2199b54502bSHong Zhang } 2209b54502bSHong Zhang 2219b54502bSHong Zhang #undef __FUNCT__ 2229b54502bSHong Zhang #define __FUNCT__ "PCSetUp_Cholesky" 2239b54502bSHong Zhang static PetscErrorCode PCSetUp_Cholesky(PC pc) 2249b54502bSHong Zhang { 2259b54502bSHong Zhang PetscErrorCode ierr; 2269b54502bSHong Zhang PetscTruth flg; 2279b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 2289b54502bSHong Zhang 2299b54502bSHong Zhang PetscFunctionBegin; 2309b54502bSHong Zhang if (dir->reusefill && pc->setupcalled) dir->info.fill = dir->actualfill; 2319b54502bSHong Zhang 2329b54502bSHong Zhang if (dir->inplace) { 2339b54502bSHong Zhang if (dir->row && dir->col && (dir->row != dir->col)) { 2349b54502bSHong Zhang ierr = ISDestroy(dir->row);CHKERRQ(ierr); 2359b54502bSHong Zhang dir->row = 0; 2369b54502bSHong Zhang } 2379b54502bSHong Zhang if (dir->col) { 2389b54502bSHong Zhang ierr = ISDestroy(dir->col);CHKERRQ(ierr); 2399b54502bSHong Zhang dir->col = 0; 2409b54502bSHong Zhang } 2419b54502bSHong Zhang ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 2429b54502bSHong Zhang if (dir->col && (dir->row != dir->col)) { /* only use row ordering for SBAIJ */ 2439b54502bSHong Zhang ierr = ISDestroy(dir->col);CHKERRQ(ierr); 2449b54502bSHong Zhang dir->col=0; 2459b54502bSHong Zhang } 246efee365bSSatish Balay if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 2479b54502bSHong Zhang ierr = MatCholeskyFactor(pc->pmat,dir->row,&dir->info);CHKERRQ(ierr); 2489b54502bSHong Zhang dir->fact = pc->pmat; 2499b54502bSHong Zhang } else { 2509b54502bSHong Zhang MatInfo info; 2519b54502bSHong Zhang if (!pc->setupcalled) { 2529b54502bSHong Zhang ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 2539bfd6278SHong Zhang /* check if dir->row == dir->col */ 2549bfd6278SHong Zhang ierr = ISEqual(dir->row,dir->col,&flg);CHKERRQ(ierr); 2559bfd6278SHong Zhang if (!flg) SETERRQ(PETSC_ERR_ARG_INCOMP,"row and column permutations must equal"); 2569bfd6278SHong Zhang ierr = ISDestroy(dir->col);CHKERRQ(ierr); /* only pass one ordering into CholeskyFactor */ 2579b54502bSHong Zhang dir->col=0; 2589bfd6278SHong Zhang 2597adad957SLisandro Dalcin ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr); 2609b54502bSHong Zhang if (flg) { 2619b54502bSHong Zhang PetscReal tol = 1.e-10; 2627adad957SLisandro Dalcin ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr); 2639b54502bSHong Zhang ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr); 2649b54502bSHong Zhang } 265efee365bSSatish Balay if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 2665c9eb25fSBarry Smith ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr); 2679b54502bSHong Zhang ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr); 2689b54502bSHong Zhang ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 2699b54502bSHong Zhang dir->actualfill = info.fill_ratio_needed; 27052e6d16bSBarry Smith ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr); 2719b54502bSHong Zhang } else if (pc->flag != SAME_NONZERO_PATTERN) { 2729b54502bSHong Zhang if (!dir->reuseordering) { 2739b54502bSHong Zhang if (dir->row && dir->col && (dir->row != dir->col)) { 2749b54502bSHong Zhang ierr = ISDestroy(dir->row);CHKERRQ(ierr); 2759b54502bSHong Zhang dir->row = 0; 2769b54502bSHong Zhang } 2779b54502bSHong Zhang if (dir->col) { 2789b54502bSHong Zhang ierr = ISDestroy(dir->col);CHKERRQ(ierr); 2799b54502bSHong Zhang dir->col =0; 2809b54502bSHong Zhang } 2819b54502bSHong Zhang ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 2829b54502bSHong Zhang if (dir->col && (dir->row != dir->col)) { /* only use row ordering for SBAIJ */ 2839b54502bSHong Zhang ierr = ISDestroy(dir->col);CHKERRQ(ierr); 2849b54502bSHong Zhang dir->col=0; 2859b54502bSHong Zhang } 2867adad957SLisandro Dalcin ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr); 2879b54502bSHong Zhang if (flg) { 2889b54502bSHong Zhang PetscReal tol = 1.e-10; 2897adad957SLisandro Dalcin ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr); 2909b54502bSHong Zhang ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr); 2919b54502bSHong Zhang } 292efee365bSSatish Balay if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 2939b54502bSHong Zhang } 2949b54502bSHong Zhang ierr = MatDestroy(dir->fact);CHKERRQ(ierr); 2955c9eb25fSBarry Smith ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr); 2969b54502bSHong Zhang ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr); 2979b54502bSHong Zhang ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 2989b54502bSHong Zhang dir->actualfill = info.fill_ratio_needed; 29952e6d16bSBarry Smith ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr); 3009b54502bSHong Zhang } 3019b54502bSHong Zhang ierr = MatCholeskyFactorNumeric(pc->pmat,&dir->info,&dir->fact);CHKERRQ(ierr); 3029b54502bSHong Zhang } 3039b54502bSHong Zhang PetscFunctionReturn(0); 3049b54502bSHong Zhang } 3059b54502bSHong Zhang 3069b54502bSHong Zhang #undef __FUNCT__ 3079b54502bSHong Zhang #define __FUNCT__ "PCDestroy_Cholesky" 3089b54502bSHong Zhang static PetscErrorCode PCDestroy_Cholesky(PC pc) 3099b54502bSHong Zhang { 3109b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 3119b54502bSHong Zhang PetscErrorCode ierr; 3129b54502bSHong Zhang 3139b54502bSHong Zhang PetscFunctionBegin; 3149b54502bSHong Zhang if (!dir->inplace && dir->fact) {ierr = MatDestroy(dir->fact);CHKERRQ(ierr);} 3159b54502bSHong Zhang if (dir->row) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);} 3169b54502bSHong Zhang if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);} 3179b54502bSHong Zhang ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 3185c9eb25fSBarry Smith ierr = PetscStrfree(dir->solvertype);CHKERRQ(ierr); 3199b54502bSHong Zhang ierr = PetscFree(dir);CHKERRQ(ierr); 3209b54502bSHong Zhang PetscFunctionReturn(0); 3219b54502bSHong Zhang } 3229b54502bSHong Zhang 3239b54502bSHong Zhang #undef __FUNCT__ 3249b54502bSHong Zhang #define __FUNCT__ "PCApply_Cholesky" 3259b54502bSHong Zhang static PetscErrorCode PCApply_Cholesky(PC pc,Vec x,Vec y) 3269b54502bSHong Zhang { 3279b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 3289b54502bSHong Zhang PetscErrorCode ierr; 3299b54502bSHong Zhang 3309b54502bSHong Zhang PetscFunctionBegin; 3319b54502bSHong Zhang if (dir->inplace) {ierr = MatSolve(pc->pmat,x,y);CHKERRQ(ierr);} 3329b54502bSHong Zhang else {ierr = MatSolve(dir->fact,x,y);CHKERRQ(ierr);} 3339b54502bSHong Zhang PetscFunctionReturn(0); 3349b54502bSHong Zhang } 3359b54502bSHong Zhang 3369b54502bSHong Zhang #undef __FUNCT__ 3379b54502bSHong Zhang #define __FUNCT__ "PCApplyTranspose_Cholesky" 3389b54502bSHong Zhang static PetscErrorCode PCApplyTranspose_Cholesky(PC pc,Vec x,Vec y) 3399b54502bSHong Zhang { 3409b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 3419b54502bSHong Zhang PetscErrorCode ierr; 3429b54502bSHong Zhang 3439b54502bSHong Zhang PetscFunctionBegin; 3449b54502bSHong Zhang if (dir->inplace) {ierr = MatSolveTranspose(pc->pmat,x,y);CHKERRQ(ierr);} 3459b54502bSHong Zhang else {ierr = MatSolveTranspose(dir->fact,x,y);CHKERRQ(ierr);} 3469b54502bSHong Zhang PetscFunctionReturn(0); 3479b54502bSHong Zhang } 3489b54502bSHong Zhang 3499b54502bSHong Zhang /* -----------------------------------------------------------------------------------*/ 3509b54502bSHong Zhang 3519b54502bSHong Zhang EXTERN_C_BEGIN 3529b54502bSHong Zhang #undef __FUNCT__ 35355ba2a51SBarry Smith #define __FUNCT__ "PCFactorSetFill_Cholesky" 35455ba2a51SBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Cholesky(PC pc,PetscReal fill) 3559b54502bSHong Zhang { 3569b54502bSHong Zhang PC_Cholesky *dir; 3579b54502bSHong Zhang 3589b54502bSHong Zhang PetscFunctionBegin; 3599b54502bSHong Zhang dir = (PC_Cholesky*)pc->data; 3609b54502bSHong Zhang dir->info.fill = fill; 3619b54502bSHong Zhang PetscFunctionReturn(0); 3629b54502bSHong Zhang } 3639b54502bSHong Zhang EXTERN_C_END 3649b54502bSHong Zhang 3659b54502bSHong Zhang EXTERN_C_BEGIN 3669b54502bSHong Zhang #undef __FUNCT__ 3672401956bSBarry Smith #define __FUNCT__ "PCFactorSetUseInPlace_Cholesky" 3682401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_Cholesky(PC pc) 3699b54502bSHong Zhang { 3709b54502bSHong Zhang PC_Cholesky *dir; 3719b54502bSHong Zhang 3729b54502bSHong Zhang PetscFunctionBegin; 3739b54502bSHong Zhang dir = (PC_Cholesky*)pc->data; 3749b54502bSHong Zhang dir->inplace = PETSC_TRUE; 3759b54502bSHong Zhang PetscFunctionReturn(0); 3769b54502bSHong Zhang } 3779b54502bSHong Zhang EXTERN_C_END 3789b54502bSHong Zhang 3799b54502bSHong Zhang EXTERN_C_BEGIN 3809b54502bSHong Zhang #undef __FUNCT__ 381e5a9bf91SBarry Smith #define __FUNCT__ "PCFactorSetMatOrderingType_Cholesky" 382e36ff184SSatish Balay PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Cholesky(PC pc,const MatOrderingType ordering) 3839b54502bSHong Zhang { 3849b54502bSHong Zhang PC_Cholesky *dir = (PC_Cholesky*)pc->data; 3859b54502bSHong Zhang PetscErrorCode ierr; 3869b54502bSHong Zhang 3879b54502bSHong Zhang PetscFunctionBegin; 3889b54502bSHong Zhang ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 3899b54502bSHong Zhang ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr); 3909b54502bSHong Zhang PetscFunctionReturn(0); 3919b54502bSHong Zhang } 3929b54502bSHong Zhang EXTERN_C_END 3939b54502bSHong Zhang 3949b54502bSHong Zhang /* -----------------------------------------------------------------------------------*/ 3959b54502bSHong Zhang 3969b54502bSHong Zhang #undef __FUNCT__ 3972401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseOrdering" 3989b54502bSHong Zhang /*@ 3992401956bSBarry Smith PCFactorSetReuseOrdering - When similar matrices are factored, this 4009b54502bSHong Zhang causes the ordering computed in the first factor to be used for all 4019b54502bSHong Zhang following factors. 4029b54502bSHong Zhang 4039b54502bSHong Zhang Collective on PC 4049b54502bSHong Zhang 4059b54502bSHong Zhang Input Parameters: 4069b54502bSHong Zhang + pc - the preconditioner context 4079b54502bSHong Zhang - flag - PETSC_TRUE to reuse else PETSC_FALSE 4089b54502bSHong Zhang 4099b54502bSHong Zhang Options Database Key: 4102401956bSBarry Smith . -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering() 4119b54502bSHong Zhang 4129b54502bSHong Zhang Level: intermediate 4139b54502bSHong Zhang 4149b54502bSHong Zhang .keywords: PC, levels, reordering, factorization, incomplete, LU 4159b54502bSHong Zhang 4162401956bSBarry Smith .seealso: PCFactorSetReuseFill() 4179b54502bSHong Zhang @*/ 4182401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering(PC pc,PetscTruth flag) 4199b54502bSHong Zhang { 4209b54502bSHong Zhang PetscErrorCode ierr,(*f)(PC,PetscTruth); 4219b54502bSHong Zhang 4229b54502bSHong Zhang PetscFunctionBegin; 4239b54502bSHong Zhang PetscValidHeaderSpecific(pc,PC_COOKIE,1); 4242401956bSBarry Smith ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",(void (**)(void))&f);CHKERRQ(ierr); 4259b54502bSHong Zhang if (f) { 4269b54502bSHong Zhang ierr = (*f)(pc,flag);CHKERRQ(ierr); 4279b54502bSHong Zhang } 4289b54502bSHong Zhang PetscFunctionReturn(0); 4299b54502bSHong Zhang } 4309b54502bSHong Zhang 4319b54502bSHong Zhang #undef __FUNCT__ 4322401956bSBarry Smith #define __FUNCT__ "PCFactorSetReuseFill" 4339b54502bSHong Zhang /*@ 4342401956bSBarry Smith PCFactorSetReuseFill - When matrices with same different nonzero structure are factored, 4352401956bSBarry Smith this causes later ones to use the fill ratio computed in the initial factorization. 4369b54502bSHong Zhang 4379b54502bSHong Zhang Collective on PC 4389b54502bSHong Zhang 4399b54502bSHong Zhang Input Parameters: 4409b54502bSHong Zhang + pc - the preconditioner context 4419b54502bSHong Zhang - flag - PETSC_TRUE to reuse else PETSC_FALSE 4429b54502bSHong Zhang 4439b54502bSHong Zhang Options Database Key: 4442401956bSBarry Smith . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 4459b54502bSHong Zhang 4469b54502bSHong Zhang Level: intermediate 4479b54502bSHong Zhang 4489b54502bSHong Zhang .keywords: PC, levels, reordering, factorization, incomplete, Cholesky 4499b54502bSHong Zhang 4502401956bSBarry Smith .seealso: PCFactorSetReuseOrdering() 4519b54502bSHong Zhang @*/ 4522401956bSBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill(PC pc,PetscTruth flag) 4539b54502bSHong Zhang { 4549b54502bSHong Zhang PetscErrorCode ierr,(*f)(PC,PetscTruth); 4559b54502bSHong Zhang 4569b54502bSHong Zhang PetscFunctionBegin; 4579b54502bSHong Zhang PetscValidHeaderSpecific(pc,PC_COOKIE,2); 4582401956bSBarry Smith ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseFill_C",(void (**)(void))&f);CHKERRQ(ierr); 4599b54502bSHong Zhang if (f) { 4609b54502bSHong Zhang ierr = (*f)(pc,flag);CHKERRQ(ierr); 4619b54502bSHong Zhang } 4629b54502bSHong Zhang PetscFunctionReturn(0); 4639b54502bSHong Zhang } 4649b54502bSHong Zhang 4659b54502bSHong Zhang /*MC 46696fc60bcSBarry Smith PCCHOLESKY - Uses a direct solver, based on Cholesky factorization, as a preconditioner 4679b54502bSHong Zhang 4689b54502bSHong Zhang Options Database Keys: 4692401956bSBarry Smith + -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering() 470*c7393fdbSBarry Smith . -pc_factor_mat_solver_package - Actives PCFactorSetMatSolverPackage() to choose the direct solver, like spooles 4712401956bSBarry Smith . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 47255ba2a51SBarry Smith . -pc_factor_fill <fill> - Sets fill amount 4732401956bSBarry Smith . -pc_factor_in_place - Activates in-place factorization 4742401956bSBarry Smith . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine 475f251bdbdSHong Zhang . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 476f251bdbdSHong Zhang - -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 477f251bdbdSHong Zhang is optional with PETSC_TRUE being the default 4789b54502bSHong Zhang 4799b54502bSHong Zhang Notes: Not all options work for all matrix formats 4809b54502bSHong Zhang 4819b54502bSHong Zhang Level: beginner 4829b54502bSHong Zhang 4839b54502bSHong Zhang Concepts: Cholesky factorization, direct solver 4849b54502bSHong Zhang 4859b54502bSHong Zhang Notes: Usually this will compute an "exact" solution in one iteration and does 4869b54502bSHong Zhang not need a Krylov method (i.e. you can use -ksp_type preonly, or 4879b54502bSHong Zhang KSPSetType(ksp,KSPPREONLY) for the Krylov method 4889b54502bSHong Zhang 4899b54502bSHong Zhang .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 490a4fd02acSBarry Smith PCILU, PCLU, PCICC, PCFactorSetReuseOrdering(), PCFactorSetReuseFill(), PCFactorGetMatrix(), 49155ba2a51SBarry Smith PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(), 492e5a9bf91SBarry Smith PCFactorSetUseInPlace(), PCFactorSetMatOrderingType(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd() 4939b54502bSHong Zhang 4949b54502bSHong Zhang M*/ 4959b54502bSHong Zhang 4969b54502bSHong Zhang EXTERN_C_BEGIN 4979b54502bSHong Zhang #undef __FUNCT__ 4989b54502bSHong Zhang #define __FUNCT__ "PCCreate_Cholesky" 499dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Cholesky(PC pc) 5009b54502bSHong Zhang { 5019b54502bSHong Zhang PetscErrorCode ierr; 5029b54502bSHong Zhang PC_Cholesky *dir; 5039b54502bSHong Zhang 5049b54502bSHong Zhang PetscFunctionBegin; 50538f2d2fdSLisandro Dalcin ierr = PetscNewLog(pc,PC_Cholesky,&dir);CHKERRQ(ierr); 5069b54502bSHong Zhang 5079b54502bSHong Zhang dir->fact = 0; 5089b54502bSHong Zhang dir->inplace = PETSC_FALSE; 5099b54502bSHong Zhang ierr = MatFactorInfoInitialize(&dir->info);CHKERRQ(ierr); 5109b54502bSHong Zhang dir->info.fill = 5.0; 5110a29876aSHong Zhang dir->info.shiftnz = 0.0; 512fbf22428SSatish Balay dir->info.shiftpd = 0.0; /* false */ 5139b54502bSHong Zhang dir->info.shift_fraction = 0.0; 5149b54502bSHong Zhang dir->info.pivotinblocks = 1.0; 5159b54502bSHong Zhang dir->col = 0; 5169b54502bSHong Zhang dir->row = 0; 5179b54502bSHong Zhang ierr = PetscStrallocpy(MATORDERING_NATURAL,&dir->ordering);CHKERRQ(ierr); 5185c9eb25fSBarry Smith ierr = PetscStrallocpy("petsc",&dir->solvertype);CHKERRQ(ierr); 5199b54502bSHong Zhang dir->reusefill = PETSC_FALSE; 5209b54502bSHong Zhang dir->reuseordering = PETSC_FALSE; 5219b54502bSHong Zhang pc->data = (void*)dir; 5229b54502bSHong Zhang 5239b54502bSHong Zhang pc->ops->destroy = PCDestroy_Cholesky; 5249b54502bSHong Zhang pc->ops->apply = PCApply_Cholesky; 5259b54502bSHong Zhang pc->ops->applytranspose = PCApplyTranspose_Cholesky; 5269b54502bSHong Zhang pc->ops->setup = PCSetUp_Cholesky; 5279b54502bSHong Zhang pc->ops->setfromoptions = PCSetFromOptions_Cholesky; 5289b54502bSHong Zhang pc->ops->view = PCView_Cholesky; 5299b54502bSHong Zhang pc->ops->applyrichardson = 0; 530a4fd02acSBarry Smith pc->ops->getfactoredmatrix = PCFactorGetMatrix_Cholesky; 5319b54502bSHong Zhang 532*c7393fdbSBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatSolverPackage_C","PCFactorSetMatSolverPackage_Cholesky", 533*c7393fdbSBarry Smith PCFactorSetMatSolverPackage_Cholesky);CHKERRQ(ierr); 534afaefe49SHong Zhang ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Cholesky", 535afaefe49SHong Zhang PCFactorSetZeroPivot_Cholesky);CHKERRQ(ierr); 536afaefe49SHong Zhang ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Cholesky", 537afaefe49SHong Zhang PCFactorSetShiftNonzero_Cholesky);CHKERRQ(ierr); 538afaefe49SHong Zhang ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Cholesky", 539afaefe49SHong Zhang PCFactorSetShiftPd_Cholesky);CHKERRQ(ierr); 540afaefe49SHong Zhang 54155ba2a51SBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Cholesky", 54255ba2a51SBarry Smith PCFactorSetFill_Cholesky);CHKERRQ(ierr); 5432401956bSBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_Cholesky", 5442401956bSBarry Smith PCFactorSetUseInPlace_Cholesky);CHKERRQ(ierr); 545e5a9bf91SBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Cholesky", 546e5a9bf91SBarry Smith PCFactorSetMatOrderingType_Cholesky);CHKERRQ(ierr); 5472401956bSBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_Cholesky", 5482401956bSBarry Smith PCFactorSetReuseOrdering_Cholesky);CHKERRQ(ierr); 5492401956bSBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_Cholesky", 5502401956bSBarry Smith PCFactorSetReuseFill_Cholesky);CHKERRQ(ierr); 5519b54502bSHong Zhang PetscFunctionReturn(0); 5529b54502bSHong Zhang } 5539b54502bSHong Zhang EXTERN_C_END 554