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 } else { 109 SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetObjective() before SNESComputeObjective()."); 110 } 111 PetscFunctionReturn(0); 112 } 113 114 115 #undef __FUNCT__ 116 #define __FUNCT__ "SNESDefaultObjectiveComputeFunctionFD" 117 PetscErrorCode SNESDefaultObjectiveComputeFunctionFD(SNES snes,Vec X,Vec F,void *ctx) { 118 /* Quadratically interpolates the change in the objective based upon a change in a particular direction */ 119 120 Vec Xh; 121 PetscErrorCode ierr; 122 PetscInt i,N,start,end; 123 PetscReal ob,ob1,ob2,ob3,fob,dx,eps=1e-6; 124 PetscScalar fv,xv; 125 126 PetscFunctionBegin; 127 ierr = VecDuplicate(X,&Xh);CHKERRQ(ierr); 128 ierr = PetscOptionsReal("-snes_fd_function_eps","Tolerance for nonzero entries in fd function","None",eps,&eps,PETSC_NULL);CHKERRQ(ierr); 129 ierr = VecSet(F,0.);CHKERRQ(ierr); 130 131 ierr = VecNorm(X,NORM_2,&fob);CHKERRQ(ierr); 132 133 ierr = VecGetSize(X,&N);CHKERRQ(ierr); 134 ierr = VecGetOwnershipRange(X,&start,&end);CHKERRQ(ierr); 135 ierr = SNESComputeObjective(snes,X,&ob);CHKERRQ(ierr); 136 137 if (fob > 0.) { 138 dx =1e-6*fob; 139 } else { 140 dx = 1e-6; 141 } 142 143 for (i=0; i<N; i++) { 144 /* compute the 1st value */ 145 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 146 if (i>= start && i<end) { 147 xv = dx; 148 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 149 } 150 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 151 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 152 ierr = SNESComputeObjective(snes,Xh,&ob1);CHKERRQ(ierr); 153 154 /* compute the 2nd value */ 155 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 156 if (i>= start && i<end) { 157 xv = 2.*dx; 158 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 159 } 160 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 161 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 162 ierr = SNESComputeObjective(snes,Xh,&ob2);CHKERRQ(ierr); 163 164 /* compute the 3rd value */ 165 ierr = VecCopy(X,Xh);CHKERRQ(ierr); 166 if (i>= start && i<end) { 167 xv = -dx; 168 ierr = VecSetValues(Xh,1,&i,&xv,ADD_VALUES);CHKERRQ(ierr); 169 } 170 ierr = VecAssemblyBegin(Xh);CHKERRQ(ierr); 171 ierr = VecAssemblyEnd(Xh);CHKERRQ(ierr); 172 ierr = SNESComputeObjective(snes,Xh,&ob3);CHKERRQ(ierr); 173 174 if (i >= start && i<end) { 175 /* set this entry to be the gradient of the objective */ 176 fv = (-ob2 + 6.*ob1 - 3.*ob -2.*ob3) / (6.*dx); 177 if (PetscAbsScalar(fv) > eps) { 178 ierr = VecSetValues(F,1,&i,&fv,INSERT_VALUES);CHKERRQ(ierr); 179 } else { 180 fv = 0.; 181 ierr = VecSetValues(F,1,&i,&fv,INSERT_VALUES);CHKERRQ(ierr); 182 } 183 } 184 } 185 186 ierr = VecDestroy(&Xh);CHKERRQ(ierr); 187 188 ierr = VecAssemblyBegin(F);CHKERRQ(ierr); 189 ierr = VecAssemblyEnd(F);CHKERRQ(ierr); 190 PetscFunctionReturn(0); 191 } 192