1 2 #include <petsc/private/snesimpl.h> /*I "petscsnes.h" I*/ 3 4 /*@ 5 SNESApplyNPC - Calls SNESSolve() on preconditioner for the SNES 6 7 Collective on snes 8 9 Input Parameters: 10 + snes - the SNES context 11 . x - input vector 12 - f - optional; the function evaluation on x 13 14 Output Parameter: 15 . y - function vector, as set by SNESSetFunction() 16 17 Note: 18 SNESComputeFunction() should be called on x before SNESApplyNPC() is called, as it is 19 with SNESComuteJacobian(). 20 21 Level: developer 22 23 .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()` 24 @*/ 25 PetscErrorCode SNESApplyNPC(SNES snes, Vec x, Vec f, Vec y) 26 { 27 PetscFunctionBegin; 28 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1); 29 PetscValidHeaderSpecific(x, VEC_CLASSID, 2); 30 PetscValidHeaderSpecific(y, VEC_CLASSID, 4); 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 PetscFunctionReturn(0); 42 } 43 PetscFunctionReturn(0); 44 } 45 46 PetscErrorCode SNESComputeFunctionDefaultNPC(SNES snes, Vec X, Vec F) 47 { 48 /* This is to be used as an argument to SNESMF -- NOT as a "function" */ 49 SNESConvergedReason reason; 50 51 PetscFunctionBegin; 52 if (snes->npc) { 53 PetscCall(SNESApplyNPC(snes, X, NULL, F)); 54 PetscCall(SNESGetConvergedReason(snes->npc, &reason)); 55 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) PetscCall(SNESSetFunctionDomainError(snes)); 56 } else { 57 PetscCall(SNESComputeFunction(snes, X, F)); 58 } 59 PetscFunctionReturn(0); 60 } 61 62 /*@ 63 SNESGetNPCFunction - Gets the current function value and its norm from a nonlinear preconditioner after `SNESSolve()` has been called on that `SNES` 64 65 Collective on snes 66 67 Input Parameter: 68 . snes - the `SNES` context 69 70 Output Parameters: 71 + F - function vector 72 - fnorm - the norm of F 73 74 Level: developer 75 76 .seealso: `SNESGetNPC()`, `SNESSetNPC()`, `SNESComputeFunction()`, `SNESApplyNPC()`, `SNESSolve()` 77 @*/ 78 PetscErrorCode SNESGetNPCFunction(SNES snes, Vec F, PetscReal *fnorm) 79 { 80 PCSide npcside; 81 SNESFunctionType functype; 82 SNESNormSchedule normschedule; 83 Vec FPC, XPC; 84 85 PetscFunctionBegin; 86 if (snes->npc) { 87 PetscCall(SNESGetNPCSide(snes->npc, &npcside)); 88 PetscCall(SNESGetFunctionType(snes->npc, &functype)); 89 PetscCall(SNESGetNormSchedule(snes->npc, &normschedule)); 90 91 /* check if the function is valid based upon how the inner solver is preconditioned */ 92 if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) { 93 PetscCall(SNESGetFunction(snes->npc, &FPC, NULL, NULL)); 94 if (FPC) { 95 if (fnorm) PetscCall(VecNorm(FPC, NORM_2, fnorm)); 96 PetscCall(VecCopy(FPC, F)); 97 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no function"); 98 } else { 99 PetscCall(SNESGetSolution(snes->npc, &XPC)); 100 if (XPC) { 101 PetscCall(SNESComputeFunction(snes->npc, XPC, F)); 102 if (fnorm) PetscCall(VecNorm(F, NORM_2, fnorm)); 103 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Nonlinear preconditioner has no solution"); 104 } 105 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "No nonlinear preconditioner set"); 106 PetscFunctionReturn(0); 107 } 108