Actual source code: ksponly.c

petsc-3.13.6 2020-09-29
Report Typos and Errors
  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;
 11:   PetscInt       lits;
 12:   Vec            Y,X,F;

 15:   if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(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) {
 28:     SNESComputeFunction(snes,X,F);
 29:   } else snes->vec_func_init_set = PETSC_FALSE;

 31:   if (snes->numbermonitors) {
 32:     PetscReal fnorm;
 33:     VecNorm(F,NORM_2,&fnorm);
 34:     SNESMonitor(snes,0,fnorm);
 35:   }

 37:   /* Call general purpose update function */
 38:   if (snes->ops->update) {
 39:     (*snes->ops->update)(snes, 0);
 40:   }

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

 45:   SNESCheckJacobianDomainerror(snes);

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

 56:   KSPGetIterationNumber(snes->ksp,&lits);
 57:   PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);
 58:   snes->iter++;

 60:   /* Take the computed step. */
 61:   VecAXPY(X,-1.0,Y);
 62:   if (snes->numbermonitors) {
 63:     PetscReal fnorm;
 64:     SNESComputeFunction(snes,X,F);
 65:     VecNorm(F,NORM_2,&fnorm);
 66:     SNESMonitor(snes,1,fnorm);
 67:   }
 68:   return(0);
 69: }

 71: static PetscErrorCode SNESSetUp_KSPONLY(SNES snes)
 72: {

 76:   SNESSetUpMatrices(snes);
 77:   return(0);
 78: }

 80: static PetscErrorCode SNESDestroy_KSPONLY(SNES snes)
 81: {

 85:   PetscFree(snes->data);
 86:   return(0);
 87: }

 89: /* -------------------------------------------------------------------------- */
 90: /*MC
 91:       SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms.
 92:       The main purpose of this solver is to solve linear problems using the SNES interface, without
 93:       any additional overhead in the form of vector operations.

 95:    Level: beginner

 97: .seealso:  SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR
 98: M*/
 99: PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes)
100: {
101:   SNES_KSPONLY   *ksponly;

105:   snes->ops->setup          = SNESSetUp_KSPONLY;
106:   snes->ops->solve          = SNESSolve_KSPONLY;
107:   snes->ops->destroy        = SNESDestroy_KSPONLY;
108:   snes->ops->setfromoptions = 0;
109:   snes->ops->view           = 0;
110:   snes->ops->reset          = 0;

112:   snes->usesksp = PETSC_TRUE;
113:   snes->usesnpc = PETSC_FALSE;

115:   snes->alwayscomputesfinalresidual = PETSC_FALSE;

117:   PetscNewLog(snes,&ksponly);
118:   snes->data = (void*)ksponly;
119:   return(0);
120: }

122: /*MC
123:       SNESKSPTRANSPOSEONLY - Nonlinear solver that only performs one Newton step and does not compute any norms.
124:       The main purpose of this solver is to solve transposed linear problems using the SNES interface, without
125:       any additional overhead in the form of vector operations within adjoint solvers.

127:    Level: beginner

129: .seealso:  SNESCreate(), SNES, SNESSetType(), SNESKSPTRANSPOSEONLY, SNESNEWTONLS, SNESNEWTONTR
130: M*/
131: PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes)
132: {
133:   SNES_KSPONLY   *kspo;

137:   SNESCreate_KSPONLY(snes);
138:   PetscObjectChangeTypeName((PetscObject)snes,SNESKSPTRANSPOSEONLY);
139:   kspo = (SNES_KSPONLY*)snes->data;
140:   kspo->transpose_solve = PETSC_TRUE;
141:   return(0);
142: }