Actual source code: pcsles.c

  1: /*$Id: pcsles.c,v 1.39 2001/04/10 19:36:17 bsmith Exp $*/
  2: /*
  3:       Defines a preconditioner that can consist of any SLES solver.
  4:     This allows embedding a Krylov method inside a preconditioner.
  5: */
 6:  #include src/sles/pc/pcimpl.h
 7:  #include petscsles.h

  9: typedef struct {
 10:   PetscTruth use_true_matrix;       /* use mat rather than pmat in inner linear solve */
 11:   SLES       sles;
 12:   int        its;                   /* total number of iterations SLES uses */
 13: } PC_SLES;

 15: static int PCApply_SLES(PC pc,Vec x,Vec y)
 16: {
 17:   int     ierr,its;
 18:   PC_SLES *jac = (PC_SLES*)pc->data;

 21:   ierr      = SLESSolve(jac->sles,x,y,&its);
 22:   jac->its += its;
 23:   return(0);
 24: }

 26: static int PCApplyTranspose_SLES(PC pc,Vec x,Vec y)
 27: {
 28:   int     ierr,its;
 29:   PC_SLES *jac = (PC_SLES*)pc->data;

 32:   ierr      = SLESSolveTranspose(jac->sles,x,y,&its);
 33:   jac->its += its;
 34:   return(0);
 35: }

 37: static int PCSetUp_SLES(PC pc)
 38: {
 39:   int     ierr;
 40:   PC_SLES *jac = (PC_SLES*)pc->data;
 41:   Mat     mat;

 44:   SLESSetFromOptions(jac->sles);
 45:   if (jac->use_true_matrix) mat = pc->mat;
 46:   else                      mat = pc->pmat;

 48:   SLESSetOperators(jac->sles,mat,pc->pmat,pc->flag);
 49:   SLESSetUp(jac->sles,pc->vec,pc->vec);
 50:   return(0);
 51: }

 53: /* Default destroy, if it has never been setup */
 54: static int PCDestroy_SLES(PC pc)
 55: {
 56:   PC_SLES *jac = (PC_SLES*)pc->data;
 57:   int     ierr;

 60:   SLESDestroy(jac->sles);
 61:   PetscFree(jac);
 62:   return(0);
 63: }

 65: static int PCView_SLES(PC pc,PetscViewer viewer)
 66: {
 67:   PC_SLES    *jac = (PC_SLES*)pc->data;
 68:   int        ierr;
 69:   PetscTruth isascii;

 72:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
 73:   if (isascii) {
 74:     if (jac->use_true_matrix) {
 75:       PetscViewerASCIIPrintf(viewer,"Using true matrix (not preconditioner matrix) on inner solven");
 76:     }
 77:     PetscViewerASCIIPrintf(viewer,"KSP and PC on SLES preconditioner follown");
 78:     PetscViewerASCIIPrintf(viewer,"---------------------------------n");
 79:   } else {
 80:     SETERRQ1(1,"Viewer type %s not supported for this object",((PetscObject)viewer)->type_name);
 81:   }
 82:   PetscViewerASCIIPushTab(viewer);
 83:   SLESView(jac->sles,viewer);
 84:   PetscViewerASCIIPopTab(viewer);
 85:   if (isascii) {
 86:     PetscViewerASCIIPrintf(viewer,"---------------------------------n");
 87:   }
 88:   return(0);
 89: }

 91: static int PCSetFromOptions_SLES(PC pc){
 92:   int        ierr;
 93:   PetscTruth flg;

 96:   PetscOptionsHead("SLES preconditioner options");
 97:     PetscOptionsName("-pc_sles_true","Use true matrix to define inner linear system, not preconditioner matrix","PCSLESSetUseTrue",&flg);
 98:     if (flg) {
 99:       PCSLESSetUseTrue(pc);
100:     }
101:   PetscOptionsTail();
102:   return(0);
103: }

105: /* ----------------------------------------------------------------------------------*/

