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;
10: /*
11: TSSetTypeFromOptions - Sets the type of ts from user options.
13: Collective on TS
15: Input Parameter:
16: . ts - The ts
18: Level: intermediate
20: .keywords: TS, set, options, database, type
21: .seealso: TSSetFromOptions(), TSSetType()
22: */
23: static int TSSetTypeFromOptions(TS ts)
24: {
25: PetscTruth opt;
26: char *defaultType;
27: char typeName[256];
28: int ierr;
31: if (ts->type_name != PETSC_NULL) {
32: defaultType = ts->type_name;
33: } else {
34: defaultType = TS_EULER;
35: }
37: if (!TSRegisterAllCalled) {
38: TSRegisterAll(PETSC_NULL);
39: }
40: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
41: if (opt == PETSC_TRUE) {
42: TSSetType(ts, typeName);
43: } else {
44: TSSetType(ts, defaultType);
45: }
46: return(0);
47: }
51: /*@
52: TSSetFromOptions - Sets various TS parameters from user options.
54: Collective on TS
56: Input Parameter:
57: . ts - the TS context obtained from TSCreate()
59: Options Database Keys:
60: + -ts_type <type> - TS_EULER, TS_BEULER, TS_PVODE, TS_PSEUDO, TS_CRANK_NICHOLSON
61: . -ts_max_steps maxsteps - maximum number of time-steps to take
62: . -ts_max_time time - maximum time to compute to
63: . -ts_dt dt - initial time step
64: . -ts_monitor - print information at each timestep
65: - -ts_xmonitor - plot information at each timestep
67: Level: beginner
69: .keywords: TS, timestep, set, options, database
71: .seealso: TSGetType
72: @*/
73: int TSSetFromOptions(TS ts)
74: {
75: PetscReal dt;
76: PetscTruth opt;
77: int ierr;
81: PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");
83: /* Handle generic TS options */
84: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
85: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
86: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
87: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
88: if (opt == PETSC_TRUE) {
89: ts->initial_time_step = ts->time_step = dt;
90: }
92: /* Monitor options */
93: PetscOptionsName("-ts_monitor","Monitor timestep size","TSDefaultMonitor",&opt);
94: if (opt == PETSC_TRUE) {
95: TSSetMonitor(ts,TSDefaultMonitor,PETSC_NULL,PETSC_NULL);
96: }
97: PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
98: if (opt == PETSC_TRUE) {
99: TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
100: }
101: PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
102: if (opt == PETSC_TRUE) {
103: TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
104: }
106: /* Handle TS type options */
107: TSSetTypeFromOptions(ts);
109: /* Handle specific TS options */
110: if (ts->ops->setfromoptions != PETSC_NULL) {
111: (*ts->ops->setfromoptions)(ts);
112: }
113: PetscOptionsEnd();
115: /* Handle subobject options */
116: switch(ts->problem_type) {
117: /* Should check for implicit/explicit */
118: case TS_LINEAR:
119: if (ts->sles != PETSC_NULL) {
120: SLESSetFromOptions(ts->sles);
121: }
122: break;
123: case TS_NONLINEAR:
124: if (ts->snes != PETSC_NULL) {
125: SNESSetFromOptions(ts->snes);
126: }
127: break;
128: default:
129: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", ts->problem_type);
130: }
132: return(0);
133: }
135: #undef __FUNCT__
137: /*@
138: TSViewFromOptions - This function visualizes the ts based upon user options.
140: Collective on TS
142: Input Parameter:
143: . ts - The ts
145: Level: intermediate
147: .keywords: TS, view, options, database
148: .seealso: TSSetFromOptions(), TSView()
149: @*/
150: int TSViewFromOptions(TS ts,const char title[])
151: {
152: PetscViewer viewer;
153: PetscDraw draw;
154: PetscTruth opt;
155: char typeName[1024];
156: char fileName[PETSC_MAX_PATH_LEN];
157: int len;
158: int ierr;
161: PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
162: if (opt == PETSC_TRUE) {
163: PetscOptionsGetString(ts->prefix, "-ts_view", typeName, 1024, &opt);
164: PetscStrlen(typeName, &len);
165: if (len > 0) {
166: PetscViewerCreate(ts->comm, &viewer);
167: PetscViewerSetType(viewer, typeName);
168: PetscOptionsGetString(ts->prefix, "-ts_view_file", fileName, 1024, &opt);
169: if (opt == PETSC_TRUE) {
170: PetscViewerSetFilename(viewer, fileName);
171: } else {
172: PetscViewerSetFilename(viewer, ts->name);
173: }
174: TSView(ts, viewer);
175: PetscViewerFlush(viewer);
176: PetscViewerDestroy(viewer);
177: } else {
178: TSView(ts, PETSC_NULL);
179: }
180: }
181: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
182: if (opt == PETSC_TRUE) {
183: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
184: PetscViewerDrawGetDraw(viewer, 0, &draw);
185: if (title != PETSC_NULL) {
186: PetscDrawSetTitle(draw, (char *)title);
187: } else {
188: PetscObjectName((PetscObject) ts); CHKERRQ(ierr) ;
189: PetscDrawSetTitle(draw, ts->name);
190: }
191: TSView(ts, viewer);
192: PetscViewerFlush(viewer);
193: PetscDrawPause(draw);
194: PetscViewerDestroy(viewer);
195: }
196: return(0);
197: }
201: /*@
202: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
203: set with TSSetRHSJacobian().
205: Collective on TS and Vec
207: Input Parameters:
208: + ts - the SNES context
209: . t - current timestep
210: - x - input vector
212: Output Parameters:
213: + A - Jacobian matrix
214: . B - optional preconditioning matrix
215: - flag - flag indicating matrix structure
217: Notes:
218: Most users should not need to explicitly call this routine, as it
219: is used internally within the nonlinear solvers.
221: See SLESSetOperators() for important information about setting the
222: flag parameter.
224: TSComputeJacobian() is valid only for TS_NONLINEAR
226: Level: developer
228: .keywords: SNES, compute, Jacobian, matrix
230: .seealso: TSSetRHSJacobian(), SLESSetOperators()
231: @*/
232: int TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
233: {
240: if (ts->problem_type != TS_NONLINEAR) {
241: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
242: }
243: if (ts->ops->rhsjacobian) {
244: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
245: *flg = DIFFERENT_NONZERO_PATTERN;
246: PetscStackPush("TS user Jacobian function");
247: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
248: PetscStackPop;
249: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
250: /* make sure user returned a correct Jacobian and preconditioner */
253: } else {
254: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
255: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
256: if (*A != *B) {
257: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
258: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
259: }
260: }
261: return(0);
262: }
266: /*
267: TSComputeRHSFunction - Evaluates the right-hand-side function.
269: Note: If the user did not provide a function but merely a matrix,
270: this routine applies the matrix.
271: */
272: int TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
273: {
281: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
282: if (ts->ops->rhsfunction) {
283: PetscStackPush("TS user right-hand-side function");
284: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
285: PetscStackPop;
286: } else {
287: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
288: MatStructure flg;
289: PetscStackPush("TS user right-hand-side matrix function");
290: (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
291: PetscStackPop;
292: }
293: MatMult(ts->A,x,y);
294: }
296: /* apply user-provided boundary conditions (only needed if these are time dependent) */
297: TSComputeRHSBoundaryConditions(ts,t,y);
298: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
300: return(0);
301: }
305: /*@C
306: TSSetRHSFunction - Sets the routine for evaluating the function,
307: F(t,u), where U_t = F(t,u).
309: Collective on TS
311: Input Parameters:
312: + ts - the TS context obtained from TSCreate()
313: . f - routine for evaluating the right-hand-side function
314: - ctx - [optional] user-defined context for private data for the
315: function evaluation routine (may be PETSC_NULL)
317: Calling sequence of func:
318: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
320: + t - current timestep
321: . u - input vector
322: . F - function vector
323: - ctx - [optional] user-defined function context
325: Important:
326: The user MUST call either this routine or TSSetRHSMatrix().
328: Level: beginner
330: .keywords: TS, timestep, set, right-hand-side, function
332: .seealso: TSSetRHSMatrix()
333: @*/
334: int TSSetRHSFunction(TS ts,int (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
335: {
339: if (ts->problem_type == TS_LINEAR) {
340: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
341: }
342: ts->ops->rhsfunction = f;
343: ts->funP = ctx;
344: return(0);
345: }
349: /*@C
350: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
351: Also sets the location to store A.
353: Collective on TS
355: Input Parameters:
356: + ts - the TS context obtained from TSCreate()
357: . A - matrix
358: . B - preconditioner matrix (usually same as A)
359: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
360: if A is not a function of t.
361: - ctx - [optional] user-defined context for private data for the
362: matrix evaluation routine (may be PETSC_NULL)
364: Calling sequence of func:
365: $ func (TS ts,PetscReal t,Mat *A,Mat *B,int *flag,void *ctx);
367: + t - current timestep
368: . A - matrix A, where U_t = A(t) U
369: . B - preconditioner matrix, usually the same as A
370: . flag - flag indicating information about the preconditioner matrix
371: structure (same as flag in SLESSetOperators())
372: - ctx - [optional] user-defined context for matrix evaluation routine
374: Notes:
375: See SLESSetOperators() for important information about setting the flag
376: output parameter in the routine func(). Be sure to read this information!
378: The routine func() takes Mat * as the matrix arguments rather than Mat.
379: This allows the matrix evaluation routine to replace A and/or B with a
380: completely new new matrix structure (not just different matrix elements)
381: when appropriate, for instance, if the nonzero structure is changing
382: throughout the global iterations.
384: Important:
385: The user MUST call either this routine or TSSetRHSFunction().
387: Level: beginner
389: .keywords: TS, timestep, set, right-hand-side, matrix
391: .seealso: TSSetRHSFunction()
392: @*/
393: int TSSetRHSMatrix(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
394: {
401: if (ts->problem_type == TS_NONLINEAR) {
402: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
403: }
405: ts->ops->rhsmatrix = f;
406: ts->jacP = ctx;
407: ts->A = A;
408: ts->B = B;
410: return(0);
411: }
415: /*@C
416: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
417: where U_t = F(U,t), as well as the location to store the matrix.
419: Collective on TS
421: Input Parameters:
422: + ts - the TS context obtained from TSCreate()
423: . A - Jacobian matrix
424: . B - preconditioner matrix (usually same as A)
425: . f - the Jacobian evaluation routine
426: - ctx - [optional] user-defined context for private data for the
427: Jacobian evaluation routine (may be PETSC_NULL)
429: Calling sequence of func:
430: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
432: + t - current timestep
433: . u - input vector
434: . A - matrix A, where U_t = A(t)u
435: . B - preconditioner matrix, usually the same as A
436: . flag - flag indicating information about the preconditioner matrix
437: structure (same as flag in SLESSetOperators())
438: - ctx - [optional] user-defined context for matrix evaluation routine
440: Notes:
441: See SLESSetOperators() for important information about setting the flag
442: output parameter in the routine func(). Be sure to read this information!
444: The routine func() takes Mat * as the matrix arguments rather than Mat.
445: This allows the matrix evaluation routine to replace A and/or B with a
446: completely new new matrix structure (not just different matrix elements)
447: when appropriate, for instance, if the nonzero structure is changing
448: throughout the global iterations.
450: Level: beginner
451:
452: .keywords: TS, timestep, set, right-hand-side, Jacobian
454: .seealso: TSDefaultComputeJacobianColor(),
455: SNESDefaultComputeJacobianColor()
457: @*/
458: int TSSetRHSJacobian(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
459: {
466: if (ts->problem_type != TS_NONLINEAR) {
467: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
468: }
470: ts->ops->rhsjacobian = f;
471: ts->jacP = ctx;
472: ts->A = A;
473: ts->B = B;
474: return(0);
475: }
479: /*
480: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
482: Note: If the user did not provide a function but merely a matrix,
483: this routine applies the matrix.
484: */
485: int TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
486: {
494: if (ts->ops->rhsbc) {
495: PetscStackPush("TS user boundary condition function");
496: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
497: PetscStackPop;
498: return(0);
499: }
501: return(0);
502: }
506: /*@C
507: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
508: boundary conditions for the function F.
510: Collective on TS
512: Input Parameters:
513: + ts - the TS context obtained from TSCreate()
514: . f - routine for evaluating the boundary condition function
515: - ctx - [optional] user-defined context for private data for the
516: function evaluation routine (may be PETSC_NULL)
518: Calling sequence of func:
519: $ func (TS ts,PetscReal t,Vec F,void *ctx);
521: + t - current timestep
522: . F - function vector
523: - ctx - [optional] user-defined function context
525: Level: intermediate
527: .keywords: TS, timestep, set, boundary conditions, function
528: @*/
529: int TSSetRHSBoundaryConditions(TS ts,int (*f)(TS,PetscReal,Vec,void*),void *ctx)
530: {
534: if (ts->problem_type != TS_LINEAR) {
535: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
536: }
537: ts->ops->rhsbc = f;
538: ts->bcP = ctx;
539: return(0);
540: }
544: /*@C
545: TSView - Prints the TS data structure.
547: Collective on TS
549: Input Parameters:
550: + ts - the TS context obtained from TSCreate()
551: - viewer - visualization context
553: Options Database Key:
554: . -ts_view - calls TSView() at end of TSStep()
556: Notes:
557: The available visualization contexts include
558: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
559: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
560: output where only the first processor opens
561: the file. All other processors send their
562: data to the first processor to print.
564: The user can open an alternative visualization context with
565: PetscViewerASCIIOpen() - output to a specified file.
567: Level: beginner
569: .keywords: TS, timestep, view
571: .seealso: PetscViewerASCIIOpen()
572: @*/
573: int TSView(TS ts,PetscViewer viewer)
574: {
575: int ierr;
576: char *type;
577: PetscTruth isascii,isstring;
581: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
585: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
586: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
587: if (isascii) {
588: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
589: TSGetType(ts,(TSType *)&type);
590: if (type) {
591: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
592: } else {
593: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
594: }
595: if (ts->ops->view) {
596: PetscViewerASCIIPushTab(viewer);
597: (*ts->ops->view)(ts,viewer);
598: PetscViewerASCIIPopTab(viewer);
599: }
600: PetscViewerASCIIPrintf(viewer," maximum steps=%d\n",ts->max_steps);
601: PetscViewerASCIIPrintf(viewer," maximum time=%g\n",ts->max_time);
602: if (ts->problem_type == TS_NONLINEAR) {
603: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%d\n",ts->nonlinear_its);
604: }
605: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%d\n",ts->linear_its);
606: } else if (isstring) {
607: TSGetType(ts,(TSType *)&type);
608: PetscViewerStringSPrintf(viewer," %-7.7s",type);
609: }
610: PetscViewerASCIIPushTab(viewer);
611: if (ts->sles) {SLESView(ts->sles,viewer);}
612: if (ts->snes) {SNESView(ts->snes,viewer);}
613: PetscViewerASCIIPopTab(viewer);
614: return(0);
615: }
620: /*@C
621: TSSetApplicationContext - Sets an optional user-defined context for
622: the timesteppers.
624: Collective on TS
626: Input Parameters:
627: + ts - the TS context obtained from TSCreate()
628: - usrP - optional user context
630: Level: intermediate
632: .keywords: TS, timestep, set, application, context
634: .seealso: TSGetApplicationContext()
635: @*/
636: int TSSetApplicationContext(TS ts,void *usrP)
637: {
640: ts->user = usrP;
641: return(0);
642: }
646: /*@C
647: TSGetApplicationContext - Gets the user-defined context for the
648: timestepper.
650: Not Collective
652: Input Parameter:
653: . ts - the TS context obtained from TSCreate()
655: Output Parameter:
656: . usrP - user context
658: Level: intermediate
660: .keywords: TS, timestep, get, application, context
662: .seealso: TSSetApplicationContext()
663: @*/
664: int TSGetApplicationContext(TS ts,void **usrP)
665: {
668: *usrP = ts->user;
669: return(0);
670: }
674: /*@
675: TSGetTimeStepNumber - Gets the current number of timesteps.
677: Not Collective
679: Input Parameter:
680: . ts - the TS context obtained from TSCreate()
682: Output Parameter:
683: . iter - number steps so far
685: Level: intermediate
687: .keywords: TS, timestep, get, iteration, number
688: @*/
689: int TSGetTimeStepNumber(TS ts,int* iter)
690: {
693: *iter = ts->steps;
694: return(0);
695: }
699: /*@
700: TSSetInitialTimeStep - Sets the initial timestep to be used,
701: as well as the initial time.
703: Collective on TS
705: Input Parameters:
706: + ts - the TS context obtained from TSCreate()
707: . initial_time - the initial time
708: - time_step - the size of the timestep
710: Level: intermediate
712: .seealso: TSSetTimeStep(), TSGetTimeStep()
714: .keywords: TS, set, initial, timestep
715: @*/
716: int TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
717: {
720: ts->time_step = time_step;
721: ts->initial_time_step = time_step;
722: ts->ptime = initial_time;
723: return(0);
724: }
728: /*@
729: TSSetTimeStep - Allows one to reset the timestep at any time,
730: useful for simple pseudo-timestepping codes.
732: Collective on TS
734: Input Parameters:
735: + ts - the TS context obtained from TSCreate()
736: - time_step - the size of the timestep
738: Level: intermediate
740: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
742: .keywords: TS, set, timestep
743: @*/
744: int TSSetTimeStep(TS ts,PetscReal time_step)
745: {
748: ts->time_step = time_step;
749: return(0);
750: }
754: /*@
755: TSGetTimeStep - Gets the current timestep size.
757: Not Collective
759: Input Parameter:
760: . ts - the TS context obtained from TSCreate()
762: Output Parameter:
763: . dt - the current timestep size
765: Level: intermediate
767: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
769: .keywords: TS, get, timestep
770: @*/
771: int TSGetTimeStep(TS ts,PetscReal* dt)
772: {
775: *dt = ts->time_step;
776: return(0);
777: }
781: /*@C
782: TSGetSolution - Returns the solution at the present timestep. It
783: is valid to call this routine inside the function that you are evaluating
784: in order to move to the new timestep. This vector not changed until
785: the solution at the next timestep has been calculated.
787: Not Collective, but Vec returned is parallel if TS is parallel
789: Input Parameter:
790: . ts - the TS context obtained from TSCreate()
792: Output Parameter:
793: . v - the vector containing the solution
795: Level: intermediate
797: .seealso: TSGetTimeStep()
799: .keywords: TS, timestep, get, solution
800: @*/
801: int TSGetSolution(TS ts,Vec *v)
802: {
805: *v = ts->vec_sol_always;
806: return(0);
807: }
809: /* ----- Routines to initialize and destroy a timestepper ---- */
812: /*@C
813: TSSetProblemType - Sets the type of problem to be solved.
815: Not collective
817: Input Parameters:
818: + ts - The TS
819: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
820: .vb
821: U_t = A U
822: U_t = A(t) U
823: U_t = F(t,U)
824: .ve
826: Level: beginner
828: .keywords: TS, problem type
829: .seealso: TSSetUp(), TSProblemType, TS
830: @*/
831: int TSSetProblemType(TS ts, TSProblemType type) {
834: ts->problem_type = type;
835: return(0);
836: }
840: /*@C
841: TSGetProblemType - Gets the type of problem to be solved.
843: Not collective
845: Input Parameter:
846: . ts - The TS
848: Output Parameter:
849: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
850: .vb
851: U_t = A U
852: U_t = A(t) U
853: U_t = F(t,U)
854: .ve
856: Level: beginner
858: .keywords: TS, problem type
859: .seealso: TSSetUp(), TSProblemType, TS
860: @*/
861: int TSGetProblemType(TS ts, TSProblemType *type) {
865: *type = ts->problem_type;
866: return(0);
867: }
871: /*@
872: TSSetUp - Sets up the internal data structures for the later use
873: of a timestepper.
875: Collective on TS
877: Input Parameter:
878: . ts - the TS context obtained from TSCreate()
880: Notes:
881: For basic use of the TS solvers the user need not explicitly call
882: TSSetUp(), since these actions will automatically occur during
883: the call to TSStep(). However, if one wishes to control this
884: phase separately, TSSetUp() should be called after TSCreate()
885: and optional routines of the form TSSetXXX(), but before TSStep().
887: Level: advanced
889: .keywords: TS, timestep, setup
891: .seealso: TSCreate(), TSStep(), TSDestroy()
892: @*/
893: int TSSetUp(TS ts)
894: {
899: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
900: if (!ts->type_name) {
901: TSSetType(ts,TS_EULER);
902: }
903: (*ts->ops->setup)(ts);
904: ts->setupcalled = 1;
905: return(0);
906: }
910: /*@C
911: TSDestroy - Destroys the timestepper context that was created
912: with TSCreate().
914: Collective on TS
916: Input Parameter:
917: . ts - the TS context obtained from TSCreate()
919: Level: beginner
921: .keywords: TS, timestepper, destroy
923: .seealso: TSCreate(), TSSetUp(), TSSolve()
924: @*/
925: int TSDestroy(TS ts)
926: {
927: int ierr,i;
931: if (--ts->refct > 0) return(0);
933: /* if memory was published with AMS then destroy it */
934: PetscObjectDepublish(ts);
936: if (ts->sles) {SLESDestroy(ts->sles);}
937: if (ts->snes) {SNESDestroy(ts->snes);}
938: (*(ts)->ops->destroy)(ts);
939: for (i=0; i<ts->numbermonitors; i++) {
940: if (ts->mdestroy[i]) {
941: (*ts->mdestroy[i])(ts->monitorcontext[i]);
942: }
943: }
944: PetscLogObjectDestroy((PetscObject)ts);
945: PetscHeaderDestroy((PetscObject)ts);
946: return(0);
947: }
951: /*@C
952: TSGetSNES - Returns the SNES (nonlinear solver) associated with
953: a TS (timestepper) context. Valid only for nonlinear problems.
955: Not Collective, but SNES is parallel if TS is parallel
957: Input Parameter:
958: . ts - the TS context obtained from TSCreate()
960: Output Parameter:
961: . snes - the nonlinear solver context
963: Notes:
964: The user can then directly manipulate the SNES context to set various
965: options, etc. Likewise, the user can then extract and manipulate the
966: SLES, KSP, and PC contexts as well.
968: TSGetSNES() does not work for integrators that do not use SNES; in
969: this case TSGetSNES() returns PETSC_NULL in snes.
971: Level: beginner
973: .keywords: timestep, get, SNES
974: @*/
975: int TSGetSNES(TS ts,SNES *snes)
976: {
979: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetSLES()");
980: *snes = ts->snes;
981: return(0);
982: }
986: /*@C
987: TSGetSLES - Returns the SLES (linear solver) associated with
988: a TS (timestepper) context.
990: Not Collective, but SLES is parallel if TS is parallel
992: Input Parameter:
993: . ts - the TS context obtained from TSCreate()
995: Output Parameter:
996: . sles - the nonlinear solver context
998: Notes:
999: The user can then directly manipulate the SLES context to set various
1000: options, etc. Likewise, the user can then extract and manipulate the
1001: KSP and PC contexts as well.
1003: TSGetSLES() does not work for integrators that do not use SLES;
1004: in this case TSGetSLES() returns PETSC_NULL in sles.
1006: Level: beginner
1008: .keywords: timestep, get, SLES
1009: @*/
1010: int TSGetSLES(TS ts,SLES *sles)
1011: {
1014: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1015: *sles = ts->sles;
1016: return(0);
1017: }
1019: /* ----------- Routines to set solver parameters ---------- */
1023: /*@
1024: TSGetDuration - Gets the maximum number of timesteps to use and
1025: maximum time for iteration.
1027: Collective on TS
1029: Input Parameters:
1030: + ts - the TS context obtained from TSCreate()
1031: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1032: - maxtime - final time to iterate to, or PETSC_NULL
1034: Level: intermediate
1036: .keywords: TS, timestep, get, maximum, iterations, time
1037: @*/
1038: int TSGetDuration(TS ts, int *maxsteps, PetscReal *maxtime)
1039: {
1042: if (maxsteps != PETSC_NULL) {
1044: *maxsteps = ts->max_steps;
1045: }
1046: if (maxtime != PETSC_NULL) {
1048: *maxtime = ts->max_time;
1049: }
1050: return(0);
1051: }
1055: /*@
1056: TSSetDuration - Sets the maximum number of timesteps to use and
1057: maximum time for iteration.
1059: Collective on TS
1061: Input Parameters:
1062: + ts - the TS context obtained from TSCreate()
1063: . maxsteps - maximum number of iterations to use
1064: - maxtime - final time to iterate to
1066: Options Database Keys:
1067: . -ts_max_steps <maxsteps> - Sets maxsteps
1068: . -ts_max_time <maxtime> - Sets maxtime
1070: Notes:
1071: The default maximum number of iterations is 5000. Default time is 5.0
1073: Level: intermediate
1075: .keywords: TS, timestep, set, maximum, iterations
1076: @*/
1077: int TSSetDuration(TS ts,int maxsteps,PetscReal maxtime)
1078: {
1081: ts->max_steps = maxsteps;
1082: ts->max_time = maxtime;
1083: return(0);
1084: }
1088: /*@
1089: TSSetSolution - Sets the initial solution vector
1090: for use by the TS routines.
1092: Collective on TS and Vec
1094: Input Parameters:
1095: + ts - the TS context obtained from TSCreate()
1096: - x - the solution vector
1098: Level: beginner
1100: .keywords: TS, timestep, set, solution, initial conditions
1101: @*/
1102: int TSSetSolution(TS ts,Vec x)
1103: {
1106: ts->vec_sol = ts->vec_sol_always = x;
1107: return(0);
1108: }
1112: /*@
1113: TSSetRhsBC - Sets the function which applies boundary conditions
1114: to the Rhs of each system.
1116: Collective on TS
1118: Input Parameters:
1119: + ts - The TS context obtained from TSCreate()
1120: - func - The function
1122: Calling sequence of func:
1123: . func (TS ts, Vec rhs, void *ctx);
1125: + rhs - The current rhs vector
1126: - ctx - The user-context
1128: Level: intermediate
1130: .keywords: TS, Rhs, boundary conditions
1131: @*/
1132: int TSSetRhsBC(TS ts, int (*func)(TS, Vec, void *))
1133: {
1136: ts->ops->applyrhsbc = func;
1137: return(0);
1138: }
1142: /*@
1143: TSDefaultRhsBC - The default boundary condition function which does nothing.
1145: Collective on TS
1147: Input Parameters:
1148: + ts - The TS context obtained from TSCreate()
1149: . rhs - The Rhs
1150: - ctx - The user-context
1152: Level: developer
1154: .keywords: TS, Rhs, boundary conditions
1155: @*/
1156: int TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1157: {
1159: return(0);
1160: }
1164: /*@
1165: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1166: to the system matrix and preconditioner of each system.
1168: Collective on TS
1170: Input Parameters:
1171: + ts - The TS context obtained from TSCreate()
1172: - func - The function
1174: Calling sequence of func:
1175: . func (TS ts, Mat A, Mat B, void *ctx);
1177: + A - The current system matrix
1178: . B - The current preconditioner
1179: - ctx - The user-context
1181: Level: intermediate
1183: .keywords: TS, System matrix, boundary conditions
1184: @*/
1185: int TSSetSystemMatrixBC(TS ts, int (*func)(TS, Mat, Mat, void *))
1186: {
1189: ts->ops->applymatrixbc = func;
1190: return(0);
1191: }
1195: /*@
1196: TSDefaultSystemMatrixBC - The default boundary condition function which
1197: does nothing.
1199: Collective on TS
1201: Input Parameters:
1202: + ts - The TS context obtained from TSCreate()
1203: . A - The system matrix
1204: . B - The preconditioner
1205: - ctx - The user-context
1207: Level: developer
1209: .keywords: TS, System matrix, boundary conditions
1210: @*/
1211: int TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1212: {
1214: return(0);
1215: }
1219: /*@
1220: TSSetSolutionBC - Sets the function which applies boundary conditions
1221: to the solution of each system. This is necessary in nonlinear systems
1222: which time dependent boundary conditions.
1224: Collective on TS
1226: Input Parameters:
1227: + ts - The TS context obtained from TSCreate()
1228: - func - The function
1230: Calling sequence of func:
1231: . func (TS ts, Vec rsol, void *ctx);
1233: + sol - The current solution vector
1234: - ctx - The user-context
1236: Level: intermediate
1238: .keywords: TS, solution, boundary conditions
1239: @*/
1240: int TSSetSolutionBC(TS ts, int (*func)(TS, Vec, void *))
1241: {
1244: ts->ops->applysolbc = func;
1245: return(0);
1246: }
1250: /*@
1251: TSDefaultSolutionBC - The default boundary condition function which
1252: does nothing.
1254: Collective on TS
1256: Input Parameters:
1257: + ts - The TS context obtained from TSCreate()
1258: . sol - The solution
1259: - ctx - The user-context
1261: Level: developer
1263: .keywords: TS, solution, boundary conditions
1264: @*/
1265: int TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1266: {
1268: return(0);
1269: }
1273: /*@
1274: TSSetPreStep - Sets the general-purpose function
1275: called once at the beginning of time stepping.
1277: Collective on TS
1279: Input Parameters:
1280: + ts - The TS context obtained from TSCreate()
1281: - func - The function
1283: Calling sequence of func:
1284: . func (TS ts);
1286: Level: intermediate
1288: .keywords: TS, timestep
1289: @*/
1290: int TSSetPreStep(TS ts, int (*func)(TS))
1291: {
1294: ts->ops->prestep = func;
1295: return(0);
1296: }
1300: /*@
1301: TSDefaultPreStep - The default pre-stepping function which does nothing.
1303: Collective on TS
1305: Input Parameters:
1306: . ts - The TS context obtained from TSCreate()
1308: Level: developer
1310: .keywords: TS, timestep
1311: @*/
1312: int TSDefaultPreStep(TS ts)
1313: {
1315: return(0);
1316: }
1320: /*@
1321: TSSetUpdate - Sets the general-purpose update function called
1322: at the beginning of every time step. This function can change
1323: the time step.
1325: Collective on TS
1327: Input Parameters:
1328: + ts - The TS context obtained from TSCreate()
1329: - func - The function
1331: Calling sequence of func:
1332: . func (TS ts, double t, double *dt);
1334: + t - The current time
1335: - dt - The current time step
1337: Level: intermediate
1339: .keywords: TS, update, timestep
1340: @*/
1341: int TSSetUpdate(TS ts, int (*func)(TS, PetscReal, PetscReal *))
1342: {
1345: ts->ops->update = func;
1346: return(0);
1347: }
1351: /*@
1352: TSDefaultUpdate - The default update function which does nothing.
1354: Collective on TS
1356: Input Parameters:
1357: + ts - The TS context obtained from TSCreate()
1358: - t - The current time
1360: Output Parameters:
1361: . dt - The current time step
1363: Level: developer
1365: .keywords: TS, update, timestep
1366: @*/
1367: int TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1368: {
1370: return(0);
1371: }
1375: /*@
1376: TSSetPostStep - Sets the general-purpose function
1377: called once at the end of time stepping.
1379: Collective on TS
1381: Input Parameters:
1382: + ts - The TS context obtained from TSCreate()
1383: - func - The function
1385: Calling sequence of func:
1386: . func (TS ts);
1388: Level: intermediate
1390: .keywords: TS, timestep
1391: @*/
1392: int TSSetPostStep(TS ts, int (*func)(TS))
1393: {
1396: ts->ops->poststep = func;
1397: return(0);
1398: }
1402: /*@
1403: TSDefaultPostStep - The default post-stepping function which does nothing.
1405: Collective on TS
1407: Input Parameters:
1408: . ts - The TS context obtained from TSCreate()
1410: Level: developer
1412: .keywords: TS, timestep
1413: @*/
1414: int TSDefaultPostStep(TS ts)
1415: {
1417: return(0);
1418: }
1420: /* ------------ Routines to set performance monitoring options ----------- */
1424: /*@C
1425: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1426: timestep to display the iteration's progress.
1428: Collective on TS
1430: Input Parameters:
1431: + ts - the TS context obtained from TSCreate()
1432: . func - monitoring routine
1433: . mctx - [optional] user-defined context for private data for the
1434: monitor routine (use PETSC_NULL if no context is desired)
1435: - monitordestroy - [optional] routine that frees monitor context
1436: (may be PETSC_NULL)
1438: Calling sequence of func:
1439: $ int func(TS ts,int steps,PetscReal time,Vec x,void *mctx)
1441: + ts - the TS context
1442: . steps - iteration number
1443: . time - current time
1444: . x - current iterate
1445: - mctx - [optional] monitoring context
1447: Notes:
1448: This routine adds an additional monitor to the list of monitors that
1449: already has been loaded.
1451: Level: intermediate
1453: .keywords: TS, timestep, set, monitor
1455: .seealso: TSDefaultMonitor(), TSClearMonitor()
1456: @*/
1457: int TSSetMonitor(TS ts,int (*monitor)(TS,int,PetscReal,Vec,void*),void *mctx,int (*mdestroy)(void*))
1458: {
1461: if (ts->numbermonitors >= MAXTSMONITORS) {
1462: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1463: }
1464: ts->monitor[ts->numbermonitors] = monitor;
1465: ts->mdestroy[ts->numbermonitors] = mdestroy;
1466: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1467: return(0);
1468: }
1472: /*@C
1473: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1475: Collective on TS
1477: Input Parameters:
1478: . ts - the TS context obtained from TSCreate()
1480: Notes:
1481: There is no way to remove a single, specific monitor.
1483: Level: intermediate
1485: .keywords: TS, timestep, set, monitor
1487: .seealso: TSDefaultMonitor(), TSSetMonitor()
1488: @*/
1489: int TSClearMonitor(TS ts)
1490: {
1493: ts->numbermonitors = 0;
1494: return(0);
1495: }
1499: int TSDefaultMonitor(TS ts,int step,PetscReal ptime,Vec v,void *ctx)
1500: {
1504: PetscPrintf(ts->comm,"timestep %d dt %g time %g\n",step,ts->time_step,ptime);
1505: return(0);
1506: }
1510: /*@
1511: TSStep - Steps the requested number of timesteps.
1513: Collective on TS
1515: Input Parameter:
1516: . ts - the TS context obtained from TSCreate()
1518: Output Parameters:
1519: + steps - number of iterations until termination
1520: - ptime - time until termination
1522: Level: beginner
1524: .keywords: TS, timestep, solve
1526: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1527: @*/
1528: int TSStep(TS ts,int *steps,PetscReal *ptime)
1529: {
1530: int ierr;
1534: if (!ts->setupcalled) {
1535: TSSetUp(ts);
1536: }
1538: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1539: (*ts->ops->prestep)(ts);
1540: (*ts->ops->step)(ts, steps, ptime);
1541: (*ts->ops->poststep)(ts);
1542: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1544: if (!PetscPreLoadingOn) {
1545: TSViewFromOptions(ts,ts->name);
1546: }
1547: return(0);
1548: }
1552: /*
1553: Runs the user provided monitor routines, if they exists.
1554: */
1555: int TSMonitor(TS ts,int step,PetscReal ptime,Vec x)
1556: {
1557: int i,ierr,n = ts->numbermonitors;
1560: for (i=0; i<n; i++) {
1561: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1562: }
1563: return(0);
1564: }
1566: /* ------------------------------------------------------------------------*/
1570: /*@C
1571: TSLGMonitorCreate - Creates a line graph context for use with
1572: TS to monitor convergence of preconditioned residual norms.
1574: Collective on TS
1576: Input Parameters:
1577: + host - the X display to open, or null for the local machine
1578: . label - the title to put in the title bar
1579: . x, y - the screen coordinates of the upper left coordinate of the window
1580: - m, n - the screen width and height in pixels
1582: Output Parameter:
1583: . draw - the drawing context
1585: Options Database Key:
1586: . -ts_xmonitor - automatically sets line graph monitor
1588: Notes:
1589: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1591: Level: intermediate
1593: .keywords: TS, monitor, line graph, residual, seealso
1595: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1597: @*/
1598: int TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1599: {
1600: PetscDraw win;
1601: int ierr;
1604: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1605: PetscDrawSetType(win,PETSC_DRAW_X);
1606: PetscDrawLGCreate(win,1,draw);
1607: PetscDrawLGIndicateDataPoints(*draw);
1609: PetscLogObjectParent(*draw,win);
1610: return(0);
1611: }
1615: int TSLGMonitor(TS ts,int n,PetscReal ptime,Vec v,void *monctx)
1616: {
1617: PetscDrawLG lg = (PetscDrawLG) monctx;
1618: PetscReal x,y = ptime;
1619: int ierr;
1622: if (!monctx) {
1623: MPI_Comm comm;
1624: PetscViewer viewer;
1626: PetscObjectGetComm((PetscObject)ts,&comm);
1627: viewer = PETSC_VIEWER_DRAW_(comm);
1628: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1629: }
1631: if (!n) {PetscDrawLGReset(lg);}
1632: x = (PetscReal)n;
1633: PetscDrawLGAddPoint(lg,&x,&y);
1634: if (n < 20 || (n % 5)) {
1635: PetscDrawLGDraw(lg);
1636: }
1637: return(0);
1638: }
1642: /*@C
1643: TSLGMonitorDestroy - Destroys a line graph context that was created
1644: with TSLGMonitorCreate().
1646: Collective on PetscDrawLG
1648: Input Parameter:
1649: . draw - the drawing context
1651: Level: intermediate
1653: .keywords: TS, monitor, line graph, destroy
1655: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1656: @*/
1657: int TSLGMonitorDestroy(PetscDrawLG drawlg)
1658: {
1659: PetscDraw draw;
1660: int ierr;
1663: PetscDrawLGGetDraw(drawlg,&draw);
1664: PetscDrawDestroy(draw);
1665: PetscDrawLGDestroy(drawlg);
1666: return(0);
1667: }
1671: /*@
1672: TSGetTime - Gets the current time.
1674: Not Collective
1676: Input Parameter:
1677: . ts - the TS context obtained from TSCreate()
1679: Output Parameter:
1680: . t - the current time
1682: Contributed by: Matthew Knepley
1684: Level: beginner
1686: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1688: .keywords: TS, get, time
1689: @*/
1690: int TSGetTime(TS ts,PetscReal* t)
1691: {
1694: *t = ts->ptime;
1695: return(0);
1696: }
1700: /*@C
1701: TSSetOptionsPrefix - Sets the prefix used for searching for all
1702: TS options in the database.
1704: Collective on TS
1706: Input Parameter:
1707: + ts - The TS context
1708: - prefix - The prefix to prepend to all option names
1710: Notes:
1711: A hyphen (-) must NOT be given at the beginning of the prefix name.
1712: The first character of all runtime options is AUTOMATICALLY the
1713: hyphen.
1715: Contributed by: Matthew Knepley
1717: Level: advanced
1719: .keywords: TS, set, options, prefix, database
1721: .seealso: TSSetFromOptions()
1723: @*/
1724: int TSSetOptionsPrefix(TS ts,const char prefix[])
1725: {
1730: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1731: switch(ts->problem_type) {
1732: case TS_NONLINEAR:
1733: SNESSetOptionsPrefix(ts->snes,prefix);
1734: break;
1735: case TS_LINEAR:
1736: SLESSetOptionsPrefix(ts->sles,prefix);
1737: break;
1738: }
1739: return(0);
1740: }
1745: /*@C
1746: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1747: TS options in the database.
1749: Collective on TS
1751: Input Parameter:
1752: + ts - The TS context
1753: - prefix - The prefix to prepend to all option names
1755: Notes:
1756: A hyphen (-) must NOT be given at the beginning of the prefix name.
1757: The first character of all runtime options is AUTOMATICALLY the
1758: hyphen.
1760: Contributed by: Matthew Knepley
1762: Level: advanced
1764: .keywords: TS, append, options, prefix, database
1766: .seealso: TSGetOptionsPrefix()
1768: @*/
1769: int TSAppendOptionsPrefix(TS ts,const char prefix[])
1770: {
1775: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1776: switch(ts->problem_type) {
1777: case TS_NONLINEAR:
1778: SNESAppendOptionsPrefix(ts->snes,prefix);
1779: break;
1780: case TS_LINEAR:
1781: SLESAppendOptionsPrefix(ts->sles,prefix);
1782: break;
1783: }
1784: return(0);
1785: }
1789: /*@C
1790: TSGetOptionsPrefix - Sets the prefix used for searching for all
1791: TS options in the database.
1793: Not Collective
1795: Input Parameter:
1796: . ts - The TS context
1798: Output Parameter:
1799: . prefix - A pointer to the prefix string used
1801: Contributed by: Matthew Knepley
1803: Notes: On the fortran side, the user should pass in a string 'prifix' of
1804: sufficient length to hold the prefix.
1806: Level: intermediate
1808: .keywords: TS, get, options, prefix, database
1810: .seealso: TSAppendOptionsPrefix()
1811: @*/
1812: int TSGetOptionsPrefix(TS ts,char *prefix[])
1813: {
1818: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1819: return(0);
1820: }
1824: /*@C
1825: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1827: Not Collective, but parallel objects are returned if TS is parallel
1829: Input Parameter:
1830: . ts - The TS context obtained from TSCreate()
1832: Output Parameters:
1833: + A - The matrix A, where U_t = A(t) U
1834: . M - The preconditioner matrix, usually the same as A
1835: - ctx - User-defined context for matrix evaluation routine
1837: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1839: Contributed by: Matthew Knepley
1841: Level: intermediate
1843: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1845: .keywords: TS, timestep, get, matrix
1847: @*/
1848: int TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1849: {
1852: if (A) *A = ts->A;
1853: if (M) *M = ts->B;
1854: if (ctx) *ctx = ts->jacP;
1855: return(0);
1856: }
1860: /*@C
1861: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1863: Not Collective, but parallel objects are returned if TS is parallel
1865: Input Parameter:
1866: . ts - The TS context obtained from TSCreate()
1868: Output Parameters:
1869: + J - The Jacobian J of F, where U_t = F(U,t)
1870: . M - The preconditioner matrix, usually the same as J
1871: - ctx - User-defined context for Jacobian evaluation routine
1873: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1875: Contributed by: Matthew Knepley
1877: Level: intermediate
1879: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1881: .keywords: TS, timestep, get, matrix, Jacobian
1882: @*/
1883: int TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1884: {
1888: TSGetRHSMatrix(ts,J,M,ctx);
1889: return(0);
1890: }
1894: /*@C
1895: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1896: VecView() for the solution at each timestep
1898: Collective on TS
1900: Input Parameters:
1901: + ts - the TS context
1902: . step - current time-step
1903: . ptime - current time
1904: - dummy - either a viewer or PETSC_NULL
1906: Level: intermediate
1908: .keywords: TS, vector, monitor, view
1910: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1911: @*/
1912: int TSVecViewMonitor(TS ts,int step,PetscReal ptime,Vec x,void *dummy)
1913: {
1914: int ierr;
1915: PetscViewer viewer = (PetscViewer) dummy;
1918: if (!viewer) {
1919: MPI_Comm comm;
1920: PetscObjectGetComm((PetscObject)ts,&comm);
1921: viewer = PETSC_VIEWER_DRAW_(comm);
1922: }
1923: VecView(x,viewer);
1924: return(0);
1925: }