1 #define PETSCKSP_DLL 2 3 #include "../src/ksp/pc/impls/factor/factor.h" /*I "petscpc.h" I*/ 4 5 /* ------------------------------------------------------------------------------------------*/ 6 7 EXTERN_C_BEGIN 8 #undef __FUNCT__ 9 #define __FUNCT__ "PCFactorSetZeroPivot_Factor" 10 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_Factor(PC pc,PetscReal z) 11 { 12 PC_Factor *ilu = (PC_Factor*)pc->data; 13 14 PetscFunctionBegin; 15 ilu->info.zeropivot = z; 16 PetscFunctionReturn(0); 17 } 18 EXTERN_C_END 19 20 EXTERN_C_BEGIN 21 #undef __FUNCT__ 22 #define __FUNCT__ "PCFactorSetShiftType_Factor" 23 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftType_Factor(PC pc,MatFactorShiftType shifttype) 24 { 25 PC_Factor *dir = (PC_Factor*)pc->data; 26 27 PetscFunctionBegin; 28 if (shifttype == (MatFactorShiftType)PETSC_DECIDE){ 29 dir->info.shifttype = (PetscReal) MAT_SHIFT_NONE; 30 } else { 31 dir->info.shifttype = (PetscReal) shifttype; 32 if (shifttype == MAT_SHIFT_NONZERO && dir->info.shiftamount == 0.0){ 33 dir->info.shiftamount = 1.e-12; /* set default amount if user has not called PCFactorSetShiftAmount() yet */ 34 } 35 } 36 PetscFunctionReturn(0); 37 } 38 39 #undef __FUNCT__ 40 #define __FUNCT__ "PCFactorSetShiftAmount_Factor" 41 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftAmount_Factor(PC pc,PetscReal shiftamount) 42 { 43 PC_Factor *dir = (PC_Factor*)pc->data; 44 45 PetscFunctionBegin; 46 if (shiftamount == (PetscReal) PETSC_DECIDE){ 47 dir->info.shiftamount = 1.e-12; 48 } else { 49 dir->info.shiftamount = shiftamount; 50 } 51 PetscFunctionReturn(0); 52 } 53 EXTERN_C_END 54 55 EXTERN_C_BEGIN 56 #undef __FUNCT__ 57 #define __FUNCT__ "PCFactorSetDropTolerance_Factor" 58 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetDropTolerance_Factor(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount) 59 { 60 PC_Factor *ilu = (PC_Factor*)pc->data; 61 62 PetscFunctionBegin; 63 if (pc->setupcalled && (!ilu->info.usedt || ((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) { 64 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Cannot change tolerance after use"); 65 } 66 ilu->info.usedt = PETSC_TRUE; 67 ilu->info.dt = dt; 68 ilu->info.dtcol = dtcol; 69 ilu->info.dtcount = dtcount; 70 ilu->info.fill = PETSC_DEFAULT; 71 PetscFunctionReturn(0); 72 } 73 EXTERN_C_END 74 75 EXTERN_C_BEGIN 76 #undef __FUNCT__ 77 #define __FUNCT__ "PCFactorSetFill_Factor" 78 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Factor(PC pc,PetscReal fill) 79 { 80 PC_Factor *dir = (PC_Factor*)pc->data; 81 82 PetscFunctionBegin; 83 dir->info.fill = fill; 84 PetscFunctionReturn(0); 85 } 86 EXTERN_C_END 87 88 EXTERN_C_BEGIN 89 #undef __FUNCT__ 90 #define __FUNCT__ "PCFactorSetMatOrderingType_Factor" 91 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Factor(PC pc,const MatOrderingType ordering) 92 { 93 PC_Factor *dir = (PC_Factor*)pc->data; 94 PetscErrorCode ierr; 95 PetscTruth flg; 96 97 PetscFunctionBegin; 98 if (!pc->setupcalled) { 99 ierr = PetscFree(dir->ordering);CHKERRQ(ierr); 100 ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr); 101 } else { 102 ierr = PetscStrcmp(dir->ordering,ordering,&flg);CHKERRQ(ierr); 103 if (!flg) { 104 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Cannot change ordering after use"); 105 } 106 } 107 PetscFunctionReturn(0); 108 } 109 EXTERN_C_END 110 111 EXTERN_C_BEGIN 112 #undef __FUNCT__ 113 #define __FUNCT__ "PCFactorSetLevels_Factor" 114 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetLevels_Factor(PC pc,PetscInt levels) 115 { 116 PC_Factor *ilu = (PC_Factor*)pc->data; 117 118 PetscFunctionBegin; 119 if (!pc->setupcalled) { 120 ilu->info.levels = levels; 121 } else if (ilu->info.usedt || ilu->info.levels != levels) { 122 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Cannot change levels after use"); 123 } 124 PetscFunctionReturn(0); 125 } 126 EXTERN_C_END 127 128 EXTERN_C_BEGIN 129 #undef __FUNCT__ 130 #define __FUNCT__ "PCFactorSetAllowDiagonalFill_Factor" 131 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetAllowDiagonalFill_Factor(PC pc) 132 { 133 PC_Factor *dir = (PC_Factor*)pc->data; 134 135 PetscFunctionBegin; 136 dir->info.diagonal_fill = 1; 137 PetscFunctionReturn(0); 138 } 139 EXTERN_C_END 140 141 142 /* ------------------------------------------------------------------------------------------*/ 143 144 EXTERN_C_BEGIN 145 #undef __FUNCT__ 146 #define __FUNCT__ "PCFactorSetPivotInBlocks_Factor" 147 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetPivotInBlocks_Factor(PC pc,PetscTruth pivot) 148 { 149 PC_Factor *dir = (PC_Factor*)pc->data; 150 151 PetscFunctionBegin; 152 dir->info.pivotinblocks = pivot ? 1.0 : 0.0; 153 PetscFunctionReturn(0); 154 } 155 EXTERN_C_END 156 157 EXTERN_C_BEGIN 158 #undef __FUNCT__ 159 #define __FUNCT__ "PCFactorGetMatrix_Factor" 160 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorGetMatrix_Factor(PC pc,Mat *mat) 161 { 162 PC_Factor *ilu = (PC_Factor*)pc->data; 163 164 PetscFunctionBegin; 165 if (!ilu->fact) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()"); 166 *mat = ilu->fact; 167 PetscFunctionReturn(0); 168 } 169 EXTERN_C_END 170 171 EXTERN_C_BEGIN 172 #undef __FUNCT__ 173 #define __FUNCT__ "PCFactorSetMatSolverPackage_Factor" 174 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatSolverPackage_Factor(PC pc,const MatSolverPackage stype) 175 { 176 PetscErrorCode ierr; 177 PC_Factor *lu = (PC_Factor*)pc->data; 178 179 PetscFunctionBegin; 180 if (lu->fact) { 181 const MatSolverPackage ltype; 182 PetscTruth flg; 183 ierr = MatFactorGetSolverPackage(lu->fact,<ype);CHKERRQ(ierr); 184 ierr = PetscStrcmp(stype,ltype,&flg);CHKERRQ(ierr); 185 if (!flg) { 186 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Cannot change solver matrix package after PC has been setup or used"); 187 } 188 } 189 ierr = PetscFree(lu->solvertype);CHKERRQ(ierr); 190 ierr = PetscStrallocpy(stype,&lu->solvertype);CHKERRQ(ierr); 191 PetscFunctionReturn(0); 192 } 193 EXTERN_C_END 194 195 EXTERN_C_BEGIN 196 #undef __FUNCT__ 197 #define __FUNCT__ "PCFactorGetMatSolverPackage_Factor" 198 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorGetMatSolverPackage_Factor(PC pc,const MatSolverPackage *stype) 199 { 200 PC_Factor *lu = (PC_Factor*)pc->data; 201 202 PetscFunctionBegin; 203 *stype = lu->solvertype; 204 PetscFunctionReturn(0); 205 } 206 EXTERN_C_END 207 208 EXTERN_C_BEGIN 209 #undef __FUNCT__ 210 #define __FUNCT__ "PCFactorSetColumnPivot_Factor" 211 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetColumnPivot_Factor(PC pc,PetscReal dtcol) 212 { 213 PC_Factor *dir = (PC_Factor*)pc->data; 214 215 PetscFunctionBegin; 216 if (dtcol < 0.0 || dtcol > 1.0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Column pivot tolerance is %G must be between 0 and 1",dtcol); 217 dir->info.dtcol = dtcol; 218 PetscFunctionReturn(0); 219 } 220 EXTERN_C_END 221 222 EXTERN_C_BEGIN 223 #undef __FUNCT__ 224 #define __FUNCT__ "PCSetFromOptions_Factor" 225 PetscErrorCode PETSCKSP_DLLEXPORT PCSetFromOptions_Factor(PC pc) 226 { 227 PC_Factor *factor = (PC_Factor*)pc->data; 228 PetscErrorCode ierr; 229 PetscTruth flg = PETSC_FALSE,set; 230 char tname[256], solvertype[64]; 231 PetscFList ordlist; 232 PetscEnum etmp; 233 234 PetscFunctionBegin; 235 if (!MatOrderingRegisterAllCalled) {ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 236 ierr = PetscOptionsTruth("-pc_factor_in_place","Form factored matrix in the same memory as the matrix","PCFactorSetUseInPlace",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 237 if (flg) { 238 ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr); 239 } 240 ierr = PetscOptionsReal("-pc_factor_fill","Expected non-zeros in factored matrix","PCFactorSetFill",((PC_Factor*)factor)->info.fill,&((PC_Factor*)factor)->info.fill,0);CHKERRQ(ierr); 241 242 ierr = PetscOptionsEnum("-pc_factor_shift_type","Shift added to diagonal","PCFactorSetShiftType", 243 MatFactorShiftTypes,(PetscEnum)(int)((PC_Factor*)factor)->info.shifttype,&etmp,&flg);CHKERRQ(ierr); 244 if (flg) { 245 ierr = PCFactorSetShiftType(pc,(MatFactorShiftType)etmp);CHKERRQ(ierr); 246 } 247 ierr = PetscOptionsReal("-pc_factor_shift_amount","Shift added to diagonal","PCFactorSetShiftAmount",((PC_Factor*)factor)->info.shiftamount,&((PC_Factor*)factor)->info.shiftamount,0);CHKERRQ(ierr); 248 249 ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",((PC_Factor*)factor)->info.zeropivot,&((PC_Factor*)factor)->info.zeropivot,0);CHKERRQ(ierr); 250 ierr = PetscOptionsReal("-pc_factor_column_pivot","Column pivot tolerance (used only for some factorization)","PCFactorSetColumnPivot",((PC_Factor*)factor)->info.dtcol,&((PC_Factor*)factor)->info.dtcol,&flg);CHKERRQ(ierr); 251 252 flg = ((PC_Factor*)factor)->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE; 253 ierr = PetscOptionsTruth("-pc_factor_pivot_in_blocks","Pivot inside matrix dense blocks for BAIJ and SBAIJ","PCFactorSetPivotInBlocks",flg,&flg,&set);CHKERRQ(ierr); 254 if (set) { 255 ierr = PCFactorSetPivotInBlocks(pc,flg);CHKERRQ(ierr); 256 } 257 258 flg = PETSC_FALSE; 259 ierr = PetscOptionsTruth("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 260 if (flg) { 261 ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr); 262 } 263 flg = PETSC_FALSE; 264 ierr = PetscOptionsTruth("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 265 if (flg) { 266 ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr); 267 } 268 269 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 270 ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in factored matrix","PCFactorSetMatOrderingType",ordlist,((PC_Factor*)factor)->ordering,tname,256,&flg);CHKERRQ(ierr); 271 if (flg) { 272 ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr); 273 } 274 275 /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */ 276 ierr = PetscOptionsString("-pc_factor_mat_solver_package","Specific direct solver to use","MatGetFactor",((PC_Factor*)factor)->solvertype,solvertype,64,&flg);CHKERRQ(ierr); 277 if (flg) { 278 ierr = PCFactorSetMatSolverPackage(pc,solvertype);CHKERRQ(ierr); 279 } 280 PetscFunctionReturn(0); 281 } 282 EXTERN_C_END 283 284 EXTERN_C_BEGIN 285 #undef __FUNCT__ 286 #define __FUNCT__ "PCView_Factor" 287 PetscErrorCode PCView_Factor(PC pc,PetscViewer viewer) 288 { 289 PC_Factor *factor = (PC_Factor*)pc->data; 290 PetscErrorCode ierr; 291 PetscTruth isstring,iascii; 292 293 PetscFunctionBegin; 294 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 295 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 296 if (iascii) { 297 if (factor->factortype == MAT_FACTOR_ILU || factor->factortype == MAT_FACTOR_ICC){ 298 if (factor->info.dt > 0) { 299 ierr = PetscViewerASCIIPrintf(viewer," drop tolerance %G\n",factor->info.dt);CHKERRQ(ierr); 300 ierr = PetscViewerASCIIPrintf(viewer," max nonzeros per row %D\n",factor->info.dtcount);CHKERRQ(ierr); 301 ierr = PetscViewerASCIIPrintf(viewer," column permutation tolerance %G\n",(PetscInt)factor->info.dtcol);CHKERRQ(ierr); 302 } else if (factor->info.levels == 1) { 303 ierr = PetscViewerASCIIPrintf(viewer," %D level of fill\n",(PetscInt)factor->info.levels);CHKERRQ(ierr); 304 } else { 305 ierr = PetscViewerASCIIPrintf(viewer," %D levels of fill\n",(PetscInt)factor->info.levels);CHKERRQ(ierr); 306 } 307 } 308 309 ierr = PetscViewerASCIIPrintf(viewer," tolerance for zero pivot %G\n",factor->info.zeropivot);CHKERRQ(ierr); 310 if (factor->info.shifttype==(PetscReal)MAT_SHIFT_POSITIVE_DEFINITE) { 311 ierr = PetscViewerASCIIPrintf(viewer," using Manteuffel shift\n");CHKERRQ(ierr); 312 } 313 if (factor->info.shifttype==(PetscReal)MAT_SHIFT_NONZERO) { 314 ierr = PetscViewerASCIIPrintf(viewer," using diagonal shift to prevent zero pivot\n");CHKERRQ(ierr); 315 } 316 if (factor->info.shifttype==(PetscReal)MAT_SHIFT_INBLOCKS) { 317 ierr = PetscViewerASCIIPrintf(viewer," using diagonal shift on blocks to prevent zero pivot\n");CHKERRQ(ierr); 318 } 319 320 ierr = PetscViewerASCIIPrintf(viewer," matrix ordering: %s\n",factor->ordering);CHKERRQ(ierr); 321 322 if (factor->fact) { 323 MatInfo info; 324 ierr = MatGetInfo(factor->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 325 ierr = PetscViewerASCIIPrintf(viewer," factor fill ratio given %G, needed %G\n",info.fill_ratio_given,info.fill_ratio_needed);CHKERRQ(ierr); 326 ierr = PetscViewerASCIIPrintf(viewer," Factored matrix follows:\n");CHKERRQ(ierr); 327 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 328 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 329 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 330 ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr); 331 ierr = MatView(factor->fact,viewer);CHKERRQ(ierr); 332 ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 333 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 334 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 335 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 336 } 337 338 } else if (isstring) { 339 MatFactorType t; 340 ierr = MatGetFactorType(factor->fact,&t);CHKERRQ(ierr); 341 if (t == MAT_FACTOR_ILU || t == MAT_FACTOR_ICC){ 342 ierr = PetscViewerStringSPrintf(viewer," lvls=%D,order=%s",(PetscInt)factor->info.levels,factor->ordering);CHKERRQ(ierr);CHKERRQ(ierr); 343 } 344 } else { 345 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported for PC_Factor",((PetscObject)viewer)->type_name); 346 } 347 PetscFunctionReturn(0); 348 } 349 EXTERN_C_END 350