1 2 #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/ 3 #include <petscdmshell.h> 4 5 6 #undef __FUNCT__ 7 #define __FUNCT__ "SNESApplyPC" 8 /*@ 9 SNESApplyPC - Runs the nonlinear preconditioner and returns the step. 10 11 Collective on SNES 12 13 Input Parameters: 14 + snes - the SNES context 15 . x - input vector 16 - f - optional; the function evaluation on x 17 18 Output Parameter: 19 . y - function vector, as set by SNESSetFunction() 20 21 Notes: 22 SNESComputeFunction() should be called on x before SNESApplyPC() is called, as it is 23 with SNESComuteJacobian(). 24 25 Level: developer 26 27 .keywords: SNES, nonlinear, compute, function 28 29 .seealso: SNESGetPC(),SNESSetPC(),SNESComputeFunction() 30 @*/ 31 PetscErrorCode SNESApplyPC(SNES snes,Vec x,Vec f,Vec y) 32 { 33 PetscErrorCode ierr; 34 35 PetscFunctionBegin; 36 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 37 PetscValidHeaderSpecific(x,VEC_CLASSID,2); 38 PetscValidHeaderSpecific(y,VEC_CLASSID,3); 39 PetscCheckSameComm(snes,1,x,2); 40 PetscCheckSameComm(snes,1,y,3); 41 ierr = VecValidValues(x,2,PETSC_TRUE);CHKERRQ(ierr); 42 if (snes->pc) { 43 if (f) { 44 ierr = SNESSetInitialFunction(snes->pc,f);CHKERRQ(ierr); 45 } 46 ierr = VecCopy(x,y);CHKERRQ(ierr); 47 ierr = PetscLogEventBegin(SNES_NPCSolve,snes->pc,x,y,0);CHKERRQ(ierr); 48 ierr = SNESSolve(snes->pc,snes->vec_rhs,y);CHKERRQ(ierr); 49 ierr = PetscLogEventEnd(SNES_NPCSolve,snes->pc,x,y,0);CHKERRQ(ierr); 50 ierr = VecAYPX(y,-1.0,x);CHKERRQ(ierr); 51 ierr = VecValidValues(y,3,PETSC_FALSE);CHKERRQ(ierr); 52 PetscFunctionReturn(0); 53 } 54 ierr = VecValidValues(y,3,PETSC_FALSE);CHKERRQ(ierr); 55 PetscFunctionReturn(0); 56 } 57 58 #undef __FUNCT__ 59 #define __FUNCT__ "SNESComputeFunctionDefaultPC" 60 PetscErrorCode SNESComputeFunctionDefaultPC(SNES snes,Vec X,Vec F) { 61 /* This is to be used as an argument to SNESMF -- NOT as a "function" */ 62 SNESConvergedReason reason; 63 PetscErrorCode ierr; 64 65 PetscFunctionBegin; 66 if (snes->pc) { 67 ierr = SNESApplyPC(snes,X,NULL,F);CHKERRQ(ierr); 68 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 69 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 70 ierr = SNESSetFunctionDomainError(snes);CHKERRQ(ierr); 71 } 72 } else { 73 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 74 } 75 PetscFunctionReturn(0); 76 } 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "SNESGetPCFunction" 80 /*@ 81 SNESGetPCFunction - Gets the function from a preconditioner after SNESSolve() has been called. 82 83 Collective on SNES 84 85 Input Parameters: 86 . snes - the SNES context 87 88 Output Parameter: 89 . F - function vector 90 . fnorm - the norm of F 91 92 Level: developer 93 94 .keywords: SNES, nonlinear, function 95 96 .seealso: SNESGetPC(),SNESSetPC(),SNESComputeFunction(),SNESApplyPC(),SNESSolve() 97 @*/ 98 PetscErrorCode SNESGetPCFunction(SNES snes,Vec F,PetscReal *fnorm) 99 { 100 PetscErrorCode ierr; 101 PCSide npcside; 102 SNESFunctionType functype; 103 SNESNormSchedule normschedule; 104 Vec FPC,XPC; 105 106 PetscFunctionBegin; 107 if (snes->pc) { 108 ierr = SNESGetPCSide(snes->pc,&npcside);CHKERRQ(ierr); 109 ierr = SNESGetFunctionType(snes->pc,&functype);CHKERRQ(ierr); 110 ierr = SNESGetNormSchedule(snes->pc,&normschedule);CHKERRQ(ierr); 111 112 /* check if the function is valid based upon how the inner solver is preconditioned */ 113 if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) { 114 ierr = SNESGetFunction(snes->pc,&FPC,NULL,NULL);CHKERRQ(ierr); 115 if (FPC) { 116 if (fnorm) {ierr = SNESGetFunctionNorm(snes->pc,fnorm);CHKERRQ(ierr);} 117 ierr = VecCopy(FPC,F);CHKERRQ(ierr); 118 } else { 119 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no function"); 120 } 121 } else { 122 ierr = SNESGetSolution(snes->pc,&XPC);CHKERRQ(ierr); 123 if (XPC) { 124 ierr = SNESComputeFunction(snes->pc,XPC,F);CHKERRQ(ierr); 125 if (fnorm) {ierr = VecNorm(F,NORM_2,fnorm);CHKERRQ(ierr);} 126 } else { 127 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no solution"); 128 } 129 } 130 } else { 131 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No preconditioner set"); 132 } 133 134 PetscFunctionReturn(0); 135 } 136