Actual source code: snespc.c
petsc-dev 2014-02-02
2: #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/
3: #include <petscdmshell.h>
8: /*@
9: SNESApplyPC - Calls the function that has been set with SNESSetFunction().
11: Collective on SNES
13: Input Parameters:
14: + snes - the SNES context
15: - x - input vector
17: Output Parameter:
18: . y - function vector, as set by SNESSetFunction()
20: Notes:
21: SNESComputeFunction() should be called on X before SNESApplyPC() is called, as it is
22: with SNESComuteJacobian().
24: Level: developer
26: .keywords: SNES, nonlinear, compute, function
28: .seealso: SNESGetPC(),SNESSetPC(),SNESComputeFunction()
29: @*/
30: PetscErrorCode SNESApplyPC(SNES snes,Vec x,Vec f,PetscReal *fnorm,Vec y)
31: {
40: VecValidValues(x,2,PETSC_TRUE);
41: if (snes->pc) {
42: if (f) {
43: SNESSetInitialFunction(snes->pc,f);
44: }
45: VecCopy(x,y);
46: PetscLogEventBegin(SNES_NPCSolve,snes->pc,x,y,0);
47: SNESSolve(snes->pc,snes->vec_rhs,y);
48: PetscLogEventEnd(SNES_NPCSolve,snes->pc,x,y,0);
49: VecAYPX(y,-1.0,x);
50: VecValidValues(y,3,PETSC_FALSE);
51: return(0);
52: }
53: VecValidValues(y,3,PETSC_FALSE);
54: return(0);
55: }
59: PetscErrorCode SNESComputeFunctionDefaultPC(SNES snes,Vec X,Vec F) {
60: /* This is to be used as an argument to SNESMF -- NOT as a "function" */
61: SNESConvergedReason reason;
62: PetscErrorCode ierr;
65: if (snes->pc) {
66: SNESApplyPC(snes,X,NULL,NULL,F);
67: SNESGetConvergedReason(snes->pc,&reason);
68: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
69: SNESSetFunctionDomainError(snes);
70: }
71: } else {
72: SNESComputeFunction(snes,X,F);
73: }
74: return(0);
75: }
79: /*@
80: SNESGetPCFunction - Gets the function from a preconditioner after SNESSolve() has been called.
82: Collective on SNES
84: Input Parameters:
85: . snes - the SNES context
87: Output Parameter:
88: . F - function vector
89: . fnorm - the norm of F
91: Level: developer
93: .keywords: SNES, nonlinear, function
95: .seealso: SNESGetPC(),SNESSetPC(),SNESComputeFunction(),SNESApplyPC(),SNESSolve()
96: @*/
97: PetscErrorCode SNESGetPCFunction(SNES snes,Vec F,PetscReal *fnorm)
98: {
99: PetscErrorCode ierr;
100: PCSide npcside;
101: SNESFunctionType functype;
102: SNESNormSchedule normschedule;
103: Vec FPC,XPC;
106: if (snes->pc) {
107: SNESGetPCSide(snes->pc,&npcside);
108: SNESGetFunctionType(snes->pc,&functype);
109: SNESGetNormSchedule(snes->pc,&normschedule);
111: /* check if the function is valid based upon how the inner solver is preconditioned */
112: if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) {
113: SNESGetFunction(snes->pc,&FPC,NULL,NULL);
114: if (FPC) {
115: if (fnorm) {SNESGetFunctionNorm(snes->pc,fnorm);}
116: VecCopy(FPC,F);
117: } else {
118: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no function");
119: }
120: } else {
121: SNESGetSolution(snes->pc,&XPC);
122: if (XPC) {
123: SNESComputeFunction(snes->pc,XPC,F);
124: if (fnorm) {VecNorm(F,NORM_2,fnorm);}
125: } else {
126: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no solution");
127: }
128: }
129: } else {
130: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No preconditioner set");
131: }
133: return(0);
134: }