Actual source code: preonly.c

  1: /*$Id: preonly.c,v 1.42 2001/05/14 21:07:54 bsmith Exp $*/

 3:  #include src/ksp/ksp/kspimpl.h

  7: static int KSPSetUp_PREONLY(KSP ksp)
  8: {
 10:   return(0);
 11: }

 15: static int  KSPSolve_PREONLY(KSP ksp)
 16: {
 17:   int        ierr;
 18:   Vec        X,B;
 19:   PetscTruth diagonalscale;

 22:   PCDiagonalScale(ksp->B,&diagonalscale);
 23:   if (diagonalscale) SETERRQ1(1,"Krylov method %s does not support diagonal scaling",ksp->type_name);
 24:   if (!ksp->guess_zero) {
 25:     SETERRQ(1,"Running KSP of preonly doesn't make sense with nonzero initial guess\n\
 26:                you probably want a KSP type of Richardson");
 27:   }

 29:   ksp->its    = 0;
 30:   X           = ksp->vec_sol;
 31:   B           = ksp->vec_rhs;
 32:   KSP_PCApply(ksp,ksp->B,B,X);
 33:   ksp->its    = 1;
 34:   ksp->reason = KSP_CONVERGED_ITS;
 35:   return(0);
 36: }

 38: /*MC
 39:      KSPPREONLY - This implements a stub method that applies ONLY the preconditioner.
 40:                   This may be used in inner iterations, where it is desired to 
 41:                   allow multiple iterations as well as the "0-iteration" case. It is 
 42:                   commonly used with the direct solver preconditioners like PCLU and PCCHOLESKY

 44:    Options Database Keys:
 45: .   see KSPSolve()

 47:    Level: beginner

 49: .seealso:  KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP

 51: M*/

 53: EXTERN_C_BEGIN
 56: int KSPCreate_PREONLY(KSP ksp)
 57: {
 59:   ksp->data                      = (void*)0;
 60:   ksp->ops->setup                = KSPSetUp_PREONLY;
 61:   ksp->ops->solve                = KSPSolve_PREONLY;
 62:   ksp->ops->destroy              = KSPDefaultDestroy;
 63:   ksp->ops->buildsolution        = KSPDefaultBuildSolution;
 64:   ksp->ops->buildresidual        = KSPDefaultBuildResidual;
 65:   ksp->ops->setfromoptions       = 0;
 66:   ksp->ops->view                 = 0;
 67:   return(0);
 68: }
 69: EXTERN_C_END