Actual source code: ksponly.c

  1: #include <petsc/private/snesimpl.h>

  3: typedef struct {
  4:   PetscBool transpose_solve;
  5: } SNES_KSPONLY;

  7: static PetscErrorCode SNESSolve_KSPONLY(SNES snes)
  8: {
  9:   SNES_KSPONLY    *ksponly = (SNES_KSPONLY *)snes->data;
 10:   PetscInt         lits;
 11:   Vec              Y, X, F;
 12:   SNESNormSchedule normschedule;

 14:   PetscFunctionBegin;
 15:   PetscCheck(!snes->xl && !snes->xu && !snes->ops->computevariablebounds, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);

 17:   snes->numFailures            = 0;
 18:   snes->numLinearSolveFailures = 0;
 19:   snes->reason                 = SNES_CONVERGED_ITERATING;
 20:   snes->iter                   = 0;
 21:   snes->norm                   = 0.0;

 23:   X = snes->vec_sol;
 24:   F = snes->vec_func;
 25:   Y = snes->vec_sol_update;

 27:   if (!snes->vec_func_init_set) PetscCall(SNESComputeFunction(snes, X, F));
 28:   else snes->vec_func_init_set = PETSC_FALSE;

 30:   PetscCall(SNESGetNormSchedule(snes, &normschedule));
 31:   if (snes->numbermonitors && (normschedule == SNES_NORM_ALWAYS || normschedule == SNES_NORM_INITIAL_ONLY || normschedule == SNES_NORM_INITIAL_FINAL_ONLY)) {
 32:     PetscReal fnorm;
 33:     PetscCall(VecNorm(F, NORM_2, &fnorm));
 34:     SNESCheckFunctionDomainError(snes, fnorm);
 35:     PetscCall(SNESMonitor(snes, 0, fnorm));
 36:   }

 38:   /* Call general purpose update function */
 39:   PetscTryTypeMethod(snes, update, 0);

 41:   /* Solve J Y = F, where J is Jacobian matrix */
 42:   PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
 43:   SNESCheckJacobianDomainError(snes);

 45:   PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
 46:   if (ksponly->transpose_solve) {
 47:     PetscCall(KSPSolveTranspose(snes->ksp, F, Y));
 48:   } else {
 49:     PetscCall(KSPSolve(snes->ksp, F, Y));
 50:   }
 51:   snes->reason = SNES_CONVERGED_ITS;
 52:   SNESCheckKSPSolve(snes);

 54:   PetscCall(KSPGetIterationNumber(snes->ksp, &lits));
 55:   PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n", snes->iter, lits));
 56:   snes->iter++;

 58:   /* Take the computed step. */
 59:   PetscCall(VecAXPY(X, -1.0, Y));

 61:   if (snes->numbermonitors && (normschedule == SNES_NORM_ALWAYS || normschedule == SNES_NORM_FINAL_ONLY || normschedule == SNES_NORM_INITIAL_FINAL_ONLY)) {
 62:     PetscReal fnorm;
 63:     PetscCall(SNESComputeFunction(snes, X, F));
 64:     PetscCall(VecNorm(F, NORM_2, &fnorm));
 65:     SNESCheckFunctionDomainError(snes, fnorm);
 66:     PetscCall(SNESMonitor(snes, 1, fnorm));
 67:   }
 68:   PetscFunctionReturn(PETSC_SUCCESS);
 69: }

 71: static PetscErrorCode SNESSetUp_KSPONLY(SNES snes)
 72: {
 73:   PetscFunctionBegin;
 74:   PetscCall(SNESSetUpMatrices(snes));
 75:   PetscFunctionReturn(PETSC_SUCCESS);
 76: }

 78: static PetscErrorCode SNESDestroy_KSPONLY(SNES snes)
 79: {
 80:   PetscFunctionBegin;
 81:   PetscCall(PetscFree(snes->data));
 82:   PetscFunctionReturn(PETSC_SUCCESS);
 83: }

 85: /*MC
 86:    SNESKSPONLY - Nonlinear solver that performs one Newton step with `KSPSolve()` and does not compute any norms.

 88:    Level: beginner

 90:    Note:
 91:    The main purpose of this solver is to solve linear problems using the `SNES` interface, without
 92:    any additional overhead in the form of vector norm operations.

 94: .seealso: [](ch_snes), `SNES`, `SNESType`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`, `SNESKSPTRANSPOSEONLY`
 95: M*/
 96: PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes)
 97: {
 98:   SNES_KSPONLY *ksponly;

100:   PetscFunctionBegin;
101:   snes->ops->setup          = SNESSetUp_KSPONLY;
102:   snes->ops->solve          = SNESSolve_KSPONLY;
103:   snes->ops->destroy        = SNESDestroy_KSPONLY;
104:   snes->ops->setfromoptions = NULL;
105:   snes->ops->view           = NULL;
106:   snes->ops->reset          = NULL;

108:   snes->usesksp = PETSC_TRUE;
109:   snes->usesnpc = PETSC_FALSE;

111:   snes->alwayscomputesfinalresidual = PETSC_FALSE;

113:   PetscCall(SNESParametersInitialize(snes));

115:   PetscCall(PetscNew(&ksponly));
116:   snes->data = (void *)ksponly;
117:   PetscFunctionReturn(PETSC_SUCCESS);
118: }

120: /*MC
121:    SNESKSPTRANSPOSEONLY - Nonlinear solver that performs one Newton step with `KSPSolveTranspose()` and does not compute any norms.

123:    Level: beginner

125:    Note:
126:    The main purpose of this solver is to solve transposed linear problems using the `SNES` interface, without
127:    any additional overhead in the form of vector operations within adjoint solvers.

129: .seealso: [](ch_snes), `SNES`, `SNESType`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESKS`, `SNESNEWTONLS`, `SNESNEWTONTR`
130: M*/
131: PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes)
132: {
133:   SNES_KSPONLY *kspo;

135:   PetscFunctionBegin;
136:   PetscCall(SNESCreate_KSPONLY(snes));
137:   PetscCall(PetscObjectChangeTypeName((PetscObject)snes, SNESKSPTRANSPOSEONLY));
138:   kspo                  = (SNES_KSPONLY *)snes->data;
139:   kspo->transpose_solve = PETSC_TRUE;
140:   PetscFunctionReturn(PETSC_SUCCESS);
141: }