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: Note:
18: SNESComputeFunction() should be called on x before SNESApplyNPC() is called, as it is
19: with SNESComuteJacobian().
21: Level: developer
23: .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()`
24: @*/
25: PetscErrorCode SNESApplyNPC(SNES snes, Vec x, Vec f, Vec y)
26: {
32: VecValidValues_Internal(x, 2, PETSC_TRUE);
33: if (snes->npc) {
34: if (f) SNESSetInitialFunction(snes->npc, f);
35: VecCopy(x, y);
36: PetscLogEventBegin(SNES_NPCSolve, snes->npc, x, y, 0);
37: SNESSolve(snes->npc, snes->vec_rhs, y);
38: PetscLogEventEnd(SNES_NPCSolve, snes->npc, x, y, 0);
39: VecAYPX(y, -1.0, x);
40: return 0;
41: }
42: return 0;
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: if (snes->npc) {
51: SNESApplyNPC(snes, X, NULL, F);
52: SNESGetConvergedReason(snes->npc, &reason);
53: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) SNESSetFunctionDomainError(snes);
54: } else {
55: SNESComputeFunction(snes, X, F);
56: }
57: return 0;
58: }
60: /*@
61: SNESGetNPCFunction - Gets the current function value and its norm from a nonlinear preconditioner after `SNESSolve()` has been called on that `SNES`
63: Collective
65: Input Parameter:
66: . snes - the `SNES` context
68: Output Parameters:
69: + F - function vector
70: - fnorm - the norm of F
72: Level: developer
74: .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()`, `SNESApplyNPC()`, `SNESSolve()`
75: @*/
76: PetscErrorCode SNESGetNPCFunction(SNES snes, Vec F, PetscReal *fnorm)
77: {
78: PCSide npcside;
79: SNESFunctionType functype;
80: SNESNormSchedule normschedule;
81: Vec FPC, XPC;
83: if (snes->npc) {
84: SNESGetNPCSide(snes->npc, &npcside);
85: SNESGetFunctionType(snes->npc, &functype);
86: SNESGetNormSchedule(snes->npc, &normschedule);
88: /* check if the function is valid based upon how the inner solver is preconditioned */
89: if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) {
90: SNESGetFunction(snes->npc, &FPC, NULL, NULL);
91: if (FPC) {
92: if (fnorm) VecNorm(FPC, NORM_2, fnorm);
93: VecCopy(FPC, F);
94: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no function");
95: } else {
96: SNESGetSolution(snes->npc, &XPC);
97: if (XPC) {
98: SNESComputeFunction(snes->npc, XPC, F);
99: if (fnorm) VecNorm(F, NORM_2, fnorm);
100: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no solution");
101: }
102: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "No nonlinear preconditioner set");
103: return 0;
104: }