xref: /petsc/src/snes/impls/ksponly/ksponly.c (revision 58d68138c660dfb4e9f5b03334792cd4f2ffd7cc)
1 #include <petsc/private/snesimpl.h> /*I   "petscsnes.h"   I*/
2 
3 typedef struct {
4   PetscBool transpose_solve;
5 } SNES_KSPONLY;
6 
7 static PetscErrorCode SNESSolve_KSPONLY(SNES snes) {
8   SNES_KSPONLY *ksponly = (SNES_KSPONLY *)snes->data;
9   PetscInt      lits;
10   Vec           Y, X, F;
11 
12   PetscFunctionBegin;
13   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);
14 
15   snes->numFailures            = 0;
16   snes->numLinearSolveFailures = 0;
17   snes->reason                 = SNES_CONVERGED_ITERATING;
18   snes->iter                   = 0;
19   snes->norm                   = 0.0;
20 
21   X = snes->vec_sol;
22   F = snes->vec_func;
23   Y = snes->vec_sol_update;
24 
25   if (!snes->vec_func_init_set) {
26     PetscCall(SNESComputeFunction(snes, X, F));
27   } else snes->vec_func_init_set = PETSC_FALSE;
28 
29   if (snes->numbermonitors) {
30     PetscReal fnorm;
31     PetscCall(VecNorm(F, NORM_2, &fnorm));
32     PetscCall(SNESMonitor(snes, 0, fnorm));
33   }
34 
35   /* Call general purpose update function */
36   PetscTryTypeMethod(snes, update, 0);
37 
38   /* Solve J Y = F, where J is Jacobian matrix */
39   PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
40 
41   SNESCheckJacobianDomainerror(snes);
42 
43   PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
44   if (ksponly->transpose_solve) {
45     PetscCall(KSPSolveTranspose(snes->ksp, F, Y));
46   } else {
47     PetscCall(KSPSolve(snes->ksp, F, Y));
48   }
49   snes->reason = SNES_CONVERGED_ITS;
50   SNESCheckKSPSolve(snes);
51 
52   PetscCall(KSPGetIterationNumber(snes->ksp, &lits));
53   PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n", snes->iter, lits));
54   snes->iter++;
55 
56   /* Take the computed step. */
57   PetscCall(VecAXPY(X, -1.0, Y));
58   if (snes->numbermonitors) {
59     PetscReal fnorm;
60     PetscCall(SNESComputeFunction(snes, X, F));
61     PetscCall(VecNorm(F, NORM_2, &fnorm));
62     PetscCall(SNESMonitor(snes, 1, fnorm));
63   }
64   PetscFunctionReturn(0);
65 }
66 
67 static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) {
68   PetscFunctionBegin;
69   PetscCall(SNESSetUpMatrices(snes));
70   PetscFunctionReturn(0);
71 }
72 
73 static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) {
74   PetscFunctionBegin;
75   PetscCall(PetscFree(snes->data));
76   PetscFunctionReturn(0);
77 }
78 
79 /* -------------------------------------------------------------------------- */
80 /*MC
81       SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms.
82       The main purpose of this solver is to solve linear problems using the SNES interface, without
83       any additional overhead in the form of vector operations.
84 
85    Level: beginner
86 
87 .seealso: `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`
88 M*/
89 PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) {
90   SNES_KSPONLY *ksponly;
91 
92   PetscFunctionBegin;
93   snes->ops->setup          = SNESSetUp_KSPONLY;
94   snes->ops->solve          = SNESSolve_KSPONLY;
95   snes->ops->destroy        = SNESDestroy_KSPONLY;
96   snes->ops->setfromoptions = NULL;
97   snes->ops->view           = NULL;
98   snes->ops->reset          = NULL;
99 
100   snes->usesksp = PETSC_TRUE;
101   snes->usesnpc = PETSC_FALSE;
102 
103   snes->alwayscomputesfinalresidual = PETSC_FALSE;
104 
105   PetscCall(PetscNewLog(snes, &ksponly));
106   snes->data = (void *)ksponly;
107   PetscFunctionReturn(0);
108 }
109 
110 /*MC
111       SNESKSPTRANSPOSEONLY - Nonlinear solver that only performs one Newton step and does not compute any norms.
112       The main purpose of this solver is to solve transposed linear problems using the SNES interface, without
113       any additional overhead in the form of vector operations within adjoint solvers.
114 
115    Level: beginner
116 
117 .seealso: `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESKSPTRANSPOSEONLY`, `SNESNEWTONLS`, `SNESNEWTONTR`
118 M*/
119 PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes) {
120   SNES_KSPONLY *kspo;
121 
122   PetscFunctionBegin;
123   PetscCall(SNESCreate_KSPONLY(snes));
124   PetscCall(PetscObjectChangeTypeName((PetscObject)snes, SNESKSPTRANSPOSEONLY));
125   kspo                  = (SNES_KSPONLY *)snes->data;
126   kspo->transpose_solve = PETSC_TRUE;
127   PetscFunctionReturn(0);
128 }
129