1 #include <petsc-private/snesimpl.h> 2 3 #undef __FUNCT__ 4 #define __FUNCT__ "SNESSetObjective" 5 /*@C 6 SNESSetObjective - Sets the objective function minimized by 7 the SNES methods. 8 9 Logically Collective on SNES 10 11 Input Parameters: 12 + snes - the SNES context 13 . func - objective evaluation routine 14 - ctx - [optional] user-defined context for private data for the 15 function evaluation routine (may be PETSC_NULL) 16 17 Calling sequence of func: 18 $ func (SNES snes,Vec x,PetscReal *obj,void *ctx); 19 20 + snes - the SNES context 21 . X - solution 22 . F - current function/gradient 23 . obj - real to hold the objective value 24 - ctx - optional user-defined objective context 25 26 Level: beginner 27 28 .keywords: SNES, nonlinear, set, objective 29 30 .seealso: SNESGetObjective(), SNESComputeObjective(), SNESSetFunction(), SNESSetJacobian() 31 @*/ 32 PetscErrorCode SNESSetObjective(SNES snes,SNESObjective func,void *ctx) 33 { 34 PetscErrorCode ierr; 35 DM dm; 36 37 PetscFunctionBegin; 38 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 39 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 40 ierr = DMSNESSetObjective(dm,func,ctx);CHKERRQ(ierr); 41 PetscFunctionReturn(0); 42 } 43 44 #undef __FUNCT__ 45 #define __FUNCT__ "SNESGetObjective" 46 /*@C 47 SNESGetObjective - Returns the objective function. 48 49 Not Collective 50 51 Input Parameter: 52 . snes - the SNES context 53 54 Output Parameter: 55 + func - the function (or PETSC_NULL) 56 - ctx - the function context (or PETSC_NULL) 57 58 Level: advanced 59 60 .keywords: SNES, nonlinear, get, objective 61 62 .seealso: SNESSetObjective(), SNESGetSolution() 63 @*/ 64 PetscErrorCode SNESGetObjective(SNES snes,SNESObjective *func,void **ctx) 65 { 66 PetscErrorCode ierr; 67 DM dm; 68 PetscFunctionBegin; 69 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 70 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 71 ierr = DMSNESGetObjective(dm,func,ctx);CHKERRQ(ierr); 72 PetscFunctionReturn(0); 73 } 74 75 #undef __FUNCT__ 76 #define __FUNCT__ "SNESComputeObjective" 77 /*@C 78 SNESComputeObjective - Computes the objective. 79 80 Collective on SNES 81 82 Input Parameter: 83 + snes - the SNES context 84 - X - the state vector 85 86 Output Parameter: 87 . ob - the objective value 88 89 Level: advanced 90 91 .keywords: SNES, nonlinear, compute, objective 92 93 .seealso: SNESSetObjective(), SNESGetSolution() 94 @*/ 95 PetscErrorCode SNESComputeObjective(SNES snes,Vec X,PetscReal *ob) 96 { 97 PetscErrorCode ierr; 98 DM dm; 99 DMSNES sdm; 100 PetscFunctionBegin; 101 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 102 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 103 PetscValidPointer(ob,3); 104 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 105 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 106 if (sdm->ops->computeobjective) { 107 ierr = (sdm->ops->computeobjective)(snes,X,ob,sdm->objectivectx);CHKERRQ(ierr); 108 } SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetObjective() before SNESComputeObjective()."); 109 PetscFunctionReturn(0); 110 } 111 112 113 #undef __FUNCT__ 114 #define __FUNCT__ "SNESDefaultObjectiveComputeFunctionFD" 115 PetscErrorCode SNESDefaultObjectiveComputeFunctionFD(SNES snes,Vec X,Vec F,void *ctx) { 116 /* Quadratically interpolates the change in the objective based upon a change in a particular direction */ 117 118 Vec Xh; 119 PetscErrorCode ierr; 120 PetscInt i,N,start,end; 121 PetscReal ob,ob1,ob2,ob3,fob,dx,eps=1e-6; 122 PetscScalar fv,xv; 123 124 PetscFunctionBegin; 125 ierr = VecDuplicate(X,&Xh);CHKERRQ(ierr); 126 ierr = PetscOptionsReal("-snes_fd_function_eps","Tolerance for nonzero entries in fd function","None",eps,&eps,PETSC_NULL);CHKERRQ(ierr); 127 ierr = VecSet(F,0.);CHKERRQ(ierr); 128 129 ierr = VecNorm(X,NORM_2,&fob);CHKERRQ(ierr); 130 131 ierr = VecGetSize(X,&N);CHKERRQ(ierr); 132 ierr = VecGetOwnershipRange(X,&start,&end);CHKERRQ(ierr); 133 ierr = SNESComputeObjective(snes,X,&ob);CHKERRQ(ierr); 134 135 if (fob > 0.) { 136 dx =1e-6*fob; 137 } else { 138 dx = 1e-6; 139 } 140 141 for (i=0; i<N; i++) { 142 /* compute the 1st value */ 143 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 144 if (i>= start && i<end) { 145 xv = dx; 146 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 147 } 148 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 149 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 150 ierr = SNESComputeObjective(snes,Xh,&ob1);CHKERRQ(ierr); 151 152 /* compute the 2nd value */ 153 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 154 if (i>= start && i<end) { 155 xv = 2.*dx; 156 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 157 } 158 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 159 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 160 ierr = SNESComputeObjective(snes,Xh,&ob2);CHKERRQ(ierr); 161 162 /* compute the 3rd value */ 163 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 164 if (i>= start && i<end) { 165 xv = -dx; 166 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 167 } 168 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 169 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 170 ierr = SNESComputeObjective(snes,Xh,&ob3);CHKERRQ(ierr); 171 172 if (i >= start && i<end) { 173 /* set this entry to be the gradient of the objective */ 174 fv = (-ob2 + 6.*ob1 - 3.*ob -2.*ob3) / (6.*dx); 175 if (PetscAbsScalar(fv) > eps) { 176 ierr = VecSetValues(F,1,&i,&fv,INSERT_VALUES);CHKERRQ(ierr); 177 } else { 178 fv = 0.; 179 ierr = VecSetValues(F,1,&i,&fv,INSERT_VALUES);CHKERRQ(ierr); 180 } 181 } 182 } 183 184 ierr = VecDestroy(&Xh);CHKERRQ(ierr); 185 186 ierr = VecAssemblyBegin(F);CHKERRQ(ierr); 187 ierr = VecAssemblyEnd(F);CHKERRQ(ierr); 188 PetscFunctionReturn(0); 189 } 190