Actual source code: snespc.c
2: #include <petsc/private/snesimpl.h>
4: /*@
5: SNESApplyNPC - Calls `SNESSolve()` on preconditioner for the `SNES`
7: Collective
9: Input Parameters:
10: + snes - the `SNES` context
11: . x - input vector
12: - f - optional; the function evaluation on `x`
14: Output Parameter:
15: . y - function vector, as set by `SNESSetFunction()`
17: Level: developer
19: Note:
20: `SNESComputeFunction()` should be called on `x` before `SNESApplyNPC()` is called, as it is
21: with `SNESComuteJacobian()`.
23: .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()`
24: @*/
25: PetscErrorCode SNESApplyNPC(SNES snes, Vec x, Vec f, Vec y)
26: {
27: PetscFunctionBegin;
31: PetscCheckSameComm(snes, 1, x, 2);
32: PetscCheckSameComm(snes, 1, y, 4);
33: PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
34: if (snes->npc) {
35: if (f) PetscCall(SNESSetInitialFunction(snes->npc, f));
36: PetscCall(VecCopy(x, y));
37: PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, x, y, 0));
38: PetscCall(SNESSolve(snes->npc, snes->vec_rhs, y));
39: PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, x, y, 0));
40: PetscCall(VecAYPX(y, -1.0, x));
41: }
42: PetscFunctionReturn(PETSC_SUCCESS);
43: }
45: PetscErrorCode SNESComputeFunctionDefaultNPC(SNES snes, Vec X, Vec F)
46: {
47: /* This is to be used as an argument to SNESMF -- NOT as a "function" */
48: SNESConvergedReason reason;
50: PetscFunctionBegin;
51: if (snes->npc) {
52: PetscCall(SNESApplyNPC(snes, X, NULL, F));
53: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
54: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) PetscCall(SNESSetFunctionDomainError(snes));
55: } else {
56: PetscCall(SNESComputeFunction(snes, X, F));
57: }
58: PetscFunctionReturn(PETSC_SUCCESS);
59: }
61: /*@
62: SNESGetNPCFunction - Gets the current function value and its norm from a nonlinear preconditioner after `SNESSolve()` has been called on that `SNES`
64: Collective
66: Input Parameter:
67: . snes - the `SNES` context
69: Output Parameters:
70: + F - function vector
71: - fnorm - the norm of `F`
73: Level: developer
75: .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()`, `SNESApplyNPC()`, `SNESSolve()`
76: @*/
77: PetscErrorCode SNESGetNPCFunction(SNES snes, Vec F, PetscReal *fnorm)
78: {
79: PCSide npcside;
80: SNESFunctionType functype;
81: SNESNormSchedule normschedule;
82: Vec FPC, XPC;
84: PetscFunctionBegin;
87: PetscCheck(snes->npc, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "No nonlinear preconditioner set");
88: PetscCall(SNESGetNPCSide(snes->npc, &npcside));
89: PetscCall(SNESGetFunctionType(snes->npc, &functype));
90: PetscCall(SNESGetNormSchedule(snes->npc, &normschedule));
92: /* check if the function is valid based upon how the inner solver is preconditioned */
93: if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) {
94: PetscCall(SNESGetFunction(snes->npc, &FPC, NULL, NULL));
95: PetscCheck(FPC, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no function");
96: if (fnorm) PetscCall(VecNorm(FPC, NORM_2, fnorm));
97: PetscCall(VecCopy(FPC, F));
98: } else {
99: PetscCall(SNESGetSolution(snes->npc, &XPC));
100: PetscCheck(XPC, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no solution");
101: PetscCall(SNESComputeFunction(snes->npc, XPC, F));
102: if (fnorm) PetscCall(VecNorm(F, NORM_2, fnorm));
103: }
104: PetscFunctionReturn(PETSC_SUCCESS);
105: }