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;


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

 21:   X = snes->vec_sol;
 22:   F = snes->vec_func;
 23:   Y = snes->vec_sol_update;

 25:   if (!snes->vec_func_init_set) {
 26:     SNESComputeFunction(snes,X,F);
 27:   } else snes->vec_func_init_set = PETSC_FALSE;

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

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

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

 43:   SNESCheckJacobianDomainerror(snes);

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

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

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

 69: static PetscErrorCode SNESSetUp_KSPONLY(SNES snes)
 70: {
 71:   SNESSetUpMatrices(snes);
 72:   return 0;
 73: }

 75: static PetscErrorCode SNESDestroy_KSPONLY(SNES snes)
 76: {
 77:   PetscFree(snes->data);
 78:   return 0;
 79: }

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

 87:    Level: beginner

 89: .seealso:  SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR
 90: M*/
 91: PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes)
 92: {
 93:   SNES_KSPONLY   *ksponly;

 95:   snes->ops->setup          = SNESSetUp_KSPONLY;
 96:   snes->ops->solve          = SNESSolve_KSPONLY;
 97:   snes->ops->destroy        = SNESDestroy_KSPONLY;
 98:   snes->ops->setfromoptions = NULL;
 99:   snes->ops->view           = NULL;
100:   snes->ops->reset          = NULL;

102:   snes->usesksp = PETSC_TRUE;
103:   snes->usesnpc = PETSC_FALSE;

105:   snes->alwayscomputesfinalresidual = PETSC_FALSE;

107:   PetscNewLog(snes,&ksponly);
108:   snes->data = (void*)ksponly;
109:   return 0;
110: }

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

117:    Level: beginner

119: .seealso:  SNESCreate(), SNES, SNESSetType(), SNESKSPTRANSPOSEONLY, SNESNEWTONLS, SNESNEWTONTR
120: M*/
121: PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes)
122: {
123:   SNES_KSPONLY   *kspo;

125:   SNESCreate_KSPONLY(snes);
126:   PetscObjectChangeTypeName((PetscObject)snes,SNESKSPTRANSPOSEONLY);
127:   kspo = (SNES_KSPONLY*)snes->data;
128:   kspo->transpose_solve = PETSC_TRUE;
129:   return 0;
130: }