Actual source code: snesob.c
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: . obj - real to hold the objective value
14: - ctx - optional user-defined objective context
16: Level: advanced
18: .seealso: SNESSetFunction(), SNESGetFunction(), SNESSetObjective(), SNESGetObjective(), SNESJacobianFunction, SNESFunction
19: M*/
21: /*@C
22: SNESSetObjective - Sets the objective function minimized by some of the SNES linesearch methods.
24: Logically Collective on SNES
26: Input Parameters:
27: + snes - the SNES context
28: . obj - objective evaluation routine; see SNESObjectiveFunction for details
29: - ctx - [optional] user-defined context for private data for the
30: function evaluation routine (may be NULL)
32: Level: intermediate
34: Note: This is not used in the SNESLINESEARCHCP line search.
36: If not provided then this defaults to the two norm of the function evaluation (set with SNESSetFunction())
38: .seealso: SNESGetObjective(), SNESComputeObjective(), SNESSetFunction(), SNESSetJacobian(), SNESObjectiveFunction
39: @*/
40: PetscErrorCode SNESSetObjective(SNES snes,PetscErrorCode (*obj)(SNES,Vec,PetscReal*,void*),void *ctx)
41: {
42: DM dm;
45: SNESGetDM(snes,&dm);
46: DMSNESSetObjective(dm,obj,ctx);
47: return 0;
48: }
50: /*@C
51: SNESGetObjective - Returns the objective function.
53: Not Collective
55: Input Parameter:
56: . snes - the SNES context
58: Output Parameters:
59: + obj - objective evaluation routine (or NULL); see SNESObjectFunction for details
60: - ctx - the function context (or NULL)
62: Level: advanced
64: .seealso: SNESSetObjective(), SNESGetSolution()
65: @*/
66: PetscErrorCode SNESGetObjective(SNES snes,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
67: {
68: DM dm;
71: SNESGetDM(snes,&dm);
72: DMSNESGetObjective(dm,obj,ctx);
73: return 0;
74: }
76: /*@C
77: SNESComputeObjective - Computes the objective.
79: Collective on SNES
81: Input Parameters:
82: + snes - the SNES context
83: - X - the state vector
85: Output Parameter:
86: . ob - the objective value
88: Level: advanced
90: .seealso: SNESSetObjective(), SNESGetSolution()
91: @*/
92: PetscErrorCode SNESComputeObjective(SNES snes,Vec X,PetscReal *ob)
93: {
94: DM dm;
95: DMSNES sdm;
100: SNESGetDM(snes,&dm);
101: DMGetDMSNES(dm,&sdm);
102: if (sdm->ops->computeobjective) {
103: PetscLogEventBegin(SNES_ObjectiveEval,snes,X,0,0);
104: (sdm->ops->computeobjective)(snes,X,ob,sdm->objectivectx);
105: PetscLogEventEnd(SNES_ObjectiveEval,snes,X,0,0);
106: } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetObjective() before SNESComputeObjective().");
107: return 0;
108: }
110: /*@C
111: SNESObjectiveComputeFunctionDefaultFD - Computes the gradient of a user provided objective
113: Collective on SNES
115: Input Parameters:
116: + snes - the SNES context
117: . X - the state vector
118: - ctx - the (ignored) function context
120: Output Parameter:
121: . F - the function value
123: Options Database Key:
124: + -snes_fd_function_eps - The differencing parameter
125: - -snes_fd_function - Compute function from user provided objective with finite difference
127: Notes:
128: SNESObjectiveComputeFunctionDefaultFD is similar in character to SNESComputeJacobianDefault.
129: Therefore, it should be used for debugging purposes only. Using it in conjunction with
130: SNESComputeJacobianDefault is excessively costly and produces a Jacobian that is quite
131: noisy. This is often necessary, but should be done with a grain of salt, even when debugging
132: small problems.
134: Note that this uses quadratic interpolation of the objective to form each value in the function.
136: Level: advanced
138: .seealso: SNESSetFunction(), SNESComputeObjective(), SNESComputeJacobianDefault()
139: @*/
140: PetscErrorCode SNESObjectiveComputeFunctionDefaultFD(SNES snes,Vec X,Vec F,void *ctx)
141: {
142: Vec Xh;
144: PetscInt i,N,start,end;
145: PetscReal ob,ob1,ob2,ob3,fob,dx,eps=1e-6;
146: PetscScalar fv,xv;
148: VecDuplicate(X,&Xh);
149: PetscOptionsBegin(PetscObjectComm((PetscObject)snes),((PetscObject)snes)->prefix,"Differencing parameters","SNES");
150: PetscOptionsReal("-snes_fd_function_eps","Tolerance for nonzero entries in fd function","None",eps,&eps,NULL);
151: PetscOptionsEnd();
152: VecSet(F,0.);
154: VecNorm(X,NORM_2,&fob);
156: VecGetSize(X,&N);
157: VecGetOwnershipRange(X,&start,&end);
158: SNESComputeObjective(snes,X,&ob);
160: if (fob > 0.) dx =1e-6*fob;
161: else dx = 1e-6;
163: for (i=0; i<N; i++) {
164: /* compute the 1st value */
165: VecCopy(X,Xh);
166: if (i>= start && i<end) {
167: xv = dx;
168: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
169: }
170: VecAssemblyBegin(Xh);
171: VecAssemblyEnd(Xh);
172: SNESComputeObjective(snes,Xh,&ob1);
174: /* compute the 2nd value */
175: VecCopy(X,Xh);
176: if (i>= start && i<end) {
177: xv = 2.*dx;
178: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
179: }
180: VecAssemblyBegin(Xh);
181: VecAssemblyEnd(Xh);
182: SNESComputeObjective(snes,Xh,&ob2);
184: /* compute the 3rd value */
185: VecCopy(X,Xh);
186: if (i>= start && i<end) {
187: xv = -dx;
188: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
189: }
190: VecAssemblyBegin(Xh);
191: VecAssemblyEnd(Xh);
192: SNESComputeObjective(snes,Xh,&ob3);
194: if (i >= start && i<end) {
195: /* set this entry to be the gradient of the objective */
196: fv = (-ob2 + 6.*ob1 - 3.*ob -2.*ob3) / (6.*dx);
197: if (PetscAbsScalar(fv) > eps) {
198: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
199: } else {
200: fv = 0.;
201: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
202: }
203: }
204: }
206: VecDestroy(&Xh);
208: VecAssemblyBegin(F);
209: VecAssemblyEnd(F);
210: return 0;
211: }