Actual source code: snespc.c
2: #include <petsc/private/snesimpl.h>
4: /*@
5: SNESApplyNPC - Calls SNESSolve() on preconditioner for the SNES
7: Collective on SNES
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: Notes:
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(x,2,PETSC_TRUE);
33: if (snes->npc) {
34: if (f) {
35: SNESSetInitialFunction(snes->npc,f);
36: }
37: VecCopy(x,y);
38: PetscLogEventBegin(SNES_NPCSolve,snes->npc,x,y,0);
39: SNESSolve(snes->npc,snes->vec_rhs,y);
40: PetscLogEventEnd(SNES_NPCSolve,snes->npc,x,y,0);
41: VecAYPX(y,-1.0,x);
42: return 0;
43: }
44: return 0;
45: }
47: PetscErrorCode SNESComputeFunctionDefaultNPC(SNES snes,Vec X,Vec F)
48: {
49: /* This is to be used as an argument to SNESMF -- NOT as a "function" */
50: SNESConvergedReason reason;
52: if (snes->npc) {
53: SNESApplyNPC(snes,X,NULL,F);
54: SNESGetConvergedReason(snes->npc,&reason);
55: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
56: SNESSetFunctionDomainError(snes);
57: }
58: } else {
59: SNESComputeFunction(snes,X,F);
60: }
61: return 0;
62: }
64: /*@
65: SNESGetNPCFunction - Gets the function from a preconditioner after SNESSolve() has been called.
67: Collective on SNES
69: Input Parameter:
70: . snes - the SNES context
72: Output Parameters:
73: + F - function vector
74: - fnorm - the norm of F
76: Level: developer
78: .seealso: SNESGetNPC(),SNESSetNPC(),SNESComputeFunction(),SNESApplyNPC(),SNESSolve()
79: @*/
80: PetscErrorCode SNESGetNPCFunction(SNES snes,Vec F,PetscReal *fnorm)
81: {
82: PCSide npcside;
83: SNESFunctionType functype;
84: SNESNormSchedule normschedule;
85: Vec FPC,XPC;
87: if (snes->npc) {
88: SNESGetNPCSide(snes->npc,&npcside);
89: SNESGetFunctionType(snes->npc,&functype);
90: 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: SNESGetFunction(snes->npc,&FPC,NULL,NULL);
95: if (FPC) {
96: if (fnorm) VecNorm(FPC,NORM_2,fnorm);
97: VecCopy(FPC,F);
98: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no function");
99: } else {
100: SNESGetSolution(snes->npc,&XPC);
101: if (XPC) {
102: SNESComputeFunction(snes->npc,XPC,F);
103: if (fnorm) VecNorm(F,NORM_2,fnorm);
104: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no solution");
105: }
106: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No preconditioner set");
107: return 0;
108: }