Actual source code: ts.c
1: #define PETSCTS_DLL
3: #include include/private/tsimpl.h
5: /* Logging support */
6: PetscCookie TS_COOKIE = 0;
7: PetscEvent TS_Step = 0, TS_PseudoComputeTimeStep = 0, TS_FunctionEval = 0, TS_JacobianEval = 0;
11: /*
12: TSSetTypeFromOptions - Sets the type of ts from user options.
14: Collective on TS
16: Input Parameter:
17: . ts - The ts
19: Level: intermediate
21: .keywords: TS, set, options, database, type
22: .seealso: TSSetFromOptions(), TSSetType()
23: */
24: static PetscErrorCode TSSetTypeFromOptions(TS ts)
25: {
26: PetscTruth opt;
27: const char *defaultType;
28: char typeName[256];
32: if (ts->type_name) {
33: defaultType = ts->type_name;
34: } else {
35: defaultType = TS_EULER;
36: }
38: if (!TSRegisterAllCalled) {
39: TSRegisterAll(PETSC_NULL);
40: }
41: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
42: if (opt) {
43: TSSetType(ts, typeName);
44: } else {
45: TSSetType(ts, defaultType);
46: }
47: return(0);
48: }
52: /*@
53: TSSetFromOptions - Sets various TS parameters from user options.
55: Collective on TS
57: Input Parameter:
58: . ts - the TS context obtained from TSCreate()
60: Options Database Keys:
61: + -ts_type <type> - TS_EULER, TS_BEULER, TS_SUNDIALS, TS_PSEUDO, TS_CRANK_NICHOLSON
62: . -ts_max_steps maxsteps - maximum number of time-steps to take
63: . -ts_max_time time - maximum time to compute to
64: . -ts_dt dt - initial time step
65: . -ts_monitor - print information at each timestep
66: - -ts_monitor_draw - plot information at each timestep
68: Level: beginner
70: .keywords: TS, timestep, set, options, database
72: .seealso: TSGetType
73: @*/
74: PetscErrorCode TSSetFromOptions(TS ts)
75: {
76: PetscReal dt;
77: PetscTruth opt,flg;
78: PetscErrorCode ierr;
79: PetscViewerASCIIMonitor monviewer;
80: char monfilename[PETSC_MAX_PATH_LEN];
84: PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");
86: /* Handle generic TS options */
87: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
88: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
89: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
90: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
91: if (opt) {
92: ts->initial_time_step = ts->time_step = dt;
93: }
95: /* Monitor options */
96: PetscOptionsString("-ts_monitor","Monitor timestep size","TSMonitorDefault","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
97: if (flg) {
98: PetscViewerASCIIMonitorCreate(ts->comm,monfilename,0,&monviewer);
99: TSMonitorSet(ts,TSMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
100: }
101: PetscOptionsName("-ts_monitor_draw","Monitor timestep size graphically","TSMonitorLG",&opt);
102: if (opt) {
103: TSMonitorSet(ts,TSMonitorLG,PETSC_NULL,PETSC_NULL);
104: }
105: PetscOptionsName("-ts_monitor_solution","Monitor solution graphically","TSMonitorSolution",&opt);
106: if (opt) {
107: TSMonitorSet(ts,TSMonitorSolution,PETSC_NULL,PETSC_NULL);
108: }
110: /* Handle TS type options */
111: TSSetTypeFromOptions(ts);
113: /* Handle specific TS options */
114: if (ts->ops->setfromoptions) {
115: (*ts->ops->setfromoptions)(ts);
116: }
117: PetscOptionsEnd();
119: /* Handle subobject options */
120: switch(ts->problem_type) {
121: /* Should check for implicit/explicit */
122: case TS_LINEAR:
123: if (ts->ksp) {
124: KSPSetOperators(ts->ksp,ts->Arhs,ts->B,DIFFERENT_NONZERO_PATTERN);
125: KSPSetFromOptions(ts->ksp);
126: }
127: break;
128: case TS_NONLINEAR:
129: if (ts->snes) {
130: /* this is a bit of a hack, but it gets the matrix information into SNES earlier
131: so that SNES and KSP have more information to pick reasonable defaults
132: before they allow users to set options */
133: SNESSetJacobian(ts->snes,ts->Arhs,ts->B,0,ts);
134: SNESSetFromOptions(ts->snes);
135: }
136: break;
137: default:
138: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
139: }
141: return(0);
142: }
144: #undef __FUNCT__
146: /*@
147: TSViewFromOptions - This function visualizes the ts based upon user options.
149: Collective on TS
151: Input Parameter:
152: . ts - The ts
154: Level: intermediate
156: .keywords: TS, view, options, database
157: .seealso: TSSetFromOptions(), TSView()
158: @*/
159: PetscErrorCode TSViewFromOptions(TS ts,const char title[])
160: {
161: PetscViewer viewer;
162: PetscDraw draw;
163: PetscTruth opt;
164: char fileName[PETSC_MAX_PATH_LEN];
168: PetscOptionsGetString(ts->prefix, "-ts_view", fileName, PETSC_MAX_PATH_LEN, &opt);
169: if (opt && !PetscPreLoadingOn) {
170: PetscViewerASCIIOpen(ts->comm,fileName,&viewer);
171: TSView(ts, viewer);
172: PetscViewerDestroy(viewer);
173: }
174: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
175: if (opt) {
176: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
177: PetscViewerDrawGetDraw(viewer, 0, &draw);
178: if (title) {
179: PetscDrawSetTitle(draw, (char *)title);
180: } else {
181: PetscObjectName((PetscObject) ts);
182: PetscDrawSetTitle(draw, ts->name);
183: }
184: TSView(ts, viewer);
185: PetscViewerFlush(viewer);
186: PetscDrawPause(draw);
187: PetscViewerDestroy(viewer);
188: }
189: return(0);
190: }
194: /*@
195: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
196: set with TSSetRHSJacobian().
198: Collective on TS and Vec
200: Input Parameters:
201: + ts - the SNES context
202: . t - current timestep
203: - x - input vector
205: Output Parameters:
206: + A - Jacobian matrix
207: . B - optional preconditioning matrix
208: - flag - flag indicating matrix structure
210: Notes:
211: Most users should not need to explicitly call this routine, as it
212: is used internally within the nonlinear solvers.
214: See KSPSetOperators() for important information about setting the
215: flag parameter.
217: TSComputeJacobian() is valid only for TS_NONLINEAR
219: Level: developer
221: .keywords: SNES, compute, Jacobian, matrix
223: .seealso: TSSetRHSJacobian(), KSPSetOperators()
224: @*/
225: PetscErrorCode TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
226: {
233: if (ts->problem_type != TS_NONLINEAR) {
234: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
235: }
236: if (ts->ops->rhsjacobian) {
238: *flg = DIFFERENT_NONZERO_PATTERN;
239: PetscStackPush("TS user Jacobian function");
240: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
241: PetscStackPop;
243: /* make sure user returned a correct Jacobian and preconditioner */
246: } else {
247: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
248: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
249: if (*A != *B) {
250: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
251: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
252: }
253: }
254: return(0);
255: }
259: /*
260: TSComputeRHSFunction - Evaluates the right-hand-side function.
262: Note: If the user did not provide a function but merely a matrix,
263: this routine applies the matrix.
264: */
265: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
266: {
275: if (ts->ops->rhsfunction) {
276: PetscStackPush("TS user right-hand-side function");
277: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
278: PetscStackPop;
279: } else {
280: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
281: MatStructure flg;
282: PetscStackPush("TS user right-hand-side matrix function");
283: (*ts->ops->rhsmatrix)(ts,t,&ts->Arhs,&ts->B,&flg,ts->jacP);
284: PetscStackPop;
285: }
286: MatMult(ts->Arhs,x,y);
287: }
291: return(0);
292: }
296: /*@C
297: TSSetRHSFunction - Sets the routine for evaluating the function,
298: F(t,u), where U_t = F(t,u).
300: Collective on TS
302: Input Parameters:
303: + ts - the TS context obtained from TSCreate()
304: . f - routine for evaluating the right-hand-side function
305: - ctx - [optional] user-defined context for private data for the
306: function evaluation routine (may be PETSC_NULL)
308: Calling sequence of func:
309: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
311: + t - current timestep
312: . u - input vector
313: . F - function vector
314: - ctx - [optional] user-defined function context
316: Important:
317: The user MUST call either this routine or TSSetMatrices().
319: Level: beginner
321: .keywords: TS, timestep, set, right-hand-side, function
323: .seealso: TSSetMatrices()
324: @*/
325: PetscErrorCode TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
326: {
330: if (ts->problem_type == TS_LINEAR) {
331: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
332: }
333: ts->ops->rhsfunction = f;
334: ts->funP = ctx;
335: return(0);
336: }
340: /*@C
341: TSSetMatrices - Sets the functions to compute the matrices Alhs and Arhs,
342: where Alhs(t) U_t = Arhs(t) U.
344: Collective on TS
346: Input Parameters:
347: + ts - the TS context obtained from TSCreate()
348: . Arhs - matrix
349: . frhs - the matrix evaluation routine for Arhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
350: if Arhs is not a function of t.
351: . Alhs - matrix or PETSC_NULL if Alhs is an indentity matrix.
352: . flhs - the matrix evaluation routine for Alhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
353: if Alhs is not a function of t.
354: . flag - flag indicating information about the matrix structure of Arhs and Alhs.
355: The available options are
356: SAME_NONZERO_PATTERN - Alhs has the same nonzero structure as Arhs
357: DIFFERENT_NONZERO_PATTERN - Alhs has different nonzero structure as Arhs
358: - ctx - [optional] user-defined context for private data for the
359: matrix evaluation routine (may be PETSC_NULL)
361: Calling sequence of func:
362: $ func(TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
364: + t - current timestep
365: . A - matrix A, where U_t = A(t) U
366: . B - preconditioner matrix, usually the same as A
367: . flag - flag indicating information about the preconditioner matrix
368: structure (same as flag in KSPSetOperators())
369: - ctx - [optional] user-defined context for matrix evaluation routine
371: Notes:
372: The routine func() takes Mat* as the matrix arguments rather than Mat.
373: This allows the matrix evaluation routine to replace Arhs or Alhs with a
374: completely new new matrix structure (not just different matrix elements)
375: when appropriate, for instance, if the nonzero structure is changing
376: throughout the global iterations.
378: Important:
379: The user MUST call either this routine or TSSetRHSFunction().
381: Level: beginner
383: .keywords: TS, timestep, set, matrix
385: .seealso: TSSetRHSFunction()
386: @*/
387: PetscErrorCode TSSetMatrices(TS ts,Mat Arhs,PetscErrorCode (*frhs)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),Mat Alhs,PetscErrorCode (*flhs)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),MatStructure flag,void *ctx)
388: {
391: if (Arhs){
394: ts->Arhs = Arhs;
395: ts->ops->rhsmatrix = frhs;
396: }
397: if (Alhs){
400: ts->Alhs = Alhs;
401: ts->ops->lhsmatrix = flhs;
402: }
403:
404: ts->jacP = ctx;
405: ts->matflg = flag;
406: return(0);
407: }
411: /*@C
412: TSGetMatrices - Returns the matrices Arhs and Alhs at the present timestep,
413: where Alhs(t) U_t = Arhs(t) U.
415: Not Collective, but parallel objects are returned if TS is parallel
417: Input Parameter:
418: . ts - The TS context obtained from TSCreate()
420: Output Parameters:
421: + Arhs - The right-hand side matrix
422: . Alhs - The left-hand side matrix
423: - ctx - User-defined context for matrix evaluation routine
425: Notes: You can pass in PETSC_NULL for any return argument you do not need.
427: Level: intermediate
429: .seealso: TSSetMatrices(), TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
431: .keywords: TS, timestep, get, matrix
433: @*/
434: PetscErrorCode TSGetMatrices(TS ts,Mat *Arhs,Mat *Alhs,void **ctx)
435: {
438: if (Arhs) *Arhs = ts->Arhs;
439: if (Alhs) *Alhs = ts->Alhs;
440: if (ctx) *ctx = ts->jacP;
441: return(0);
442: }
446: /*@C
447: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
448: where U_t = F(U,t), as well as the location to store the matrix.
449: Use TSSetMatrices() for linear problems.
451: Collective on TS
453: Input Parameters:
454: + ts - the TS context obtained from TSCreate()
455: . A - Jacobian matrix
456: . B - preconditioner matrix (usually same as A)
457: . f - the Jacobian evaluation routine
458: - ctx - [optional] user-defined context for private data for the
459: Jacobian evaluation routine (may be PETSC_NULL)
461: Calling sequence of func:
462: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
464: + t - current timestep
465: . u - input vector
466: . A - matrix A, where U_t = A(t)u
467: . B - preconditioner matrix, usually the same as A
468: . flag - flag indicating information about the preconditioner matrix
469: structure (same as flag in KSPSetOperators())
470: - ctx - [optional] user-defined context for matrix evaluation routine
472: Notes:
473: See KSPSetOperators() for important information about setting the flag
474: output parameter in the routine func(). Be sure to read this information!
476: The routine func() takes Mat * as the matrix arguments rather than Mat.
477: This allows the matrix evaluation routine to replace A and/or B with a
478: completely new new matrix structure (not just different matrix elements)
479: when appropriate, for instance, if the nonzero structure is changing
480: throughout the global iterations.
482: Level: beginner
483:
484: .keywords: TS, timestep, set, right-hand-side, Jacobian
486: .seealso: TSDefaultComputeJacobianColor(),
487: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetMatrices()
489: @*/
490: PetscErrorCode TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
491: {
498: if (ts->problem_type != TS_NONLINEAR) {
499: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetMatrices()");
500: }
502: ts->ops->rhsjacobian = f;
503: ts->jacP = ctx;
504: ts->Arhs = A;
505: ts->B = B;
506: return(0);
507: }
511: /*@C
512: TSView - Prints the TS data structure.
514: Collective on TS
516: Input Parameters:
517: + ts - the TS context obtained from TSCreate()
518: - viewer - visualization context
520: Options Database Key:
521: . -ts_view - calls TSView() at end of TSStep()
523: Notes:
524: The available visualization contexts include
525: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
526: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
527: output where only the first processor opens
528: the file. All other processors send their
529: data to the first processor to print.
531: The user can open an alternative visualization context with
532: PetscViewerASCIIOpen() - output to a specified file.
534: Level: beginner
536: .keywords: TS, timestep, view
538: .seealso: PetscViewerASCIIOpen()
539: @*/
540: PetscErrorCode TSView(TS ts,PetscViewer viewer)
541: {
543: char *type;
544: PetscTruth iascii,isstring;
548: if (!viewer) {
549: PetscViewerASCIIGetStdout(ts->comm,&viewer);
550: }
554: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
555: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
556: if (iascii) {
557: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
558: TSGetType(ts,(TSType *)&type);
559: if (type) {
560: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
561: } else {
562: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
563: }
564: if (ts->ops->view) {
565: PetscViewerASCIIPushTab(viewer);
566: (*ts->ops->view)(ts,viewer);
567: PetscViewerASCIIPopTab(viewer);
568: }
569: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
570: PetscViewerASCIIPrintf(viewer," maximum time=%G\n",ts->max_time);
571: if (ts->problem_type == TS_NONLINEAR) {
572: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
573: }
574: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
575: } else if (isstring) {
576: TSGetType(ts,(TSType *)&type);
577: PetscViewerStringSPrintf(viewer," %-7.7s",type);
578: }
579: PetscViewerASCIIPushTab(viewer);
580: if (ts->ksp) {KSPView(ts->ksp,viewer);}
581: if (ts->snes) {SNESView(ts->snes,viewer);}
582: PetscViewerASCIIPopTab(viewer);
583: return(0);
584: }
589: /*@C
590: TSSetApplicationContext - Sets an optional user-defined context for
591: the timesteppers.
593: Collective on TS
595: Input Parameters:
596: + ts - the TS context obtained from TSCreate()
597: - usrP - optional user context
599: Level: intermediate
601: .keywords: TS, timestep, set, application, context
603: .seealso: TSGetApplicationContext()
604: @*/
605: PetscErrorCode TSSetApplicationContext(TS ts,void *usrP)
606: {
609: ts->user = usrP;
610: return(0);
611: }
615: /*@C
616: TSGetApplicationContext - Gets the user-defined context for the
617: timestepper.
619: Not Collective
621: Input Parameter:
622: . ts - the TS context obtained from TSCreate()
624: Output Parameter:
625: . usrP - user context
627: Level: intermediate
629: .keywords: TS, timestep, get, application, context
631: .seealso: TSSetApplicationContext()
632: @*/
633: PetscErrorCode TSGetApplicationContext(TS ts,void **usrP)
634: {
637: *usrP = ts->user;
638: return(0);
639: }
643: /*@
644: TSGetTimeStepNumber - Gets the current number of timesteps.
646: Not Collective
648: Input Parameter:
649: . ts - the TS context obtained from TSCreate()
651: Output Parameter:
652: . iter - number steps so far
654: Level: intermediate
656: .keywords: TS, timestep, get, iteration, number
657: @*/
658: PetscErrorCode TSGetTimeStepNumber(TS ts,PetscInt* iter)
659: {
663: *iter = ts->steps;
664: return(0);
665: }
669: /*@
670: TSSetInitialTimeStep - Sets the initial timestep to be used,
671: as well as the initial time.
673: Collective on TS
675: Input Parameters:
676: + ts - the TS context obtained from TSCreate()
677: . initial_time - the initial time
678: - time_step - the size of the timestep
680: Level: intermediate
682: .seealso: TSSetTimeStep(), TSGetTimeStep()
684: .keywords: TS, set, initial, timestep
685: @*/
686: PetscErrorCode TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
687: {
690: ts->time_step = time_step;
691: ts->initial_time_step = time_step;
692: ts->ptime = initial_time;
693: return(0);
694: }
698: /*@
699: TSSetTimeStep - Allows one to reset the timestep at any time,
700: useful for simple pseudo-timestepping codes.
702: Collective on TS
704: Input Parameters:
705: + ts - the TS context obtained from TSCreate()
706: - time_step - the size of the timestep
708: Level: intermediate
710: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
712: .keywords: TS, set, timestep
713: @*/
714: PetscErrorCode TSSetTimeStep(TS ts,PetscReal time_step)
715: {
718: ts->time_step = time_step;
719: return(0);
720: }
724: /*@
725: TSGetTimeStep - Gets the current timestep size.
727: Not Collective
729: Input Parameter:
730: . ts - the TS context obtained from TSCreate()
732: Output Parameter:
733: . dt - the current timestep size
735: Level: intermediate
737: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
739: .keywords: TS, get, timestep
740: @*/
741: PetscErrorCode TSGetTimeStep(TS ts,PetscReal* dt)
742: {
746: *dt = ts->time_step;
747: return(0);
748: }
752: /*@
753: TSGetSolution - Returns the solution at the present timestep. It
754: is valid to call this routine inside the function that you are evaluating
755: in order to move to the new timestep. This vector not changed until
756: the solution at the next timestep has been calculated.
758: Not Collective, but Vec returned is parallel if TS is parallel
760: Input Parameter:
761: . ts - the TS context obtained from TSCreate()
763: Output Parameter:
764: . v - the vector containing the solution
766: Level: intermediate
768: .seealso: TSGetTimeStep()
770: .keywords: TS, timestep, get, solution
771: @*/
772: PetscErrorCode TSGetSolution(TS ts,Vec *v)
773: {
777: *v = ts->vec_sol_always;
778: return(0);
779: }
781: /* ----- Routines to initialize and destroy a timestepper ---- */
784: /*@
785: TSSetProblemType - Sets the type of problem to be solved.
787: Not collective
789: Input Parameters:
790: + ts - The TS
791: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
792: .vb
793: U_t = A U
794: U_t = A(t) U
795: U_t = F(t,U)
796: .ve
798: Level: beginner
800: .keywords: TS, problem type
801: .seealso: TSSetUp(), TSProblemType, TS
802: @*/
803: PetscErrorCode TSSetProblemType(TS ts, TSProblemType type)
804: {
807: ts->problem_type = type;
808: return(0);
809: }
813: /*@C
814: TSGetProblemType - Gets the type of problem to be solved.
816: Not collective
818: Input Parameter:
819: . ts - The TS
821: Output Parameter:
822: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
823: .vb
824: U_t = A U
825: U_t = A(t) U
826: U_t = F(t,U)
827: .ve
829: Level: beginner
831: .keywords: TS, problem type
832: .seealso: TSSetUp(), TSProblemType, TS
833: @*/
834: PetscErrorCode TSGetProblemType(TS ts, TSProblemType *type)
835: {
839: *type = ts->problem_type;
840: return(0);
841: }
845: /*@
846: TSSetUp - Sets up the internal data structures for the later use
847: of a timestepper.
849: Collective on TS
851: Input Parameter:
852: . ts - the TS context obtained from TSCreate()
854: Notes:
855: For basic use of the TS solvers the user need not explicitly call
856: TSSetUp(), since these actions will automatically occur during
857: the call to TSStep(). However, if one wishes to control this
858: phase separately, TSSetUp() should be called after TSCreate()
859: and optional routines of the form TSSetXXX(), but before TSStep().
861: Level: advanced
863: .keywords: TS, timestep, setup
865: .seealso: TSCreate(), TSStep(), TSDestroy()
866: @*/
867: PetscErrorCode TSSetUp(TS ts)
868: {
873: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
874: if (!ts->type_name) {
875: TSSetType(ts,TS_EULER);
876: }
877: (*ts->ops->setup)(ts);
878: ts->setupcalled = 1;
879: return(0);
880: }
884: /*@
885: TSDestroy - Destroys the timestepper context that was created
886: with TSCreate().
888: Collective on TS
890: Input Parameter:
891: . ts - the TS context obtained from TSCreate()
893: Level: beginner
895: .keywords: TS, timestepper, destroy
897: .seealso: TSCreate(), TSSetUp(), TSSolve()
898: @*/
899: PetscErrorCode TSDestroy(TS ts)
900: {
905: if (--ts->refct > 0) return(0);
907: /* if memory was published with AMS then destroy it */
908: PetscObjectDepublish(ts);
909: if (ts->A) {MatDestroy(ts->A);CHKERRQ(ierr)}
910: if (ts->ksp) {KSPDestroy(ts->ksp);}
911: if (ts->snes) {SNESDestroy(ts->snes);}
912: if (ts->ops->destroy) {(*(ts)->ops->destroy)(ts);}
913: TSMonitorCancel(ts);
914: PetscHeaderDestroy(ts);
915: return(0);
916: }
920: /*@
921: TSGetSNES - Returns the SNES (nonlinear solver) associated with
922: a TS (timestepper) context. Valid only for nonlinear problems.
924: Not Collective, but SNES is parallel if TS is parallel
926: Input Parameter:
927: . ts - the TS context obtained from TSCreate()
929: Output Parameter:
930: . snes - the nonlinear solver context
932: Notes:
933: The user can then directly manipulate the SNES context to set various
934: options, etc. Likewise, the user can then extract and manipulate the
935: KSP, KSP, and PC contexts as well.
937: TSGetSNES() does not work for integrators that do not use SNES; in
938: this case TSGetSNES() returns PETSC_NULL in snes.
940: Level: beginner
942: .keywords: timestep, get, SNES
943: @*/
944: PetscErrorCode TSGetSNES(TS ts,SNES *snes)
945: {
949: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
950: *snes = ts->snes;
951: return(0);
952: }
956: /*@
957: TSGetKSP - Returns the KSP (linear solver) associated with
958: a TS (timestepper) context.
960: Not Collective, but KSP is parallel if TS is parallel
962: Input Parameter:
963: . ts - the TS context obtained from TSCreate()
965: Output Parameter:
966: . ksp - the nonlinear solver context
968: Notes:
969: The user can then directly manipulate the KSP context to set various
970: options, etc. Likewise, the user can then extract and manipulate the
971: KSP and PC contexts as well.
973: TSGetKSP() does not work for integrators that do not use KSP;
974: in this case TSGetKSP() returns PETSC_NULL in ksp.
976: Level: beginner
978: .keywords: timestep, get, KSP
979: @*/
980: PetscErrorCode TSGetKSP(TS ts,KSP *ksp)
981: {
985: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
986: *ksp = ts->ksp;
987: return(0);
988: }
990: /* ----------- Routines to set solver parameters ---------- */
994: /*@
995: TSGetDuration - Gets the maximum number of timesteps to use and
996: maximum time for iteration.
998: Collective on TS
1000: Input Parameters:
1001: + ts - the TS context obtained from TSCreate()
1002: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1003: - maxtime - final time to iterate to, or PETSC_NULL
1005: Level: intermediate
1007: .keywords: TS, timestep, get, maximum, iterations, time
1008: @*/
1009: PetscErrorCode TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1010: {
1013: if (maxsteps) {
1015: *maxsteps = ts->max_steps;
1016: }
1017: if (maxtime ) {
1019: *maxtime = ts->max_time;
1020: }
1021: return(0);
1022: }
1026: /*@
1027: TSSetDuration - Sets the maximum number of timesteps to use and
1028: maximum time for iteration.
1030: Collective on TS
1032: Input Parameters:
1033: + ts - the TS context obtained from TSCreate()
1034: . maxsteps - maximum number of iterations to use
1035: - maxtime - final time to iterate to
1037: Options Database Keys:
1038: . -ts_max_steps <maxsteps> - Sets maxsteps
1039: . -ts_max_time <maxtime> - Sets maxtime
1041: Notes:
1042: The default maximum number of iterations is 5000. Default time is 5.0
1044: Level: intermediate
1046: .keywords: TS, timestep, set, maximum, iterations
1047: @*/
1048: PetscErrorCode TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1049: {
1052: ts->max_steps = maxsteps;
1053: ts->max_time = maxtime;
1054: return(0);
1055: }
1059: /*@
1060: TSSetSolution - Sets the initial solution vector
1061: for use by the TS routines.
1063: Collective on TS and Vec
1065: Input Parameters:
1066: + ts - the TS context obtained from TSCreate()
1067: - x - the solution vector
1069: Level: beginner
1071: .keywords: TS, timestep, set, solution, initial conditions
1072: @*/
1073: PetscErrorCode TSSetSolution(TS ts,Vec x)
1074: {
1078: ts->vec_sol = ts->vec_sol_always = x;
1079: return(0);
1080: }
1084: /*@C
1085: TSSetPreStep - Sets the general-purpose function
1086: called once at the beginning of time stepping.
1088: Collective on TS
1090: Input Parameters:
1091: + ts - The TS context obtained from TSCreate()
1092: - func - The function
1094: Calling sequence of func:
1095: . func (TS ts);
1097: Level: intermediate
1099: .keywords: TS, timestep
1100: @*/
1101: PetscErrorCode TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1102: {
1105: ts->ops->prestep = func;
1106: return(0);
1107: }
1111: /*@
1112: TSDefaultPreStep - The default pre-stepping function which does nothing.
1114: Collective on TS
1116: Input Parameters:
1117: . ts - The TS context obtained from TSCreate()
1119: Level: developer
1121: .keywords: TS, timestep
1122: @*/
1123: PetscErrorCode TSDefaultPreStep(TS ts)
1124: {
1126: return(0);
1127: }
1131: /*@C
1132: TSSetUpdate - Sets the general-purpose update function called
1133: at the beginning of every time step. This function can change
1134: the time step.
1136: Collective on TS
1138: Input Parameters:
1139: + ts - The TS context obtained from TSCreate()
1140: - func - The function
1142: Calling sequence of func:
1143: . func (TS ts, double t, double *dt);
1145: + t - The current time
1146: - dt - The current time step
1148: Level: intermediate
1150: .keywords: TS, update, timestep
1151: @*/
1152: PetscErrorCode TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1153: {
1156: ts->ops->update = func;
1157: return(0);
1158: }
1162: /*@
1163: TSDefaultUpdate - The default update function which does nothing.
1165: Collective on TS
1167: Input Parameters:
1168: + ts - The TS context obtained from TSCreate()
1169: - t - The current time
1171: Output Parameters:
1172: . dt - The current time step
1174: Level: developer
1176: .keywords: TS, update, timestep
1177: @*/
1178: PetscErrorCode TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1179: {
1181: return(0);
1182: }
1186: /*@C
1187: TSSetPostStep - Sets the general-purpose function
1188: called once at the end of time stepping.
1190: Collective on TS
1192: Input Parameters:
1193: + ts - The TS context obtained from TSCreate()
1194: - func - The function
1196: Calling sequence of func:
1197: . func (TS ts);
1199: Level: intermediate
1201: .keywords: TS, timestep
1202: @*/
1203: PetscErrorCode TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1204: {
1207: ts->ops->poststep = func;
1208: return(0);
1209: }
1213: /*@
1214: TSDefaultPostStep - The default post-stepping function which does nothing.
1216: Collective on TS
1218: Input Parameters:
1219: . ts - The TS context obtained from TSCreate()
1221: Level: developer
1223: .keywords: TS, timestep
1224: @*/
1225: PetscErrorCode TSDefaultPostStep(TS ts)
1226: {
1228: return(0);
1229: }
1231: /* ------------ Routines to set performance monitoring options ----------- */
1235: /*@C
1236: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
1237: timestep to display the iteration's progress.
1239: Collective on TS
1241: Input Parameters:
1242: + ts - the TS context obtained from TSCreate()
1243: . func - monitoring routine
1244: . mctx - [optional] user-defined context for private data for the
1245: monitor routine (use PETSC_NULL if no context is desired)
1246: - monitordestroy - [optional] routine that frees monitor context
1247: (may be PETSC_NULL)
1249: Calling sequence of func:
1250: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1252: + ts - the TS context
1253: . steps - iteration number
1254: . time - current time
1255: . x - current iterate
1256: - mctx - [optional] monitoring context
1258: Notes:
1259: This routine adds an additional monitor to the list of monitors that
1260: already has been loaded.
1262: Level: intermediate
1264: .keywords: TS, timestep, set, monitor
1266: .seealso: TSMonitorDefault(), TSMonitorCancel()
1267: @*/
1268: PetscErrorCode TSMonitorSet(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1269: {
1272: if (ts->numbermonitors >= MAXTSMONITORS) {
1273: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1274: }
1275: ts->monitor[ts->numbermonitors] = monitor;
1276: ts->mdestroy[ts->numbermonitors] = mdestroy;
1277: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1278: return(0);
1279: }
1283: /*@C
1284: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
1286: Collective on TS
1288: Input Parameters:
1289: . ts - the TS context obtained from TSCreate()
1291: Notes:
1292: There is no way to remove a single, specific monitor.
1294: Level: intermediate
1296: .keywords: TS, timestep, set, monitor
1298: .seealso: TSMonitorDefault(), TSMonitorSet()
1299: @*/
1300: PetscErrorCode TSMonitorCancel(TS ts)
1301: {
1303: PetscInt i;
1307: for (i=0; i<ts->numbermonitors; i++) {
1308: if (ts->mdestroy[i]) {
1309: (*ts->mdestroy[i])(ts->monitorcontext[i]);
1310: }
1311: }
1312: ts->numbermonitors = 0;
1313: return(0);
1314: }
1318: /*@
1319: TSMonitorDefault - Sets the Default monitor
1321: Level: intermediate
1323: .keywords: TS, set, monitor
1325: .seealso: TSMonitorDefault(), TSMonitorSet()
1326: @*/
1327: PetscErrorCode TSMonitorDefault(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1328: {
1329: PetscErrorCode ierr;
1330: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor)ctx;
1333: if (!ctx) {
1334: PetscViewerASCIIMonitorCreate(ts->comm,"stdout",0,&viewer);
1335: }
1336: PetscViewerASCIIMonitorPrintf(viewer,"timestep %D dt %G time %G\n",step,ts->time_step,ptime);
1337: if (!ctx) {
1338: PetscViewerASCIIMonitorDestroy(viewer);
1339: }
1340: return(0);
1341: }
1345: /*@
1346: TSStep - Steps the requested number of timesteps.
1348: Collective on TS
1350: Input Parameter:
1351: . ts - the TS context obtained from TSCreate()
1353: Output Parameters:
1354: + steps - number of iterations until termination
1355: - ptime - time until termination
1357: Level: beginner
1359: .keywords: TS, timestep, solve
1361: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1362: @*/
1363: PetscErrorCode TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1364: {
1369: if (!ts->setupcalled) {
1370: TSSetUp(ts);
1371: }
1374: (*ts->ops->prestep)(ts);
1375: (*ts->ops->step)(ts, steps, ptime);
1376: (*ts->ops->poststep)(ts);
1379: if (!PetscPreLoadingOn) {
1380: TSViewFromOptions(ts,ts->name);
1381: }
1382: return(0);
1383: }
1387: /*@
1388: TSSolve - Steps the requested number of timesteps.
1390: Collective on TS
1392: Input Parameter:
1393: + ts - the TS context obtained from TSCreate()
1394: - x - the solution vector, or PETSC_NULL if it was set with TSSetSolution()
1396: Level: beginner
1398: .keywords: TS, timestep, solve
1400: .seealso: TSCreate(), TSSetSolution(), TSStep()
1401: @*/
1402: PetscErrorCode TSSolve(TS ts, Vec x)
1403: {
1404: PetscInt steps;
1405: PetscReal ptime;
1409: /* set solution vector if provided */
1410: if (x) { TSSetSolution(ts, x); }
1411: /* reset time step and iteration counters */
1412: ts->steps = 0; ts->linear_its = 0; ts->nonlinear_its = 0;
1413: /* steps the requested number of timesteps. */
1414: TSStep(ts, &steps, &ptime);
1415: return(0);
1416: }
1420: /*
1421: Runs the user provided monitor routines, if they exists.
1422: */
1423: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1424: {
1426: PetscInt i,n = ts->numbermonitors;
1429: for (i=0; i<n; i++) {
1430: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1431: }
1432: return(0);
1433: }
1435: /* ------------------------------------------------------------------------*/
1439: /*@C
1440: TSMonitorLGCreate - Creates a line graph context for use with
1441: TS to monitor convergence of preconditioned residual norms.
1443: Collective on TS
1445: Input Parameters:
1446: + host - the X display to open, or null for the local machine
1447: . label - the title to put in the title bar
1448: . x, y - the screen coordinates of the upper left coordinate of the window
1449: - m, n - the screen width and height in pixels
1451: Output Parameter:
1452: . draw - the drawing context
1454: Options Database Key:
1455: . -ts_monitor_draw - automatically sets line graph monitor
1457: Notes:
1458: Use TSMonitorLGDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1460: Level: intermediate
1462: .keywords: TS, monitor, line graph, residual, seealso
1464: .seealso: TSMonitorLGDestroy(), TSMonitorSet()
1466: @*/
1467: PetscErrorCode TSMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1468: {
1469: PetscDraw win;
1473: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1474: PetscDrawSetType(win,PETSC_DRAW_X);
1475: PetscDrawLGCreate(win,1,draw);
1476: PetscDrawLGIndicateDataPoints(*draw);
1478: PetscLogObjectParent(*draw,win);
1479: return(0);
1480: }
1484: PetscErrorCode TSMonitorLG(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1485: {
1486: PetscDrawLG lg = (PetscDrawLG) monctx;
1487: PetscReal x,y = ptime;
1491: if (!monctx) {
1492: MPI_Comm comm;
1493: PetscViewer viewer;
1495: PetscObjectGetComm((PetscObject)ts,&comm);
1496: viewer = PETSC_VIEWER_DRAW_(comm);
1497: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1498: }
1500: if (!n) {PetscDrawLGReset(lg);}
1501: x = (PetscReal)n;
1502: PetscDrawLGAddPoint(lg,&x,&y);
1503: if (n < 20 || (n % 5)) {
1504: PetscDrawLGDraw(lg);
1505: }
1506: return(0);
1507: }
1511: /*@C
1512: TSMonitorLGDestroy - Destroys a line graph context that was created
1513: with TSMonitorLGCreate().
1515: Collective on PetscDrawLG
1517: Input Parameter:
1518: . draw - the drawing context
1520: Level: intermediate
1522: .keywords: TS, monitor, line graph, destroy
1524: .seealso: TSMonitorLGCreate(), TSMonitorSet(), TSMonitorLG();
1525: @*/
1526: PetscErrorCode TSMonitorLGDestroy(PetscDrawLG drawlg)
1527: {
1528: PetscDraw draw;
1532: PetscDrawLGGetDraw(drawlg,&draw);
1533: PetscDrawDestroy(draw);
1534: PetscDrawLGDestroy(drawlg);
1535: return(0);
1536: }
1540: /*@
1541: TSGetTime - Gets the current time.
1543: Not Collective
1545: Input Parameter:
1546: . ts - the TS context obtained from TSCreate()
1548: Output Parameter:
1549: . t - the current time
1551: Contributed by: Matthew Knepley
1553: Level: beginner
1555: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1557: .keywords: TS, get, time
1558: @*/
1559: PetscErrorCode TSGetTime(TS ts,PetscReal* t)
1560: {
1564: *t = ts->ptime;
1565: return(0);
1566: }
1570: /*@
1571: TSSetTime - Allows one to reset the time.
1573: Collective on TS
1575: Input Parameters:
1576: + ts - the TS context obtained from TSCreate()
1577: - time - the time
1579: Level: intermediate
1581: .seealso: TSGetTime(), TSSetDuration()
1583: .keywords: TS, set, time
1584: @*/
1585: PetscErrorCode TSSetTime(TS ts, PetscReal t)
1586: {
1589: ts->ptime = t;
1590: return(0);
1591: }
1595: /*@C
1596: TSSetOptionsPrefix - Sets the prefix used for searching for all
1597: TS options in the database.
1599: Collective on TS
1601: Input Parameter:
1602: + ts - The TS context
1603: - prefix - The prefix to prepend to all option names
1605: Notes:
1606: A hyphen (-) must NOT be given at the beginning of the prefix name.
1607: The first character of all runtime options is AUTOMATICALLY the
1608: hyphen.
1610: Contributed by: Matthew Knepley
1612: Level: advanced
1614: .keywords: TS, set, options, prefix, database
1616: .seealso: TSSetFromOptions()
1618: @*/
1619: PetscErrorCode TSSetOptionsPrefix(TS ts,const char prefix[])
1620: {
1625: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1626: switch(ts->problem_type) {
1627: case TS_NONLINEAR:
1628: if (ts->snes) {
1629: SNESSetOptionsPrefix(ts->snes,prefix);
1630: }
1631: break;
1632: case TS_LINEAR:
1633: if (ts->ksp) {
1634: KSPSetOptionsPrefix(ts->ksp,prefix);
1635: }
1636: break;
1637: }
1638: return(0);
1639: }
1644: /*@C
1645: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1646: TS options in the database.
1648: Collective on TS
1650: Input Parameter:
1651: + ts - The TS context
1652: - prefix - The prefix to prepend to all option names
1654: Notes:
1655: A hyphen (-) must NOT be given at the beginning of the prefix name.
1656: The first character of all runtime options is AUTOMATICALLY the
1657: hyphen.
1659: Contributed by: Matthew Knepley
1661: Level: advanced
1663: .keywords: TS, append, options, prefix, database
1665: .seealso: TSGetOptionsPrefix()
1667: @*/
1668: PetscErrorCode TSAppendOptionsPrefix(TS ts,const char prefix[])
1669: {
1674: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1675: switch(ts->problem_type) {
1676: case TS_NONLINEAR:
1677: if (ts->snes) {
1678: SNESAppendOptionsPrefix(ts->snes,prefix);
1679: }
1680: break;
1681: case TS_LINEAR:
1682: if (ts->ksp) {
1683: KSPAppendOptionsPrefix(ts->ksp,prefix);
1684: }
1685: break;
1686: }
1687: return(0);
1688: }
1692: /*@C
1693: TSGetOptionsPrefix - Sets the prefix used for searching for all
1694: TS options in the database.
1696: Not Collective
1698: Input Parameter:
1699: . ts - The TS context
1701: Output Parameter:
1702: . prefix - A pointer to the prefix string used
1704: Contributed by: Matthew Knepley
1706: Notes: On the fortran side, the user should pass in a string 'prifix' of
1707: sufficient length to hold the prefix.
1709: Level: intermediate
1711: .keywords: TS, get, options, prefix, database
1713: .seealso: TSAppendOptionsPrefix()
1714: @*/
1715: PetscErrorCode TSGetOptionsPrefix(TS ts,const char *prefix[])
1716: {
1722: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1723: return(0);
1724: }
1728: /*@C
1729: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1731: Not Collective, but parallel objects are returned if TS is parallel
1733: Input Parameter:
1734: . ts - The TS context obtained from TSCreate()
1736: Output Parameters:
1737: + J - The Jacobian J of F, where U_t = F(U,t)
1738: . M - The preconditioner matrix, usually the same as J
1739: - ctx - User-defined context for Jacobian evaluation routine
1741: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1743: Level: intermediate
1745: .seealso: TSGetTimeStep(), TSGetMatrices(), TSGetTime(), TSGetTimeStepNumber()
1747: .keywords: TS, timestep, get, matrix, Jacobian
1748: @*/
1749: PetscErrorCode TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1750: {
1752: if (J) *J = ts->Arhs;
1753: if (M) *M = ts->B;
1754: if (ctx) *ctx = ts->jacP;
1755: return(0);
1756: }
1760: /*@C
1761: TSMonitorSolution - Monitors progress of the TS solvers by calling
1762: VecView() for the solution at each timestep
1764: Collective on TS
1766: Input Parameters:
1767: + ts - the TS context
1768: . step - current time-step
1769: . ptime - current time
1770: - dummy - either a viewer or PETSC_NULL
1772: Level: intermediate
1774: .keywords: TS, vector, monitor, view
1776: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
1777: @*/
1778: PetscErrorCode TSMonitorSolution(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1779: {
1781: PetscViewer viewer = (PetscViewer) dummy;
1784: if (!dummy) {
1785: viewer = PETSC_VIEWER_DRAW_(ts->comm);
1786: }
1787: VecView(x,viewer);
1788: return(0);
1789: }