Actual source code: snesob.c
petsc-3.12.5 2020-03-29
1: #include <petsc/private/snesimpl.h>
3: /*MC
4: SNESObjectiveFunction - functional form used to convey the objective function to the nonlinear solver
6: Synopsis:
7: #include <petscsnes.h>
8: SNESObjectiveFunction(SNES snes,Vec x,PetscReal *obj,void *ctx);
10: Input Parameters:
11: + snes - the SNES context
12: . X - solution
13: . F - current function/gradient
14: . obj - real to hold the objective value
15: - ctx - optional user-defined objective context
17: Level: advanced
19: .seealso: SNESSetFunction(), SNESGetFunction(), SNESSetObjective(), SNESGetObjective(), SNESJacobianFunction, SNESFunction
20: M*/
23: /*@C
24: SNESSetObjective - Sets the objective function minimized by some of the SNES linesearch methods.
26: Logically Collective on SNES
28: Input Parameters:
29: + snes - the SNES context
30: . obj - objective evaluation routine; see SNESObjectiveFunction for details
31: - ctx - [optional] user-defined context for private data for the
32: function evaluation routine (may be NULL)
34: Level: intermediate
36: Note: This is not used in the SNESLINESEARCHCP line search.
38: If not provided then this defaults to the two norm of the function evaluation (set with SNESSetFunction())
40: .seealso: SNESGetObjective(), SNESComputeObjective(), SNESSetFunction(), SNESSetJacobian(), SNESObjectiveFunction
41: @*/
42: PetscErrorCode SNESSetObjective(SNES snes,PetscErrorCode (*obj)(SNES,Vec,PetscReal*,void*),void *ctx)
43: {
45: DM dm;
49: SNESGetDM(snes,&dm);
50: DMSNESSetObjective(dm,obj,ctx);
51: return(0);
52: }
54: /*@C
55: SNESGetObjective - Returns the objective function.
57: Not Collective
59: Input Parameter:
60: . snes - the SNES context
62: Output Parameter:
63: + obj - objective evaluation routine (or NULL); see SNESObjectFunction for details
64: - ctx - the function context (or NULL)
66: Level: advanced
68: .seealso: SNESSetObjective(), SNESGetSolution()
69: @*/
70: PetscErrorCode SNESGetObjective(SNES snes,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
71: {
73: DM dm;
77: SNESGetDM(snes,&dm);
78: DMSNESGetObjective(dm,obj,ctx);
79: return(0);
80: }
82: /*@C
83: SNESComputeObjective - Computes the objective.
85: Collective on SNES
87: Input Parameter:
88: + snes - the SNES context
89: - X - the state vector
91: Output Parameter:
92: . ob - the objective value
94: Level: advanced
96: .seealso: SNESSetObjective(), SNESGetSolution()
97: @*/
98: PetscErrorCode SNESComputeObjective(SNES snes,Vec X,PetscReal *ob)
99: {
101: DM dm;
102: DMSNES sdm;
108: SNESGetDM(snes,&dm);
109: DMGetDMSNES(dm,&sdm);
110: if (sdm->ops->computeobjective) {
111: PetscLogEventBegin(SNES_ObjectiveEval,snes,X,0,0);
112: (sdm->ops->computeobjective)(snes,X,ob,sdm->objectivectx);
113: PetscLogEventEnd(SNES_ObjectiveEval,snes,X,0,0);
114: } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetObjective() before SNESComputeObjective().");
115: return(0);
116: }
119: /*@C
120: SNESObjectiveComputeFunctionDefaultFD - Computes the gradient of a user provided objective
122: Collective on SNES
124: Input Parameter:
125: + snes - the SNES context
126: . X - the state vector
127: - ctx - the (ignored) function context
129: Output Parameter:
130: . F - the function value
132: Options Database Key:
133: + -snes_fd_function_eps - The differencing parameter
134: - -snes_fd_function - Compute function from user provided objective with finite difference
136: Notes:
137: SNESObjectiveComputeFunctionDefaultFD is similar in character to SNESComputeJacobianDefault.
138: Therefore, it should be used for debugging purposes only. Using it in conjunction with
139: SNESComputeJacobianDefault is excessively costly and produces a Jacobian that is quite
140: noisy. This is often necessary, but should be done with a grain of salt, even when debugging
141: small problems.
143: Note that this uses quadratic interpolation of the objective to form each value in the function.
145: Level: advanced
147: .seealso: SNESSetFunction(), SNESComputeObjective(), SNESComputeJacobianDefault()
148: @*/
149: PetscErrorCode SNESObjectiveComputeFunctionDefaultFD(SNES snes,Vec X,Vec F,void *ctx)
150: {
151: Vec Xh;
153: PetscInt i,N,start,end;
154: PetscReal ob,ob1,ob2,ob3,fob,dx,eps=1e-6;
155: PetscScalar fv,xv;
158: VecDuplicate(X,&Xh);
159: PetscOptionsBegin(PetscObjectComm((PetscObject)snes),((PetscObject)snes)->prefix,"Differencing parameters","SNES");
160: PetscOptionsReal("-snes_fd_function_eps","Tolerance for nonzero entries in fd function","None",eps,&eps,NULL);
161: PetscOptionsEnd();
162: VecSet(F,0.);
164: VecNorm(X,NORM_2,&fob);
166: VecGetSize(X,&N);
167: VecGetOwnershipRange(X,&start,&end);
168: SNESComputeObjective(snes,X,&ob);
170: if (fob > 0.) dx =1e-6*fob;
171: else dx = 1e-6;
173: for (i=0; i<N; i++) {
174: /* compute the 1st value */
175: VecCopy(X,Xh);
176: if (i>= start && i<end) {
177: xv = dx;
178: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
179: }
180: VecAssemblyBegin(Xh);
181: VecAssemblyEnd(Xh);
182: SNESComputeObjective(snes,Xh,&ob1);
184: /* compute the 2nd value */
185: VecCopy(X,Xh);
186: if (i>= start && i<end) {
187: xv = 2.*dx;
188: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
189: }
190: VecAssemblyBegin(Xh);
191: VecAssemblyEnd(Xh);
192: SNESComputeObjective(snes,Xh,&ob2);
194: /* compute the 3rd value */
195: VecCopy(X,Xh);
196: if (i>= start && i<end) {
197: xv = -dx;
198: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
199: }
200: VecAssemblyBegin(Xh);
201: VecAssemblyEnd(Xh);
202: SNESComputeObjective(snes,Xh,&ob3);
204: if (i >= start && i<end) {
205: /* set this entry to be the gradient of the objective */
206: fv = (-ob2 + 6.*ob1 - 3.*ob -2.*ob3) / (6.*dx);
207: if (PetscAbsScalar(fv) > eps) {
208: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
209: } else {
210: fv = 0.;
211: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
212: }
213: }
214: }
216: VecDestroy(&Xh);
218: VecAssemblyBegin(F);
219: VecAssemblyEnd(F);
220: return(0);
221: }