Actual source code: preonly.c

petsc-master 2020-07-04
Report Typos and Errors

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

  4: static PetscErrorCode KSPSetUp_PREONLY(KSP ksp)
  5: {
  7:   return(0);
  8: }

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

 17:   PCGetDiagonalScale(ksp->pc,&diagonalscale);
 18:   if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Krylov method %s does not support diagonal scaling",((PetscObject)ksp)->type_name);
 19:   if (!ksp->guess_zero) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_USER,"Running KSP of preonly doesn't make sense with nonzero initial guess\n\
 20:                you probably want a KSP type of Richardson");
 21:   ksp->its = 0;
 22:   KSP_PCApply(ksp,ksp->vec_rhs,ksp->vec_sol);
 23:   PCGetFailedReason(ksp->pc,&pcreason);
 24:   if (pcreason) {
 25:     ksp->reason = KSP_DIVERGED_PC_FAILED;
 26:   } else {
 27:     ksp->its    = 1;
 28:     ksp->reason = KSP_CONVERGED_ITS;
 29:     if (PetscDefined(USE_DEBUG)) {
 30:       PetscReal norm;
 31:       VecNorm(ksp->vec_sol,NORM_2,&norm);
 32:       if (PetscIsInfOrNanReal(norm)) {
 33:         ksp->reason = KSP_DIVERGED_NANORINF;
 34:       }
 35:     }
 36:   }
 37:   return(0);
 38: }

 40: static PetscErrorCode KSPMatSolve_PREONLY(KSP ksp, Mat B, Mat X)
 41: {
 43:   PetscBool      diagonalscale;
 44:   PCFailedReason pcreason;

 47:   PCGetDiagonalScale(ksp->pc,&diagonalscale);
 48:   if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Krylov method %s does not support diagonal scaling",((PetscObject)ksp)->type_name);
 49:   if (!ksp->guess_zero) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_USER,"Running KSP of preonly doesn't make sense with nonzero initial guess\n\
 50:                you probably want a KSP type of Richardson");
 51:   ksp->its = 0;
 52:   PCMatApply(ksp->pc,B,X);
 53:   PCGetFailedReason(ksp->pc,&pcreason);
 54:   if (pcreason) {
 55:     ksp->reason = KSP_DIVERGED_PC_FAILED;
 56:   } else {
 57:     ksp->its    = 1;
 58:     ksp->reason = KSP_CONVERGED_ITS;
 59:     if (PetscDefined(USE_DEBUG)) {
 60:       PetscReal norm;
 61:       MatNorm(X,NORM_FROBENIUS,&norm);
 62:       if (PetscIsInfOrNanReal(norm)) {
 63:         ksp->reason = KSP_DIVERGED_NANORINF;
 64:       }
 65:     }
 66:   }
 67:   return(0);
 68: }

 70: /*MC
 71:      KSPPREONLY - This implements a method that applies ONLY the preconditioner exactly once.
 72:                   This may be used in inner iterations, where it is desired to
 73:                   allow multiple iterations as well as the "0-iteration" case. It is
 74:                   commonly used with the direct solver preconditioners like PCLU and PCCHOLESKY

 76:    Options Database Keys:
 77: .   -ksp_type preonly

 79:    Level: beginner

 81:    Notes:
 82:     Since this does not involve an iteration the basic KSP parameters such as tolerances and iteration counts
 83:           do not apply

 85:     To apply multiple preconditioners in a simple iteration use KSPRICHARDSON

 87:    Developer Notes:
 88:     Even though this method does not use any norms, the user is allowed to set the KSPNormType to any value.
 89:     This is so the users does not have to change KSPNormType options when they switch from other KSP methods to this one.

 91: .seealso:  KSPCreate(), KSPSetType(), KSPType, KSP, KSPRICHARDSON, KSPCHEBYSHEV

 93: M*/

 95: PETSC_EXTERN PetscErrorCode KSPCreate_PREONLY(KSP ksp)
 96: {

100:   KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_LEFT,3);
101:   KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_RIGHT,2);
102:   KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,2);
103:   KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_RIGHT,2);
104:   KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,2);
105:   KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_RIGHT,2);
106:   KSPSetSupportedNorm(ksp,KSP_NORM_NATURAL,PC_LEFT,2);

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