Actual source code: preonly.c


  2: #include <petsc/private/kspimpl.h>

  4: static PetscErrorCode KSPSetUp_PREONLY(KSP ksp)
  5: {
  6:   PetscFunctionBegin;
  7:   PetscFunctionReturn(PETSC_SUCCESS);
  8: }

 10: static PetscErrorCode KSPSolve_PREONLY(KSP ksp)
 11: {
 12:   PetscBool      diagonalscale;
 13:   PCFailedReason pcreason;

 15:   PetscFunctionBegin;
 16:   PetscCall(PCGetDiagonalScale(ksp->pc, &diagonalscale));
 17:   PetscCheck(!diagonalscale, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Krylov method %s does not support diagonal scaling", ((PetscObject)ksp)->type_name);
 18:   PetscCheck(ksp->guess_zero, PetscObjectComm((PetscObject)ksp), PETSC_ERR_USER, "Running KSP of preonly doesn't make sense with nonzero initial guess\n\
 19:                you probably want a KSP type of Richardson");
 20:   ksp->its = 0;
 21:   PetscCall(KSP_PCApply(ksp, ksp->vec_rhs, ksp->vec_sol));
 22:   PetscCall(PCGetFailedReasonRank(ksp->pc, &pcreason));
 23:   /* Note: only some ranks may have this set; this may lead to problems if the caller assumes ksp->reason is set on all processes or just uses the result */
 24:   if (pcreason) {
 25:     PetscCall(VecSetInf(ksp->vec_sol));
 26:     ksp->reason = KSP_DIVERGED_PC_FAILED;
 27:   } else {
 28:     ksp->its    = 1;
 29:     ksp->reason = KSP_CONVERGED_ITS;
 30:   }
 31:   if (ksp->numbermonitors) {
 32:     Vec       v;
 33:     PetscReal norm;
 34:     Mat       A;

 36:     PetscCall(VecNorm(ksp->vec_rhs, NORM_2, &norm));
 37:     PetscCall(KSPMonitor(ksp, 0, norm));
 38:     PetscCall(VecDuplicate(ksp->vec_rhs, &v));
 39:     PetscCall(PCGetOperators(ksp->pc, &A, NULL));
 40:     PetscCall(KSP_MatMult(ksp, A, ksp->vec_sol, v));
 41:     PetscCall(VecAYPX(v, -1.0, ksp->vec_rhs));
 42:     PetscCall(VecNorm(v, NORM_2, &norm));
 43:     PetscCall(VecDestroy(&v));
 44:     PetscCall(KSPMonitor(ksp, 1, norm));
 45:   }
 46:   PetscFunctionReturn(PETSC_SUCCESS);
 47: }

 49: static PetscErrorCode KSPMatSolve_PREONLY(KSP ksp, Mat B, Mat X)
 50: {
 51:   PetscBool      diagonalscale;
 52:   PCFailedReason pcreason;

 54:   PetscFunctionBegin;
 55:   PetscCall(PCGetDiagonalScale(ksp->pc, &diagonalscale));
 56:   PetscCheck(!diagonalscale, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Krylov method %s does not support diagonal scaling", ((PetscObject)ksp)->type_name);
 57:   PetscCheck(ksp->guess_zero, PetscObjectComm((PetscObject)ksp), PETSC_ERR_USER, "Running KSP of preonly doesn't make sense with nonzero initial guess\n\
 58:                you probably want a KSP type of Richardson");
 59:   ksp->its = 0;
 60:   PetscCall(KSP_PCMatApply(ksp, B, X));
 61:   PetscCall(PCGetFailedReason(ksp->pc, &pcreason));
 62:   /* Note: only some ranks may have this set; this may lead to problems if the caller assumes ksp->reason is set on all processes or just uses the result */
 63:   if (pcreason) {
 64:     PetscCall(MatSetInf(X));
 65:     ksp->reason = KSP_DIVERGED_PC_FAILED;
 66:   } else {
 67:     ksp->its    = 1;
 68:     ksp->reason = KSP_CONVERGED_ITS;
 69:   }
 70:   PetscFunctionReturn(PETSC_SUCCESS);
 71: }

 73: /*MC
 74:      KSPPREONLY - This implements a method that applies ONLY the preconditioner exactly once.
 75:                   This may be used in inner iterations, where it is desired to
 76:                   allow multiple iterations as well as the "0-iteration" case. It is
 77:                   commonly used with the direct solver preconditioners like `PCLU` and `PCCHOLESKY`.
 78:                   There is an alias of `KSPNONE`.

 80:    Options Database Key:
 81: .   -ksp_type preonly - use preconditioner only

 83:    Level: beginner

 85:    Notes:
 86:    Since this does not involve an iteration the basic `KSP` parameters such as tolerances and iteration counts
 87:    do not apply

 89:    To apply multiple preconditioners in a simple iteration use `KSPRICHARDSON`

 91:    Developer Note:
 92:    Even though this method does not use any norms, the user is allowed to set the `KSPNormType` to any value.
 93:    This is so the users does not have to change `KSPNormType` options when they switch from other `KSP` methods to this one.

 95: .seealso: [](ch_ksp), `KSPCreate()`, `KSPSetType()`, `KSPType`, `KSP`, `KSPRICHARDSON`, `KSPCHEBYSHEV`
 96: M*/

 98: PETSC_EXTERN PetscErrorCode KSPCreate_PREONLY(KSP ksp)
 99: {
100:   PetscFunctionBegin;
101:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_NONE, PC_LEFT, 3));
102:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_NONE, PC_RIGHT, 2));
103:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_PRECONDITIONED, PC_LEFT, 2));
104:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_PRECONDITIONED, PC_RIGHT, 2));
105:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_UNPRECONDITIONED, PC_LEFT, 2));
106:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_UNPRECONDITIONED, PC_RIGHT, 2));
107:   PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_NATURAL, PC_LEFT, 2));

109:   ksp->data                = NULL;
110:   ksp->ops->setup          = KSPSetUp_PREONLY;
111:   ksp->ops->solve          = KSPSolve_PREONLY;
112:   ksp->ops->matsolve       = KSPMatSolve_PREONLY;
113:   ksp->ops->destroy        = KSPDestroyDefault;
114:   ksp->ops->buildsolution  = KSPBuildSolutionDefault;
115:   ksp->ops->buildresidual  = KSPBuildResidualDefault;
116:   ksp->ops->setfromoptions = NULL;
117:   ksp->ops->view           = NULL;
118:   PetscFunctionReturn(PETSC_SUCCESS);
119: }