Actual source code: preonly.c
1: /*$Id: preonly.c,v 1.42 2001/05/14 21:07:54 bsmith Exp $*/
3: /*
4: This implements a stub method that applies ONLY the preconditioner.
5: This may be used in inner iterations, where it is desired to
6: allow multiple iterations as well as the "0-iteration" case
7: */
8: #include src/sles/ksp/kspimpl.h
10: static int KSPSetUp_PREONLY(KSP ksp)
11: {
13: return(0);
14: }
16: static int KSPSolve_PREONLY(KSP ksp,int *its)
17: {
18: int ierr;
19: Vec X,B;
20: PetscTruth diagonalscale;
23: ierr = PCDiagonalScale(ksp->B,&diagonalscale);
24: if (diagonalscale) SETERRQ1(1,"Krylov method %s does not support diagonal scaling",ksp->type_name);
25: if (!ksp->guess_zero) {
26: SETERRQ(1,"Running KSP of preonly doesn't make sense with nonzero initial guessn
27: you probably want a KSP type of Richardson");
28: }
30: ksp->its = 0;
31: X = ksp->vec_sol;
32: B = ksp->vec_rhs;
33: ierr = KSP_PCApply(ksp,ksp->B,B,X);
34: *its = 1;
35: ksp->its = 1;
36: ksp->reason = KSP_CONVERGED_ITS;
37: return(0);
38: }
40: EXTERN_C_BEGIN
41: int KSPCreate_PREONLY(KSP ksp)
42: {
44: ksp->data = (void*)0;
45: ksp->ops->setup = KSPSetUp_PREONLY;
46: ksp->ops->solve = KSPSolve_PREONLY;
47: ksp->ops->destroy = KSPDefaultDestroy;
48: ksp->ops->buildsolution = KSPDefaultBuildSolution;
49: ksp->ops->buildresidual = KSPDefaultBuildResidual;
50: ksp->ops->setfromoptions = 0;
51: ksp->ops->view = 0;
52: return(0);
53: }
54: EXTERN_C_END