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