Actual source code: ts.c
1: #include "src/ts/tsimpl.h" /*I "petscts.h" I*/
3: /* Logging support */
4: PetscCookie TS_COOKIE = 0;
5: PetscEvent TS_Step = 0, TS_PseudoComputeTimeStep = 0, TS_FunctionEval = 0, TS_JacobianEval = 0;
9: /*
10: TSSetTypeFromOptions - Sets the type of ts from user options.
12: Collective on TS
14: Input Parameter:
15: . ts - The ts
17: Level: intermediate
19: .keywords: TS, set, options, database, type
20: .seealso: TSSetFromOptions(), TSSetType()
21: */
22: static PetscErrorCode TSSetTypeFromOptions(TS ts)
23: {
24: PetscTruth opt;
25: const char *defaultType;
26: char typeName[256];
30: if (ts->type_name != PETSC_NULL) {
31: defaultType = ts->type_name;
32: } else {
33: defaultType = TS_EULER;
34: }
36: if (!TSRegisterAllCalled) {
37: TSRegisterAll(PETSC_NULL);
38: }
39: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
40: if (opt) {
41: TSSetType(ts, typeName);
42: } else {
43: TSSetType(ts, defaultType);
44: }
45: return(0);
46: }
50: /*@
51: TSSetFromOptions - Sets various TS parameters from user options.
53: Collective on TS
55: Input Parameter:
56: . ts - the TS context obtained from TSCreate()
58: Options Database Keys:
59: + -ts_type <type> - TS_EULER, TS_BEULER, TS_PVODE, TS_PSEUDO, TS_CRANK_NICHOLSON
60: . -ts_max_steps maxsteps - maximum number of time-steps to take
61: . -ts_max_time time - maximum time to compute to
62: . -ts_dt dt - initial time step
63: . -ts_monitor - print information at each timestep
64: - -ts_xmonitor - plot information at each timestep
66: Level: beginner
68: .keywords: TS, timestep, set, options, database
70: .seealso: TSGetType
71: @*/
72: PetscErrorCode TSSetFromOptions(TS ts)
73: {
74: PetscReal dt;
75: PetscTruth opt;
80: PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");
82: /* Handle generic TS options */
83: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
84: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
85: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
86: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
87: if (opt) {
88: ts->initial_time_step = ts->time_step = dt;
89: }
91: /* Monitor options */
92: PetscOptionsName("-ts_monitor","Monitor timestep size","TSDefaultMonitor",&opt);
93: if (opt) {
94: TSSetMonitor(ts,TSDefaultMonitor,PETSC_NULL,PETSC_NULL);
95: }
96: PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
97: if (opt) {
98: TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
99: }
100: PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
101: if (opt) {
102: TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
103: }
105: /* Handle TS type options */
106: TSSetTypeFromOptions(ts);
108: /* Handle specific TS options */
109: if (ts->ops->setfromoptions != PETSC_NULL) {
110: (*ts->ops->setfromoptions)(ts);
111: }
112: PetscOptionsEnd();
114: /* Handle subobject options */
115: switch(ts->problem_type) {
116: /* Should check for implicit/explicit */
117: case TS_LINEAR:
118: if (ts->ksp != PETSC_NULL) {
119: KSPSetOperators(ts->ksp,ts->A,ts->B,DIFFERENT_NONZERO_PATTERN);
120: KSPSetFromOptions(ts->ksp);
121: }
122: break;
123: case TS_NONLINEAR:
124: if (ts->snes != PETSC_NULL) {
125: /* this is a bit of a hack, but it gets the matrix information into SNES earlier
126: so that SNES and KSP have more information to pick reasonable defaults
127: before they allow users to set options */
128: SNESSetJacobian(ts->snes,ts->A,ts->B,0,ts);
129: SNESSetFromOptions(ts->snes);
130: }
131: break;
132: default:
133: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
134: }
136: return(0);
137: }
139: #undef __FUNCT__
141: /*@
142: TSViewFromOptions - This function visualizes the ts based upon user options.
144: Collective on TS
146: Input Parameter:
147: . ts - The ts
149: Level: intermediate
151: .keywords: TS, view, options, database
152: .seealso: TSSetFromOptions(), TSView()
153: @*/
154: PetscErrorCode TSViewFromOptions(TS ts,const char title[])
155: {
156: PetscViewer viewer;
157: PetscDraw draw;
158: PetscTruth opt;
159: char typeName[1024];
160: char fileName[PETSC_MAX_PATH_LEN];
161: size_t len;
165: PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
166: if (opt) {
167: PetscOptionsGetString(ts->prefix, "-ts_view", typeName, 1024, &opt);
168: PetscStrlen(typeName, &len);
169: if (len > 0) {
170: PetscViewerCreate(ts->comm, &viewer);
171: PetscViewerSetType(viewer, typeName);
172: PetscOptionsGetString(ts->prefix, "-ts_view_file", fileName, 1024, &opt);
173: if (opt) {
174: PetscViewerSetFilename(viewer, fileName);
175: } else {
176: PetscViewerSetFilename(viewer, ts->name);
177: }
178: TSView(ts, viewer);
179: PetscViewerFlush(viewer);
180: PetscViewerDestroy(viewer);
181: } else {
182: TSView(ts, PETSC_NULL);
183: }
184: }
185: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
186: if (opt) {
187: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
188: PetscViewerDrawGetDraw(viewer, 0, &draw);
189: if (title) {
190: PetscDrawSetTitle(draw, (char *)title);
191: } else {
192: PetscObjectName((PetscObject) ts);
193: PetscDrawSetTitle(draw, ts->name);
194: }
195: TSView(ts, viewer);
196: PetscViewerFlush(viewer);
197: PetscDrawPause(draw);
198: PetscViewerDestroy(viewer);
199: }
200: return(0);
201: }
205: /*@
206: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
207: set with TSSetRHSJacobian().
209: Collective on TS and Vec
211: Input Parameters:
212: + ts - the SNES context
213: . t - current timestep
214: - x - input vector
216: Output Parameters:
217: + A - Jacobian matrix
218: . B - optional preconditioning matrix
219: - flag - flag indicating matrix structure
221: Notes:
222: Most users should not need to explicitly call this routine, as it
223: is used internally within the nonlinear solvers.
225: See KSPSetOperators() for important information about setting the
226: flag parameter.
228: TSComputeJacobian() is valid only for TS_NONLINEAR
230: Level: developer
232: .keywords: SNES, compute, Jacobian, matrix
234: .seealso: TSSetRHSJacobian(), KSPSetOperators()
235: @*/
236: PetscErrorCode TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
237: {
244: if (ts->problem_type != TS_NONLINEAR) {
245: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
246: }
247: if (ts->ops->rhsjacobian) {
248: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
249: *flg = DIFFERENT_NONZERO_PATTERN;
250: PetscStackPush("TS user Jacobian function");
251: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
252: PetscStackPop;
253: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
254: /* make sure user returned a correct Jacobian and preconditioner */
257: } else {
258: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
259: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
260: if (*A != *B) {
261: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
262: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
263: }
264: }
265: return(0);
266: }
270: /*
271: TSComputeRHSFunction - Evaluates the right-hand-side function.
273: Note: If the user did not provide a function but merely a matrix,
274: this routine applies the matrix.
275: */
276: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
277: {
285: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
286: if (ts->ops->rhsfunction) {
287: PetscStackPush("TS user right-hand-side function");
288: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
289: PetscStackPop;
290: } else {
291: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
292: MatStructure flg;
293: PetscStackPush("TS user right-hand-side matrix function");
294: (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
295: PetscStackPop;
296: }
297: MatMult(ts->A,x,y);
298: }
300: /* apply user-provided boundary conditions (only needed if these are time dependent) */
301: TSComputeRHSBoundaryConditions(ts,t,y);
302: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
304: return(0);
305: }
309: /*@C
310: TSSetRHSFunction - Sets the routine for evaluating the function,
311: F(t,u), where U_t = F(t,u).
313: Collective on TS
315: Input Parameters:
316: + ts - the TS context obtained from TSCreate()
317: . f - routine for evaluating the right-hand-side function
318: - ctx - [optional] user-defined context for private data for the
319: function evaluation routine (may be PETSC_NULL)
321: Calling sequence of func:
322: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
324: + t - current timestep
325: . u - input vector
326: . F - function vector
327: - ctx - [optional] user-defined function context
329: Important:
330: The user MUST call either this routine or TSSetRHSMatrix().
332: Level: beginner
334: .keywords: TS, timestep, set, right-hand-side, function
336: .seealso: TSSetRHSMatrix()
337: @*/
338: PetscErrorCode TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
339: {
343: if (ts->problem_type == TS_LINEAR) {
344: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
345: }
346: ts->ops->rhsfunction = f;
347: ts->funP = ctx;
348: return(0);
349: }
353: /*@C
354: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
355: Also sets the location to store A.
357: Collective on TS
359: Input Parameters:
360: + ts - the TS context obtained from TSCreate()
361: . A - matrix
362: . B - preconditioner matrix (usually same as A)
363: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
364: if A is not a function of t.
365: - ctx - [optional] user-defined context for private data for the
366: matrix evaluation routine (may be PETSC_NULL)
368: Calling sequence of func:
369: $ func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
371: + t - current timestep
372: . A - matrix A, where U_t = A(t) U
373: . B - preconditioner matrix, usually the same as A
374: . flag - flag indicating information about the preconditioner matrix
375: structure (same as flag in KSPSetOperators())
376: - ctx - [optional] user-defined context for matrix evaluation routine
378: Notes:
379: See KSPSetOperators() for important information about setting the flag
380: output parameter in the routine func(). Be sure to read this information!
382: The routine func() takes Mat * as the matrix arguments rather than Mat.
383: This allows the matrix evaluation routine to replace A and/or B with a
384: completely new new matrix structure (not just different matrix elements)
385: when appropriate, for instance, if the nonzero structure is changing
386: throughout the global iterations.
388: Important:
389: The user MUST call either this routine or TSSetRHSFunction().
391: Level: beginner
393: .keywords: TS, timestep, set, right-hand-side, matrix
395: .seealso: TSSetRHSFunction()
396: @*/
397: PetscErrorCode TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
398: {
405: if (ts->problem_type == TS_NONLINEAR) {
406: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
407: }
409: ts->ops->rhsmatrix = f;
410: ts->jacP = ctx;
411: ts->A = A;
412: ts->B = B;
414: return(0);
415: }
419: /*@C
420: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
421: where U_t = F(U,t), as well as the location to store the matrix.
422: Use TSSetRHSMatrix() for linear problems.
424: Collective on TS
426: Input Parameters:
427: + ts - the TS context obtained from TSCreate()
428: . A - Jacobian matrix
429: . B - preconditioner matrix (usually same as A)
430: . f - the Jacobian evaluation routine
431: - ctx - [optional] user-defined context for private data for the
432: Jacobian evaluation routine (may be PETSC_NULL)
434: Calling sequence of func:
435: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
437: + t - current timestep
438: . u - input vector
439: . A - matrix A, where U_t = A(t)u
440: . B - preconditioner matrix, usually the same as A
441: . flag - flag indicating information about the preconditioner matrix
442: structure (same as flag in KSPSetOperators())
443: - ctx - [optional] user-defined context for matrix evaluation routine
445: Notes:
446: See KSPSetOperators() for important information about setting the flag
447: output parameter in the routine func(). Be sure to read this information!
449: The routine func() takes Mat * as the matrix arguments rather than Mat.
450: This allows the matrix evaluation routine to replace A and/or B with a
451: completely new new matrix structure (not just different matrix elements)
452: when appropriate, for instance, if the nonzero structure is changing
453: throughout the global iterations.
455: Level: beginner
456:
457: .keywords: TS, timestep, set, right-hand-side, Jacobian
459: .seealso: TSDefaultComputeJacobianColor(),
460: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()
462: @*/
463: PetscErrorCode TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
464: {
471: if (ts->problem_type != TS_NONLINEAR) {
472: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
473: }
475: ts->ops->rhsjacobian = f;
476: ts->jacP = ctx;
477: ts->A = A;
478: ts->B = B;
479: return(0);
480: }
484: /*
485: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
487: Note: If the user did not provide a function but merely a matrix,
488: this routine applies the matrix.
489: */
490: PetscErrorCode TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
491: {
499: if (ts->ops->rhsbc) {
500: PetscStackPush("TS user boundary condition function");
501: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
502: PetscStackPop;
503: return(0);
504: }
506: return(0);
507: }
511: /*@C
512: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
513: boundary conditions for the function F.
515: Collective on TS
517: Input Parameters:
518: + ts - the TS context obtained from TSCreate()
519: . f - routine for evaluating the boundary condition function
520: - ctx - [optional] user-defined context for private data for the
521: function evaluation routine (may be PETSC_NULL)
523: Calling sequence of func:
524: $ func (TS ts,PetscReal t,Vec F,void *ctx);
526: + t - current timestep
527: . F - function vector
528: - ctx - [optional] user-defined function context
530: Level: intermediate
532: .keywords: TS, timestep, set, boundary conditions, function
533: @*/
534: PetscErrorCode TSSetRHSBoundaryConditions(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
535: {
539: if (ts->problem_type != TS_LINEAR) {
540: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
541: }
542: ts->ops->rhsbc = f;
543: ts->bcP = ctx;
544: return(0);
545: }
549: /*@C
550: TSView - Prints the TS data structure.
552: Collective on TS
554: Input Parameters:
555: + ts - the TS context obtained from TSCreate()
556: - viewer - visualization context
558: Options Database Key:
559: . -ts_view - calls TSView() at end of TSStep()
561: Notes:
562: The available visualization contexts include
563: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
564: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
565: output where only the first processor opens
566: the file. All other processors send their
567: data to the first processor to print.
569: The user can open an alternative visualization context with
570: PetscViewerASCIIOpen() - output to a specified file.
572: Level: beginner
574: .keywords: TS, timestep, view
576: .seealso: PetscViewerASCIIOpen()
577: @*/
578: PetscErrorCode TSView(TS ts,PetscViewer viewer)
579: {
581: char *type;
582: PetscTruth iascii,isstring;
586: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
590: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
591: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
592: if (iascii) {
593: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
594: TSGetType(ts,(TSType *)&type);
595: if (type) {
596: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
597: } else {
598: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
599: }
600: if (ts->ops->view) {
601: PetscViewerASCIIPushTab(viewer);
602: (*ts->ops->view)(ts,viewer);
603: PetscViewerASCIIPopTab(viewer);
604: }
605: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
606: PetscViewerASCIIPrintf(viewer," maximum time=%g\n",ts->max_time);
607: if (ts->problem_type == TS_NONLINEAR) {
608: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
609: }
610: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
611: } else if (isstring) {
612: TSGetType(ts,(TSType *)&type);
613: PetscViewerStringSPrintf(viewer," %-7.7s",type);
614: }
615: PetscViewerASCIIPushTab(viewer);
616: if (ts->ksp) {KSPView(ts->ksp,viewer);}
617: if (ts->snes) {SNESView(ts->snes,viewer);}
618: PetscViewerASCIIPopTab(viewer);
619: return(0);
620: }
625: /*@C
626: TSSetApplicationContext - Sets an optional user-defined context for
627: the timesteppers.
629: Collective on TS
631: Input Parameters:
632: + ts - the TS context obtained from TSCreate()
633: - usrP - optional user context
635: Level: intermediate
637: .keywords: TS, timestep, set, application, context
639: .seealso: TSGetApplicationContext()
640: @*/
641: PetscErrorCode TSSetApplicationContext(TS ts,void *usrP)
642: {
645: ts->user = usrP;
646: return(0);
647: }
651: /*@C
652: TSGetApplicationContext - Gets the user-defined context for the
653: timestepper.
655: Not Collective
657: Input Parameter:
658: . ts - the TS context obtained from TSCreate()
660: Output Parameter:
661: . usrP - user context
663: Level: intermediate
665: .keywords: TS, timestep, get, application, context
667: .seealso: TSSetApplicationContext()
668: @*/
669: PetscErrorCode TSGetApplicationContext(TS ts,void **usrP)
670: {
673: *usrP = ts->user;
674: return(0);
675: }
679: /*@
680: TSGetTimeStepNumber - Gets the current number of timesteps.
682: Not Collective
684: Input Parameter:
685: . ts - the TS context obtained from TSCreate()
687: Output Parameter:
688: . iter - number steps so far
690: Level: intermediate
692: .keywords: TS, timestep, get, iteration, number
693: @*/
694: PetscErrorCode TSGetTimeStepNumber(TS ts,PetscInt* iter)
695: {
699: *iter = ts->steps;
700: return(0);
701: }
705: /*@
706: TSSetInitialTimeStep - Sets the initial timestep to be used,
707: as well as the initial time.
709: Collective on TS
711: Input Parameters:
712: + ts - the TS context obtained from TSCreate()
713: . initial_time - the initial time
714: - time_step - the size of the timestep
716: Level: intermediate
718: .seealso: TSSetTimeStep(), TSGetTimeStep()
720: .keywords: TS, set, initial, timestep
721: @*/
722: PetscErrorCode TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
723: {
726: ts->time_step = time_step;
727: ts->initial_time_step = time_step;
728: ts->ptime = initial_time;
729: return(0);
730: }
734: /*@
735: TSSetTimeStep - Allows one to reset the timestep at any time,
736: useful for simple pseudo-timestepping codes.
738: Collective on TS
740: Input Parameters:
741: + ts - the TS context obtained from TSCreate()
742: - time_step - the size of the timestep
744: Level: intermediate
746: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
748: .keywords: TS, set, timestep
749: @*/
750: PetscErrorCode TSSetTimeStep(TS ts,PetscReal time_step)
751: {
754: ts->time_step = time_step;
755: return(0);
756: }
760: /*@
761: TSGetTimeStep - Gets the current timestep size.
763: Not Collective
765: Input Parameter:
766: . ts - the TS context obtained from TSCreate()
768: Output Parameter:
769: . dt - the current timestep size
771: Level: intermediate
773: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
775: .keywords: TS, get, timestep
776: @*/
777: PetscErrorCode TSGetTimeStep(TS ts,PetscReal* dt)
778: {
782: *dt = ts->time_step;
783: return(0);
784: }
788: /*@C
789: TSGetSolution - Returns the solution at the present timestep. It
790: is valid to call this routine inside the function that you are evaluating
791: in order to move to the new timestep. This vector not changed until
792: the solution at the next timestep has been calculated.
794: Not Collective, but Vec returned is parallel if TS is parallel
796: Input Parameter:
797: . ts - the TS context obtained from TSCreate()
799: Output Parameter:
800: . v - the vector containing the solution
802: Level: intermediate
804: .seealso: TSGetTimeStep()
806: .keywords: TS, timestep, get, solution
807: @*/
808: PetscErrorCode TSGetSolution(TS ts,Vec *v)
809: {
813: *v = ts->vec_sol_always;
814: return(0);
815: }
817: /* ----- Routines to initialize and destroy a timestepper ---- */
820: /*@C
821: TSSetProblemType - Sets the type of problem to be solved.
823: Not collective
825: Input Parameters:
826: + ts - The TS
827: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
828: .vb
829: U_t = A U
830: U_t = A(t) U
831: U_t = F(t,U)
832: .ve
834: Level: beginner
836: .keywords: TS, problem type
837: .seealso: TSSetUp(), TSProblemType, TS
838: @*/
839: PetscErrorCode TSSetProblemType(TS ts, TSProblemType type)
840: {
843: ts->problem_type = type;
844: return(0);
845: }
849: /*@C
850: TSGetProblemType - Gets the type of problem to be solved.
852: Not collective
854: Input Parameter:
855: . ts - The TS
857: Output Parameter:
858: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
859: .vb
860: U_t = A U
861: U_t = A(t) U
862: U_t = F(t,U)
863: .ve
865: Level: beginner
867: .keywords: TS, problem type
868: .seealso: TSSetUp(), TSProblemType, TS
869: @*/
870: PetscErrorCode TSGetProblemType(TS ts, TSProblemType *type)
871: {
875: *type = ts->problem_type;
876: return(0);
877: }
881: /*@
882: TSSetUp - Sets up the internal data structures for the later use
883: of a timestepper.
885: Collective on TS
887: Input Parameter:
888: . ts - the TS context obtained from TSCreate()
890: Notes:
891: For basic use of the TS solvers the user need not explicitly call
892: TSSetUp(), since these actions will automatically occur during
893: the call to TSStep(). However, if one wishes to control this
894: phase separately, TSSetUp() should be called after TSCreate()
895: and optional routines of the form TSSetXXX(), but before TSStep().
897: Level: advanced
899: .keywords: TS, timestep, setup
901: .seealso: TSCreate(), TSStep(), TSDestroy()
902: @*/
903: PetscErrorCode TSSetUp(TS ts)
904: {
909: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
910: if (!ts->type_name) {
911: TSSetType(ts,TS_EULER);
912: }
913: (*ts->ops->setup)(ts);
914: ts->setupcalled = 1;
915: return(0);
916: }
920: /*@C
921: TSDestroy - Destroys the timestepper context that was created
922: with TSCreate().
924: Collective on TS
926: Input Parameter:
927: . ts - the TS context obtained from TSCreate()
929: Level: beginner
931: .keywords: TS, timestepper, destroy
933: .seealso: TSCreate(), TSSetUp(), TSSolve()
934: @*/
935: PetscErrorCode TSDestroy(TS ts)
936: {
938: PetscInt i;
942: if (--ts->refct > 0) return(0);
944: /* if memory was published with AMS then destroy it */
945: PetscObjectDepublish(ts);
947: if (ts->ksp) {KSPDestroy(ts->ksp);}
948: if (ts->snes) {SNESDestroy(ts->snes);}
949: (*(ts)->ops->destroy)(ts);
950: for (i=0; i<ts->numbermonitors; i++) {
951: if (ts->mdestroy[i]) {
952: (*ts->mdestroy[i])(ts->monitorcontext[i]);
953: }
954: }
955: PetscLogObjectDestroy((PetscObject)ts);
956: PetscHeaderDestroy((PetscObject)ts);
957: return(0);
958: }
962: /*@C
963: TSGetSNES - Returns the SNES (nonlinear solver) associated with
964: a TS (timestepper) context. Valid only for nonlinear problems.
966: Not Collective, but SNES is parallel if TS is parallel
968: Input Parameter:
969: . ts - the TS context obtained from TSCreate()
971: Output Parameter:
972: . snes - the nonlinear solver context
974: Notes:
975: The user can then directly manipulate the SNES context to set various
976: options, etc. Likewise, the user can then extract and manipulate the
977: KSP, KSP, and PC contexts as well.
979: TSGetSNES() does not work for integrators that do not use SNES; in
980: this case TSGetSNES() returns PETSC_NULL in snes.
982: Level: beginner
984: .keywords: timestep, get, SNES
985: @*/
986: PetscErrorCode TSGetSNES(TS ts,SNES *snes)
987: {
991: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
992: *snes = ts->snes;
993: return(0);
994: }
998: /*@C
999: TSGetKSP - Returns the KSP (linear solver) associated with
1000: a TS (timestepper) context.
1002: Not Collective, but KSP is parallel if TS is parallel
1004: Input Parameter:
1005: . ts - the TS context obtained from TSCreate()
1007: Output Parameter:
1008: . ksp - the nonlinear solver context
1010: Notes:
1011: The user can then directly manipulate the KSP context to set various
1012: options, etc. Likewise, the user can then extract and manipulate the
1013: KSP and PC contexts as well.
1015: TSGetKSP() does not work for integrators that do not use KSP;
1016: in this case TSGetKSP() returns PETSC_NULL in ksp.
1018: Level: beginner
1020: .keywords: timestep, get, KSP
1021: @*/
1022: PetscErrorCode TSGetKSP(TS ts,KSP *ksp)
1023: {
1027: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1028: *ksp = ts->ksp;
1029: return(0);
1030: }
1032: /* ----------- Routines to set solver parameters ---------- */
1036: /*@
1037: TSGetDuration - Gets the maximum number of timesteps to use and
1038: maximum time for iteration.
1040: Collective on TS
1042: Input Parameters:
1043: + ts - the TS context obtained from TSCreate()
1044: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1045: - maxtime - final time to iterate to, or PETSC_NULL
1047: Level: intermediate
1049: .keywords: TS, timestep, get, maximum, iterations, time
1050: @*/
1051: PetscErrorCode TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1052: {
1055: if (maxsteps != PETSC_NULL) {
1057: *maxsteps = ts->max_steps;
1058: }
1059: if (maxtime != PETSC_NULL) {
1061: *maxtime = ts->max_time;
1062: }
1063: return(0);
1064: }
1068: /*@
1069: TSSetDuration - Sets the maximum number of timesteps to use and
1070: maximum time for iteration.
1072: Collective on TS
1074: Input Parameters:
1075: + ts - the TS context obtained from TSCreate()
1076: . maxsteps - maximum number of iterations to use
1077: - maxtime - final time to iterate to
1079: Options Database Keys:
1080: . -ts_max_steps <maxsteps> - Sets maxsteps
1081: . -ts_max_time <maxtime> - Sets maxtime
1083: Notes:
1084: The default maximum number of iterations is 5000. Default time is 5.0
1086: Level: intermediate
1088: .keywords: TS, timestep, set, maximum, iterations
1089: @*/
1090: PetscErrorCode TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1091: {
1094: ts->max_steps = maxsteps;
1095: ts->max_time = maxtime;
1096: return(0);
1097: }
1101: /*@
1102: TSSetSolution - Sets the initial solution vector
1103: for use by the TS routines.
1105: Collective on TS and Vec
1107: Input Parameters:
1108: + ts - the TS context obtained from TSCreate()
1109: - x - the solution vector
1111: Level: beginner
1113: .keywords: TS, timestep, set, solution, initial conditions
1114: @*/
1115: PetscErrorCode TSSetSolution(TS ts,Vec x)
1116: {
1120: ts->vec_sol = ts->vec_sol_always = x;
1121: return(0);
1122: }
1126: /*@C
1127: TSSetRhsBC - Sets the function which applies boundary conditions
1128: to the Rhs of each system.
1130: Collective on TS
1132: Input Parameters:
1133: + ts - The TS context obtained from TSCreate()
1134: - func - The function
1136: Calling sequence of func:
1137: . func (TS ts, Vec rhs, void *ctx);
1139: + rhs - The current rhs vector
1140: - ctx - The user-context
1142: Level: intermediate
1144: .keywords: TS, Rhs, boundary conditions
1145: @*/
1146: PetscErrorCode TSSetRhsBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1147: {
1150: ts->ops->applyrhsbc = func;
1151: return(0);
1152: }
1156: /*@
1157: TSDefaultRhsBC - The default boundary condition function which does nothing.
1159: Collective on TS
1161: Input Parameters:
1162: + ts - The TS context obtained from TSCreate()
1163: . rhs - The Rhs
1164: - ctx - The user-context
1166: Level: developer
1168: .keywords: TS, Rhs, boundary conditions
1169: @*/
1170: PetscErrorCode TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1171: {
1173: return(0);
1174: }
1178: /*@C
1179: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1180: to the system matrix and preconditioner of each system.
1182: Collective on TS
1184: Input Parameters:
1185: + ts - The TS context obtained from TSCreate()
1186: - func - The function
1188: Calling sequence of func:
1189: . func (TS ts, Mat A, Mat B, void *ctx);
1191: + A - The current system matrix
1192: . B - The current preconditioner
1193: - ctx - The user-context
1195: Level: intermediate
1197: .keywords: TS, System matrix, boundary conditions
1198: @*/
1199: PetscErrorCode TSSetSystemMatrixBC(TS ts, PetscErrorCode (*func)(TS, Mat, Mat, void *))
1200: {
1203: ts->ops->applymatrixbc = func;
1204: return(0);
1205: }
1209: /*@
1210: TSDefaultSystemMatrixBC - The default boundary condition function which
1211: does nothing.
1213: Collective on TS
1215: Input Parameters:
1216: + ts - The TS context obtained from TSCreate()
1217: . A - The system matrix
1218: . B - The preconditioner
1219: - ctx - The user-context
1221: Level: developer
1223: .keywords: TS, System matrix, boundary conditions
1224: @*/
1225: PetscErrorCode TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1226: {
1228: return(0);
1229: }
1233: /*@C
1234: TSSetSolutionBC - Sets the function which applies boundary conditions
1235: to the solution of each system. This is necessary in nonlinear systems
1236: which time dependent boundary conditions.
1238: Collective on TS
1240: Input Parameters:
1241: + ts - The TS context obtained from TSCreate()
1242: - func - The function
1244: Calling sequence of func:
1245: . func (TS ts, Vec rsol, void *ctx);
1247: + sol - The current solution vector
1248: - ctx - The user-context
1250: Level: intermediate
1252: .keywords: TS, solution, boundary conditions
1253: @*/
1254: PetscErrorCode TSSetSolutionBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1255: {
1258: ts->ops->applysolbc = func;
1259: return(0);
1260: }
1264: /*@
1265: TSDefaultSolutionBC - The default boundary condition function which
1266: does nothing.
1268: Collective on TS
1270: Input Parameters:
1271: + ts - The TS context obtained from TSCreate()
1272: . sol - The solution
1273: - ctx - The user-context
1275: Level: developer
1277: .keywords: TS, solution, boundary conditions
1278: @*/
1279: PetscErrorCode TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1280: {
1282: return(0);
1283: }
1287: /*@C
1288: TSSetPreStep - Sets the general-purpose function
1289: called once at the beginning of time stepping.
1291: Collective on TS
1293: Input Parameters:
1294: + ts - The TS context obtained from TSCreate()
1295: - func - The function
1297: Calling sequence of func:
1298: . func (TS ts);
1300: Level: intermediate
1302: .keywords: TS, timestep
1303: @*/
1304: PetscErrorCode TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1305: {
1308: ts->ops->prestep = func;
1309: return(0);
1310: }
1314: /*@
1315: TSDefaultPreStep - The default pre-stepping function which does nothing.
1317: Collective on TS
1319: Input Parameters:
1320: . ts - The TS context obtained from TSCreate()
1322: Level: developer
1324: .keywords: TS, timestep
1325: @*/
1326: PetscErrorCode TSDefaultPreStep(TS ts)
1327: {
1329: return(0);
1330: }
1334: /*@C
1335: TSSetUpdate - Sets the general-purpose update function called
1336: at the beginning of every time step. This function can change
1337: the time step.
1339: Collective on TS
1341: Input Parameters:
1342: + ts - The TS context obtained from TSCreate()
1343: - func - The function
1345: Calling sequence of func:
1346: . func (TS ts, double t, double *dt);
1348: + t - The current time
1349: - dt - The current time step
1351: Level: intermediate
1353: .keywords: TS, update, timestep
1354: @*/
1355: PetscErrorCode TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1356: {
1359: ts->ops->update = func;
1360: return(0);
1361: }
1365: /*@
1366: TSDefaultUpdate - The default update function which does nothing.
1368: Collective on TS
1370: Input Parameters:
1371: + ts - The TS context obtained from TSCreate()
1372: - t - The current time
1374: Output Parameters:
1375: . dt - The current time step
1377: Level: developer
1379: .keywords: TS, update, timestep
1380: @*/
1381: PetscErrorCode TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1382: {
1384: return(0);
1385: }
1389: /*@C
1390: TSSetPostStep - Sets the general-purpose function
1391: called once at the end of time stepping.
1393: Collective on TS
1395: Input Parameters:
1396: + ts - The TS context obtained from TSCreate()
1397: - func - The function
1399: Calling sequence of func:
1400: . func (TS ts);
1402: Level: intermediate
1404: .keywords: TS, timestep
1405: @*/
1406: PetscErrorCode TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1407: {
1410: ts->ops->poststep = func;
1411: return(0);
1412: }
1416: /*@
1417: TSDefaultPostStep - The default post-stepping function which does nothing.
1419: Collective on TS
1421: Input Parameters:
1422: . ts - The TS context obtained from TSCreate()
1424: Level: developer
1426: .keywords: TS, timestep
1427: @*/
1428: PetscErrorCode TSDefaultPostStep(TS ts)
1429: {
1431: return(0);
1432: }
1434: /* ------------ Routines to set performance monitoring options ----------- */
1438: /*@C
1439: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1440: timestep to display the iteration's progress.
1442: Collective on TS
1444: Input Parameters:
1445: + ts - the TS context obtained from TSCreate()
1446: . func - monitoring routine
1447: . mctx - [optional] user-defined context for private data for the
1448: monitor routine (use PETSC_NULL if no context is desired)
1449: - monitordestroy - [optional] routine that frees monitor context
1450: (may be PETSC_NULL)
1452: Calling sequence of func:
1453: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1455: + ts - the TS context
1456: . steps - iteration number
1457: . time - current time
1458: . x - current iterate
1459: - mctx - [optional] monitoring context
1461: Notes:
1462: This routine adds an additional monitor to the list of monitors that
1463: already has been loaded.
1465: Level: intermediate
1467: .keywords: TS, timestep, set, monitor
1469: .seealso: TSDefaultMonitor(), TSClearMonitor()
1470: @*/
1471: PetscErrorCode TSSetMonitor(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1472: {
1475: if (ts->numbermonitors >= MAXTSMONITORS) {
1476: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1477: }
1478: ts->monitor[ts->numbermonitors] = monitor;
1479: ts->mdestroy[ts->numbermonitors] = mdestroy;
1480: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1481: return(0);
1482: }
1486: /*@C
1487: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1489: Collective on TS
1491: Input Parameters:
1492: . ts - the TS context obtained from TSCreate()
1494: Notes:
1495: There is no way to remove a single, specific monitor.
1497: Level: intermediate
1499: .keywords: TS, timestep, set, monitor
1501: .seealso: TSDefaultMonitor(), TSSetMonitor()
1502: @*/
1503: PetscErrorCode TSClearMonitor(TS ts)
1504: {
1507: ts->numbermonitors = 0;
1508: return(0);
1509: }
1513: PetscErrorCode TSDefaultMonitor(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1514: {
1518: PetscPrintf(ts->comm,"timestep %D dt %g time %g\n",step,ts->time_step,ptime);
1519: return(0);
1520: }
1524: /*@
1525: TSStep - Steps the requested number of timesteps.
1527: Collective on TS
1529: Input Parameter:
1530: . ts - the TS context obtained from TSCreate()
1532: Output Parameters:
1533: + steps - number of iterations until termination
1534: - ptime - time until termination
1536: Level: beginner
1538: .keywords: TS, timestep, solve
1540: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1541: @*/
1542: PetscErrorCode TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1543: {
1548: if (!ts->setupcalled) {
1549: TSSetUp(ts);
1550: }
1552: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1553: (*ts->ops->prestep)(ts);
1554: (*ts->ops->step)(ts, steps, ptime);
1555: (*ts->ops->poststep)(ts);
1556: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1558: if (!PetscPreLoadingOn) {
1559: TSViewFromOptions(ts,ts->name);
1560: }
1561: return(0);
1562: }
1566: /*
1567: Runs the user provided monitor routines, if they exists.
1568: */
1569: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1570: {
1572: PetscInt i,n = ts->numbermonitors;
1575: for (i=0; i<n; i++) {
1576: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1577: }
1578: return(0);
1579: }
1581: /* ------------------------------------------------------------------------*/
1585: /*@C
1586: TSLGMonitorCreate - Creates a line graph context for use with
1587: TS to monitor convergence of preconditioned residual norms.
1589: Collective on TS
1591: Input Parameters:
1592: + host - the X display to open, or null for the local machine
1593: . label - the title to put in the title bar
1594: . x, y - the screen coordinates of the upper left coordinate of the window
1595: - m, n - the screen width and height in pixels
1597: Output Parameter:
1598: . draw - the drawing context
1600: Options Database Key:
1601: . -ts_xmonitor - automatically sets line graph monitor
1603: Notes:
1604: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1606: Level: intermediate
1608: .keywords: TS, monitor, line graph, residual, seealso
1610: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1612: @*/
1613: PetscErrorCode TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1614: {
1615: PetscDraw win;
1619: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1620: PetscDrawSetType(win,PETSC_DRAW_X);
1621: PetscDrawLGCreate(win,1,draw);
1622: PetscDrawLGIndicateDataPoints(*draw);
1624: PetscLogObjectParent(*draw,win);
1625: return(0);
1626: }
1630: PetscErrorCode TSLGMonitor(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1631: {
1632: PetscDrawLG lg = (PetscDrawLG) monctx;
1633: PetscReal x,y = ptime;
1637: if (!monctx) {
1638: MPI_Comm comm;
1639: PetscViewer viewer;
1641: PetscObjectGetComm((PetscObject)ts,&comm);
1642: viewer = PETSC_VIEWER_DRAW_(comm);
1643: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1644: }
1646: if (!n) {PetscDrawLGReset(lg);}
1647: x = (PetscReal)n;
1648: PetscDrawLGAddPoint(lg,&x,&y);
1649: if (n < 20 || (n % 5)) {
1650: PetscDrawLGDraw(lg);
1651: }
1652: return(0);
1653: }
1657: /*@C
1658: TSLGMonitorDestroy - Destroys a line graph context that was created
1659: with TSLGMonitorCreate().
1661: Collective on PetscDrawLG
1663: Input Parameter:
1664: . draw - the drawing context
1666: Level: intermediate
1668: .keywords: TS, monitor, line graph, destroy
1670: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1671: @*/
1672: PetscErrorCode TSLGMonitorDestroy(PetscDrawLG drawlg)
1673: {
1674: PetscDraw draw;
1678: PetscDrawLGGetDraw(drawlg,&draw);
1679: PetscDrawDestroy(draw);
1680: PetscDrawLGDestroy(drawlg);
1681: return(0);
1682: }
1686: /*@
1687: TSGetTime - Gets the current time.
1689: Not Collective
1691: Input Parameter:
1692: . ts - the TS context obtained from TSCreate()
1694: Output Parameter:
1695: . t - the current time
1697: Contributed by: Matthew Knepley
1699: Level: beginner
1701: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1703: .keywords: TS, get, time
1704: @*/
1705: PetscErrorCode TSGetTime(TS ts,PetscReal* t)
1706: {
1710: *t = ts->ptime;
1711: return(0);
1712: }
1716: /*@C
1717: TSSetOptionsPrefix - Sets the prefix used for searching for all
1718: TS options in the database.
1720: Collective on TS
1722: Input Parameter:
1723: + ts - The TS context
1724: - prefix - The prefix to prepend to all option names
1726: Notes:
1727: A hyphen (-) must NOT be given at the beginning of the prefix name.
1728: The first character of all runtime options is AUTOMATICALLY the
1729: hyphen.
1731: Contributed by: Matthew Knepley
1733: Level: advanced
1735: .keywords: TS, set, options, prefix, database
1737: .seealso: TSSetFromOptions()
1739: @*/
1740: PetscErrorCode TSSetOptionsPrefix(TS ts,const char prefix[])
1741: {
1746: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1747: switch(ts->problem_type) {
1748: case TS_NONLINEAR:
1749: SNESSetOptionsPrefix(ts->snes,prefix);
1750: break;
1751: case TS_LINEAR:
1752: KSPSetOptionsPrefix(ts->ksp,prefix);
1753: break;
1754: }
1755: return(0);
1756: }
1761: /*@C
1762: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1763: TS options in the database.
1765: Collective on TS
1767: Input Parameter:
1768: + ts - The TS context
1769: - prefix - The prefix to prepend to all option names
1771: Notes:
1772: A hyphen (-) must NOT be given at the beginning of the prefix name.
1773: The first character of all runtime options is AUTOMATICALLY the
1774: hyphen.
1776: Contributed by: Matthew Knepley
1778: Level: advanced
1780: .keywords: TS, append, options, prefix, database
1782: .seealso: TSGetOptionsPrefix()
1784: @*/
1785: PetscErrorCode TSAppendOptionsPrefix(TS ts,const char prefix[])
1786: {
1791: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1792: switch(ts->problem_type) {
1793: case TS_NONLINEAR:
1794: SNESAppendOptionsPrefix(ts->snes,prefix);
1795: break;
1796: case TS_LINEAR:
1797: KSPAppendOptionsPrefix(ts->ksp,prefix);
1798: break;
1799: }
1800: return(0);
1801: }
1805: /*@C
1806: TSGetOptionsPrefix - Sets the prefix used for searching for all
1807: TS options in the database.
1809: Not Collective
1811: Input Parameter:
1812: . ts - The TS context
1814: Output Parameter:
1815: . prefix - A pointer to the prefix string used
1817: Contributed by: Matthew Knepley
1819: Notes: On the fortran side, the user should pass in a string 'prifix' of
1820: sufficient length to hold the prefix.
1822: Level: intermediate
1824: .keywords: TS, get, options, prefix, database
1826: .seealso: TSAppendOptionsPrefix()
1827: @*/
1828: PetscErrorCode TSGetOptionsPrefix(TS ts,char *prefix[])
1829: {
1835: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1836: return(0);
1837: }
1841: /*@C
1842: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1844: Not Collective, but parallel objects are returned if TS is parallel
1846: Input Parameter:
1847: . ts - The TS context obtained from TSCreate()
1849: Output Parameters:
1850: + A - The matrix A, where U_t = A(t) U
1851: . M - The preconditioner matrix, usually the same as A
1852: - ctx - User-defined context for matrix evaluation routine
1854: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1856: Contributed by: Matthew Knepley
1858: Level: intermediate
1860: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1862: .keywords: TS, timestep, get, matrix
1864: @*/
1865: PetscErrorCode TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1866: {
1869: if (A) *A = ts->A;
1870: if (M) *M = ts->B;
1871: if (ctx) *ctx = ts->jacP;
1872: return(0);
1873: }
1877: /*@C
1878: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1880: Not Collective, but parallel objects are returned if TS is parallel
1882: Input Parameter:
1883: . ts - The TS context obtained from TSCreate()
1885: Output Parameters:
1886: + J - The Jacobian J of F, where U_t = F(U,t)
1887: . M - The preconditioner matrix, usually the same as J
1888: - ctx - User-defined context for Jacobian evaluation routine
1890: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1892: Contributed by: Matthew Knepley
1894: Level: intermediate
1896: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1898: .keywords: TS, timestep, get, matrix, Jacobian
1899: @*/
1900: PetscErrorCode TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1901: {
1905: TSGetRHSMatrix(ts,J,M,ctx);
1906: return(0);
1907: }
1911: /*@C
1912: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1913: VecView() for the solution at each timestep
1915: Collective on TS
1917: Input Parameters:
1918: + ts - the TS context
1919: . step - current time-step
1920: . ptime - current time
1921: - dummy - either a viewer or PETSC_NULL
1923: Level: intermediate
1925: .keywords: TS, vector, monitor, view
1927: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1928: @*/
1929: PetscErrorCode TSVecViewMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1930: {
1932: PetscViewer viewer = (PetscViewer) dummy;
1935: if (!viewer) {
1936: MPI_Comm comm;
1937: PetscObjectGetComm((PetscObject)ts,&comm);
1938: viewer = PETSC_VIEWER_DRAW_(comm);
1939: }
1940: VecView(x,viewer);
1941: return(0);
1942: }