Actual source code: tsdiscgrad.c
1: /*
2: Code for timestepping with discrete gradient integrators
3: */
4: #include <petsc/private/tsimpl.h>
5: #include <petscdm.h>
7: PetscBool DGCite = PETSC_FALSE;
8: const char DGCitation[] = "@article{Gonzalez1996,\n"
9: " title = {Time integration and discrete Hamiltonian systems},\n"
10: " author = {Oscar Gonzalez},\n"
11: " journal = {Journal of Nonlinear Science},\n"
12: " volume = {6},\n"
13: " pages = {449--467},\n"
14: " doi = {10.1007/978-1-4612-1246-1_10},\n"
15: " year = {1996}\n}\n";
17: typedef struct {
18: PetscReal stage_time;
19: Vec X0, X, Xdot;
20: void *funcCtx;
21: PetscBool gonzalez;
22: PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *);
23: PetscErrorCode (*Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *);
24: PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *);
25: } TS_DiscGrad;
27: static PetscErrorCode TSDiscGradGetX0AndXdot(TS ts, DM dm, Vec *X0, Vec *Xdot)
28: {
29: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
31: if (X0) {
32: if (dm && dm != ts->dm) DMGetNamedGlobalVector(dm, "TSDiscGrad_X0", X0);
33: else {*X0 = ts->vec_sol;}
34: }
35: if (Xdot) {
36: if (dm && dm != ts->dm) DMGetNamedGlobalVector(dm, "TSDiscGrad_Xdot", Xdot);
37: else {*Xdot = dg->Xdot;}
38: }
39: return 0;
40: }
42: static PetscErrorCode TSDiscGradRestoreX0AndXdot(TS ts, DM dm, Vec *X0, Vec *Xdot)
43: {
44: if (X0) {
45: if (dm && dm != ts->dm) DMRestoreNamedGlobalVector(dm, "TSDiscGrad_X0", X0);
46: }
47: if (Xdot) {
48: if (dm && dm != ts->dm) DMRestoreNamedGlobalVector(dm, "TSDiscGrad_Xdot", Xdot);
49: }
50: return 0;
51: }
53: static PetscErrorCode DMCoarsenHook_TSDiscGrad(DM fine, DM coarse, void *ctx)
54: {
55: return 0;
56: }
58: static PetscErrorCode DMRestrictHook_TSDiscGrad(DM fine, Mat restrct, Vec rscale, Mat inject, DM coarse, void *ctx)
59: {
60: TS ts = (TS) ctx;
61: Vec X0, Xdot, X0_c, Xdot_c;
63: TSDiscGradGetX0AndXdot(ts, fine, &X0, &Xdot);
64: TSDiscGradGetX0AndXdot(ts, coarse, &X0_c, &Xdot_c);
65: MatRestrict(restrct, X0, X0_c);
66: MatRestrict(restrct, Xdot, Xdot_c);
67: VecPointwiseMult(X0_c, rscale, X0_c);
68: VecPointwiseMult(Xdot_c, rscale, Xdot_c);
69: TSDiscGradRestoreX0AndXdot(ts, fine, &X0, &Xdot);
70: TSDiscGradRestoreX0AndXdot(ts, coarse, &X0_c, &Xdot_c);
71: return 0;
72: }
74: static PetscErrorCode DMSubDomainHook_TSDiscGrad(DM dm, DM subdm, void *ctx)
75: {
76: return 0;
77: }
79: static PetscErrorCode DMSubDomainRestrictHook_TSDiscGrad(DM dm, VecScatter gscat, VecScatter lscat, DM subdm, void *ctx)
80: {
81: TS ts = (TS) ctx;
82: Vec X0, Xdot, X0_sub, Xdot_sub;
84: TSDiscGradGetX0AndXdot(ts, dm, &X0, &Xdot);
85: TSDiscGradGetX0AndXdot(ts, subdm, &X0_sub, &Xdot_sub);
87: VecScatterBegin(gscat, X0, X0_sub, INSERT_VALUES, SCATTER_FORWARD);
88: VecScatterEnd(gscat, X0, X0_sub, INSERT_VALUES, SCATTER_FORWARD);
90: VecScatterBegin(gscat, Xdot, Xdot_sub, INSERT_VALUES, SCATTER_FORWARD);
91: VecScatterEnd(gscat, Xdot, Xdot_sub, INSERT_VALUES, SCATTER_FORWARD);
93: TSDiscGradRestoreX0AndXdot(ts, dm, &X0, &Xdot);
94: TSDiscGradRestoreX0AndXdot(ts, subdm, &X0_sub, &Xdot_sub);
95: return 0;
96: }
98: static PetscErrorCode TSSetUp_DiscGrad(TS ts)
99: {
100: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
101: DM dm;
103: if (!dg->X) VecDuplicate(ts->vec_sol, &dg->X);
104: if (!dg->X0) VecDuplicate(ts->vec_sol, &dg->X0);
105: if (!dg->Xdot) VecDuplicate(ts->vec_sol, &dg->Xdot);
107: TSGetDM(ts, &dm);
108: DMCoarsenHookAdd(dm, DMCoarsenHook_TSDiscGrad, DMRestrictHook_TSDiscGrad, ts);
109: DMSubDomainHookAdd(dm, DMSubDomainHook_TSDiscGrad, DMSubDomainRestrictHook_TSDiscGrad, ts);
110: return 0;
111: }
113: static PetscErrorCode TSSetFromOptions_DiscGrad(PetscOptionItems *PetscOptionsObject, TS ts)
114: {
115: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
117: PetscOptionsHead(PetscOptionsObject, "Discrete Gradients ODE solver options");
118: {
119: PetscOptionsBool("-ts_discgrad_gonzalez","Use Gonzalez term in discrete gradients formulation","TSDiscGradUseGonzalez",dg->gonzalez,&dg->gonzalez,NULL);
120: }
121: PetscOptionsTail();
122: return 0;
123: }
125: static PetscErrorCode TSView_DiscGrad(TS ts,PetscViewer viewer)
126: {
127: PetscBool iascii;
129: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
130: if (iascii) {
131: PetscViewerASCIIPrintf(viewer," Discrete Gradients\n");
132: }
133: return 0;
134: }
136: static PetscErrorCode TSDiscGradIsGonzalez_DiscGrad(TS ts,PetscBool *gonzalez)
137: {
138: TS_DiscGrad *dg = (TS_DiscGrad*)ts->data;
140: *gonzalez = dg->gonzalez;
141: return 0;
142: }
144: static PetscErrorCode TSDiscGradUseGonzalez_DiscGrad(TS ts,PetscBool flg)
145: {
146: TS_DiscGrad *dg = (TS_DiscGrad*)ts->data;
148: dg->gonzalez = flg;
149: return 0;
150: }
152: static PetscErrorCode TSReset_DiscGrad(TS ts)
153: {
154: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
156: VecDestroy(&dg->X);
157: VecDestroy(&dg->X0);
158: VecDestroy(&dg->Xdot);
159: return 0;
160: }
162: static PetscErrorCode TSDestroy_DiscGrad(TS ts)
163: {
164: DM dm;
166: TSReset_DiscGrad(ts);
167: TSGetDM(ts, &dm);
168: if (dm) {
169: DMCoarsenHookRemove(dm, DMCoarsenHook_TSDiscGrad, DMRestrictHook_TSDiscGrad, ts);
170: DMSubDomainHookRemove(dm, DMSubDomainHook_TSDiscGrad, DMSubDomainRestrictHook_TSDiscGrad, ts);
171: }
172: PetscFree(ts->data);
173: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradGetFormulation_C",NULL);
174: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradSetFormulation_C",NULL);
175: return 0;
176: }
178: static PetscErrorCode TSInterpolate_DiscGrad(TS ts, PetscReal t, Vec X)
179: {
180: TS_DiscGrad *dg = (TS_DiscGrad*)ts->data;
181: PetscReal dt = t - ts->ptime;
183: VecCopy(ts->vec_sol, dg->X);
184: VecWAXPY(X, dt, dg->Xdot, dg->X);
185: return 0;
186: }
188: static PetscErrorCode TSDiscGrad_SNESSolve(TS ts, Vec b, Vec x)
189: {
190: SNES snes;
191: PetscInt nits, lits;
193: TSGetSNES(ts, &snes);
194: SNESSolve(snes, b, x);
195: SNESGetIterationNumber(snes, &nits);
196: SNESGetLinearSolveIterations(snes, &lits);
197: ts->snes_its += nits;
198: ts->ksp_its += lits;
199: return 0;
200: }
202: static PetscErrorCode TSStep_DiscGrad(TS ts)
203: {
204: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
205: TSAdapt adapt;
206: TSStepStatus status = TS_STEP_INCOMPLETE;
207: PetscInt rejections = 0;
208: PetscBool stageok, accept = PETSC_TRUE;
209: PetscReal next_time_step = ts->time_step;
211: TSGetAdapt(ts, &adapt);
212: if (!ts->steprollback) VecCopy(ts->vec_sol, dg->X0);
214: while (!ts->reason && status != TS_STEP_COMPLETE) {
215: PetscReal shift = 1/(0.5*ts->time_step);
217: dg->stage_time = ts->ptime + 0.5*ts->time_step;
219: VecCopy(dg->X0, dg->X);
220: TSPreStage(ts, dg->stage_time);
221: TSDiscGrad_SNESSolve(ts, NULL, dg->X);
222: TSPostStage(ts, dg->stage_time, 0, &dg->X);
223: TSAdaptCheckStage(adapt, ts, dg->stage_time, dg->X, &stageok);
224: if (!stageok) goto reject_step;
226: status = TS_STEP_PENDING;
227: VecAXPBYPCZ(dg->Xdot, -shift, shift, 0, dg->X0, dg->X);
228: VecAXPY(ts->vec_sol, ts->time_step, dg->Xdot);
229: TSAdaptChoose(adapt, ts, ts->time_step, NULL, &next_time_step, &accept);
230: status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE;
231: if (!accept) {
232: VecCopy(dg->X0, ts->vec_sol);
233: ts->time_step = next_time_step;
234: goto reject_step;
235: }
236: ts->ptime += ts->time_step;
237: ts->time_step = next_time_step;
238: break;
240: reject_step:
241: ts->reject++; accept = PETSC_FALSE;
242: if (!ts->reason && ts->max_reject >= 0 && ++rejections > ts->max_reject) {
243: ts->reason = TS_DIVERGED_STEP_REJECTED;
244: PetscInfo(ts, "Step=%D, step rejections %D greater than current TS allowed, stopping solve\n", ts->steps, rejections);
245: }
246: }
247: return 0;
248: }
250: static PetscErrorCode TSGetStages_DiscGrad(TS ts, PetscInt *ns, Vec **Y)
251: {
252: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
254: if (ns) *ns = 1;
255: if (Y) *Y = &(dg->X);
256: return 0;
257: }
259: /*
260: This defines the nonlinear equation that is to be solved with SNES
261: G(U) = F[t0 + 0.5*dt, U, (U-U0)/dt] = 0
262: */
264: /* x = (x+x')/2 */
265: /* NEED TO CALCULATE x_{n+1} from x and x_{n}*/
266: static PetscErrorCode SNESTSFormFunction_DiscGrad(SNES snes, Vec x, Vec y, TS ts)
267: {
269: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
270: PetscReal norm, shift = 1/(0.5*ts->time_step);
271: PetscInt n;
272: Vec X0, Xdot, Xp, Xdiff;
273: Mat S;
274: PetscScalar F=0, F0=0, Gp;
275: Vec G, SgF;
276: DM dm, dmsave;
278: SNESGetDM(snes, &dm);
280: VecDuplicate(y, &Xp);
281: VecDuplicate(y, &Xdiff);
282: VecDuplicate(y, &SgF);
283: VecDuplicate(y, &G);
285: VecGetLocalSize(y, &n);
286: MatCreate(PETSC_COMM_WORLD, &S);
287: MatSetSizes(S, PETSC_DECIDE, PETSC_DECIDE, n, n);
288: MatSetFromOptions(S);
289: MatSetUp(S);
290: MatAssemblyBegin(S,MAT_FINAL_ASSEMBLY);
291: MatAssemblyEnd(S,MAT_FINAL_ASSEMBLY);
293: TSDiscGradGetX0AndXdot(ts, dm, &X0, &Xdot);
294: VecAXPBYPCZ(Xdot, -shift, shift, 0, X0, x); /* Xdot = shift (x - X0) */
296: VecAXPBYPCZ(Xp, -1, 2, 0, X0, x); /* Xp = 2*x - X0 + (0)*Xmid */
297: VecAXPBYPCZ(Xdiff, -1, 1, 0, X0, Xp); /* Xdiff = xp - X0 + (0)*Xdiff */
299: if (dg->gonzalez) {
300: (*dg->Sfunc)(ts, dg->stage_time, x , S, dg->funcCtx);
301: (*dg->Ffunc)(ts, dg->stage_time, Xp, &F, dg->funcCtx);
302: (*dg->Ffunc)(ts, dg->stage_time, X0, &F0, dg->funcCtx);
303: (*dg->Gfunc)(ts, dg->stage_time, x , G, dg->funcCtx);
305: /* Adding Extra Gonzalez Term */
306: VecDot(Xdiff, G, &Gp);
307: VecNorm(Xdiff, NORM_2, &norm);
308: if (norm < PETSC_SQRT_MACHINE_EPSILON) {
309: Gp = 0;
310: } else {
311: /* Gp = (1/|xn+1 - xn|^2) * (F(xn+1) - F(xn) - Gp) */
312: Gp = (F - F0 - Gp)/PetscSqr(norm);
313: }
314: VecAXPY(G, Gp, Xdiff);
315: MatMult(S, G , SgF); /* S*gradF */
317: } else {
318: (*dg->Sfunc)(ts, dg->stage_time, x, S, dg->funcCtx);
319: (*dg->Gfunc)(ts, dg->stage_time, x, G, dg->funcCtx);
321: MatMult(S, G , SgF); /* Xdot = S*gradF */
322: }
323: /* DM monkey-business allows user code to call TSGetDM() inside of functions evaluated on levels of FAS */
324: dmsave = ts->dm;
325: ts->dm = dm;
326: VecAXPBYPCZ(y, 1, -1, 0, Xdot, SgF);
327: ts->dm = dmsave;
328: TSDiscGradRestoreX0AndXdot(ts, dm, &X0, &Xdot);
330: VecDestroy(&Xp);
331: VecDestroy(&Xdiff);
332: VecDestroy(&SgF);
333: VecDestroy(&G);
334: MatDestroy(&S);
336: return 0;
337: }
339: static PetscErrorCode SNESTSFormJacobian_DiscGrad(SNES snes, Vec x, Mat A, Mat B, TS ts)
340: {
341: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
342: PetscReal shift = 1/(0.5*ts->time_step);
343: Vec Xdot;
344: DM dm,dmsave;
346: SNESGetDM(snes, &dm);
347: /* Xdot has already been computed in SNESTSFormFunction_DiscGrad (SNES guarantees this) */
348: TSDiscGradGetX0AndXdot(ts, dm, NULL, &Xdot);
350: dmsave = ts->dm;
351: ts->dm = dm;
352: TSComputeIJacobian(ts, dg->stage_time, x, Xdot, shift, A, B, PETSC_FALSE);
353: ts->dm = dmsave;
354: TSDiscGradRestoreX0AndXdot(ts, dm, NULL, &Xdot);
355: return 0;
356: }
358: static PetscErrorCode TSDiscGradGetFormulation_DiscGrad(TS ts, PetscErrorCode (**Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (**Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (**Gfunc)(TS, PetscReal, Vec, Vec, void *), void *ctx)
359: {
360: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
362: *Sfunc = dg->Sfunc;
363: *Ffunc = dg->Ffunc;
364: *Gfunc = dg->Gfunc;
365: return 0;
366: }
368: static PetscErrorCode TSDiscGradSetFormulation_DiscGrad(TS ts, PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (*Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *), void *ctx)
369: {
370: TS_DiscGrad *dg = (TS_DiscGrad *) ts->data;
372: dg->Sfunc = Sfunc;
373: dg->Ffunc = Ffunc;
374: dg->Gfunc = Gfunc;
375: dg->funcCtx = ctx;
376: return 0;
377: }
379: /*MC
380: TSDISCGRAD - ODE solver using the discrete gradients version of the implicit midpoint method
382: Notes:
383: This is the implicit midpoint rule, with an optional term that guarantees the discrete gradient property. This
384: timestepper applies to systems of the form
385: $ u_t = S(u) grad F(u)
386: where S(u) is a linear operator, and F is a functional of u.
388: Level: intermediate
390: .seealso: TSCreate(), TSSetType(), TS, TSDISCGRAD, TSDiscGradSetFormulation()
391: M*/
392: PETSC_EXTERN PetscErrorCode TSCreate_DiscGrad(TS ts)
393: {
394: TS_DiscGrad *th;
396: PetscCitationsRegister(DGCitation, &DGCite);
397: ts->ops->reset = TSReset_DiscGrad;
398: ts->ops->destroy = TSDestroy_DiscGrad;
399: ts->ops->view = TSView_DiscGrad;
400: ts->ops->setfromoptions = TSSetFromOptions_DiscGrad;
401: ts->ops->setup = TSSetUp_DiscGrad;
402: ts->ops->step = TSStep_DiscGrad;
403: ts->ops->interpolate = TSInterpolate_DiscGrad;
404: ts->ops->getstages = TSGetStages_DiscGrad;
405: ts->ops->snesfunction = SNESTSFormFunction_DiscGrad;
406: ts->ops->snesjacobian = SNESTSFormJacobian_DiscGrad;
407: ts->default_adapt_type = TSADAPTNONE;
409: ts->usessnes = PETSC_TRUE;
411: PetscNewLog(ts,&th);
412: ts->data = (void*)th;
414: th->gonzalez = PETSC_FALSE;
416: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradGetFormulation_C",TSDiscGradGetFormulation_DiscGrad);
417: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradSetFormulation_C",TSDiscGradSetFormulation_DiscGrad);
418: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradIsGonzalez_C",TSDiscGradIsGonzalez_DiscGrad);
419: PetscObjectComposeFunction((PetscObject)ts,"TSDiscGradUseGonzalez_C",TSDiscGradUseGonzalez_DiscGrad);
420: return 0;
421: }
423: /*@C
424: TSDiscGradGetFormulation - Get the construction method for S, F, and grad F from the formulation u_t = S grad F
426: Not Collective
428: Input Parameter:
429: . ts - timestepping context
431: Output Parameters:
432: + Sfunc - constructor for the S matrix from the formulation
433: . Ffunc - functional F from the formulation
434: . Gfunc - constructor for the gradient of F from the formulation
435: - ctx - the user context
437: Calling sequence of Sfunc:
438: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Mat S, void *)
440: Calling sequence of Ffunc:
441: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, PetscScalar *F, void *)
443: Calling sequence of Gfunc:
444: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Vec G, void *)
446: Level: intermediate
448: .seealso: TSDiscGradSetFormulation()
449: @*/
450: PetscErrorCode TSDiscGradGetFormulation(TS ts, PetscErrorCode (**Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (**Ffunc)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (**Gfunc)(TS, PetscReal, Vec, Vec, void *), void *ctx)
451: {
456: PetscUseMethod(ts,"TSDiscGradGetFormulation_C",(TS,PetscErrorCode(**Sfunc)(TS,PetscReal,Vec,Mat,void*),PetscErrorCode(**Ffunc)(TS,PetscReal,Vec,PetscScalar*,void*),PetscErrorCode(**Gfunc)(TS,PetscReal,Vec,Vec,void*), void*),(ts,Sfunc,Ffunc,Gfunc,ctx));
457: return 0;
458: }
460: /*@C
461: TSDiscGradSetFormulation - Set the construction method for S, F, and grad F from the formulation u_t = S(u) grad F(u)
463: Not Collective
465: Input Parameters:
466: + ts - timestepping context
467: . Sfunc - constructor for the S matrix from the formulation
468: . Ffunc - functional F from the formulation
469: - Gfunc - constructor for the gradient of F from the formulation
470: Calling sequence of Sfunc:
471: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Mat S, void *)
473: Calling sequence of Ffunc:
474: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, PetscScalar *F, void *)
476: Calling sequence of Gfunc:
477: $ PetscErrorCode func(TS ts, PetscReal time, Vec u, Vec G, void *)
479: Level: Intermediate
481: .seealso: TSDiscGradGetFormulation()
482: @*/
483: PetscErrorCode TSDiscGradSetFormulation(TS ts, PetscErrorCode (*Sfunc)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (*Ffunc)(TS, PetscReal, Vec , PetscScalar *, void *), PetscErrorCode (*Gfunc)(TS, PetscReal, Vec, Vec, void *), void *ctx)
484: {
489: PetscTryMethod(ts,"TSDiscGradSetFormulation_C",(TS,PetscErrorCode(*Sfunc)(TS,PetscReal,Vec,Mat,void*),PetscErrorCode(*Ffunc)(TS,PetscReal,Vec,PetscScalar*,void*),PetscErrorCode(*Gfunc)(TS,PetscReal,Vec,Vec,void*), void*),(ts,Sfunc,Ffunc,Gfunc,ctx));
490: return 0;
491: }
493: /*@
494: TSDiscGradIsGonzalez - Checks flag for whether to use additional conservative terms in discrete gradient formulation.
496: Not Collective
498: Input Parameter:
499: . ts - timestepping context
501: Output Parameter:
502: . gonzalez - PETSC_TRUE when using the Gonzalez term
504: Level: Advanced
506: .seealso: TSDiscGradUseGonzalez(), TSDISCGRAD
507: @*/
508: PetscErrorCode TSDiscGradIsGonzalez(TS ts,PetscBool *gonzalez)
509: {
512: PetscUseMethod(ts,"TSDiscGradIsGonzalez_C",(TS,PetscBool*),(ts,gonzalez));
513: return 0;
514: }
516: /*@
517: TSDiscGradUseGonzalez - Sets discrete gradient formulation with or without additional conservative terms. Without flag, the discrete gradients timestepper is just backwards euler
519: Not Collective
521: Input Parameters:
522: + ts - timestepping context
523: - flg - PETSC_TRUE to use the Gonzalez term
525: Options Database:
526: . -ts_discgrad_gonzalez <flg> - use the Gonzalez term for the discrete gradient formulation
528: Level: Intermediate
530: .seealso: TSDISCGRAD
531: @*/
532: PetscErrorCode TSDiscGradUseGonzalez(TS ts,PetscBool flg)
533: {
535: PetscTryMethod(ts,"TSDiscGradUseGonzalez_C",(TS,PetscBool),(ts,flg));
536: return 0;
537: }