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: {
43: DM dm;
47: SNESGetDM(snes,&dm);
48: DMSNESSetObjective(dm,obj,ctx);
49: return(0);
50: }
52: /*@C
53: SNESGetObjective - Returns the objective function.
55: Not Collective
57: Input Parameter:
58: . snes - the SNES context
60: Output Parameters:
61: + obj - objective evaluation routine (or NULL); see SNESObjectFunction for details
62: - ctx - the function context (or NULL)
64: Level: advanced
66: .seealso: SNESSetObjective(), SNESGetSolution()
67: @*/
68: PetscErrorCode SNESGetObjective(SNES snes,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
69: {
71: DM dm;
75: SNESGetDM(snes,&dm);
76: DMSNESGetObjective(dm,obj,ctx);
77: return(0);
78: }
80: /*@C
81: SNESComputeObjective - Computes the objective.
83: Collective on SNES
85: Input Parameters:
86: + snes - the SNES context
87: - X - the state vector
89: Output Parameter:
90: . ob - the objective value
92: Level: advanced
94: .seealso: SNESSetObjective(), SNESGetSolution()
95: @*/
96: PetscErrorCode SNESComputeObjective(SNES snes,Vec X,PetscReal *ob)
97: {
99: DM dm;
100: DMSNES sdm;
106: SNESGetDM(snes,&dm);
107: DMGetDMSNES(dm,&sdm);
108: if (sdm->ops->computeobjective) {
109: PetscLogEventBegin(SNES_ObjectiveEval,snes,X,0,0);
110: (sdm->ops->computeobjective)(snes,X,ob,sdm->objectivectx);
111: PetscLogEventEnd(SNES_ObjectiveEval,snes,X,0,0);
112: } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetObjective() before SNESComputeObjective().");
113: return(0);
114: }
116: /*@C
117: SNESObjectiveComputeFunctionDefaultFD - Computes the gradient of a user provided objective
119: Collective on SNES
121: Input Parameters:
122: + snes - the SNES context
123: . X - the state vector
124: - ctx - the (ignored) function context
126: Output Parameter:
127: . F - the function value
129: Options Database Key:
130: + -snes_fd_function_eps - The differencing parameter
131: - -snes_fd_function - Compute function from user provided objective with finite difference
133: Notes:
134: SNESObjectiveComputeFunctionDefaultFD is similar in character to SNESComputeJacobianDefault.
135: Therefore, it should be used for debugging purposes only. Using it in conjunction with
136: SNESComputeJacobianDefault is excessively costly and produces a Jacobian that is quite
137: noisy. This is often necessary, but should be done with a grain of salt, even when debugging
138: small problems.
140: Note that this uses quadratic interpolation of the objective to form each value in the function.
142: Level: advanced
144: .seealso: SNESSetFunction(), SNESComputeObjective(), SNESComputeJacobianDefault()
145: @*/
146: PetscErrorCode SNESObjectiveComputeFunctionDefaultFD(SNES snes,Vec X,Vec F,void *ctx)
147: {
148: Vec Xh;
150: PetscInt i,N,start,end;
151: PetscReal ob,ob1,ob2,ob3,fob,dx,eps=1e-6;
152: PetscScalar fv,xv;
155: VecDuplicate(X,&Xh);
156: PetscOptionsBegin(PetscObjectComm((PetscObject)snes),((PetscObject)snes)->prefix,"Differencing parameters","SNES");
157: PetscOptionsReal("-snes_fd_function_eps","Tolerance for nonzero entries in fd function","None",eps,&eps,NULL);
158: PetscOptionsEnd();
159: VecSet(F,0.);
161: VecNorm(X,NORM_2,&fob);
163: VecGetSize(X,&N);
164: VecGetOwnershipRange(X,&start,&end);
165: SNESComputeObjective(snes,X,&ob);
167: if (fob > 0.) dx =1e-6*fob;
168: else dx = 1e-6;
170: for (i=0; i<N; i++) {
171: /* compute the 1st value */
172: VecCopy(X,Xh);
173: if (i>= start && i<end) {
174: xv = dx;
175: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
176: }
177: VecAssemblyBegin(Xh);
178: VecAssemblyEnd(Xh);
179: SNESComputeObjective(snes,Xh,&ob1);
181: /* compute the 2nd value */
182: VecCopy(X,Xh);
183: if (i>= start && i<end) {
184: xv = 2.*dx;
185: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
186: }
187: VecAssemblyBegin(Xh);
188: VecAssemblyEnd(Xh);
189: SNESComputeObjective(snes,Xh,&ob2);
191: /* compute the 3rd value */
192: VecCopy(X,Xh);
193: if (i>= start && i<end) {
194: xv = -dx;
195: VecSetValues(Xh,1,&i,&xv,ADD_VALUES);
196: }
197: VecAssemblyBegin(Xh);
198: VecAssemblyEnd(Xh);
199: SNESComputeObjective(snes,Xh,&ob3);
201: if (i >= start && i<end) {
202: /* set this entry to be the gradient of the objective */
203: fv = (-ob2 + 6.*ob1 - 3.*ob -2.*ob3) / (6.*dx);
204: if (PetscAbsScalar(fv) > eps) {
205: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
206: } else {
207: fv = 0.;
208: VecSetValues(F,1,&i,&fv,INSERT_VALUES);
209: }
210: }
211: }
213: VecDestroy(&Xh);
215: VecAssemblyBegin(F);
216: VecAssemblyEnd(F);
217: return(0);
218: }