Actual source code: modpcf.c

 2:  #include petscksp.h
  5: /*@C
  6:    KSPFGMRESSetModifyPC - Sets the routine used by FGMRES to modify the preconditioner.

  8:    Collective on KSP

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

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

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

 25:    Options Database Keys:
 26:    -ksp_fgmres_modifypcnochange
 27:    -ksp_fgmres_modifypcksp

 29:    Level: intermediate

 31:    Contributed by Allison Baker

 33:    Notes:
 34:    Several modifypc routines are predefined, including
 35:     KSPFGMRESModifyPCNoChange()
 36:     KSPFGMRESModifyPCKSP()

 38: .seealso: KSPFGMRESModifyPCNoChange(), KSPFGMRESModifyPCKSP()

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

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


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

 59: /*@C

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

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

 71:    Level: intermediate

 73:    Contributed by Allison Baker

 75: You can use this as a template!

 77: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()

 79: @*/
 80: PetscErrorCode KSPFGMRESModifyPCNoChange(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void* dummy)
 81: {

 84:   return(0);
 85: }

 89: /*@C

 91:  KSPFGMRESModifyPCKSP - modifies the attributes of the
 92:      GMRES preconditioner.  It serves as an example (not as something 
 93:      useful!) 

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

102:    Level: intermediate

104:    Contributed by Allison Baker

106:  This could be used as a template!

108: .seealso: KSPFGMRESSetModifyPC(), KSPFGMRESModifyPCKSP()

110: @*/
111: PetscErrorCode KSPFGMRESModifyPCKSP(KSP ksp,PetscInt total_its,PetscInt loc_its,PetscReal res_norm,void *dummy)
112: {
113:   PC             pc;
115:   PetscInt       maxits;
116:   KSP            sub_ksp;
117:   PetscReal      rtol,abstol,dtol;
118:   PetscTruth     isksp;

121:   KSPGetPC(ksp,&pc);

123:   PetscTypeCompare((PetscObject)pc,PCKSP,&isksp);
124:   if (isksp) {
125:     PCKSPGetKSP(pc,&sub_ksp);
126: 
127:     /* note that at this point you could check the type of KSP with KSPGetType() */

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

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