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 = MAT_SHIFT_NONE; 30 } else { 31 dir->info.shifttype = shifttype; 32 } 33 PetscFunctionReturn(0); 34 } 35 36 #undef __FUNCT__ 37 #define __FUNCT__ "PCFactorSetShiftAmount_Factor" 38 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftAmount_Factor(PC pc,PetscReal shiftamount) 39 { 40 PC_Factor *dir = (PC_Factor*)pc->data; 41 42 PetscFunctionBegin; 43 if (shiftamount == (PetscReal) PETSC_DECIDE){ 44 dir->info.shiftamount = 1.e-12; 45 } else { 46 dir->info.shiftamount = shiftamount; 47 } 48 PetscFunctionReturn(0); 49 } 50 EXTERN_C_END 51 52 EXTERN_C_BEGIN 53 #undef __FUNCT__ 54 #define __FUNCT__ "PCFactorSetShiftNonzero_Factor" 55 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_Factor(PC pc,PetscReal shift) 56 { 57 PC_Factor *dir = (PC_Factor*)pc->data; 58 59 PetscFunctionBegin; 60 if (shift == (PetscReal) PETSC_DECIDE) { 61 dir->info.shiftnz = 1.e-12; 62 } else { 63 dir->info.shiftnz = shift; 64 } 65 PetscFunctionReturn(0); 66 } 67 EXTERN_C_END 68 69 EXTERN_C_BEGIN 70 #undef __FUNCT__ 71 #define __FUNCT__ "PCFactorSetShiftPd_Factor" 72 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_Factor(PC pc,PetscTruth shift) 73 { 74 PC_Factor *dir = (PC_Factor*)pc->data; 75 76 PetscFunctionBegin; 77 if (shift) { 78 dir->info.shiftpd = 1.0; 79 } else { 80 dir->info.shiftpd = 0.0; 81 } 82 PetscFunctionReturn(0); 83 } 84 EXTERN_C_END 85 86 EXTERN_C_BEGIN 87 #undef __FUNCT__ 88 #define __FUNCT__ "PCFactorSetDropTolerance_Factor" 89 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetDropTolerance_Factor(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount) 90 { 91 PC_Factor *ilu = (PC_Factor*)pc->data; 92 93 PetscFunctionBegin; 94 if (pc->setupcalled && (!ilu->info.usedt || ((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) { 95 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change tolerance after use"); 96 } 97 ilu->info.usedt = PETSC_TRUE; 98 ilu->info.dt = dt; 99 ilu->info.dtcol = dtcol; 100 ilu->info.dtcount = dtcount; 101 ilu->info.fill = PETSC_DEFAULT; 102 PetscFunctionReturn(0); 103 } 104 EXTERN_C_END 105 106 EXTERN_C_BEGIN 107 #undef __FUNCT__ 108 #define __FUNCT__ "PCFactorSetFill_Factor" 109 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Factor(PC pc,PetscReal fill) 110 { 111 PC_Factor *dir = (PC_Factor*)pc->data; 112 113 PetscFunctionBegin; 114 dir->info.fill = fill; 115 PetscFunctionReturn(0); 116 } 117 EXTERN_C_END 118 119 EXTERN_C_BEGIN 120 #undef __FUNCT__ 121 #define __FUNCT__ "PCFactorSetMatOrderingType_Factor" 122 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Factor(PC pc,const MatOrderingType ordering) 123 { 124 PC_Factor *dir = (PC_Factor*)pc->data; 125 PetscErrorCode ierr; 126 PetscTruth flg; 127 128 PetscFunctionBegin; 129 if (!pc->setupcalled) { 130 ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 131 ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr); 132 } else { 133 ierr = PetscStrcmp(dir->ordering,ordering,&flg);CHKERRQ(ierr); 134 if (!flg) { 135 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change ordering after use"); 136 } 137 } 138 PetscFunctionReturn(0); 139 } 140 EXTERN_C_END 141 142 EXTERN_C_BEGIN 143 #undef __FUNCT__ 144 #define __FUNCT__ "PCFactorSetLevels_Factor" 145 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetLevels_Factor(PC pc,PetscInt levels) 146 { 147 PC_Factor *ilu = (PC_Factor*)pc->data; 148 149 PetscFunctionBegin; 150 if (!pc->setupcalled) { 151 ilu->info.levels = levels; 152 } else if (ilu->info.usedt || ilu->info.levels != levels) { 153 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change levels after use"); 154 } 155 PetscFunctionReturn(0); 156 } 157 EXTERN_C_END 158 159 EXTERN_C_BEGIN 160 #undef __FUNCT__ 161 #define __FUNCT__ "PCFactorSetAllowDiagonalFill_Factor" 162 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetAllowDiagonalFill_Factor(PC pc) 163 { 164 PC_Factor *dir = (PC_Factor*)pc->data; 165 166 PetscFunctionBegin; 167 dir->info.diagonal_fill = 1; 168 PetscFunctionReturn(0); 169 } 170 EXTERN_C_END 171 172 173 /* ------------------------------------------------------------------------------------------*/ 174 175 EXTERN_C_BEGIN 176 #undef __FUNCT__ 177 #define __FUNCT__ "PCFactorSetPivotInBlocks_Factor" 178 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetPivotInBlocks_Factor(PC pc,PetscTruth pivot) 179 { 180 PC_Factor *dir = (PC_Factor*)pc->data; 181 182 PetscFunctionBegin; 183 dir->info.pivotinblocks = pivot ? 1.0 : 0.0; 184 PetscFunctionReturn(0); 185 } 186 EXTERN_C_END 187 188 EXTERN_C_BEGIN 189 #undef __FUNCT__ 190 #define __FUNCT__ "PCFactorSetShiftInBlocks_Factor" 191 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftInBlocks_Factor(PC pc,PetscReal shift) 192 { 193 PC_Factor *dir = (PC_Factor*)pc->data; 194 195 PetscFunctionBegin; 196 if (shift == PETSC_DEFAULT) { 197 dir->info.shiftinblocks = 1.e-12; 198 } else { 199 dir->info.shiftinblocks = shift; 200 } 201 PetscFunctionReturn(0); 202 } 203 EXTERN_C_END 204 205 206 EXTERN_C_BEGIN 207 #undef __FUNCT__ 208 #define __FUNCT__ "PCFactorGetMatrix_Factor" 209 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorGetMatrix_Factor(PC pc,Mat *mat) 210 { 211 PC_Factor *ilu = (PC_Factor*)pc->data; 212 213 PetscFunctionBegin; 214 if (!ilu->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()"); 215 *mat = ilu->fact; 216 PetscFunctionReturn(0); 217 } 218 EXTERN_C_END 219 220 EXTERN_C_BEGIN 221 #undef __FUNCT__ 222 #define __FUNCT__ "PCFactorSetMatSolverPackage_Factor" 223 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatSolverPackage_Factor(PC pc,const MatSolverPackage stype) 224 { 225 PetscErrorCode ierr; 226 PC_Factor *lu = (PC_Factor*)pc->data; 227 228 PetscFunctionBegin; 229 if (lu->fact) { 230 const MatSolverPackage ltype; 231 PetscTruth flg; 232 ierr = MatFactorGetSolverPackage(lu->fact,<ype);CHKERRQ(ierr); 233 ierr = PetscStrcmp(stype,ltype,&flg);CHKERRQ(ierr); 234 if (!flg) { 235 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change solver matrix package after PC has been setup or used"); 236 } 237 } 238 ierr = PetscStrfree(lu->solvertype);CHKERRQ(ierr); 239 ierr = PetscStrallocpy(stype,&lu->solvertype);CHKERRQ(ierr); 240 PetscFunctionReturn(0); 241 } 242 EXTERN_C_END 243 244 EXTERN_C_BEGIN 245 #undef __FUNCT__ 246 #define __FUNCT__ "PCFactorGetMatSolverPackage_Factor" 247 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorGetMatSolverPackage_Factor(PC pc,const MatSolverPackage *stype) 248 { 249 PC_Factor *lu = (PC_Factor*)pc->data; 250 251 PetscFunctionBegin; 252 *stype = lu->solvertype; 253 PetscFunctionReturn(0); 254 } 255 EXTERN_C_END 256 257 EXTERN_C_BEGIN 258 #undef __FUNCT__ 259 #define __FUNCT__ "PCFactorSetColumnPivot_Factor" 260 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetColumnPivot_Factor(PC pc,PetscReal dtcol) 261 { 262 PC_Factor *dir = (PC_Factor*)pc->data; 263 264 PetscFunctionBegin; 265 if (dtcol < 0.0 || dtcol > 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Column pivot tolerance is %G must be between 0 and 1",dtcol); 266 dir->info.dtcol = dtcol; 267 PetscFunctionReturn(0); 268 } 269 EXTERN_C_END 270 271 EXTERN_C_BEGIN 272 #undef __FUNCT__ 273 #define __FUNCT__ "PCSetFromOptions_Factor" 274 PetscErrorCode PETSCKSP_DLLEXPORT PCSetFromOptions_Factor(PC pc) 275 { 276 PC_Factor *factor = (PC_Factor*)pc->data; 277 PetscErrorCode ierr; 278 PetscTruth flg = PETSC_FALSE,set; 279 char tname[256], solvertype[64]; 280 PetscFList ordlist; 281 282 PetscFunctionBegin; 283 if (!MatOrderingRegisterAllCalled) {ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 284 ierr = PetscOptionsTruth("-pc_factor_in_place","Form factored matrix in the same memory as the matrix","PCFactorSetUseInPlace",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 285 if (flg) { 286 ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr); 287 } 288 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); 289 290 ierr = PetscOptionsEnum("-pc_factor_shift_type","Shift added to diagonal","PCFactorSetShiftType", 291 MatFactorShiftTypes,(PetscEnum)((PC_Factor*)factor)->info.shifttype,(PetscEnum*)&((PC_Factor*)factor)->info.shifttype,&flg);CHKERRQ(ierr); 292 293 ierr = PetscOptionsReal("-pc_factor_shift_amount","Shift added to diagonal","PCFactorSetShiftAmount",((PC_Factor*)factor)->info.shiftamount,&((PC_Factor*)factor)->info.shiftamount,0);CHKERRQ(ierr); 294 295 #if defined(MV) 296 flg = PETSC_FALSE; 297 ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 298 if (flg) { 299 ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr); 300 } 301 ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",((PC_Factor*)factor)->info.shiftnz,&((PC_Factor*)factor)->info.shiftnz,0);CHKERRQ(ierr); 302 flg = PETSC_FALSE; 303 ierr = PetscOptionsTruth("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 304 if (flg) { 305 ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 306 } 307 flg = PETSC_FALSE; 308 ierr = PetscOptionsName("-pc_factor_shift_in_blocks","Shift added to diagonal dense blocks","PCFactorSetShiftInBlocks",&flg);CHKERRQ(ierr); 309 if (flg) { 310 ierr = PCFactorSetShiftInBlocks(pc,(PetscReal) PETSC_DEFAULT);CHKERRQ(ierr); 311 } 312 ierr = PetscOptionsReal("-pc_factor_shift_in_blocks","Shift added to diagonal dense blocks","PCFactorSetShiftInBlocks",((PC_Factor*)factor)->info.shiftinblocks,&((PC_Factor*)factor)->info.shiftinblocks,0);CHKERRQ(ierr); 313 #endif 314 315 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); 316 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); 317 318 flg = ((PC_Factor*)factor)->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE; 319 ierr = PetscOptionsTruth("-pc_factor_pivot_in_blocks","Pivot inside matrix dense blocks for BAIJ and SBAIJ","PCFactorSetPivotInBlocks",flg,&flg,&set);CHKERRQ(ierr); 320 if (set) { 321 ierr = PCFactorSetPivotInBlocks(pc,flg);CHKERRQ(ierr); 322 } 323 324 flg = PETSC_FALSE; 325 ierr = PetscOptionsTruth("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 326 if (flg) { 327 ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr); 328 } 329 flg = PETSC_FALSE; 330 ierr = PetscOptionsTruth("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 331 if (flg) { 332 ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr); 333 } 334 335 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 336 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); 337 if (flg) { 338 ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr); 339 } 340 341 /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */ 342 ierr = PetscOptionsString("-pc_factor_mat_solver_package","Specific direct solver to use","MatGetFactor",((PC_Factor*)factor)->solvertype,solvertype,64,&flg);CHKERRQ(ierr); 343 if (flg) { 344 ierr = PCFactorSetMatSolverPackage(pc,solvertype);CHKERRQ(ierr); 345 } 346 347 /* 348 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); 349 350 flg = ((PC_Factor*)factor)->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE; 351 ierr = PetscOptionsTruth("-pc_factor_pivot_in_blocks","Pivot inside matrix dense blocks for BAIJ and SBAIJ","PCFactorSetPivotInBlocks",flg,&flg,&set);CHKERRQ(ierr); 352 if (set) { 353 ierr = PCFactorSetPivotInBlocks(pc,flg);CHKERRQ(ierr); 354 } 355 */ 356 PetscFunctionReturn(0); 357 } 358 EXTERN_C_END 359