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 Notes: 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 CHKERRQ(VecValidValues(x,2,PETSC_TRUE)); 34 if (snes->npc) { 35 if (f) { 36 CHKERRQ(SNESSetInitialFunction(snes->npc,f)); 37 } 38 CHKERRQ(VecCopy(x,y)); 39 CHKERRQ(PetscLogEventBegin(SNES_NPCSolve,snes->npc,x,y,0)); 40 CHKERRQ(SNESSolve(snes->npc,snes->vec_rhs,y)); 41 CHKERRQ(PetscLogEventEnd(SNES_NPCSolve,snes->npc,x,y,0)); 42 CHKERRQ(VecAYPX(y,-1.0,x)); 43 PetscFunctionReturn(0); 44 } 45 PetscFunctionReturn(0); 46 } 47 48 PetscErrorCode SNESComputeFunctionDefaultNPC(SNES snes,Vec X,Vec F) 49 { 50 /* This is to be used as an argument to SNESMF -- NOT as a "function" */ 51 SNESConvergedReason reason; 52 53 PetscFunctionBegin; 54 if (snes->npc) { 55 CHKERRQ(SNESApplyNPC(snes,X,NULL,F)); 56 CHKERRQ(SNESGetConvergedReason(snes->npc,&reason)); 57 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 58 CHKERRQ(SNESSetFunctionDomainError(snes)); 59 } 60 } else { 61 CHKERRQ(SNESComputeFunction(snes,X,F)); 62 } 63 PetscFunctionReturn(0); 64 } 65 66 /*@ 67 SNESGetNPCFunction - Gets the function from a preconditioner after SNESSolve() has been called. 68 69 Collective on SNES 70 71 Input Parameter: 72 . snes - the SNES context 73 74 Output Parameters: 75 + F - function vector 76 - fnorm - the norm of F 77 78 Level: developer 79 80 .seealso: SNESGetNPC(),SNESSetNPC(),SNESComputeFunction(),SNESApplyNPC(),SNESSolve() 81 @*/ 82 PetscErrorCode SNESGetNPCFunction(SNES snes,Vec F,PetscReal *fnorm) 83 { 84 PCSide npcside; 85 SNESFunctionType functype; 86 SNESNormSchedule normschedule; 87 Vec FPC,XPC; 88 89 PetscFunctionBegin; 90 if (snes->npc) { 91 CHKERRQ(SNESGetNPCSide(snes->npc,&npcside)); 92 CHKERRQ(SNESGetFunctionType(snes->npc,&functype)); 93 CHKERRQ(SNESGetNormSchedule(snes->npc,&normschedule)); 94 95 /* check if the function is valid based upon how the inner solver is preconditioned */ 96 if (normschedule != SNES_NORM_NONE && normschedule != SNES_NORM_INITIAL_ONLY && (npcside == PC_RIGHT || functype == SNES_FUNCTION_UNPRECONDITIONED)) { 97 CHKERRQ(SNESGetFunction(snes->npc,&FPC,NULL,NULL)); 98 if (FPC) { 99 if (fnorm) CHKERRQ(VecNorm(FPC,NORM_2,fnorm)); 100 CHKERRQ(VecCopy(FPC,F)); 101 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no function"); 102 } else { 103 CHKERRQ(SNESGetSolution(snes->npc,&XPC)); 104 if (XPC) { 105 CHKERRQ(SNESComputeFunction(snes->npc,XPC,F)); 106 if (fnorm) CHKERRQ(VecNorm(F,NORM_2,fnorm)); 107 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Preconditioner has no solution"); 108 } 109 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"No preconditioner set"); 110 PetscFunctionReturn(0); 111 } 112