Actual source code: snesgs.c
petsc-3.4.5 2014-06-29
1: #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/
3: typedef struct {
4: PetscInt sweeps; /* number of sweeps through the local subdomain before neighbor communication */
5: PetscInt max_its; /* maximum iterations of the inner pointblock solver */
6: PetscReal rtol; /* relative tolerance of the inner pointblock solver */
7: PetscReal abstol; /* absolute tolerance of the inner pointblock solver */
8: PetscReal stol; /* step tolerance of the inner pointblock solver */
9: } SNES_GS;
14: /*@
15: SNESGSSetTolerances - Sets various parameters used in convergence tests.
17: Logically Collective on SNES
19: Input Parameters:
20: + snes - the SNES context
21: . abstol - absolute convergence tolerance
22: . rtol - relative convergence tolerance
23: . stol - convergence tolerance in terms of the norm of the change in the solution between steps, || delta x || < stol*|| x ||
24: - maxit - maximum number of iterations
27: Options Database Keys:
28: + -snes_gs_atol <abstol> - Sets abstol
29: . -snes_gs_rtol <rtol> - Sets rtol
30: . -snes_gs_stol <stol> - Sets stol
31: - -snes_max_it <maxit> - Sets maxit
33: Level: intermediate
35: .keywords: SNES, nonlinear, gauss-seidel, set, convergence, tolerances
37: .seealso: SNESSetTrustRegionTolerance()
38: @*/
39: PetscErrorCode SNESGSSetTolerances(SNES snes,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit)
40: {
41: SNES_GS *gs = (SNES_GS*)snes->data;
46: if (abstol != PETSC_DEFAULT) {
47: if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %G must be non-negative",abstol);
48: gs->abstol = abstol;
49: }
50: if (rtol != PETSC_DEFAULT) {
51: if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %G must be non-negative and less than 1.0",rtol);
52: gs->rtol = rtol;
53: }
54: if (stol != PETSC_DEFAULT) {
55: if (stol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %G must be non-negative",stol);
56: gs->stol = stol;
57: }
58: if (maxit != PETSC_DEFAULT) {
59: if (maxit < 0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit);
60: gs->max_its = maxit;
61: }
62: return(0);
63: }
67: /*@
68: SNESGSGetTolerances - Gets various parameters used in convergence tests.
70: Not Collective
72: Input Parameters:
73: + snes - the SNES context
74: . atol - absolute convergence tolerance
75: . rtol - relative convergence tolerance
76: . stol - convergence tolerance in terms of the norm
77: of the change in the solution between steps
78: - maxit - maximum number of iterations
80: Notes:
81: The user can specify NULL for any parameter that is not needed.
83: Level: intermediate
85: .keywords: SNES, nonlinear, get, convergence, tolerances
87: .seealso: SNESSetTolerances()
88: @*/
89: PetscErrorCode SNESGSGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit)
90: {
91: SNES_GS *gs = (SNES_GS*)snes->data;
95: if (atol) *atol = gs->abstol;
96: if (rtol) *rtol = gs->rtol;
97: if (stol) *stol = gs->stol;
98: if (maxit) *maxit = gs->max_its;
99: return(0);
100: }
104: /*@
105: SNESGSSetSweeps - Sets the number of sweeps of GS to use.
107: Input Parameters:
108: + snes - the SNES context
109: - sweeps - the number of sweeps of GS to perform.
111: Level: intermediate
113: .keywords: SNES, nonlinear, set, Gauss-Siedel
115: .seealso: SNESSetGS(), SNESGetGS(), SNESSetPC(), SNESGSGetSweeps()
116: @*/
118: PetscErrorCode SNESGSSetSweeps(SNES snes, PetscInt sweeps)
119: {
120: SNES_GS *gs = (SNES_GS*)snes->data;
124: gs->sweeps = sweeps;
125: return(0);
126: }
130: /*@
131: SNESGSGetSweeps - Gets the number of sweeps GS will use.
133: Input Parameters:
134: . snes - the SNES context
136: Output Parameters:
137: . sweeps - the number of sweeps of GS to perform.
139: Level: intermediate
141: .keywords: SNES, nonlinear, set, Gauss-Siedel
143: .seealso: SNESSetGS(), SNESGetGS(), SNESSetPC(), SNESGSSetSweeps()
144: @*/
145: PetscErrorCode SNESGSGetSweeps(SNES snes, PetscInt * sweeps)
146: {
147: SNES_GS *gs = (SNES_GS*)snes->data;
151: *sweeps = gs->sweeps;
152: return(0);
153: }
158: PetscErrorCode SNESDefaultApplyGS(SNES snes,Vec X,Vec F,void *ctx)
159: {
161: /* see if there's a coloring on the DM */
162: return(0);
163: }
167: PetscErrorCode SNESReset_GS(SNES snes)
168: {
170: return(0);
171: }
175: PetscErrorCode SNESDestroy_GS(SNES snes)
176: {
180: SNESReset_GS(snes);
181: PetscFree(snes->data);
182: return(0);
183: }
187: PetscErrorCode SNESSetUp_GS(SNES snes)
188: {
190: return(0);
191: }
195: PetscErrorCode SNESSetFromOptions_GS(SNES snes)
196: {
197: SNES_GS *gs = (SNES_GS*)snes->data;
199: PetscInt sweeps,max_its=PETSC_DEFAULT;
200: PetscReal rtol=PETSC_DEFAULT,atol=PETSC_DEFAULT,stol=PETSC_DEFAULT;
201: PetscBool flg,flg1,flg2,flg3;
204: PetscOptionsHead("SNES GS options");
205: /* GS Options */
206: PetscOptionsInt("-snes_gs_sweeps","Number of sweeps of GS to apply","SNESComputeGS",gs->sweeps,&sweeps,&flg);
207: if (flg) {
208: SNESGSSetSweeps(snes,sweeps);
209: }
210: PetscOptionsReal("-snes_gs_atol","Number of sweeps of GS to apply","SNESComputeGS",gs->abstol,&atol,&flg);
211: PetscOptionsReal("-snes_gs_rtol","Number of sweeps of GS to apply","SNESComputeGS",gs->rtol,&rtol,&flg1);
212: PetscOptionsReal("-snes_gs_stol","Number of sweeps of GS to apply","SNESComputeGS",gs->stol,&stol,&flg2);
213: PetscOptionsInt("-snes_gs_max_it","Number of sweeps of GS to apply","SNESComputeGS",gs->max_its,&max_its,&flg3);
214: if (flg || flg1 || flg2 || flg3) {
215: SNESGSSetTolerances(snes,atol,rtol,stol,max_its);
216: }
217: PetscOptionsTail();
218: return(0);
219: }
223: PetscErrorCode SNESView_GS(SNES snes, PetscViewer viewer)
224: {
226: return(0);
227: }
231: PetscErrorCode SNESSolve_GS(SNES snes)
232: {
233: Vec F;
234: Vec X;
235: Vec B;
236: PetscInt i;
237: PetscReal fnorm;
239: SNESNormType normtype;
242: X = snes->vec_sol;
243: F = snes->vec_func;
244: B = snes->vec_rhs;
246: PetscObjectAMSTakeAccess((PetscObject)snes);
247: snes->iter = 0;
248: snes->norm = 0.;
249: PetscObjectAMSGrantAccess((PetscObject)snes);
250: snes->reason = SNES_CONVERGED_ITERATING;
252: SNESGetNormType(snes, &normtype);
253: if (normtype == SNES_NORM_FUNCTION || normtype == SNES_NORM_INITIAL_ONLY || normtype == SNES_NORM_INITIAL_FINAL_ONLY) {
254: /* compute the initial function and preconditioned update delX */
255: if (!snes->vec_func_init_set) {
256: SNESComputeFunction(snes,X,F);
257: if (snes->domainerror) {
258: snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
259: return(0);
260: }
261: } else snes->vec_func_init_set = PETSC_FALSE;
263: /* convergence test */
264: if (!snes->norm_init_set) {
265: VecNorm(F, NORM_2, &fnorm); /* fnorm <- ||F|| */
266: if (PetscIsInfOrNanReal(fnorm)) {
267: snes->reason = SNES_DIVERGED_FNORM_NAN;
268: return(0);
269: }
270: } else {
271: fnorm = snes->norm_init;
272: snes->norm_init_set = PETSC_FALSE;
273: }
274: PetscObjectAMSTakeAccess((PetscObject)snes);
275: snes->iter = 0;
276: snes->norm = fnorm;
277: PetscObjectAMSGrantAccess((PetscObject)snes);
278: SNESLogConvergenceHistory(snes,snes->norm,0);
279: SNESMonitor(snes,0,snes->norm);
281: /* set parameter for default relative tolerance convergence test */
282: snes->ttol = fnorm*snes->rtol;
284: /* test convergence */
285: (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
286: if (snes->reason) return(0);
287: } else {
288: PetscObjectAMSGrantAccess((PetscObject)snes);
289: SNESLogConvergenceHistory(snes,snes->norm,0);
290: SNESMonitor(snes,0,snes->norm);
291: }
294: /* Call general purpose update function */
295: if (snes->ops->update) {
296: (*snes->ops->update)(snes, snes->iter);
297: }
299: for (i = 0; i < snes->max_its; i++) {
300: SNESComputeGS(snes, B, X);
301: /* only compute norms if requested or about to exit due to maximum iterations */
302: if (normtype == SNES_NORM_FUNCTION || ((i == snes->max_its - 1) && (normtype == SNES_NORM_INITIAL_FINAL_ONLY || normtype == SNES_NORM_FINAL_ONLY))) {
303: SNESComputeFunction(snes,X,F);
304: if (snes->domainerror) {
305: snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
306: return(0);
307: }
308: VecNorm(F, NORM_2, &fnorm); /* fnorm <- ||F|| */
309: if (PetscIsInfOrNanReal(fnorm)) {
310: snes->reason = SNES_DIVERGED_FNORM_NAN;
311: return(0);
312: }
313: }
314: /* Monitor convergence */
315: PetscObjectAMSTakeAccess((PetscObject)snes);
316: snes->iter = i+1;
317: snes->norm = fnorm;
318: PetscObjectAMSGrantAccess((PetscObject)snes);
319: SNESLogConvergenceHistory(snes,snes->norm,0);
320: SNESMonitor(snes,snes->iter,snes->norm);
321: /* Test for convergence */
322: if (normtype == SNES_NORM_FUNCTION) (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
323: if (snes->reason) return(0);
324: /* Call general purpose update function */
325: if (snes->ops->update) {
326: (*snes->ops->update)(snes, snes->iter);
327: }
328: }
329: if (normtype == SNES_NORM_FUNCTION) {
330: if (i == snes->max_its) {
331: PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",snes->max_its);
332: if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT;
333: }
334: } else if (!snes->reason) snes->reason = SNES_CONVERGED_ITS; /* GS is meant to be used as a preconditioner */
335: return(0);
336: }
338: /*MC
339: SNESGS - Just calls the user-provided solution routine provided with SNESSetGS()
341: Level: advanced
343: Notes:
344: the Gauss-Seidel smoother is inherited through composition. If a solver has been created with SNESGetPC(), it will have
345: its parent's Gauss-Seidel routine associated with it.
347: .seealso: SNESCreate(), SNES, SNESSetType(), SNESSetGS(), SNESType (for list of available types)
348: M*/
352: PETSC_EXTERN PetscErrorCode SNESCreate_GS(SNES snes)
353: {
354: SNES_GS *gs;
358: snes->ops->destroy = SNESDestroy_GS;
359: snes->ops->setup = SNESSetUp_GS;
360: snes->ops->setfromoptions = SNESSetFromOptions_GS;
361: snes->ops->view = SNESView_GS;
362: snes->ops->solve = SNESSolve_GS;
363: snes->ops->reset = SNESReset_GS;
365: snes->usesksp = PETSC_FALSE;
366: snes->usespc = PETSC_FALSE;
368: if (!snes->tolerancesset) {
369: snes->max_its = 10000;
370: snes->max_funcs = 10000;
371: }
373: PetscNewLog(snes, SNES_GS, &gs);
375: gs->sweeps = 1;
376: gs->rtol = 1e-5;
377: gs->abstol = 1e-15;
378: gs->stol = 1e-12;
379: gs->max_its = 50;
381: snes->data = (void*) gs;
382: return(0);
383: }