Actual source code: ts.c
1: /* $Id: ts.c,v 1.43 2001/09/07 20:12:01 bsmith Exp $ */
2: #include src/ts/tsimpl.h
4: /* Logging support */
5: int TS_COOKIE;
6: int TS_Step, TS_PseudoComputeTimeStep, TS_FunctionEval, TS_JacobianEval;
8: /*
9: TSSetTypeFromOptions - Sets the type of ts from user options.
11: Collective on TS
13: Input Parameter:
14: . ts - The ts
16: Level: intermediate
18: .keywords: TS, set, options, database, type
19: .seealso: TSSetFromOptions(), TSSetType()
20: */
21: static int TSSetTypeFromOptions(TS ts)
22: {
23: PetscTruth opt;
24: char *defaultType;
25: char typeName[256];
26: int ierr;
29: if (ts->type_name != PETSC_NULL) {
30: defaultType = ts->type_name;
31: } else {
32: defaultType = TS_EULER;
33: }
35: if (!TSRegisterAllCalled) {
36: TSRegisterAll(PETSC_NULL);
37: }
38: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
39: if (opt == PETSC_TRUE) {
40: TSSetType(ts, typeName);
41: } else {
42: TSSetType(ts, defaultType);
43: }
44: return(0);
45: }
47: /*@
48: TSSetFromOptions - Sets various TS parameters from user options.
50: Collective on TS
52: Input Parameter:
53: . ts - the TS context obtained from TSCreate()
55: Options Database Keys:
56: + -ts_type <type> - TS_EULER, TS_BEULER, TS_PVODE, TS_PSEUDO, TS_CRANK_NICHOLSON
57: . -ts_max_steps maxsteps - maximum number of time-steps to take
58: . -ts_max_time time - maximum time to compute to
59: . -ts_dt dt - initial time step
60: . -ts_monitor - print information at each timestep
61: - -ts_xmonitor - plot information at each timestep
63: Level: beginner
65: .keywords: TS, timestep, set, options, database
67: .seealso: TSGetType
68: @*/
69: int TSSetFromOptions(TS ts)
70: {
71: PetscReal dt;
72: PetscTruth opt;
73: int ierr;
77: PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");
79: /* Handle generic TS options */
80: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
81: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
82: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
83: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
84: if (opt == PETSC_TRUE) {
85: ts->initial_time_step = ts->time_step = dt;
86: }
88: /* Monitor options */
89: PetscOptionsName("-ts_monitor","Monitor timestep size","TSDefaultMonitor",&opt);
90: if (opt == PETSC_TRUE) {
91: TSSetMonitor(ts,TSDefaultMonitor,PETSC_NULL,PETSC_NULL);
92: }
93: PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
94: if (opt == PETSC_TRUE) {
95: TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
96: }
97: PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
98: if (opt == PETSC_TRUE) {
99: TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
100: }
102: /* Handle TS type options */
103: TSSetTypeFromOptions(ts);
105: /* Handle specific TS options */
106: if (ts->ops->setfromoptions != PETSC_NULL) {
107: (*ts->ops->setfromoptions)(ts);
108: }
109: PetscOptionsEnd();
111: /* Handle subobject options */
112: switch(ts->problem_type) {
113: /* Should check for implicit/explicit */
114: case TS_LINEAR:
115: if (ts->sles != PETSC_NULL) {
116: SLESSetFromOptions(ts->sles);
117: }
118: break;
119: case TS_NONLINEAR:
120: if (ts->snes != PETSC_NULL) {
121: SNESSetFromOptions(ts->snes);
122: }
123: break;
124: default:
125: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", ts->problem_type);
126: }
128: TSViewFromOptions(ts, ts->name);
129: return(0);
130: }
132: /*@
133: TSViewFromOptions - This function visualizes the ts based upon user options.
135: Collective on TS
137: Input Parameter:
138: . ts - The ts
140: Level: intermediate
142: .keywords: TS, view, options, database
143: .seealso: TSSetFromOptions(), TSView()
144: @*/
145: int TSViewFromOptions(TS ts, char *title)
146: {
147: PetscViewer viewer;
148: PetscDraw draw;
149: PetscTruth opt;
150: char *titleStr;
151: char typeName[1024];
152: char fileName[1024];
153: int len;
154: int ierr;
157: PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
158: if (opt == PETSC_TRUE) {
159: PetscOptionsGetString(ts->prefix, "-ts_view", typeName, 1024, &opt);
160: PetscStrlen(typeName, &len);
161: if (len > 0) {
162: PetscViewerCreate(ts->comm, &viewer);
163: PetscViewerSetType(viewer, typeName);
164: PetscOptionsGetString(ts->prefix, "-ts_view_file", fileName, 1024, &opt);
165: if (opt == PETSC_TRUE) {
166: PetscViewerSetFilename(viewer, fileName);
167: } else {
168: PetscViewerSetFilename(viewer, ts->name);
169: }
170: TSView(ts, viewer);
171: PetscViewerFlush(viewer);
172: PetscViewerDestroy(viewer);
173: } else {
174: TSView(ts, PETSC_NULL);
175: }
176: }
177: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
178: if (opt == PETSC_TRUE) {
179: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
180: PetscViewerDrawGetDraw(viewer, 0, &draw);
181: if (title != PETSC_NULL) {
182: titleStr = title;
183: } else {
184: PetscObjectName((PetscObject) ts); CHKERRQ(ierr) ;
185: titleStr = ts->name;
186: }
187: PetscDrawSetTitle(draw, titleStr);
188: TSView(ts, viewer);
189: PetscViewerFlush(viewer);
190: PetscDrawPause(draw);
191: PetscViewerDestroy(viewer);
192: }
193: return(0);
194: }
196: /*@
197: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
198: set with TSSetRHSJacobian().
200: Collective on TS and Vec
202: Input Parameters:
203: + ts - the SNES context
204: . t - current timestep
205: - x - input vector
207: Output Parameters:
208: + A - Jacobian matrix
209: . B - optional preconditioning matrix
210: - flag - flag indicating matrix structure
212: Notes:
213: Most users should not need to explicitly call this routine, as it
214: is used internally within the nonlinear solvers.
216: See SLESSetOperators() for important information about setting the
217: flag parameter.
219: TSComputeJacobian() is valid only for TS_NONLINEAR
221: Level: developer
223: .keywords: SNES, compute, Jacobian, matrix
225: .seealso: TSSetRHSJacobian(), SLESSetOperators()
226: @*/
227: int TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
228: {
235: if (ts->problem_type != TS_NONLINEAR) {
236: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
237: }
238: if (ts->ops->rhsjacobian) {
239: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
240: *flg = DIFFERENT_NONZERO_PATTERN;
241: PetscStackPush("TS user Jacobian function");
242: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
243: PetscStackPop;
244: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
245: /* make sure user returned a correct Jacobian and preconditioner */
248: } else {
249: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
250: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
251: if (*A != *B) {
252: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
253: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
254: }
255: }
256: return(0);
257: }
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: int TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
266: {
274: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
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->A,&ts->B,&flg,ts->jacP);
284: PetscStackPop;
285: }
286: MatMult(ts->A,x,y);
287: }
289: /* apply user-provided boundary conditions (only needed if these are time dependent) */
290: TSComputeRHSBoundaryConditions(ts,t,y);
291: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
293: return(0);
294: }
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 TSSetRHSMatrix().
319: Level: beginner
321: .keywords: TS, timestep, set, right-hand-side, function
323: .seealso: TSSetRHSMatrix()
324: @*/
325: int TSSetRHSFunction(TS ts,int (*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: }
338: /*@C
339: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
340: Also sets the location to store A.
342: Collective on TS
344: Input Parameters:
345: + ts - the TS context obtained from TSCreate()
346: . A - matrix
347: . B - preconditioner matrix (usually same as A)
348: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
349: if A is not a function of t.
350: - ctx - [optional] user-defined context for private data for the
351: matrix evaluation routine (may be PETSC_NULL)
353: Calling sequence of func:
354: $ func (TS ts,PetscReal t,Mat *A,Mat *B,int *flag,void *ctx);
356: + t - current timestep
357: . A - matrix A, where U_t = A(t) U
358: . B - preconditioner matrix, usually the same as A
359: . flag - flag indicating information about the preconditioner matrix
360: structure (same as flag in SLESSetOperators())
361: - ctx - [optional] user-defined context for matrix evaluation routine
363: Notes:
364: See SLESSetOperators() for important information about setting the flag
365: output parameter in the routine func(). Be sure to read this information!
367: The routine func() takes Mat * as the matrix arguments rather than Mat.
368: This allows the matrix evaluation routine to replace A and/or B with a
369: completely new new matrix structure (not just different matrix elements)
370: when appropriate, for instance, if the nonzero structure is changing
371: throughout the global iterations.
373: Important:
374: The user MUST call either this routine or TSSetRHSFunction().
376: Level: beginner
378: .keywords: TS, timestep, set, right-hand-side, matrix
380: .seealso: TSSetRHSFunction()
381: @*/
382: int TSSetRHSMatrix(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
383: {
390: if (ts->problem_type == TS_NONLINEAR) {
391: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
392: }
394: ts->ops->rhsmatrix = f;
395: ts->jacP = ctx;
396: ts->A = A;
397: ts->B = B;
399: return(0);
400: }
402: /*@C
403: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
404: where U_t = F(U,t), as well as the location to store the matrix.
406: Collective on TS
408: Input Parameters:
409: + ts - the TS context obtained from TSCreate()
410: . A - Jacobian matrix
411: . B - preconditioner matrix (usually same as A)
412: . f - the Jacobian evaluation routine
413: - ctx - [optional] user-defined context for private data for the
414: Jacobian evaluation routine (may be PETSC_NULL)
416: Calling sequence of func:
417: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
419: + t - current timestep
420: . u - input vector
421: . A - matrix A, where U_t = A(t)u
422: . B - preconditioner matrix, usually the same as A
423: . flag - flag indicating information about the preconditioner matrix
424: structure (same as flag in SLESSetOperators())
425: - ctx - [optional] user-defined context for matrix evaluation routine
427: Notes:
428: See SLESSetOperators() for important information about setting the flag
429: output parameter in the routine func(). Be sure to read this information!
431: The routine func() takes Mat * as the matrix arguments rather than Mat.
432: This allows the matrix evaluation routine to replace A and/or B with a
433: completely new new matrix structure (not just different matrix elements)
434: when appropriate, for instance, if the nonzero structure is changing
435: throughout the global iterations.
437: Level: beginner
438:
439: .keywords: TS, timestep, set, right-hand-side, Jacobian
441: .seealso: TSDefaultComputeJacobianColor(),
442: SNESDefaultComputeJacobianColor()
444: @*/
445: int TSSetRHSJacobian(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
446: {
453: if (ts->problem_type != TS_NONLINEAR) {
454: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
455: }
457: ts->ops->rhsjacobian = f;
458: ts->jacP = ctx;
459: ts->A = A;
460: ts->B = B;
461: return(0);
462: }
464: /*
465: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
467: Note: If the user did not provide a function but merely a matrix,
468: this routine applies the matrix.
469: */
470: int TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
471: {
479: if (ts->ops->rhsbc) {
480: PetscStackPush("TS user boundary condition function");
481: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
482: PetscStackPop;
483: return(0);
484: }
486: return(0);
487: }
489: /*@C
490: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
491: boundary conditions for the function F.
493: Collective on TS
495: Input Parameters:
496: + ts - the TS context obtained from TSCreate()
497: . f - routine for evaluating the boundary condition function
498: - ctx - [optional] user-defined context for private data for the
499: function evaluation routine (may be PETSC_NULL)
501: Calling sequence of func:
502: $ func (TS ts,PetscReal t,Vec F,void *ctx);
504: + t - current timestep
505: . F - function vector
506: - ctx - [optional] user-defined function context
508: Level: intermediate
510: .keywords: TS, timestep, set, boundary conditions, function
511: @*/
512: int TSSetRHSBoundaryConditions(TS ts,int (*f)(TS,PetscReal,Vec,void*),void *ctx)
513: {
517: if (ts->problem_type != TS_LINEAR) {
518: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
519: }
520: ts->ops->rhsbc = f;
521: ts->bcP = ctx;
522: return(0);
523: }
525: /*@C
526: TSView - Prints the TS data structure.
528: Collective on TS
530: Input Parameters:
531: + ts - the TS context obtained from TSCreate()
532: - viewer - visualization context
534: Options Database Key:
535: . -ts_view - calls TSView() at end of TSStep()
537: Notes:
538: The available visualization contexts include
539: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
540: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
541: output where only the first processor opens
542: the file. All other processors send their
543: data to the first processor to print.
545: The user can open an alternative visualization context with
546: PetscViewerASCIIOpen() - output to a specified file.
548: Level: beginner
550: .keywords: TS, timestep, view
552: .seealso: PetscViewerASCIIOpen()
553: @*/
554: int TSView(TS ts,PetscViewer viewer)
555: {
556: int ierr;
557: char *type;
558: PetscTruth isascii,isstring;
562: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
566: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
567: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
568: if (isascii) {
569: PetscViewerASCIIPrintf(viewer,"TS Object:n");
570: TSGetType(ts,(TSType *)&type);
571: if (type) {
572: PetscViewerASCIIPrintf(viewer," type: %sn",type);
573: } else {
574: PetscViewerASCIIPrintf(viewer," type: not yet setn");
575: }
576: if (ts->ops->view) {
577: PetscViewerASCIIPushTab(viewer);
578: (*ts->ops->view)(ts,viewer);
579: PetscViewerASCIIPopTab(viewer);
580: }
581: PetscViewerASCIIPrintf(viewer," maximum steps=%dn",ts->max_steps);
582: PetscViewerASCIIPrintf(viewer," maximum time=%gn",ts->max_time);
583: if (ts->problem_type == TS_NONLINEAR) {
584: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%dn",ts->nonlinear_its);
585: }
586: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%dn",ts->linear_its);
587: } else if (isstring) {
588: TSGetType(ts,(TSType *)&type);
589: PetscViewerStringSPrintf(viewer," %-7.7s",type);
590: }
591: PetscViewerASCIIPushTab(viewer);
592: if (ts->sles) {SLESView(ts->sles,viewer);}
593: if (ts->snes) {SNESView(ts->snes,viewer);}
594: PetscViewerASCIIPopTab(viewer);
595: return(0);
596: }
599: /*@C
600: TSSetApplicationContext - Sets an optional user-defined context for
601: the timesteppers.
603: Collective on TS
605: Input Parameters:
606: + ts - the TS context obtained from TSCreate()
607: - usrP - optional user context
609: Level: intermediate
611: .keywords: TS, timestep, set, application, context
613: .seealso: TSGetApplicationContext()
614: @*/
615: int TSSetApplicationContext(TS ts,void *usrP)
616: {
619: ts->user = usrP;
620: return(0);
621: }
623: /*@C
624: TSGetApplicationContext - Gets the user-defined context for the
625: timestepper.
627: Not Collective
629: Input Parameter:
630: . ts - the TS context obtained from TSCreate()
632: Output Parameter:
633: . usrP - user context
635: Level: intermediate
637: .keywords: TS, timestep, get, application, context
639: .seealso: TSSetApplicationContext()
640: @*/
641: int TSGetApplicationContext(TS ts,void **usrP)
642: {
645: *usrP = ts->user;
646: return(0);
647: }
649: /*@
650: TSGetTimeStepNumber - Gets the current number of timesteps.
652: Not Collective
654: Input Parameter:
655: . ts - the TS context obtained from TSCreate()
657: Output Parameter:
658: . iter - number steps so far
660: Level: intermediate
662: .keywords: TS, timestep, get, iteration, number
663: @*/
664: int TSGetTimeStepNumber(TS ts,int* iter)
665: {
668: *iter = ts->steps;
669: return(0);
670: }
672: /*@
673: TSSetInitialTimeStep - Sets the initial timestep to be used,
674: as well as the initial time.
676: Collective on TS
678: Input Parameters:
679: + ts - the TS context obtained from TSCreate()
680: . initial_time - the initial time
681: - time_step - the size of the timestep
683: Level: intermediate
685: .seealso: TSSetTimeStep(), TSGetTimeStep()
687: .keywords: TS, set, initial, timestep
688: @*/
689: int TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
690: {
693: ts->time_step = time_step;
694: ts->initial_time_step = time_step;
695: ts->ptime = initial_time;
696: return(0);
697: }
699: /*@
700: TSSetTimeStep - Allows one to reset the timestep at any time,
701: useful for simple pseudo-timestepping codes.
703: Collective on TS
705: Input Parameters:
706: + ts - the TS context obtained from TSCreate()
707: - time_step - the size of the timestep
709: Level: intermediate
711: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
713: .keywords: TS, set, timestep
714: @*/
715: int TSSetTimeStep(TS ts,PetscReal time_step)
716: {
719: ts->time_step = time_step;
720: return(0);
721: }
723: /*@
724: TSGetTimeStep - Gets the current timestep size.
726: Not Collective
728: Input Parameter:
729: . ts - the TS context obtained from TSCreate()
731: Output Parameter:
732: . dt - the current timestep size
734: Level: intermediate
736: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
738: .keywords: TS, get, timestep
739: @*/
740: int TSGetTimeStep(TS ts,PetscReal* dt)
741: {
744: *dt = ts->time_step;
745: return(0);
746: }
748: /*@C
749: TSGetSolution - Returns the solution at the present timestep. It
750: is valid to call this routine inside the function that you are evaluating
751: in order to move to the new timestep. This vector not changed until
752: the solution at the next timestep has been calculated.
754: Not Collective, but Vec returned is parallel if TS is parallel
756: Input Parameter:
757: . ts - the TS context obtained from TSCreate()
759: Output Parameter:
760: . v - the vector containing the solution
762: Level: intermediate
764: .seealso: TSGetTimeStep()
766: .keywords: TS, timestep, get, solution
767: @*/
768: int TSGetSolution(TS ts,Vec *v)
769: {
772: *v = ts->vec_sol_always;
773: return(0);
774: }
776: /* ----- Routines to initialize and destroy a timestepper ---- */
777: /*@C
778: TSSetProblemType - Sets the type of problem to be solved.
780: Not collective
782: Input Parameters:
783: + ts - The TS
784: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
785: .vb
786: U_t = A U
787: U_t = A(t) U
788: U_t = F(t,U)
789: .ve
791: Level: beginner
793: .keywords: TS, problem type
794: .seealso: TSSetUp(), TSProblemType, TS
795: @*/
796: int TSSetProblemType(TS ts, TSProblemType type) {
799: ts->problem_type = type;
800: return(0);
801: }
803: /*@C
804: TSGetProblemType - Gets the type of problem to be solved.
806: Not collective
808: Input Parameter:
809: . ts - The TS
811: Output Parameter:
812: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
813: .vb
814: U_t = A U
815: U_t = A(t) U
816: U_t = F(t,U)
817: .ve
819: Level: beginner
821: .keywords: TS, problem type
822: .seealso: TSSetUp(), TSProblemType, TS
823: @*/
824: int TSGetProblemType(TS ts, TSProblemType *type) {
828: *type = ts->problem_type;
829: return(0);
830: }
832: /*@
833: TSSetUp - Sets up the internal data structures for the later use
834: of a timestepper.
836: Collective on TS
838: Input Parameter:
839: . ts - the TS context obtained from TSCreate()
841: Notes:
842: For basic use of the TS solvers the user need not explicitly call
843: TSSetUp(), since these actions will automatically occur during
844: the call to TSStep(). However, if one wishes to control this
845: phase separately, TSSetUp() should be called after TSCreate()
846: and optional routines of the form TSSetXXX(), but before TSStep().
848: Level: advanced
850: .keywords: TS, timestep, setup
852: .seealso: TSCreate(), TSStep(), TSDestroy()
853: @*/
854: int TSSetUp(TS ts)
855: {
860: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
861: if (!ts->type_name) {
862: TSSetType(ts,TS_EULER);
863: }
864: (*ts->ops->setup)(ts);
865: ts->setupcalled = 1;
866: return(0);
867: }
869: /*@C
870: TSDestroy - Destroys the timestepper context that was created
871: with TSCreate().
873: Collective on TS
875: Input Parameter:
876: . ts - the TS context obtained from TSCreate()
878: Level: beginner
880: .keywords: TS, timestepper, destroy
882: .seealso: TSCreate(), TSSetUp(), TSSolve()
883: @*/
884: int TSDestroy(TS ts)
885: {
886: int ierr,i;
890: if (--ts->refct > 0) return(0);
892: /* if memory was published with AMS then destroy it */
893: PetscObjectDepublish(ts);
895: if (ts->sles) {SLESDestroy(ts->sles);}
896: if (ts->snes) {SNESDestroy(ts->snes);}
897: (*(ts)->ops->destroy)(ts);
898: for (i=0; i<ts->numbermonitors; i++) {
899: if (ts->mdestroy[i]) {
900: (*ts->mdestroy[i])(ts->monitorcontext[i]);
901: }
902: }
903: PetscLogObjectDestroy((PetscObject)ts);
904: PetscHeaderDestroy((PetscObject)ts);
905: return(0);
906: }
908: /*@C
909: TSGetSNES - Returns the SNES (nonlinear solver) associated with
910: a TS (timestepper) context. Valid only for nonlinear problems.
912: Not Collective, but SNES is parallel if TS is parallel
914: Input Parameter:
915: . ts - the TS context obtained from TSCreate()
917: Output Parameter:
918: . snes - the nonlinear solver context
920: Notes:
921: The user can then directly manipulate the SNES context to set various
922: options, etc. Likewise, the user can then extract and manipulate the
923: SLES, KSP, and PC contexts as well.
925: TSGetSNES() does not work for integrators that do not use SNES; in
926: this case TSGetSNES() returns PETSC_NULL in snes.
928: Level: beginner
930: .keywords: timestep, get, SNES
931: @*/
932: int TSGetSNES(TS ts,SNES *snes)
933: {
936: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetSLES()");
937: *snes = ts->snes;
938: return(0);
939: }
941: /*@C
942: TSGetSLES - Returns the SLES (linear solver) associated with
943: a TS (timestepper) context.
945: Not Collective, but SLES is parallel if TS is parallel
947: Input Parameter:
948: . ts - the TS context obtained from TSCreate()
950: Output Parameter:
951: . sles - the nonlinear solver context
953: Notes:
954: The user can then directly manipulate the SLES context to set various
955: options, etc. Likewise, the user can then extract and manipulate the
956: KSP and PC contexts as well.
958: TSGetSLES() does not work for integrators that do not use SLES;
959: in this case TSGetSLES() returns PETSC_NULL in sles.
961: Level: beginner
963: .keywords: timestep, get, SLES
964: @*/
965: int TSGetSLES(TS ts,SLES *sles)
966: {
969: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
970: *sles = ts->sles;
971: return(0);
972: }
974: /* ----------- Routines to set solver parameters ---------- */
976: /*@
977: TSSetDuration - Sets the maximum number of timesteps to use and
978: maximum time for iteration.
980: Collective on TS
982: Input Parameters:
983: + ts - the TS context obtained from TSCreate()
984: . maxsteps - maximum number of iterations to use
985: - maxtime - final time to iterate to
987: Options Database Keys:
988: . -ts_max_steps <maxsteps> - Sets maxsteps
989: . -ts_max_time <maxtime> - Sets maxtime
991: Notes:
992: The default maximum number of iterations is 5000. Default time is 5.0
994: Level: intermediate
996: .keywords: TS, timestep, set, maximum, iterations
997: @*/
998: int TSSetDuration(TS ts,int maxsteps,PetscReal maxtime)
999: {
1002: ts->max_steps = maxsteps;
1003: ts->max_time = maxtime;
1004: return(0);
1005: }
1007: /*@
1008: TSSetSolution - Sets the initial solution vector
1009: for use by the TS routines.
1011: Collective on TS and Vec
1013: Input Parameters:
1014: + ts - the TS context obtained from TSCreate()
1015: - x - the solution vector
1017: Level: beginner
1019: .keywords: TS, timestep, set, solution, initial conditions
1020: @*/
1021: int TSSetSolution(TS ts,Vec x)
1022: {
1025: ts->vec_sol = ts->vec_sol_always = x;
1026: return(0);
1027: }
1029: /*@
1030: TSSetRhsBC - Sets the function which applies boundary conditions
1031: to the Rhs of each system.
1033: Collective on TS
1035: Input Parameters:
1036: + ts - The TS context obtained from TSCreate()
1037: - func - The function
1039: Calling sequence of func:
1040: . func (TS ts, Vec rhs, void *ctx);
1042: + rhs - The current rhs vector
1043: - ctx - The user-context
1045: Level: intermediate
1047: .keywords: TS, Rhs, boundary conditions
1048: @*/
1049: int TSSetRhsBC(TS ts, int (*func)(TS, Vec, void *))
1050: {
1053: ts->ops->applyrhsbc = func;
1054: return(0);
1055: }
1057: /*@
1058: TSDefaultRhsBC - The default boundary condition function which does nothing.
1060: Collective on TS
1062: Input Parameters:
1063: + ts - The TS context obtained from TSCreate()
1064: . rhs - The Rhs
1065: - ctx - The user-context
1067: Level: developer
1069: .keywords: TS, Rhs, boundary conditions
1070: @*/
1071: int TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1072: {
1074: return(0);
1075: }
1077: /*@
1078: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1079: to the system matrix and preconditioner of each system.
1081: Collective on TS
1083: Input Parameters:
1084: + ts - The TS context obtained from TSCreate()
1085: - func - The function
1087: Calling sequence of func:
1088: . func (TS ts, Mat A, Mat B, void *ctx);
1090: + A - The current system matrix
1091: . B - The current preconditioner
1092: - ctx - The user-context
1094: Level: intermediate
1096: .keywords: TS, System matrix, boundary conditions
1097: @*/
1098: int TSSetSystemMatrixBC(TS ts, int (*func)(TS, Mat, Mat, void *))
1099: {
1102: ts->ops->applymatrixbc = func;
1103: return(0);
1104: }
1106: /*@
1107: TSDefaultSystemMatrixBC - The default boundary condition function which
1108: does nothing.
1110: Collective on TS
1112: Input Parameters:
1113: + ts - The TS context obtained from TSCreate()
1114: . A - The system matrix
1115: . B - The preconditioner
1116: - ctx - The user-context
1118: Level: developer
1120: .keywords: TS, System matrix, boundary conditions
1121: @*/
1122: int TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1123: {
1125: return(0);
1126: }
1128: /*@
1129: TSSetSolutionBC - Sets the function which applies boundary conditions
1130: to the solution of each system. This is necessary in nonlinear systems
1131: which time dependent boundary conditions.
1133: Collective on TS
1135: Input Parameters:
1136: + ts - The TS context obtained from TSCreate()
1137: - func - The function
1139: Calling sequence of func:
1140: . func (TS ts, Vec rsol, void *ctx);
1142: + sol - The current solution vector
1143: - ctx - The user-context
1145: Level: intermediate
1147: .keywords: TS, solution, boundary conditions
1148: @*/
1149: int TSSetSolutionBC(TS ts, int (*func)(TS, Vec, void *))
1150: {
1153: ts->ops->applysolbc = func;
1154: return(0);
1155: }
1157: /*@
1158: TSDefaultSolutionBC - The default boundary condition function which
1159: does nothing.
1161: Collective on TS
1163: Input Parameters:
1164: + ts - The TS context obtained from TSCreate()
1165: . sol - The solution
1166: - ctx - The user-context
1168: Level: developer
1170: .keywords: TS, solution, boundary conditions
1171: @*/
1172: int TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1173: {
1175: return(0);
1176: }
1178: /*@
1179: TSSetPreStep - Sets the general-purpose function
1180: called once at the beginning of time stepping.
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);
1191: Level: intermediate
1193: .keywords: TS, timestep
1194: @*/
1195: int TSSetPreStep(TS ts, int (*func)(TS))
1196: {
1199: ts->ops->prestep = func;
1200: return(0);
1201: }
1203: /*@
1204: TSDefaultPreStep - The default pre-stepping function which does nothing.
1206: Collective on TS
1208: Input Parameters:
1209: . ts - The TS context obtained from TSCreate()
1211: Level: developer
1213: .keywords: TS, timestep
1214: @*/
1215: int TSDefaultPreStep(TS ts)
1216: {
1218: return(0);
1219: }
1221: /*@
1222: TSSetUpdate - Sets the general-purpose update function called
1223: at the beginning of every time step. This function can change
1224: the time step.
1226: Collective on TS
1228: Input Parameters:
1229: + ts - The TS context obtained from TSCreate()
1230: - func - The function
1232: Calling sequence of func:
1233: . func (TS ts, double t, double *dt);
1235: + t - The current time
1236: - dt - The current time step
1238: Level: intermediate
1240: .keywords: TS, update, timestep
1241: @*/
1242: int TSSetUpdate(TS ts, int (*func)(TS, PetscReal, PetscReal *))
1243: {
1246: ts->ops->update = func;
1247: return(0);
1248: }
1250: /*@
1251: TSDefaultUpdate - The default update function which does nothing.
1253: Collective on TS
1255: Input Parameters:
1256: + ts - The TS context obtained from TSCreate()
1257: - t - The current time
1259: Output Parameters:
1260: . dt - The current time step
1262: Level: developer
1264: .keywords: TS, update, timestep
1265: @*/
1266: int TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1267: {
1269: return(0);
1270: }
1272: /*@
1273: TSSetPostStep - Sets the general-purpose function
1274: called once at the end of time stepping.
1276: Collective on TS
1278: Input Parameters:
1279: + ts - The TS context obtained from TSCreate()
1280: - func - The function
1282: Calling sequence of func:
1283: . func (TS ts);
1285: Level: intermediate
1287: .keywords: TS, timestep
1288: @*/
1289: int TSSetPostStep(TS ts, int (*func)(TS))
1290: {
1293: ts->ops->poststep = func;
1294: return(0);
1295: }
1297: /*@
1298: TSDefaultPostStep - The default post-stepping function which does nothing.
1300: Collective on TS
1302: Input Parameters:
1303: . ts - The TS context obtained from TSCreate()
1305: Level: developer
1307: .keywords: TS, timestep
1308: @*/
1309: int TSDefaultPostStep(TS ts)
1310: {
1312: return(0);
1313: }
1315: /* ------------ Routines to set performance monitoring options ----------- */
1317: /*@C
1318: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1319: timestep to display the iteration's progress.
1321: Collective on TS
1323: Input Parameters:
1324: + ts - the TS context obtained from TSCreate()
1325: . func - monitoring routine
1326: . mctx - [optional] user-defined context for private data for the
1327: monitor routine (use PETSC_NULL if no context is desired)
1328: - monitordestroy - [optional] routine that frees monitor context
1329: (may be PETSC_NULL)
1331: Calling sequence of func:
1332: $ int func(TS ts,int steps,PetscReal time,Vec x,void *mctx)
1334: + ts - the TS context
1335: . steps - iteration number
1336: . time - current time
1337: . x - current iterate
1338: - mctx - [optional] monitoring context
1340: Notes:
1341: This routine adds an additional monitor to the list of monitors that
1342: already has been loaded.
1344: Level: intermediate
1346: .keywords: TS, timestep, set, monitor
1348: .seealso: TSDefaultMonitor(), TSClearMonitor()
1349: @*/
1350: int TSSetMonitor(TS ts,int (*monitor)(TS,int,PetscReal,Vec,void*),void *mctx,int (*mdestroy)(void*))
1351: {
1354: if (ts->numbermonitors >= MAXTSMONITORS) {
1355: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1356: }
1357: ts->monitor[ts->numbermonitors] = monitor;
1358: ts->mdestroy[ts->numbermonitors] = mdestroy;
1359: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1360: return(0);
1361: }
1363: /*@C
1364: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1366: Collective on TS
1368: Input Parameters:
1369: . ts - the TS context obtained from TSCreate()
1371: Notes:
1372: There is no way to remove a single, specific monitor.
1374: Level: intermediate
1376: .keywords: TS, timestep, set, monitor
1378: .seealso: TSDefaultMonitor(), TSSetMonitor()
1379: @*/
1380: int TSClearMonitor(TS ts)
1381: {
1384: ts->numbermonitors = 0;
1385: return(0);
1386: }
1388: int TSDefaultMonitor(TS ts,int step,PetscReal ptime,Vec v,void *ctx)
1389: {
1393: PetscPrintf(ts->comm,"timestep %d dt %g time %gn",step,ts->time_step,ptime);
1394: return(0);
1395: }
1397: /*@
1398: TSStep - Steps the requested number of timesteps.
1400: Collective on TS
1402: Input Parameter:
1403: . ts - the TS context obtained from TSCreate()
1405: Output Parameters:
1406: + steps - number of iterations until termination
1407: - ptime - time until termination
1409: Level: beginner
1411: .keywords: TS, timestep, solve
1413: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1414: @*/
1415: int TSStep(TS ts,int *steps,PetscReal *ptime)
1416: {
1417: PetscTruth opt;
1418: int ierr;
1422: if (!ts->setupcalled) {
1423: TSSetUp(ts);
1424: }
1426: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1427: (*ts->ops->prestep)(ts);
1428: (*ts->ops->step)(ts, steps, ptime);
1429: (*ts->ops->poststep)(ts);
1430: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1432: PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
1433: if ((opt == PETSC_TRUE) && !PetscPreLoadingOn) {
1434: TSView(ts, PETSC_VIEWER_STDOUT_WORLD);
1435: }
1436: return(0);
1437: }
1439: /*
1440: Runs the user provided monitor routines, if they exists.
1441: */
1442: int TSMonitor(TS ts,int step,PetscReal ptime,Vec x)
1443: {
1444: int i,ierr,n = ts->numbermonitors;
1447: for (i=0; i<n; i++) {
1448: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1449: }
1450: return(0);
1451: }
1453: /* ------------------------------------------------------------------------*/
1455: /*@C
1456: TSLGMonitorCreate - Creates a line graph context for use with
1457: TS to monitor convergence of preconditioned residual norms.
1459: Collective on TS
1461: Input Parameters:
1462: + host - the X display to open, or null for the local machine
1463: . label - the title to put in the title bar
1464: . x, y - the screen coordinates of the upper left coordinate of the window
1465: - m, n - the screen width and height in pixels
1467: Output Parameter:
1468: . draw - the drawing context
1470: Options Database Key:
1471: . -ts_xmonitor - automatically sets line graph monitor
1473: Notes:
1474: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1476: Level: intermediate
1478: .keywords: TS, monitor, line graph, residual, seealso
1480: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1482: @*/
1483: int TSLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1484: {
1485: PetscDraw win;
1486: int ierr;
1489: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1490: PetscDrawSetType(win,PETSC_DRAW_X);
1491: PetscDrawLGCreate(win,1,draw);
1492: PetscDrawLGIndicateDataPoints(*draw);
1494: PetscLogObjectParent(*draw,win);
1495: return(0);
1496: }
1498: int TSLGMonitor(TS ts,int n,PetscReal ptime,Vec v,void *monctx)
1499: {
1500: PetscDrawLG lg = (PetscDrawLG) monctx;
1501: PetscReal x,y = ptime;
1502: int ierr;
1505: if (!monctx) {
1506: MPI_Comm comm;
1507: PetscViewer viewer;
1509: ierr = PetscObjectGetComm((PetscObject)ts,&comm);
1510: viewer = PETSC_VIEWER_DRAW_(comm);
1511: ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);
1512: }
1514: if (!n) {PetscDrawLGReset(lg);}
1515: x = (PetscReal)n;
1516: PetscDrawLGAddPoint(lg,&x,&y);
1517: if (n < 20 || (n % 5)) {
1518: PetscDrawLGDraw(lg);
1519: }
1520: return(0);
1521: }
1523: /*@C
1524: TSLGMonitorDestroy - Destroys a line graph context that was created
1525: with TSLGMonitorCreate().
1527: Collective on PetscDrawLG
1529: Input Parameter:
1530: . draw - the drawing context
1532: Level: intermediate
1534: .keywords: TS, monitor, line graph, destroy
1536: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1537: @*/
1538: int TSLGMonitorDestroy(PetscDrawLG drawlg)
1539: {
1540: PetscDraw draw;
1541: int ierr;
1544: PetscDrawLGGetDraw(drawlg,&draw);
1545: PetscDrawDestroy(draw);
1546: PetscDrawLGDestroy(drawlg);
1547: return(0);
1548: }
1550: /*@
1551: TSGetTime - Gets the current time.
1553: Not Collective
1555: Input Parameter:
1556: . ts - the TS context obtained from TSCreate()
1558: Output Parameter:
1559: . t - the current time
1561: Contributed by: Matthew Knepley
1563: Level: beginner
1565: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1567: .keywords: TS, get, time
1568: @*/
1569: int TSGetTime(TS ts,PetscReal* t)
1570: {
1573: *t = ts->ptime;
1574: return(0);
1575: }
1577: /*@C
1578: TSSetOptionsPrefix - Sets the prefix used for searching for all
1579: TS options in the database.
1581: Collective on TS
1583: Input Parameter:
1584: + ts - The TS context
1585: - prefix - The prefix to prepend to all option names
1587: Notes:
1588: A hyphen (-) must NOT be given at the beginning of the prefix name.
1589: The first character of all runtime options is AUTOMATICALLY the
1590: hyphen.
1592: Contributed by: Matthew Knepley
1594: Level: advanced
1596: .keywords: TS, set, options, prefix, database
1598: .seealso: TSSetFromOptions()
1600: @*/
1601: int TSSetOptionsPrefix(TS ts,char *prefix)
1602: {
1607: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1608: switch(ts->problem_type) {
1609: case TS_NONLINEAR:
1610: SNESSetOptionsPrefix(ts->snes,prefix);
1611: break;
1612: case TS_LINEAR:
1613: SLESSetOptionsPrefix(ts->sles,prefix);
1614: break;
1615: }
1616: return(0);
1617: }
1620: /*@C
1621: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1622: TS options in the database.
1624: Collective on TS
1626: Input Parameter:
1627: + ts - The TS context
1628: - prefix - The prefix to prepend to all option names
1630: Notes:
1631: A hyphen (-) must NOT be given at the beginning of the prefix name.
1632: The first character of all runtime options is AUTOMATICALLY the
1633: hyphen.
1635: Contributed by: Matthew Knepley
1637: Level: advanced
1639: .keywords: TS, append, options, prefix, database
1641: .seealso: TSGetOptionsPrefix()
1643: @*/
1644: int TSAppendOptionsPrefix(TS ts,char *prefix)
1645: {
1650: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1651: switch(ts->problem_type) {
1652: case TS_NONLINEAR:
1653: SNESAppendOptionsPrefix(ts->snes,prefix);
1654: break;
1655: case TS_LINEAR:
1656: SLESAppendOptionsPrefix(ts->sles,prefix);
1657: break;
1658: }
1659: return(0);
1660: }
1662: /*@C
1663: TSGetOptionsPrefix - Sets the prefix used for searching for all
1664: TS options in the database.
1666: Not Collective
1668: Input Parameter:
1669: . ts - The TS context
1671: Output Parameter:
1672: . prefix - A pointer to the prefix string used
1674: Contributed by: Matthew Knepley
1676: Notes: On the fortran side, the user should pass in a string 'prifix' of
1677: sufficient length to hold the prefix.
1679: Level: intermediate
1681: .keywords: TS, get, options, prefix, database
1683: .seealso: TSAppendOptionsPrefix()
1684: @*/
1685: int TSGetOptionsPrefix(TS ts,char **prefix)
1686: {
1691: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1692: return(0);
1693: }
1695: /*@C
1696: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1698: Not Collective, but parallel objects are returned if TS is parallel
1700: Input Parameter:
1701: . ts - The TS context obtained from TSCreate()
1703: Output Parameters:
1704: + A - The matrix A, where U_t = A(t) U
1705: . M - The preconditioner matrix, usually the same as A
1706: - ctx - User-defined context for matrix evaluation routine
1708: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1710: Contributed by: Matthew Knepley
1712: Level: intermediate
1714: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1716: .keywords: TS, timestep, get, matrix
1718: @*/
1719: int TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1720: {
1723: if (A) *A = ts->A;
1724: if (M) *M = ts->B;
1725: if (ctx) *ctx = ts->jacP;
1726: return(0);
1727: }
1729: /*@C
1730: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1732: Not Collective, but parallel objects are returned if TS is parallel
1734: Input Parameter:
1735: . ts - The TS context obtained from TSCreate()
1737: Output Parameters:
1738: + J - The Jacobian J of F, where U_t = F(U,t)
1739: . M - The preconditioner matrix, usually the same as J
1740: - ctx - User-defined context for Jacobian evaluation routine
1742: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1744: Contributed by: Matthew Knepley
1746: Level: intermediate
1748: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1750: .keywords: TS, timestep, get, matrix, Jacobian
1751: @*/
1752: int TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1753: {
1757: TSGetRHSMatrix(ts,J,M,ctx);
1758: return(0);
1759: }
1761: /*@C
1762: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1763: VecView() for the solution at each timestep
1765: Collective on TS
1767: Input Parameters:
1768: + ts - the TS context
1769: . step - current time-step
1770: . ptime - current time
1771: - dummy - either a viewer or PETSC_NULL
1773: Level: intermediate
1775: .keywords: TS, vector, monitor, view
1777: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1778: @*/
1779: int TSVecViewMonitor(TS ts,int step,PetscReal ptime,Vec x,void *dummy)
1780: {
1781: int ierr;
1782: PetscViewer viewer = (PetscViewer) dummy;
1785: if (!viewer) {
1786: MPI_Comm comm;
1787: ierr = PetscObjectGetComm((PetscObject)ts,&comm);
1788: viewer = PETSC_VIEWER_DRAW_(comm);
1789: }
1790: VecView(x,viewer);
1791: return(0);
1792: }