Actual source code: modpcf.c

  1: /* $Id: modpcf.c,v 1.14 2001/04/10 19:36:35 bsmith Exp $*/

 3:  #include petscsles.h
  4: /*@C
  5:    KSPFGMRESSetModifyPC - Sets the routine used by FGMRES to modify the preconditioner.

  7:    Collective on KSP

  9:    Input Parameters:
 10: +  ksp - iterative context obtained from KSPCreate
 11: .  fcn - modifypc function
 12: .  ctx - optional contex
 13: -  d - optional context destroy routine

 15:    Calling Sequence of function:
 16:     int fcn(KSP ksp,int total_its,int loc_its,PetscReal res_norm,void*ctx);

 18:     ksp - the ksp context being used.
 19:     total_its     - the total number of FGMRES iterations that have occurred.    
 20:     loc_its       - the number of FGMRES iterations since last restart.
 21:     res_norm      - the current residual norm.
 22:     ctx           - optional context variable

 24:    Options Database Keys:
 25:    -ksp_fgmres_modifypcnochange
 26:    -ksp_fgmres_modifypcsles

 28:    Level: intermediate

 30:    Contributed by Allison Baker

 32:    Notes:
 33:    Several modifypc routines are predefined, including
 34:     KSPFGMRESModifyPCNoChange()
 35:     KSPFGMRESModifyPCSLES()

 37: .seealso: KSPFGMRESModifyPCNoChange(), KSPFGMRESModifyPCSLES()

 39: @*/
 40: int KSPFGMRESSetModifyPC(KSP ksp,int (*fcn)(KSP,int,int,PetscReal,void*),void* ctx,int (*d)(void*))
 41: {
 42:   int ierr,(*f)(KSP,int (*)(KSP,int,int,PetscReal,void*),void*,int (*)(void*));

 46:   PetscObjectQueryFunction((PetscObject)ksp,"KSPFGMRESSetModifyPC_C",(void (**)(void))&f);
 47:   if (f) {
 48:     (*f)(ksp,fcn,ctx,d);
 49:   }
 50:   return(0);
 51: }



 55: /* The following are different routines used to modify the preconditioner */

 57: /*@C

 59:   KSPFGMRESModifyPCNoChange - this is the default used by fgmres - it doesn't change the preconditioner. 

 61:   Input Parameters:
 62: +    ksp - the ksp context being used.
 63: .    total_its     - the total number of FGMRES iterations that have occurred.    
 64: .    loc_its       - the number of FGMRES iterations since last restart.
 65:                     a restart (so number of Krylov directions to be computed)
 66: .    res_norm      - the current residual norm.
 67: -    dummy         - context variable, unused in this routine

 69:    Level: intermediate

 71:    Contributed by Allison Baker

 73: You can use this as a template!

 75: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCSLES()

 77: @*/
 78: int KSPFGMRESModifyPCNoChange(KSP ksp,int total_its,int loc_its,PetscReal res_norm,void* dummy)
 79: {

 82:   return(0);
 83: }

 85: /*@C

 87:  KSPFGMRESModifyPCSLES - modifies the attributes of the
 88:      GMRES preconditioner.  It serves as an example (not as something 
 89:      useful!) 

 91:   Input Parameters:
 92: +    ksp - the ksp context being used.
 93: .    total_its     - the total number of FGMRES iterations that have occurred.    
 94: .    loc_its       - the number of FGMRES iterations since last restart.
 95: .    res_norm      - the current residual norm.
 96: -    dummy         - context, not used here

 98:    Level: intermediate

100:    Contributed by Allison Baker

102:  This could be used as a template!

104: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCSLES()

106: @*/
107: int KSPFGMRESModifyPCSLES(KSP ksp,int total_its,int loc_its,PetscReal res_norm,void *dummy)
108: {
109:   PC         pc;
110:   int        ierr,maxits;
111:   SLES       sub_sles;
112:   KSP        sub_ksp;
113:   PetscReal  rtol,atol,dtol;
114:   PetscTruth issles;


118:   KSPGetPC(ksp,&pc);

120:   PetscTypeCompare((PetscObject)pc,PCSLES,&issles);
121:   if (issles) {
122:     PCSLESGetSLES(pc,&sub_sles);
123:     SLESGetKSP(sub_sles,&sub_ksp);
124: 
125:     /* note that at this point you could check the type of KSP with KSPGetType() */

127:     /* Now we can use functions such as KSPGMRESSetRestart() or 
128:       KSPGMRESSetOrthogonalization() or KSPSetTolerances() */

130:     KSPGetTolerances(sub_ksp,&rtol,&atol,&dtol,&maxits);
131:     if (loc_its == 0) {
132:       rtol = .1;
133:     } else {
134:       rtol *= .9;
135:     }
136:     KSPSetTolerances(sub_ksp,rtol,atol,dtol,maxits);
137:   }
138:   return(0);
139: }