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 = 0;
6: int TS_Step = 0, TS_PseudoComputeTimeStep = 0, TS_FunctionEval = 0, TS_JacobianEval = 0;
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: const 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->ksp != PETSC_NULL) {
120: KSPSetFromOptions(ts->ksp);
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 KSPSetOperators() 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(), KSPSetOperators()
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 KSPSetOperators())
372: - ctx - [optional] user-defined context for matrix evaluation routine
374: Notes:
375: See KSPSetOperators() 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.
418: Use TSSetRHSMatrix() for linear problems.
420: Collective on TS
422: Input Parameters:
423: + ts - the TS context obtained from TSCreate()
424: . A - Jacobian matrix
425: . B - preconditioner matrix (usually same as A)
426: . f - the Jacobian evaluation routine
427: - ctx - [optional] user-defined context for private data for the
428: Jacobian evaluation routine (may be PETSC_NULL)
430: Calling sequence of func:
431: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
433: + t - current timestep
434: . u - input vector
435: . A - matrix A, where U_t = A(t)u
436: . B - preconditioner matrix, usually the same as A
437: . flag - flag indicating information about the preconditioner matrix
438: structure (same as flag in KSPSetOperators())
439: - ctx - [optional] user-defined context for matrix evaluation routine
441: Notes:
442: See KSPSetOperators() for important information about setting the flag
443: output parameter in the routine func(). Be sure to read this information!
445: The routine func() takes Mat * as the matrix arguments rather than Mat.
446: This allows the matrix evaluation routine to replace A and/or B with a
447: completely new new matrix structure (not just different matrix elements)
448: when appropriate, for instance, if the nonzero structure is changing
449: throughout the global iterations.
451: Level: beginner
452:
453: .keywords: TS, timestep, set, right-hand-side, Jacobian
455: .seealso: TSDefaultComputeJacobianColor(),
456: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()
458: @*/
459: int TSSetRHSJacobian(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
460: {
467: if (ts->problem_type != TS_NONLINEAR) {
468: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
469: }
471: ts->ops->rhsjacobian = f;
472: ts->jacP = ctx;
473: ts->A = A;
474: ts->B = B;
475: return(0);
476: }
480: /*
481: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
483: Note: If the user did not provide a function but merely a matrix,
484: this routine applies the matrix.
485: */
486: int TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
487: {
495: if (ts->ops->rhsbc) {
496: PetscStackPush("TS user boundary condition function");
497: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
498: PetscStackPop;
499: return(0);
500: }
502: return(0);
503: }
507: /*@C
508: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
509: boundary conditions for the function F.
511: Collective on TS
513: Input Parameters:
514: + ts - the TS context obtained from TSCreate()
515: . f - routine for evaluating the boundary condition function
516: - ctx - [optional] user-defined context for private data for the
517: function evaluation routine (may be PETSC_NULL)
519: Calling sequence of func:
520: $ func (TS ts,PetscReal t,Vec F,void *ctx);
522: + t - current timestep
523: . F - function vector
524: - ctx - [optional] user-defined function context
526: Level: intermediate
528: .keywords: TS, timestep, set, boundary conditions, function
529: @*/
530: int TSSetRHSBoundaryConditions(TS ts,int (*f)(TS,PetscReal,Vec,void*),void *ctx)
531: {
535: if (ts->problem_type != TS_LINEAR) {
536: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
537: }
538: ts->ops->rhsbc = f;
539: ts->bcP = ctx;
540: return(0);
541: }
545: /*@C
546: TSView - Prints the TS data structure.
548: Collective on TS
550: Input Parameters:
551: + ts - the TS context obtained from TSCreate()
552: - viewer - visualization context
554: Options Database Key:
555: . -ts_view - calls TSView() at end of TSStep()
557: Notes:
558: The available visualization contexts include
559: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
560: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
561: output where only the first processor opens
562: the file. All other processors send their
563: data to the first processor to print.
565: The user can open an alternative visualization context with
566: PetscViewerASCIIOpen() - output to a specified file.
568: Level: beginner
570: .keywords: TS, timestep, view
572: .seealso: PetscViewerASCIIOpen()
573: @*/
574: int TSView(TS ts,PetscViewer viewer)
575: {
576: int ierr;
577: char *type;
578: PetscTruth isascii,isstring;
582: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
586: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
587: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
588: if (isascii) {
589: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
590: TSGetType(ts,(TSType *)&type);
591: if (type) {
592: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
593: } else {
594: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
595: }
596: if (ts->ops->view) {
597: PetscViewerASCIIPushTab(viewer);
598: (*ts->ops->view)(ts,viewer);
599: PetscViewerASCIIPopTab(viewer);
600: }
601: PetscViewerASCIIPrintf(viewer," maximum steps=%d\n",ts->max_steps);
602: PetscViewerASCIIPrintf(viewer," maximum time=%g\n",ts->max_time);
603: if (ts->problem_type == TS_NONLINEAR) {
604: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%d\n",ts->nonlinear_its);
605: }
606: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%d\n",ts->linear_its);
607: } else if (isstring) {
608: TSGetType(ts,(TSType *)&type);
609: PetscViewerStringSPrintf(viewer," %-7.7s",type);
610: }
611: PetscViewerASCIIPushTab(viewer);
612: if (ts->ksp) {KSPView(ts->ksp,viewer);}
613: if (ts->snes) {SNESView(ts->snes,viewer);}
614: PetscViewerASCIIPopTab(viewer);
615: return(0);
616: }
621: /*@C
622: TSSetApplicationContext - Sets an optional user-defined context for
623: the timesteppers.
625: Collective on TS
627: Input Parameters:
628: + ts - the TS context obtained from TSCreate()
629: - usrP - optional user context
631: Level: intermediate
633: .keywords: TS, timestep, set, application, context
635: .seealso: TSGetApplicationContext()
636: @*/
637: int TSSetApplicationContext(TS ts,void *usrP)
638: {
641: ts->user = usrP;
642: return(0);
643: }
647: /*@C
648: TSGetApplicationContext - Gets the user-defined context for the
649: timestepper.
651: Not Collective
653: Input Parameter:
654: . ts - the TS context obtained from TSCreate()
656: Output Parameter:
657: . usrP - user context
659: Level: intermediate
661: .keywords: TS, timestep, get, application, context
663: .seealso: TSSetApplicationContext()
664: @*/
665: int TSGetApplicationContext(TS ts,void **usrP)
666: {
669: *usrP = ts->user;
670: return(0);
671: }
675: /*@
676: TSGetTimeStepNumber - Gets the current number of timesteps.
678: Not Collective
680: Input Parameter:
681: . ts - the TS context obtained from TSCreate()
683: Output Parameter:
684: . iter - number steps so far
686: Level: intermediate
688: .keywords: TS, timestep, get, iteration, number
689: @*/
690: int TSGetTimeStepNumber(TS ts,int* iter)
691: {
695: *iter = ts->steps;
696: return(0);
697: }
701: /*@
702: TSSetInitialTimeStep - Sets the initial timestep to be used,
703: as well as the initial time.
705: Collective on TS
707: Input Parameters:
708: + ts - the TS context obtained from TSCreate()
709: . initial_time - the initial time
710: - time_step - the size of the timestep
712: Level: intermediate
714: .seealso: TSSetTimeStep(), TSGetTimeStep()
716: .keywords: TS, set, initial, timestep
717: @*/
718: int TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
719: {
722: ts->time_step = time_step;
723: ts->initial_time_step = time_step;
724: ts->ptime = initial_time;
725: return(0);
726: }
730: /*@
731: TSSetTimeStep - Allows one to reset the timestep at any time,
732: useful for simple pseudo-timestepping codes.
734: Collective on TS
736: Input Parameters:
737: + ts - the TS context obtained from TSCreate()
738: - time_step - the size of the timestep
740: Level: intermediate
742: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
744: .keywords: TS, set, timestep
745: @*/
746: int TSSetTimeStep(TS ts,PetscReal time_step)
747: {
750: ts->time_step = time_step;
751: return(0);
752: }
756: /*@
757: TSGetTimeStep - Gets the current timestep size.
759: Not Collective
761: Input Parameter:
762: . ts - the TS context obtained from TSCreate()
764: Output Parameter:
765: . dt - the current timestep size
767: Level: intermediate
769: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
771: .keywords: TS, get, timestep
772: @*/
773: int TSGetTimeStep(TS ts,PetscReal* dt)
774: {
778: *dt = ts->time_step;
779: return(0);
780: }
784: /*@C
785: TSGetSolution - Returns the solution at the present timestep. It
786: is valid to call this routine inside the function that you are evaluating
787: in order to move to the new timestep. This vector not changed until
788: the solution at the next timestep has been calculated.
790: Not Collective, but Vec returned is parallel if TS is parallel
792: Input Parameter:
793: . ts - the TS context obtained from TSCreate()
795: Output Parameter:
796: . v - the vector containing the solution
798: Level: intermediate
800: .seealso: TSGetTimeStep()
802: .keywords: TS, timestep, get, solution
803: @*/
804: int TSGetSolution(TS ts,Vec *v)
805: {
809: *v = ts->vec_sol_always;
810: return(0);
811: }
813: /* ----- Routines to initialize and destroy a timestepper ---- */
816: /*@C
817: TSSetProblemType - Sets the type of problem to be solved.
819: Not collective
821: Input Parameters:
822: + ts - The TS
823: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
824: .vb
825: U_t = A U
826: U_t = A(t) U
827: U_t = F(t,U)
828: .ve
830: Level: beginner
832: .keywords: TS, problem type
833: .seealso: TSSetUp(), TSProblemType, TS
834: @*/
835: int TSSetProblemType(TS ts, TSProblemType type) {
838: ts->problem_type = type;
839: return(0);
840: }
844: /*@C
845: TSGetProblemType - Gets the type of problem to be solved.
847: Not collective
849: Input Parameter:
850: . ts - The TS
852: Output Parameter:
853: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
854: .vb
855: U_t = A U
856: U_t = A(t) U
857: U_t = F(t,U)
858: .ve
860: Level: beginner
862: .keywords: TS, problem type
863: .seealso: TSSetUp(), TSProblemType, TS
864: @*/
865: int TSGetProblemType(TS ts, TSProblemType *type) {
869: *type = ts->problem_type;
870: return(0);
871: }
875: /*@
876: TSSetUp - Sets up the internal data structures for the later use
877: of a timestepper.
879: Collective on TS
881: Input Parameter:
882: . ts - the TS context obtained from TSCreate()
884: Notes:
885: For basic use of the TS solvers the user need not explicitly call
886: TSSetUp(), since these actions will automatically occur during
887: the call to TSStep(). However, if one wishes to control this
888: phase separately, TSSetUp() should be called after TSCreate()
889: and optional routines of the form TSSetXXX(), but before TSStep().
891: Level: advanced
893: .keywords: TS, timestep, setup
895: .seealso: TSCreate(), TSStep(), TSDestroy()
896: @*/
897: int TSSetUp(TS ts)
898: {
903: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
904: if (!ts->type_name) {
905: TSSetType(ts,TS_EULER);
906: }
907: (*ts->ops->setup)(ts);
908: ts->setupcalled = 1;
909: return(0);
910: }
914: /*@C
915: TSDestroy - Destroys the timestepper context that was created
916: with TSCreate().
918: Collective on TS
920: Input Parameter:
921: . ts - the TS context obtained from TSCreate()
923: Level: beginner
925: .keywords: TS, timestepper, destroy
927: .seealso: TSCreate(), TSSetUp(), TSSolve()
928: @*/
929: int TSDestroy(TS ts)
930: {
931: int ierr,i;
935: if (--ts->refct > 0) return(0);
937: /* if memory was published with AMS then destroy it */
938: PetscObjectDepublish(ts);
940: if (ts->ksp) {KSPDestroy(ts->ksp);}
941: if (ts->snes) {SNESDestroy(ts->snes);}
942: (*(ts)->ops->destroy)(ts);
943: for (i=0; i<ts->numbermonitors; i++) {
944: if (ts->mdestroy[i]) {
945: (*ts->mdestroy[i])(ts->monitorcontext[i]);
946: }
947: }
948: PetscLogObjectDestroy((PetscObject)ts);
949: PetscHeaderDestroy((PetscObject)ts);
950: return(0);
951: }
955: /*@C
956: TSGetSNES - Returns the SNES (nonlinear solver) associated with
957: a TS (timestepper) context. Valid only for nonlinear problems.
959: Not Collective, but SNES is parallel if TS is parallel
961: Input Parameter:
962: . ts - the TS context obtained from TSCreate()
964: Output Parameter:
965: . snes - the nonlinear solver context
967: Notes:
968: The user can then directly manipulate the SNES context to set various
969: options, etc. Likewise, the user can then extract and manipulate the
970: KSP, KSP, and PC contexts as well.
972: TSGetSNES() does not work for integrators that do not use SNES; in
973: this case TSGetSNES() returns PETSC_NULL in snes.
975: Level: beginner
977: .keywords: timestep, get, SNES
978: @*/
979: int TSGetSNES(TS ts,SNES *snes)
980: {
984: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
985: *snes = ts->snes;
986: return(0);
987: }
991: /*@C
992: TSGetKSP - Returns the KSP (linear solver) associated with
993: a TS (timestepper) context.
995: Not Collective, but KSP is parallel if TS is parallel
997: Input Parameter:
998: . ts - the TS context obtained from TSCreate()
1000: Output Parameter:
1001: . ksp - the nonlinear solver context
1003: Notes:
1004: The user can then directly manipulate the KSP context to set various
1005: options, etc. Likewise, the user can then extract and manipulate the
1006: KSP and PC contexts as well.
1008: TSGetKSP() does not work for integrators that do not use KSP;
1009: in this case TSGetKSP() returns PETSC_NULL in ksp.
1011: Level: beginner
1013: .keywords: timestep, get, KSP
1014: @*/
1015: int TSGetKSP(TS ts,KSP *ksp)
1016: {
1020: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1021: *ksp = ts->ksp;
1022: return(0);
1023: }
1025: /* ----------- Routines to set solver parameters ---------- */
1029: /*@
1030: TSGetDuration - Gets the maximum number of timesteps to use and
1031: maximum time for iteration.
1033: Collective on TS
1035: Input Parameters:
1036: + ts - the TS context obtained from TSCreate()
1037: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1038: - maxtime - final time to iterate to, or PETSC_NULL
1040: Level: intermediate
1042: .keywords: TS, timestep, get, maximum, iterations, time
1043: @*/
1044: int TSGetDuration(TS ts, int *maxsteps, PetscReal *maxtime)
1045: {
1048: if (maxsteps != PETSC_NULL) {
1050: *maxsteps = ts->max_steps;
1051: }
1052: if (maxtime != PETSC_NULL) {
1054: *maxtime = ts->max_time;
1055: }
1056: return(0);
1057: }
1061: /*@
1062: TSSetDuration - Sets the maximum number of timesteps to use and
1063: maximum time for iteration.
1065: Collective on TS
1067: Input Parameters:
1068: + ts - the TS context obtained from TSCreate()
1069: . maxsteps - maximum number of iterations to use
1070: - maxtime - final time to iterate to
1072: Options Database Keys:
1073: . -ts_max_steps <maxsteps> - Sets maxsteps
1074: . -ts_max_time <maxtime> - Sets maxtime
1076: Notes:
1077: The default maximum number of iterations is 5000. Default time is 5.0
1079: Level: intermediate
1081: .keywords: TS, timestep, set, maximum, iterations
1082: @*/
1083: int TSSetDuration(TS ts,int maxsteps,PetscReal maxtime)
1084: {
1087: ts->max_steps = maxsteps;
1088: ts->max_time = maxtime;
1089: return(0);
1090: }
1094: /*@
1095: TSSetSolution - Sets the initial solution vector
1096: for use by the TS routines.
1098: Collective on TS and Vec
1100: Input Parameters:
1101: + ts - the TS context obtained from TSCreate()
1102: - x - the solution vector
1104: Level: beginner
1106: .keywords: TS, timestep, set, solution, initial conditions
1107: @*/
1108: int TSSetSolution(TS ts,Vec x)
1109: {
1113: ts->vec_sol = ts->vec_sol_always = x;
1114: return(0);
1115: }
1119: /*@C
1120: TSSetRhsBC - Sets the function which applies boundary conditions
1121: to the Rhs of each system.
1123: Collective on TS
1125: Input Parameters:
1126: + ts - The TS context obtained from TSCreate()
1127: - func - The function
1129: Calling sequence of func:
1130: . func (TS ts, Vec rhs, void *ctx);
1132: + rhs - The current rhs vector
1133: - ctx - The user-context
1135: Level: intermediate
1137: .keywords: TS, Rhs, boundary conditions
1138: @*/
1139: int TSSetRhsBC(TS ts, int (*func)(TS, Vec, void *))
1140: {
1143: ts->ops->applyrhsbc = func;
1144: return(0);
1145: }
1149: /*@
1150: TSDefaultRhsBC - The default boundary condition function which does nothing.
1152: Collective on TS
1154: Input Parameters:
1155: + ts - The TS context obtained from TSCreate()
1156: . rhs - The Rhs
1157: - ctx - The user-context
1159: Level: developer
1161: .keywords: TS, Rhs, boundary conditions
1162: @*/
1163: int TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1164: {
1166: return(0);
1167: }
1171: /*@C
1172: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1173: to the system matrix and preconditioner of each system.
1175: Collective on TS
1177: Input Parameters:
1178: + ts - The TS context obtained from TSCreate()
1179: - func - The function
1181: Calling sequence of func:
1182: . func (TS ts, Mat A, Mat B, void *ctx);
1184: + A - The current system matrix
1185: . B - The current preconditioner
1186: - ctx - The user-context
1188: Level: intermediate
1190: .keywords: TS, System matrix, boundary conditions
1191: @*/
1192: int TSSetSystemMatrixBC(TS ts, int (*func)(TS, Mat, Mat, void *))
1193: {
1196: ts->ops->applymatrixbc = func;
1197: return(0);
1198: }
1202: /*@
1203: TSDefaultSystemMatrixBC - The default boundary condition function which
1204: does nothing.
1206: Collective on TS
1208: Input Parameters:
1209: + ts - The TS context obtained from TSCreate()
1210: . A - The system matrix
1211: . B - The preconditioner
1212: - ctx - The user-context
1214: Level: developer
1216: .keywords: TS, System matrix, boundary conditions
1217: @*/
1218: int TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1219: {
1221: return(0);
1222: }
1226: /*@C
1227: TSSetSolutionBC - Sets the function which applies boundary conditions
1228: to the solution of each system. This is necessary in nonlinear systems
1229: which time dependent boundary conditions.
1231: Collective on TS
1233: Input Parameters:
1234: + ts - The TS context obtained from TSCreate()
1235: - func - The function
1237: Calling sequence of func:
1238: . func (TS ts, Vec rsol, void *ctx);
1240: + sol - The current solution vector
1241: - ctx - The user-context
1243: Level: intermediate
1245: .keywords: TS, solution, boundary conditions
1246: @*/
1247: int TSSetSolutionBC(TS ts, int (*func)(TS, Vec, void *))
1248: {
1251: ts->ops->applysolbc = func;
1252: return(0);
1253: }
1257: /*@
1258: TSDefaultSolutionBC - The default boundary condition function which
1259: does nothing.
1261: Collective on TS
1263: Input Parameters:
1264: + ts - The TS context obtained from TSCreate()
1265: . sol - The solution
1266: - ctx - The user-context
1268: Level: developer
1270: .keywords: TS, solution, boundary conditions
1271: @*/
1272: int TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1273: {
1275: return(0);
1276: }
1280: /*@C
1281: TSSetPreStep - Sets the general-purpose function
1282: called once at the beginning of time stepping.
1284: Collective on TS
1286: Input Parameters:
1287: + ts - The TS context obtained from TSCreate()
1288: - func - The function
1290: Calling sequence of func:
1291: . func (TS ts);
1293: Level: intermediate
1295: .keywords: TS, timestep
1296: @*/
1297: int TSSetPreStep(TS ts, int (*func)(TS))
1298: {
1301: ts->ops->prestep = func;
1302: return(0);
1303: }
1307: /*@
1308: TSDefaultPreStep - The default pre-stepping function which does nothing.
1310: Collective on TS
1312: Input Parameters:
1313: . ts - The TS context obtained from TSCreate()
1315: Level: developer
1317: .keywords: TS, timestep
1318: @*/
1319: int TSDefaultPreStep(TS ts)
1320: {
1322: return(0);
1323: }
1327: /*@C
1328: TSSetUpdate - Sets the general-purpose update function called
1329: at the beginning of every time step. This function can change
1330: the time step.
1332: Collective on TS
1334: Input Parameters:
1335: + ts - The TS context obtained from TSCreate()
1336: - func - The function
1338: Calling sequence of func:
1339: . func (TS ts, double t, double *dt);
1341: + t - The current time
1342: - dt - The current time step
1344: Level: intermediate
1346: .keywords: TS, update, timestep
1347: @*/
1348: int TSSetUpdate(TS ts, int (*func)(TS, PetscReal, PetscReal *))
1349: {
1352: ts->ops->update = func;
1353: return(0);
1354: }
1358: /*@
1359: TSDefaultUpdate - The default update function which does nothing.
1361: Collective on TS
1363: Input Parameters:
1364: + ts - The TS context obtained from TSCreate()
1365: - t - The current time
1367: Output Parameters:
1368: . dt - The current time step
1370: Level: developer
1372: .keywords: TS, update, timestep
1373: @*/
1374: int TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1375: {
1377: return(0);
1378: }
1382: /*@C
1383: TSSetPostStep - Sets the general-purpose function
1384: called once at the end of time stepping.
1386: Collective on TS
1388: Input Parameters:
1389: + ts - The TS context obtained from TSCreate()
1390: - func - The function
1392: Calling sequence of func:
1393: . func (TS ts);
1395: Level: intermediate
1397: .keywords: TS, timestep
1398: @*/
1399: int TSSetPostStep(TS ts, int (*func)(TS))
1400: {
1403: ts->ops->poststep = func;
1404: return(0);
1405: }
1409: /*@
1410: TSDefaultPostStep - The default post-stepping function which does nothing.
1412: Collective on TS
1414: Input Parameters:
1415: . ts - The TS context obtained from TSCreate()
1417: Level: developer
1419: .keywords: TS, timestep
1420: @*/
1421: int TSDefaultPostStep(TS ts)
1422: {
1424: return(0);
1425: }
1427: /* ------------ Routines to set performance monitoring options ----------- */
1431: /*@C
1432: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1433: timestep to display the iteration's progress.
1435: Collective on TS
1437: Input Parameters:
1438: + ts - the TS context obtained from TSCreate()
1439: . func - monitoring routine
1440: . mctx - [optional] user-defined context for private data for the
1441: monitor routine (use PETSC_NULL if no context is desired)
1442: - monitordestroy - [optional] routine that frees monitor context
1443: (may be PETSC_NULL)
1445: Calling sequence of func:
1446: $ int func(TS ts,int steps,PetscReal time,Vec x,void *mctx)
1448: + ts - the TS context
1449: . steps - iteration number
1450: . time - current time
1451: . x - current iterate
1452: - mctx - [optional] monitoring context
1454: Notes:
1455: This routine adds an additional monitor to the list of monitors that
1456: already has been loaded.
1458: Level: intermediate
1460: .keywords: TS, timestep, set, monitor
1462: .seealso: TSDefaultMonitor(), TSClearMonitor()
1463: @*/
1464: int TSSetMonitor(TS ts,int (*monitor)(TS,int,PetscReal,Vec,void*),void *mctx,int (*mdestroy)(void*))
1465: {
1468: if (ts->numbermonitors >= MAXTSMONITORS) {
1469: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1470: }
1471: ts->monitor[ts->numbermonitors] = monitor;
1472: ts->mdestroy[ts->numbermonitors] = mdestroy;
1473: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1474: return(0);
1475: }
1479: /*@C
1480: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1482: Collective on TS
1484: Input Parameters:
1485: . ts - the TS context obtained from TSCreate()
1487: Notes:
1488: There is no way to remove a single, specific monitor.
1490: Level: intermediate
1492: .keywords: TS, timestep, set, monitor
1494: .seealso: TSDefaultMonitor(), TSSetMonitor()
1495: @*/
1496: int TSClearMonitor(TS ts)
1497: {
1500: ts->numbermonitors = 0;
1501: return(0);
1502: }
1506: int TSDefaultMonitor(TS ts,int step,PetscReal ptime,Vec v,void *ctx)
1507: {
1511: PetscPrintf(ts->comm,"timestep %d dt %g time %g\n",step,ts->time_step,ptime);
1512: return(0);
1513: }
1517: /*@
1518: TSStep - Steps the requested number of timesteps.
1520: Collective on TS
1522: Input Parameter:
1523: . ts - the TS context obtained from TSCreate()
1525: Output Parameters:
1526: + steps - number of iterations until termination
1527: - ptime - time until termination
1529: Level: beginner
1531: .keywords: TS, timestep, solve
1533: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1534: @*/
1535: int TSStep(TS ts,int *steps,PetscReal *ptime)
1536: {
1537: int ierr;
1541: if (!ts->setupcalled) {
1542: TSSetUp(ts);
1543: }
1545: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1546: (*ts->ops->prestep)(ts);
1547: (*ts->ops->step)(ts, steps, ptime);
1548: (*ts->ops->poststep)(ts);
1549: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1551: if (!PetscPreLoadingOn) {
1552: TSViewFromOptions(ts,ts->name);
1553: }
1554: return(0);
1555: }
1559: /*
1560: Runs the user provided monitor routines, if they exists.
1561: */
1562: int TSMonitor(TS ts,int step,PetscReal ptime,Vec x)
1563: {
1564: int i,ierr,n = ts->numbermonitors;
1567: for (i=0; i<n; i++) {
1568: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1569: }
1570: return(0);
1571: }
1573: /* ------------------------------------------------------------------------*/
1577: /*@C
1578: TSLGMonitorCreate - Creates a line graph context for use with
1579: TS to monitor convergence of preconditioned residual norms.
1581: Collective on TS
1583: Input Parameters:
1584: + host - the X display to open, or null for the local machine
1585: . label - the title to put in the title bar
1586: . x, y - the screen coordinates of the upper left coordinate of the window
1587: - m, n - the screen width and height in pixels
1589: Output Parameter:
1590: . draw - the drawing context
1592: Options Database Key:
1593: . -ts_xmonitor - automatically sets line graph monitor
1595: Notes:
1596: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1598: Level: intermediate
1600: .keywords: TS, monitor, line graph, residual, seealso
1602: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1604: @*/
1605: int TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1606: {
1607: PetscDraw win;
1608: int ierr;
1611: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1612: PetscDrawSetType(win,PETSC_DRAW_X);
1613: PetscDrawLGCreate(win,1,draw);
1614: PetscDrawLGIndicateDataPoints(*draw);
1616: PetscLogObjectParent(*draw,win);
1617: return(0);
1618: }
1622: int TSLGMonitor(TS ts,int n,PetscReal ptime,Vec v,void *monctx)
1623: {
1624: PetscDrawLG lg = (PetscDrawLG) monctx;
1625: PetscReal x,y = ptime;
1626: int ierr;
1629: if (!monctx) {
1630: MPI_Comm comm;
1631: PetscViewer viewer;
1633: PetscObjectGetComm((PetscObject)ts,&comm);
1634: viewer = PETSC_VIEWER_DRAW_(comm);
1635: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1636: }
1638: if (!n) {PetscDrawLGReset(lg);}
1639: x = (PetscReal)n;
1640: PetscDrawLGAddPoint(lg,&x,&y);
1641: if (n < 20 || (n % 5)) {
1642: PetscDrawLGDraw(lg);
1643: }
1644: return(0);
1645: }
1649: /*@C
1650: TSLGMonitorDestroy - Destroys a line graph context that was created
1651: with TSLGMonitorCreate().
1653: Collective on PetscDrawLG
1655: Input Parameter:
1656: . draw - the drawing context
1658: Level: intermediate
1660: .keywords: TS, monitor, line graph, destroy
1662: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1663: @*/
1664: int TSLGMonitorDestroy(PetscDrawLG drawlg)
1665: {
1666: PetscDraw draw;
1667: int ierr;
1670: PetscDrawLGGetDraw(drawlg,&draw);
1671: PetscDrawDestroy(draw);
1672: PetscDrawLGDestroy(drawlg);
1673: return(0);
1674: }
1678: /*@
1679: TSGetTime - Gets the current time.
1681: Not Collective
1683: Input Parameter:
1684: . ts - the TS context obtained from TSCreate()
1686: Output Parameter:
1687: . t - the current time
1689: Contributed by: Matthew Knepley
1691: Level: beginner
1693: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1695: .keywords: TS, get, time
1696: @*/
1697: int TSGetTime(TS ts,PetscReal* t)
1698: {
1702: *t = ts->ptime;
1703: return(0);
1704: }
1708: /*@C
1709: TSSetOptionsPrefix - Sets the prefix used for searching for all
1710: TS options in the database.
1712: Collective on TS
1714: Input Parameter:
1715: + ts - The TS context
1716: - prefix - The prefix to prepend to all option names
1718: Notes:
1719: A hyphen (-) must NOT be given at the beginning of the prefix name.
1720: The first character of all runtime options is AUTOMATICALLY the
1721: hyphen.
1723: Contributed by: Matthew Knepley
1725: Level: advanced
1727: .keywords: TS, set, options, prefix, database
1729: .seealso: TSSetFromOptions()
1731: @*/
1732: int TSSetOptionsPrefix(TS ts,const char prefix[])
1733: {
1738: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1739: switch(ts->problem_type) {
1740: case TS_NONLINEAR:
1741: SNESSetOptionsPrefix(ts->snes,prefix);
1742: break;
1743: case TS_LINEAR:
1744: KSPSetOptionsPrefix(ts->ksp,prefix);
1745: break;
1746: }
1747: return(0);
1748: }
1753: /*@C
1754: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1755: TS options in the database.
1757: Collective on TS
1759: Input Parameter:
1760: + ts - The TS context
1761: - prefix - The prefix to prepend to all option names
1763: Notes:
1764: A hyphen (-) must NOT be given at the beginning of the prefix name.
1765: The first character of all runtime options is AUTOMATICALLY the
1766: hyphen.
1768: Contributed by: Matthew Knepley
1770: Level: advanced
1772: .keywords: TS, append, options, prefix, database
1774: .seealso: TSGetOptionsPrefix()
1776: @*/
1777: int TSAppendOptionsPrefix(TS ts,const char prefix[])
1778: {
1783: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1784: switch(ts->problem_type) {
1785: case TS_NONLINEAR:
1786: SNESAppendOptionsPrefix(ts->snes,prefix);
1787: break;
1788: case TS_LINEAR:
1789: KSPAppendOptionsPrefix(ts->ksp,prefix);
1790: break;
1791: }
1792: return(0);
1793: }
1797: /*@C
1798: TSGetOptionsPrefix - Sets the prefix used for searching for all
1799: TS options in the database.
1801: Not Collective
1803: Input Parameter:
1804: . ts - The TS context
1806: Output Parameter:
1807: . prefix - A pointer to the prefix string used
1809: Contributed by: Matthew Knepley
1811: Notes: On the fortran side, the user should pass in a string 'prifix' of
1812: sufficient length to hold the prefix.
1814: Level: intermediate
1816: .keywords: TS, get, options, prefix, database
1818: .seealso: TSAppendOptionsPrefix()
1819: @*/
1820: int TSGetOptionsPrefix(TS ts,char *prefix[])
1821: {
1827: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1828: return(0);
1829: }
1833: /*@C
1834: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1836: Not Collective, but parallel objects are returned if TS is parallel
1838: Input Parameter:
1839: . ts - The TS context obtained from TSCreate()
1841: Output Parameters:
1842: + A - The matrix A, where U_t = A(t) U
1843: . M - The preconditioner matrix, usually the same as A
1844: - ctx - User-defined context for matrix evaluation routine
1846: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1848: Contributed by: Matthew Knepley
1850: Level: intermediate
1852: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1854: .keywords: TS, timestep, get, matrix
1856: @*/
1857: int TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1858: {
1861: if (A) *A = ts->A;
1862: if (M) *M = ts->B;
1863: if (ctx) *ctx = ts->jacP;
1864: return(0);
1865: }
1869: /*@C
1870: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1872: Not Collective, but parallel objects are returned if TS is parallel
1874: Input Parameter:
1875: . ts - The TS context obtained from TSCreate()
1877: Output Parameters:
1878: + J - The Jacobian J of F, where U_t = F(U,t)
1879: . M - The preconditioner matrix, usually the same as J
1880: - ctx - User-defined context for Jacobian evaluation routine
1882: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1884: Contributed by: Matthew Knepley
1886: Level: intermediate
1888: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1890: .keywords: TS, timestep, get, matrix, Jacobian
1891: @*/
1892: int TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1893: {
1897: TSGetRHSMatrix(ts,J,M,ctx);
1898: return(0);
1899: }
1903: /*@C
1904: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1905: VecView() for the solution at each timestep
1907: Collective on TS
1909: Input Parameters:
1910: + ts - the TS context
1911: . step - current time-step
1912: . ptime - current time
1913: - dummy - either a viewer or PETSC_NULL
1915: Level: intermediate
1917: .keywords: TS, vector, monitor, view
1919: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1920: @*/
1921: int TSVecViewMonitor(TS ts,int step,PetscReal ptime,Vec x,void *dummy)
1922: {
1923: int ierr;
1924: PetscViewer viewer = (PetscViewer) dummy;
1927: if (!viewer) {
1928: MPI_Comm comm;
1929: PetscObjectGetComm((PetscObject)ts,&comm);
1930: viewer = PETSC_VIEWER_DRAW_(comm);
1931: }
1932: VecView(x,viewer);
1933: return(0);
1934: }