Actual source code: factimpl.c

  1: #define PETSCKSP_DLL

 3:  #include ../src/ksp/pc/impls/factor/factor.h

  5: /* ------------------------------------------------------------------------------------------*/

 10: PetscErrorCode  PCFactorSetZeroPivot_Factor(PC pc,PetscReal z)
 11: {
 12:   PC_Factor *ilu = (PC_Factor*)pc->data;

 15:   ilu->info.zeropivot = z;
 16:   return(0);
 17: }

 23: PetscErrorCode  PCFactorSetShiftNonzero_Factor(PC pc,PetscReal shift)
 24: {
 25:   PC_Factor *dir = (PC_Factor*)pc->data;

 28:   if (shift == (PetscReal) PETSC_DECIDE) {
 29:     dir->info.shiftnz = 1.e-12;
 30:   } else {
 31:     dir->info.shiftnz = shift;
 32:   }
 33:   return(0);
 34: }

 40: PetscErrorCode  PCFactorSetShiftPd_Factor(PC pc,PetscTruth shift)
 41: {
 42:   PC_Factor *dir = (PC_Factor*)pc->data;
 43: 
 45:   if (shift) {
 46:     dir->info.shiftpd = 1.0;
 47:   } else {
 48:     dir->info.shiftpd = 0.0;
 49:   }
 50:   return(0);
 51: }


 58: PetscErrorCode  PCFactorSetUseDropTolerance_Factor(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount)
 59: {
 60:   PC_Factor         *ilu = (PC_Factor*)pc->data;

 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_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:   return(0);
 72: }

 78: PetscErrorCode  PCFactorSetFill_Factor(PC pc,PetscReal fill)
 79: {
 80:   PC_Factor *dir = (PC_Factor*)pc->data;

 83:   dir->info.fill = fill;
 84:   return(0);
 85: }

 91: PetscErrorCode  PCFactorSetMatOrderingType_Factor(PC pc,const MatOrderingType ordering)
 92: {
 93:   PC_Factor      *dir = (PC_Factor*)pc->data;
 95:   PetscTruth     flg;
 96: 
 98:   if (!pc->setupcalled) {
 99:      PetscStrfree(dir->ordering);
100:      PetscStrallocpy(ordering,&dir->ordering);
101:   } else {
102:     PetscStrcmp(dir->ordering,ordering,&flg);
103:     if (!flg) {
104:       SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change ordering after use");
105:     }
106:   }
107:   return(0);
108: }

114: PetscErrorCode  PCFactorSetLevels_Factor(PC pc,PetscInt levels)
115: {
116:   PC_Factor      *ilu = (PC_Factor*)pc->data;

119:   if (!pc->setupcalled) {
120:     ilu->info.levels = levels;
121:   } else if (ilu->info.usedt || ilu->info.levels != levels) {
122:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change levels after use");
123:   }
124:   return(0);
125: }

131: PetscErrorCode  PCFactorSetAllowDiagonalFill_Factor(PC pc)
132: {
133:   PC_Factor *dir = (PC_Factor*)pc->data;

136:   dir->info.diagonal_fill = 1;
137:   return(0);
138: }


142: /* ------------------------------------------------------------------------------------------*/

147: PetscErrorCode  PCFactorSetPivotInBlocks_Factor(PC pc,PetscTruth pivot)
148: {
149:   PC_Factor *dir = (PC_Factor*)pc->data;

152:   dir->info.pivotinblocks = pivot ? 1.0 : 0.0;
153:   return(0);
154: }

160: PetscErrorCode  PCFactorSetShiftInBlocks_Factor(PC pc,PetscReal shift)
161: {
162:   PC_Factor *dir = (PC_Factor*)pc->data;

165:   if (shift == PETSC_DEFAULT) {
166:     dir->info.shiftinblocks = 1.e-12;
167:   } else {
168:     dir->info.shiftinblocks = shift;
169:   }
170:   return(0);
171: }


178: PetscErrorCode  PCFactorGetMatrix_Factor(PC pc,Mat *mat)
179: {
180:   PC_Factor *ilu = (PC_Factor*)pc->data;

183:   if (!ilu->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()");
184:   *mat = ilu->fact;
185:   return(0);
186: }

192: PetscErrorCode  PCFactorSetMatSolverPackage_Factor(PC pc,const MatSolverPackage stype)
193: {
195:   PC_Factor      *lu = (PC_Factor*)pc->data;

198:   if (lu->fact) {
199:     const MatSolverPackage ltype;
200:     PetscTruth             flg;
201:     MatFactorGetSolverPackage(lu->fact,&ltype);
202:     PetscStrcmp(stype,ltype,&flg);
203:     if (!flg) {
204:       SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot change solver matrix package after PC has been setup or used");
205:     }
206:   }
207:   PetscStrfree(lu->solvertype);
208:   PetscStrallocpy(stype,&lu->solvertype);
209:   return(0);
210: }

216: PetscErrorCode  PCFactorGetMatSolverPackage_Factor(PC pc,const MatSolverPackage *stype)
217: {
218:   PC_Factor      *lu = (PC_Factor*)pc->data;

221:   *stype = lu->solvertype;
222:   return(0);
223: }

229: PetscErrorCode  PCFactorSetPivoting_Factor(PC pc,PetscReal dtcol)
230: {
231:   PC_Factor *dir = (PC_Factor*)pc->data;

234:   if (dtcol < 0.0 || dtcol > 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Column pivot tolerance is %G must be between 0 and 1",dtcol);
235:   dir->info.dtcol = dtcol;
236:   return(0);
237: }