107: EXTERN_C_BEGIN
108: int PCSLESSetUseTrue_SLES(PC pc)
109: {
110:   PC_SLES   *jac;

113:   jac                  = (PC_SLES*)pc->data;
114:   jac->use_true_matrix = PETSC_TRUE;
115:   return(0);
116: }
117: EXTERN_C_END

119: EXTERN_C_BEGIN
120: int PCSLESGetSLES_SLES(PC pc,SLES *sles)
121: {
122:   PC_SLES   *jac;

125:   jac          = (PC_SLES*)pc->data;
126:   *sles        = jac->sles;
127:   return(0);
128: }
129: EXTERN_C_END

131: /*@
132:    PCSLESSetUseTrue - Sets a flag to indicate that the true matrix (rather than
133:    the matrix used to define the preconditioner) is used to compute the
134:    residual inside the inner solve.

136:    Collective on PC

138:    Input Parameters:
139: .  pc - the preconditioner context

141:    Options Database Key:
142: .  -pc_sles_true - Activates PCSLESSetUseTrue()

144:    Note:
145:    For the common case in which the preconditioning and linear 
146:    system matrices are identical, this routine is unnecessary.

148:    Level: advanced

150: .keywords:  PC, SLES, set, true, local, flag

152: .seealso: PCSetOperators(), PCBJacobiSetUseTrueLocal()
153: @*/
154: int PCSLESSetUseTrue(PC pc)
155: {
156:   int ierr,(*f)(PC);

160:   PetscObjectQueryFunction((PetscObject)pc,"PCSLESSetUseTrue_C",(void (**)(void))&f);
161:   if (f) {
162:     (*f)(pc);
163:   }
164:   return(0);
165: }

167: /*@C
168:    PCSLESGetSLES - Gets the SLES context for a SLES PC.

170:    Not Collective but SLES returned is parallel if PC was parallel

172:    Input Parameter:
173: .  pc - the preconditioner context

175:    Output Parameters:
176: .  sles - the PC solver

178:    Notes:
179:    You must call SLESSetUp() before calling PCSLESGetSLES().

181:    Level: advanced

183: .keywords:  PC, SLES, get, context
184: @*/
185: int PCSLESGetSLES(PC pc,SLES *sles)
186: {
187:   int ierr,(*f)(PC,SLES*);

191:   if (!pc->setupcalled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call SLESSetUp first");
192:   PetscObjectQueryFunction((PetscObject)pc,"PCSLESGetSLES_C",(void (**)(void))&f);
193:   if (f) {
194:     (*f)(pc,sles);
195:   }
196:   return(0);
197: }

199: /* ----------------------------------------------------------------------------------*/

201: EXTERN_C_BEGIN
202: int PCCreate_SLES(PC pc)
203: {
204:   int       ierr;
205:   char      *prefix;
206:   PC_SLES   *jac;

209:   PetscNew(PC_SLES,&jac);
210:   PetscLogObjectMemory(pc,sizeof(PC_SLES));
211:   pc->ops->apply              = PCApply_SLES;
212:   pc->ops->applytranspose     = PCApplyTranspose_SLES;
213:   pc->ops->setup              = PCSetUp_SLES;
214:   pc->ops->destroy            = PCDestroy_SLES;
215:   pc->ops->setfromoptions     = PCSetFromOptions_SLES;
216:   pc->ops->view               = PCView_SLES;
217:   pc->ops->applyrichardson    = 0;

219:   pc->data               = (void*)jac;
220:   ierr                   = SLESCreate(pc->comm,&jac->sles);

222:   PCGetOptionsPrefix(pc,&prefix);
223:   SLESSetOptionsPrefix(jac->sles,prefix);
224:   SLESAppendOptionsPrefix(jac->sles,"sles_");
225:   jac->use_true_matrix = PETSC_FALSE;
226:   jac->its             = 0;

228:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSLESSetUseTrue_C","PCSLESSetUseTrue_SLES",
229:                     PCSLESSetUseTrue_SLES);
230:   PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSLESGetSLES_C","PCSLESGetSLES_SLES",
231:                     PCSLESGetSLES_SLES);

233:   return(0);
234: }
235: EXTERN_C_END