Actual source code: ts.c
1: #define PETSCTS_DLL
3: #include private/tsimpl.h
5: /* Logging support */
6: PetscCookie TS_COOKIE;
7: PetscLogEvent TS_Step, TS_PseudoComputeTimeStep, TS_FunctionEval, TS_JacobianEval;
11: /*
12: TSSetTypeFromOptions - Sets the type of ts from user options.
14: Collective on TS
16: Input Parameter:
17: . ts - The ts
19: Level: intermediate
21: .keywords: TS, set, options, database, type
22: .seealso: TSSetFromOptions(), TSSetType()
23: */
24: static PetscErrorCode TSSetTypeFromOptions(TS ts)
25: {
26: PetscTruth opt;
27: const char *defaultType;
28: char typeName[256];
32: if (((PetscObject)ts)->type_name) {
33: defaultType = ((PetscObject)ts)->type_name;
34: } else {
35: defaultType = TS_EULER;
36: }
38: if (!TSRegisterAllCalled) {TSRegisterAll(PETSC_NULL);}
39: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
40: if (opt) {
41: TSSetType(ts, typeName);
42: } else {
43: TSSetType(ts, defaultType);
44: }
45: return(0);
46: }
50: /*@
51: TSSetFromOptions - Sets various TS parameters from user options.
53: Collective on TS
55: Input Parameter:
56: . ts - the TS context obtained from TSCreate()
58: Options Database Keys:
59: + -ts_type <type> - TS_EULER, TS_BEULER, TS_SUNDIALS, TS_PSEUDO, TS_CRANK_NICHOLSON
60: . -ts_max_steps maxsteps - maximum number of time-steps to take
61: . -ts_max_time time - maximum time to compute to
62: . -ts_dt dt - initial time step
63: . -ts_monitor - print information at each timestep
64: - -ts_monitor_draw - plot information at each timestep
66: Level: beginner
68: .keywords: TS, timestep, set, options, database
70: .seealso: TSGetType()
71: @*/
72: PetscErrorCode TSSetFromOptions(TS ts)
73: {
74: PetscReal dt;
75: PetscTruth opt,flg;
76: PetscErrorCode ierr;
77: PetscViewerASCIIMonitor monviewer;
78: char monfilename[PETSC_MAX_PATH_LEN];
82: PetscOptionsBegin(((PetscObject)ts)->comm, ((PetscObject)ts)->prefix, "Time step options", "TS");
84: /* Handle generic TS options */
85: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
86: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
87: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
88: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
89: if (opt) {
90: ts->initial_time_step = ts->time_step = dt;
91: }
93: /* Monitor options */
94: PetscOptionsString("-ts_monitor","Monitor timestep size","TSMonitorDefault","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
95: if (flg) {
96: PetscViewerASCIIMonitorCreate(((PetscObject)ts)->comm,monfilename,((PetscObject)ts)->tablevel,&monviewer);
97: TSMonitorSet(ts,TSMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
98: }
99: PetscOptionsName("-ts_monitor_draw","Monitor timestep size graphically","TSMonitorLG",&opt);
100: if (opt) {
101: TSMonitorSet(ts,TSMonitorLG,PETSC_NULL,PETSC_NULL);
102: }
103: PetscOptionsName("-ts_monitor_solution","Monitor solution graphically","TSMonitorSolution",&opt);
104: if (opt) {
105: TSMonitorSet(ts,TSMonitorSolution,PETSC_NULL,PETSC_NULL);
106: }
108: /* Handle TS type options */
109: TSSetTypeFromOptions(ts);
111: /* Handle specific TS options */
112: if (ts->ops->setfromoptions) {
113: (*ts->ops->setfromoptions)(ts);
114: }
115: PetscOptionsEnd();
117: /* Handle subobject options */
118: switch(ts->problem_type) {
119: /* Should check for implicit/explicit */
120: case TS_LINEAR:
121: if (ts->ksp) {
122: KSPSetOperators(ts->ksp,ts->Arhs,ts->B,DIFFERENT_NONZERO_PATTERN);
123: KSPSetFromOptions(ts->ksp);
124: }
125: break;
126: case TS_NONLINEAR:
127: if (ts->snes) {
128: /* this is a bit of a hack, but it gets the matrix information into SNES earlier
129: so that SNES and KSP have more information to pick reasonable defaults
130: before they allow users to set options */
131: SNESSetJacobian(ts->snes,ts->Arhs,ts->B,0,ts);
132: SNESSetFromOptions(ts->snes);
133: }
134: break;
135: default:
136: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
137: }
139: return(0);
140: }
142: #undef __FUNCT__
144: /*@
145: TSViewFromOptions - This function visualizes the ts based upon user options.
147: Collective on TS
149: Input Parameter:
150: . ts - The ts
152: Level: intermediate
154: .keywords: TS, view, options, database
155: .seealso: TSSetFromOptions(), TSView()
156: @*/
157: PetscErrorCode TSViewFromOptions(TS ts,const char title[])
158: {
159: PetscViewer viewer;
160: PetscDraw draw;
161: PetscTruth opt;
162: char fileName[PETSC_MAX_PATH_LEN];
166: PetscOptionsGetString(((PetscObject)ts)->prefix, "-ts_view", fileName, PETSC_MAX_PATH_LEN, &opt);
167: if (opt && !PetscPreLoadingOn) {
168: PetscViewerASCIIOpen(((PetscObject)ts)->comm,fileName,&viewer);
169: TSView(ts, viewer);
170: PetscViewerDestroy(viewer);
171: }
172: PetscOptionsHasName(((PetscObject)ts)->prefix, "-ts_view_draw", &opt);
173: if (opt) {
174: PetscViewerDrawOpen(((PetscObject)ts)->comm, 0, 0, 0, 0, 300, 300, &viewer);
175: PetscViewerDrawGetDraw(viewer, 0, &draw);
176: if (title) {
177: PetscDrawSetTitle(draw, (char *)title);
178: } else {
179: PetscObjectName((PetscObject)ts);
180: PetscDrawSetTitle(draw, ((PetscObject)ts)->name);
181: }
182: TSView(ts, viewer);
183: PetscViewerFlush(viewer);
184: PetscDrawPause(draw);
185: PetscViewerDestroy(viewer);
186: }
187: return(0);
188: }
192: /*@
193: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
194: set with TSSetRHSJacobian().
196: Collective on TS and Vec
198: Input Parameters:
199: + ts - the SNES context
200: . t - current timestep
201: - x - input vector
203: Output Parameters:
204: + A - Jacobian matrix
205: . B - optional preconditioning matrix
206: - flag - flag indicating matrix structure
208: Notes:
209: Most users should not need to explicitly call this routine, as it
210: is used internally within the nonlinear solvers.
212: See KSPSetOperators() for important information about setting the
213: flag parameter.
215: TSComputeJacobian() is valid only for TS_NONLINEAR
217: Level: developer
219: .keywords: SNES, compute, Jacobian, matrix
221: .seealso: TSSetRHSJacobian(), KSPSetOperators()
222: @*/
223: PetscErrorCode TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
224: {
231: if (ts->problem_type != TS_NONLINEAR) {
232: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
233: }
234: if (ts->ops->rhsjacobian) {
235: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
236: *flg = DIFFERENT_NONZERO_PATTERN;
237: PetscStackPush("TS user Jacobian function");
238: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
239: PetscStackPop;
240: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
241: /* make sure user returned a correct Jacobian and preconditioner */
244: } else {
245: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
246: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
247: if (*A != *B) {
248: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
249: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
250: }
251: }
252: return(0);
253: }
257: /*
258: TSComputeRHSFunction - Evaluates the right-hand-side function.
260: Note: If the user did not provide a function but merely a matrix,
261: this routine applies the matrix.
262: */
263: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
264: {
272: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
273: if (ts->ops->rhsfunction) {
274: PetscStackPush("TS user right-hand-side function");
275: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
276: PetscStackPop;
277: } else {
278: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
279: MatStructure flg;
280: PetscStackPush("TS user right-hand-side matrix function");
281: (*ts->ops->rhsmatrix)(ts,t,&ts->Arhs,&ts->B,&flg,ts->jacP);
282: PetscStackPop;
283: }
284: MatMult(ts->Arhs,x,y);
285: }
287: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
289: return(0);
290: }
294: /*@C
295: TSSetRHSFunction - Sets the routine for evaluating the function,
296: F(t,u), where U_t = F(t,u).
298: Collective on TS
300: Input Parameters:
301: + ts - the TS context obtained from TSCreate()
302: . f - routine for evaluating the right-hand-side function
303: - ctx - [optional] user-defined context for private data for the
304: function evaluation routine (may be PETSC_NULL)
306: Calling sequence of func:
307: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
309: + t - current timestep
310: . u - input vector
311: . F - function vector
312: - ctx - [optional] user-defined function context
314: Important:
315: The user MUST call either this routine or TSSetMatrices().
317: Level: beginner
319: .keywords: TS, timestep, set, right-hand-side, function
321: .seealso: TSSetMatrices()
322: @*/
323: PetscErrorCode TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
324: {
328: if (ts->problem_type == TS_LINEAR) {
329: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
330: }
331: ts->ops->rhsfunction = f;
332: ts->funP = ctx;
333: return(0);
334: }
338: /*@C
339: TSSetMatrices - Sets the functions to compute the matrices Alhs and Arhs,
340: where Alhs(t) U_t = Arhs(t) U.
342: Collective on TS
344: Input Parameters:
345: + ts - the TS context obtained from TSCreate()
346: . Arhs - matrix
347: . frhs - the matrix evaluation routine for Arhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
348: if Arhs is not a function of t.
349: . Alhs - matrix or PETSC_NULL if Alhs is an indentity matrix.
350: . flhs - the matrix evaluation routine for Alhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
351: if Alhs is not a function of t.
352: . flag - flag indicating information about the matrix structure of Arhs and Alhs.
353: The available options are
354: SAME_NONZERO_PATTERN - Alhs has the same nonzero structure as Arhs
355: DIFFERENT_NONZERO_PATTERN - Alhs has different nonzero structure as Arhs
356: - ctx - [optional] user-defined context for private data for the
357: matrix evaluation routine (may be PETSC_NULL)
359: Calling sequence of func:
360: $ func(TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
362: + t - current timestep
363: . A - matrix A, where U_t = A(t) U
364: . B - preconditioner matrix, usually the same as A
365: . flag - flag indicating information about the preconditioner matrix
366: structure (same as flag in KSPSetOperators())
367: - ctx - [optional] user-defined context for matrix evaluation routine
369: Notes:
370: The routine func() takes Mat* as the matrix arguments rather than Mat.
371: This allows the matrix evaluation routine to replace Arhs or Alhs with a
372: completely new new matrix structure (not just different matrix elements)
373: when appropriate, for instance, if the nonzero structure is changing
374: throughout the global iterations.
376: Important:
377: The user MUST call either this routine or TSSetRHSFunction().
379: Level: beginner
381: .keywords: TS, timestep, set, matrix
383: .seealso: TSSetRHSFunction()
384: @*/
385: PetscErrorCode TSSetMatrices(TS ts,Mat Arhs,PetscErrorCode (*frhs)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),Mat Alhs,PetscErrorCode (*flhs)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),MatStructure flag,void *ctx)
386: {
389: if (Arhs){
392: ts->Arhs = Arhs;
393: ts->ops->rhsmatrix = frhs;
394: }
395: if (Alhs){
398: ts->Alhs = Alhs;
399: ts->ops->lhsmatrix = flhs;
400: }
401:
402: ts->jacP = ctx;
403: ts->matflg = flag;
404: return(0);
405: }
409: /*@C
410: TSGetMatrices - Returns the matrices Arhs and Alhs at the present timestep,
411: where Alhs(t) U_t = Arhs(t) U.
413: Not Collective, but parallel objects are returned if TS is parallel
415: Input Parameter:
416: . ts - The TS context obtained from TSCreate()
418: Output Parameters:
419: + Arhs - The right-hand side matrix
420: . Alhs - The left-hand side matrix
421: - ctx - User-defined context for matrix evaluation routine
423: Notes: You can pass in PETSC_NULL for any return argument you do not need.
425: Level: intermediate
427: .seealso: TSSetMatrices(), TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
429: .keywords: TS, timestep, get, matrix
431: @*/
432: PetscErrorCode TSGetMatrices(TS ts,Mat *Arhs,Mat *Alhs,void **ctx)
433: {
436: if (Arhs) *Arhs = ts->Arhs;
437: if (Alhs) *Alhs = ts->Alhs;
438: if (ctx) *ctx = ts->jacP;
439: return(0);
440: }
444: /*@C
445: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
446: where U_t = F(U,t), as well as the location to store the matrix.
447: Use TSSetMatrices() for linear problems.
449: Collective on TS
451: Input Parameters:
452: + ts - the TS context obtained from TSCreate()
453: . A - Jacobian matrix
454: . B - preconditioner matrix (usually same as A)
455: . f - the Jacobian evaluation routine
456: - ctx - [optional] user-defined context for private data for the
457: Jacobian evaluation routine (may be PETSC_NULL)
459: Calling sequence of func:
460: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
462: + t - current timestep
463: . u - input vector
464: . A - matrix A, where U_t = A(t)u
465: . B - preconditioner matrix, usually the same as A
466: . flag - flag indicating information about the preconditioner matrix
467: structure (same as flag in KSPSetOperators())
468: - ctx - [optional] user-defined context for matrix evaluation routine
470: Notes:
471: See KSPSetOperators() for important information about setting the flag
472: output parameter in the routine func(). Be sure to read this information!
474: The routine func() takes Mat * as the matrix arguments rather than Mat.
475: This allows the matrix evaluation routine to replace A and/or B with a
476: completely new matrix structure (not just different matrix elements)
477: when appropriate, for instance, if the nonzero structure is changing
478: throughout the global iterations.
480: Level: beginner
481:
482: .keywords: TS, timestep, set, right-hand-side, Jacobian
484: .seealso: TSDefaultComputeJacobianColor(),
485: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetMatrices()
487: @*/
488: PetscErrorCode TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
489: {
496: if (ts->problem_type != TS_NONLINEAR) {
497: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetMatrices()");
498: }
500: ts->ops->rhsjacobian = f;
501: ts->jacP = ctx;
502: ts->Arhs = A;
503: ts->B = B;
504: return(0);
505: }
509: /*@C
510: TSView - Prints the TS data structure.
512: Collective on TS
514: Input Parameters:
515: + ts - the TS context obtained from TSCreate()
516: - viewer - visualization context
518: Options Database Key:
519: . -ts_view - calls TSView() at end of TSStep()
521: Notes:
522: The available visualization contexts include
523: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
524: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
525: output where only the first processor opens
526: the file. All other processors send their
527: data to the first processor to print.
529: The user can open an alternative visualization context with
530: PetscViewerASCIIOpen() - output to a specified file.
532: Level: beginner
534: .keywords: TS, timestep, view
536: .seealso: PetscViewerASCIIOpen()
537: @*/
538: PetscErrorCode TSView(TS ts,PetscViewer viewer)
539: {
541: const TSType type;
542: PetscTruth iascii,isstring;
546: if (!viewer) {
547: PetscViewerASCIIGetStdout(((PetscObject)ts)->comm,&viewer);
548: }
552: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
553: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
554: if (iascii) {
555: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
556: TSGetType(ts,&type);
557: if (type) {
558: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
559: } else {
560: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
561: }
562: if (ts->ops->view) {
563: PetscViewerASCIIPushTab(viewer);
564: (*ts->ops->view)(ts,viewer);
565: PetscViewerASCIIPopTab(viewer);
566: }
567: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
568: PetscViewerASCIIPrintf(viewer," maximum time=%G\n",ts->max_time);
569: if (ts->problem_type == TS_NONLINEAR) {
570: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
571: }
572: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
573: } else if (isstring) {
574: TSGetType(ts,&type);
575: PetscViewerStringSPrintf(viewer," %-7.7s",type);
576: }
577: PetscViewerASCIIPushTab(viewer);
578: if (ts->ksp) {KSPView(ts->ksp,viewer);}
579: if (ts->snes) {SNESView(ts->snes,viewer);}
580: PetscViewerASCIIPopTab(viewer);
581: return(0);
582: }
587: /*@C
588: TSSetApplicationContext - Sets an optional user-defined context for
589: the timesteppers.
591: Collective on TS
593: Input Parameters:
594: + ts - the TS context obtained from TSCreate()
595: - usrP - optional user context
597: Level: intermediate
599: .keywords: TS, timestep, set, application, context
601: .seealso: TSGetApplicationContext()
602: @*/
603: PetscErrorCode TSSetApplicationContext(TS ts,void *usrP)
604: {
607: ts->user = usrP;
608: return(0);
609: }
613: /*@C
614: TSGetApplicationContext - Gets the user-defined context for the
615: timestepper.
617: Not Collective
619: Input Parameter:
620: . ts - the TS context obtained from TSCreate()
622: Output Parameter:
623: . usrP - user context
625: Level: intermediate
627: .keywords: TS, timestep, get, application, context
629: .seealso: TSSetApplicationContext()
630: @*/
631: PetscErrorCode TSGetApplicationContext(TS ts,void **usrP)
632: {
635: *usrP = ts->user;
636: return(0);
637: }
641: /*@
642: TSGetTimeStepNumber - Gets the current number of timesteps.
644: Not Collective
646: Input Parameter:
647: . ts - the TS context obtained from TSCreate()
649: Output Parameter:
650: . iter - number steps so far
652: Level: intermediate
654: .keywords: TS, timestep, get, iteration, number
655: @*/
656: PetscErrorCode TSGetTimeStepNumber(TS ts,PetscInt* iter)
657: {
661: *iter = ts->steps;
662: return(0);
663: }
667: /*@
668: TSSetInitialTimeStep - Sets the initial timestep to be used,
669: as well as the initial time.
671: Collective on TS
673: Input Parameters:
674: + ts - the TS context obtained from TSCreate()
675: . initial_time - the initial time
676: - time_step - the size of the timestep
678: Level: intermediate
680: .seealso: TSSetTimeStep(), TSGetTimeStep()
682: .keywords: TS, set, initial, timestep
683: @*/
684: PetscErrorCode TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
685: {
688: ts->time_step = time_step;
689: ts->initial_time_step = time_step;
690: ts->ptime = initial_time;
691: return(0);
692: }
696: /*@
697: TSSetTimeStep - Allows one to reset the timestep at any time,
698: useful for simple pseudo-timestepping codes.
700: Collective on TS
702: Input Parameters:
703: + ts - the TS context obtained from TSCreate()
704: - time_step - the size of the timestep
706: Level: intermediate
708: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
710: .keywords: TS, set, timestep
711: @*/
712: PetscErrorCode TSSetTimeStep(TS ts,PetscReal time_step)
713: {
716: ts->time_step = time_step;
717: return(0);
718: }
722: /*@
723: TSGetTimeStep - Gets the current timestep size.
725: Not Collective
727: Input Parameter:
728: . ts - the TS context obtained from TSCreate()
730: Output Parameter:
731: . dt - the current timestep size
733: Level: intermediate
735: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
737: .keywords: TS, get, timestep
738: @*/
739: PetscErrorCode TSGetTimeStep(TS ts,PetscReal* dt)
740: {
744: *dt = ts->time_step;
745: return(0);
746: }
750: /*@
751: TSGetSolution - Returns the solution at the present timestep. It
752: is valid to call this routine inside the function that you are evaluating
753: in order to move to the new timestep. This vector not changed until
754: the solution at the next timestep has been calculated.
756: Not Collective, but Vec returned is parallel if TS is parallel
758: Input Parameter:
759: . ts - the TS context obtained from TSCreate()
761: Output Parameter:
762: . v - the vector containing the solution
764: Level: intermediate
766: .seealso: TSGetTimeStep()
768: .keywords: TS, timestep, get, solution
769: @*/
770: PetscErrorCode TSGetSolution(TS ts,Vec *v)
771: {
775: *v = ts->vec_sol_always;
776: return(0);
777: }
779: /* ----- Routines to initialize and destroy a timestepper ---- */
782: /*@
783: TSSetProblemType - Sets the type of problem to be solved.
785: Not collective
787: Input Parameters:
788: + ts - The TS
789: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
790: .vb
791: U_t = A U
792: U_t = A(t) U
793: U_t = F(t,U)
794: .ve
796: Level: beginner
798: .keywords: TS, problem type
799: .seealso: TSSetUp(), TSProblemType, TS
800: @*/
801: PetscErrorCode TSSetProblemType(TS ts, TSProblemType type)
802: {
805: ts->problem_type = type;
806: return(0);
807: }
811: /*@C
812: TSGetProblemType - Gets the type of problem to be solved.
814: Not collective
816: Input Parameter:
817: . ts - The TS
819: Output Parameter:
820: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
821: .vb
822: U_t = A U
823: U_t = A(t) U
824: U_t = F(t,U)
825: .ve
827: Level: beginner
829: .keywords: TS, problem type
830: .seealso: TSSetUp(), TSProblemType, TS
831: @*/
832: PetscErrorCode TSGetProblemType(TS ts, TSProblemType *type)
833: {
837: *type = ts->problem_type;
838: return(0);
839: }
843: /*@
844: TSSetUp - Sets up the internal data structures for the later use
845: of a timestepper.
847: Collective on TS
849: Input Parameter:
850: . ts - the TS context obtained from TSCreate()
852: Notes:
853: For basic use of the TS solvers the user need not explicitly call
854: TSSetUp(), since these actions will automatically occur during
855: the call to TSStep(). However, if one wishes to control this
856: phase separately, TSSetUp() should be called after TSCreate()
857: and optional routines of the form TSSetXXX(), but before TSStep().
859: Level: advanced
861: .keywords: TS, timestep, setup
863: .seealso: TSCreate(), TSStep(), TSDestroy()
864: @*/
865: PetscErrorCode TSSetUp(TS ts)
866: {
871: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
872: if (!((PetscObject)ts)->type_name) {
873: TSSetType(ts,TS_EULER);
874: }
875: (*ts->ops->setup)(ts);
876: ts->setupcalled = 1;
877: return(0);
878: }
882: /*@
883: TSDestroy - Destroys the timestepper context that was created
884: with TSCreate().
886: Collective on TS
888: Input Parameter:
889: . ts - the TS context obtained from TSCreate()
891: Level: beginner
893: .keywords: TS, timestepper, destroy
895: .seealso: TSCreate(), TSSetUp(), TSSolve()
896: @*/
897: PetscErrorCode TSDestroy(TS ts)
898: {
903: if (--((PetscObject)ts)->refct > 0) return(0);
905: /* if memory was published with AMS then destroy it */
906: PetscObjectDepublish(ts);
907: if (ts->A) {MatDestroy(ts->A);CHKERRQ(ierr)}
908: if (ts->ksp) {KSPDestroy(ts->ksp);}
909: if (ts->snes) {SNESDestroy(ts->snes);}
910: if (ts->ops->destroy) {(*(ts)->ops->destroy)(ts);}
911: TSMonitorCancel(ts);
912: PetscHeaderDestroy(ts);
913: return(0);
914: }
918: /*@
919: TSGetSNES - Returns the SNES (nonlinear solver) associated with
920: a TS (timestepper) context. Valid only for nonlinear problems.
922: Not Collective, but SNES is parallel if TS is parallel
924: Input Parameter:
925: . ts - the TS context obtained from TSCreate()
927: Output Parameter:
928: . snes - the nonlinear solver context
930: Notes:
931: The user can then directly manipulate the SNES context to set various
932: options, etc. Likewise, the user can then extract and manipulate the
933: KSP, KSP, and PC contexts as well.
935: TSGetSNES() does not work for integrators that do not use SNES; in
936: this case TSGetSNES() returns PETSC_NULL in snes.
938: Level: beginner
940: .keywords: timestep, get, SNES
941: @*/
942: PetscErrorCode TSGetSNES(TS ts,SNES *snes)
943: {
947: if (((PetscObject)ts)->type_name == PETSC_NULL)
948: SETERRQ(PETSC_ERR_ARG_NULL,"SNES is not created yet. Call TSSetType() first");
949: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
950: *snes = ts->snes;
951: return(0);
952: }
956: /*@
957: TSGetKSP - Returns the KSP (linear solver) associated with
958: a TS (timestepper) context.
960: Not Collective, but KSP is parallel if TS is parallel
962: Input Parameter:
963: . ts - the TS context obtained from TSCreate()
965: Output Parameter:
966: . ksp - the nonlinear solver context
968: Notes:
969: The user can then directly manipulate the KSP context to set various
970: options, etc. Likewise, the user can then extract and manipulate the
971: KSP and PC contexts as well.
973: TSGetKSP() does not work for integrators that do not use KSP;
974: in this case TSGetKSP() returns PETSC_NULL in ksp.
976: Level: beginner
978: .keywords: timestep, get, KSP
979: @*/
980: PetscErrorCode TSGetKSP(TS ts,KSP *ksp)
981: {
985: if (((PetscObject)ts)->type_name == PETSC_NULL)
986: SETERRQ(PETSC_ERR_ARG_NULL,"KSP is not created yet. Call TSSetType() first");
987: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
988: *ksp = ts->ksp;
989: return(0);
990: }
992: /* ----------- Routines to set solver parameters ---------- */
996: /*@
997: TSGetDuration - Gets the maximum number of timesteps to use and
998: maximum time for iteration.
1000: Collective on TS
1002: Input Parameters:
1003: + ts - the TS context obtained from TSCreate()
1004: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1005: - maxtime - final time to iterate to, or PETSC_NULL
1007: Level: intermediate
1009: .keywords: TS, timestep, get, maximum, iterations, time
1010: @*/
1011: PetscErrorCode TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1012: {
1015: if (maxsteps) {
1017: *maxsteps = ts->max_steps;
1018: }
1019: if (maxtime ) {
1021: *maxtime = ts->max_time;
1022: }
1023: return(0);
1024: }
1028: /*@
1029: TSSetDuration - Sets the maximum number of timesteps to use and
1030: maximum time for iteration.
1032: Collective on TS
1034: Input Parameters:
1035: + ts - the TS context obtained from TSCreate()
1036: . maxsteps - maximum number of iterations to use
1037: - maxtime - final time to iterate to
1039: Options Database Keys:
1040: . -ts_max_steps <maxsteps> - Sets maxsteps
1041: . -ts_max_time <maxtime> - Sets maxtime
1043: Notes:
1044: The default maximum number of iterations is 5000. Default time is 5.0
1046: Level: intermediate
1048: .keywords: TS, timestep, set, maximum, iterations
1049: @*/
1050: PetscErrorCode TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1051: {
1054: ts->max_steps = maxsteps;
1055: ts->max_time = maxtime;
1056: return(0);
1057: }
1061: /*@
1062: TSSetSolution - Sets the initial solution vector
1063: for use by the TS routines.
1065: Collective on TS and Vec
1067: Input Parameters:
1068: + ts - the TS context obtained from TSCreate()
1069: - x - the solution vector
1071: Level: beginner
1073: .keywords: TS, timestep, set, solution, initial conditions
1074: @*/
1075: PetscErrorCode TSSetSolution(TS ts,Vec x)
1076: {
1080: ts->vec_sol = ts->vec_sol_always = x;
1081: return(0);
1082: }
1086: /*@C
1087: TSSetPreStep - Sets the general-purpose function
1088: called once at the beginning of time stepping.
1090: Collective on TS
1092: Input Parameters:
1093: + ts - The TS context obtained from TSCreate()
1094: - func - The function
1096: Calling sequence of func:
1097: . func (TS ts);
1099: Level: intermediate
1101: .keywords: TS, timestep
1102: @*/
1103: PetscErrorCode TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1104: {
1107: ts->ops->prestep = func;
1108: return(0);
1109: }
1113: /*@
1114: TSDefaultPreStep - The default pre-stepping function which does nothing.
1116: Collective on TS
1118: Input Parameters:
1119: . ts - The TS context obtained from TSCreate()
1121: Level: developer
1123: .keywords: TS, timestep
1124: @*/
1125: PetscErrorCode TSDefaultPreStep(TS ts)
1126: {
1128: return(0);
1129: }
1133: /*@C
1134: TSSetUpdate - Sets the general-purpose update function called
1135: at the beginning of every time step. This function can change
1136: the time step.
1138: Collective on TS
1140: Input Parameters:
1141: + ts - The TS context obtained from TSCreate()
1142: - func - The function
1144: Calling sequence of func:
1145: . func (TS ts, double t, double *dt);
1147: + t - The current time
1148: - dt - The current time step
1150: Level: intermediate
1152: .keywords: TS, update, timestep
1153: @*/
1154: PetscErrorCode TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1155: {
1158: ts->ops->update = func;
1159: return(0);
1160: }
1164: /*@C
1165: TSSetPostUpdate - Sets the general-purpose update function called
1166: after every time step.
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, double t, double *dt);
1177: + t - The current time
1178: - dt - The current time step
1180: Level: intermediate
1182: .keywords: TS, update, timestep
1183: @*/
1184: PetscErrorCode TSSetPostUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1185: {
1188: ts->ops->postupdate = func;
1189: return(0);
1190: }
1194: /*@
1195: TSDefaultUpdate - The default update function which does nothing.
1197: Collective on TS
1199: Input Parameters:
1200: + ts - The TS context obtained from TSCreate()
1201: - t - The current time
1203: Output Parameters:
1204: . dt - The current time step
1206: Level: developer
1208: .keywords: TS, update, timestep
1209: @*/
1210: PetscErrorCode TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1211: {
1213: return(0);
1214: }
1218: /*@C
1219: TSSetPostStep - Sets the general-purpose function
1220: called once at the end of time stepping.
1222: Collective on TS
1224: Input Parameters:
1225: + ts - The TS context obtained from TSCreate()
1226: - func - The function
1228: Calling sequence of func:
1229: . func (TS ts);
1231: Level: intermediate
1233: .keywords: TS, timestep
1234: @*/
1235: PetscErrorCode TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1236: {
1239: ts->ops->poststep = func;
1240: return(0);
1241: }
1245: /*@
1246: TSDefaultPostStep - The default post-stepping function which does nothing.
1248: Collective on TS
1250: Input Parameters:
1251: . ts - The TS context obtained from TSCreate()
1253: Level: developer
1255: .keywords: TS, timestep
1256: @*/
1257: PetscErrorCode TSDefaultPostStep(TS ts)
1258: {
1260: return(0);
1261: }
1263: /* ------------ Routines to set performance monitoring options ----------- */
1267: /*@C
1268: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
1269: timestep to display the iteration's progress.
1271: Collective on TS
1273: Input Parameters:
1274: + ts - the TS context obtained from TSCreate()
1275: . func - monitoring routine
1276: . mctx - [optional] user-defined context for private data for the
1277: monitor routine (use PETSC_NULL if no context is desired)
1278: - monitordestroy - [optional] routine that frees monitor context
1279: (may be PETSC_NULL)
1281: Calling sequence of func:
1282: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1284: + ts - the TS context
1285: . steps - iteration number
1286: . time - current time
1287: . x - current iterate
1288: - mctx - [optional] monitoring context
1290: Notes:
1291: This routine adds an additional monitor to the list of monitors that
1292: already has been loaded.
1294: Fortran notes: Only a single monitor function can be set for each TS object
1296: Level: intermediate
1298: .keywords: TS, timestep, set, monitor
1300: .seealso: TSMonitorDefault(), TSMonitorCancel()
1301: @*/
1302: PetscErrorCode TSMonitorSet(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1303: {
1306: if (ts->numbermonitors >= MAXTSMONITORS) {
1307: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1308: }
1309: ts->monitor[ts->numbermonitors] = monitor;
1310: ts->mdestroy[ts->numbermonitors] = mdestroy;
1311: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1312: return(0);
1313: }
1317: /*@C
1318: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
1320: Collective on TS
1322: Input Parameters:
1323: . ts - the TS context obtained from TSCreate()
1325: Notes:
1326: There is no way to remove a single, specific monitor.
1328: Level: intermediate
1330: .keywords: TS, timestep, set, monitor
1332: .seealso: TSMonitorDefault(), TSMonitorSet()
1333: @*/
1334: PetscErrorCode TSMonitorCancel(TS ts)
1335: {
1337: PetscInt i;
1341: for (i=0; i<ts->numbermonitors; i++) {
1342: if (ts->mdestroy[i]) {
1343: (*ts->mdestroy[i])(ts->monitorcontext[i]);
1344: }
1345: }
1346: ts->numbermonitors = 0;
1347: return(0);
1348: }
1352: /*@
1353: TSMonitorDefault - Sets the Default monitor
1355: Level: intermediate
1357: .keywords: TS, set, monitor
1359: .seealso: TSMonitorDefault(), TSMonitorSet()
1360: @*/
1361: PetscErrorCode TSMonitorDefault(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1362: {
1363: PetscErrorCode ierr;
1364: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor)ctx;
1367: if (!ctx) {
1368: PetscViewerASCIIMonitorCreate(((PetscObject)ts)->comm,"stdout",0,&viewer);
1369: }
1370: PetscViewerASCIIMonitorPrintf(viewer,"timestep %D dt %G time %G\n",step,ts->time_step,ptime);
1371: if (!ctx) {
1372: PetscViewerASCIIMonitorDestroy(viewer);
1373: }
1374: return(0);
1375: }
1379: /*@
1380: TSStep - Steps the requested number of timesteps.
1382: Collective on TS
1384: Input Parameter:
1385: . ts - the TS context obtained from TSCreate()
1387: Output Parameters:
1388: + steps - number of iterations until termination
1389: - ptime - time until termination
1391: Level: beginner
1393: .keywords: TS, timestep, solve
1395: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1396: @*/
1397: PetscErrorCode TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1398: {
1403: if (!ts->setupcalled) {
1404: TSSetUp(ts);
1405: }
1407: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1408: (*ts->ops->prestep)(ts);
1409: (*ts->ops->step)(ts, steps, ptime);
1410: (*ts->ops->poststep)(ts);
1411: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1413: if (!PetscPreLoadingOn) {
1414: TSViewFromOptions(ts,((PetscObject)ts)->name);
1415: }
1416: return(0);
1417: }
1421: /*@
1422: TSSolve - Steps the requested number of timesteps.
1424: Collective on TS
1426: Input Parameter:
1427: + ts - the TS context obtained from TSCreate()
1428: - x - the solution vector, or PETSC_NULL if it was set with TSSetSolution()
1430: Level: beginner
1432: .keywords: TS, timestep, solve
1434: .seealso: TSCreate(), TSSetSolution(), TSStep()
1435: @*/
1436: PetscErrorCode TSSolve(TS ts, Vec x)
1437: {
1438: PetscInt steps;
1439: PetscReal ptime;
1443: /* set solution vector if provided */
1444: if (x) { TSSetSolution(ts, x); }
1445: /* reset time step and iteration counters */
1446: ts->steps = 0; ts->linear_its = 0; ts->nonlinear_its = 0;
1447: /* steps the requested number of timesteps. */
1448: TSStep(ts, &steps, &ptime);
1449: return(0);
1450: }
1454: /*
1455: Runs the user provided monitor routines, if they exists.
1456: */
1457: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1458: {
1460: PetscInt i,n = ts->numbermonitors;
1463: for (i=0; i<n; i++) {
1464: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1465: }
1466: return(0);
1467: }
1469: /* ------------------------------------------------------------------------*/
1473: /*@C
1474: TSMonitorLGCreate - Creates a line graph context for use with
1475: TS to monitor convergence of preconditioned residual norms.
1477: Collective on TS
1479: Input Parameters:
1480: + host - the X display to open, or null for the local machine
1481: . label - the title to put in the title bar
1482: . x, y - the screen coordinates of the upper left coordinate of the window
1483: - m, n - the screen width and height in pixels
1485: Output Parameter:
1486: . draw - the drawing context
1488: Options Database Key:
1489: . -ts_monitor_draw - automatically sets line graph monitor
1491: Notes:
1492: Use TSMonitorLGDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1494: Level: intermediate
1496: .keywords: TS, monitor, line graph, residual, seealso
1498: .seealso: TSMonitorLGDestroy(), TSMonitorSet()
1500: @*/
1501: PetscErrorCode TSMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1502: {
1503: PetscDraw win;
1507: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1508: PetscDrawSetType(win,PETSC_DRAW_X);
1509: PetscDrawLGCreate(win,1,draw);
1510: PetscDrawLGIndicateDataPoints(*draw);
1512: PetscLogObjectParent(*draw,win);
1513: return(0);
1514: }
1518: PetscErrorCode TSMonitorLG(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1519: {
1520: PetscDrawLG lg = (PetscDrawLG) monctx;
1521: PetscReal x,y = ptime;
1525: if (!monctx) {
1526: MPI_Comm comm;
1527: PetscViewer viewer;
1529: PetscObjectGetComm((PetscObject)ts,&comm);
1530: viewer = PETSC_VIEWER_DRAW_(comm);
1531: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1532: }
1534: if (!n) {PetscDrawLGReset(lg);}
1535: x = (PetscReal)n;
1536: PetscDrawLGAddPoint(lg,&x,&y);
1537: if (n < 20 || (n % 5)) {
1538: PetscDrawLGDraw(lg);
1539: }
1540: return(0);
1541: }
1545: /*@C
1546: TSMonitorLGDestroy - Destroys a line graph context that was created
1547: with TSMonitorLGCreate().
1549: Collective on PetscDrawLG
1551: Input Parameter:
1552: . draw - the drawing context
1554: Level: intermediate
1556: .keywords: TS, monitor, line graph, destroy
1558: .seealso: TSMonitorLGCreate(), TSMonitorSet(), TSMonitorLG();
1559: @*/
1560: PetscErrorCode TSMonitorLGDestroy(PetscDrawLG drawlg)
1561: {
1562: PetscDraw draw;
1566: PetscDrawLGGetDraw(drawlg,&draw);
1567: PetscDrawDestroy(draw);
1568: PetscDrawLGDestroy(drawlg);
1569: return(0);
1570: }
1574: /*@
1575: TSGetTime - Gets the current time.
1577: Not Collective
1579: Input Parameter:
1580: . ts - the TS context obtained from TSCreate()
1582: Output Parameter:
1583: . t - the current time
1585: Contributed by: Matthew Knepley
1587: Level: beginner
1589: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1591: .keywords: TS, get, time
1592: @*/
1593: PetscErrorCode TSGetTime(TS ts,PetscReal* t)
1594: {
1598: *t = ts->ptime;
1599: return(0);
1600: }
1604: /*@
1605: TSSetTime - Allows one to reset the time.
1607: Collective on TS
1609: Input Parameters:
1610: + ts - the TS context obtained from TSCreate()
1611: - time - the time
1613: Level: intermediate
1615: .seealso: TSGetTime(), TSSetDuration()
1617: .keywords: TS, set, time
1618: @*/
1619: PetscErrorCode TSSetTime(TS ts, PetscReal t)
1620: {
1623: ts->ptime = t;
1624: return(0);
1625: }
1629: /*@C
1630: TSSetOptionsPrefix - Sets the prefix used for searching for all
1631: TS options in the database.
1633: Collective on TS
1635: Input Parameter:
1636: + ts - The TS context
1637: - prefix - The prefix to prepend to all option names
1639: Notes:
1640: A hyphen (-) must NOT be given at the beginning of the prefix name.
1641: The first character of all runtime options is AUTOMATICALLY the
1642: hyphen.
1644: Contributed by: Matthew Knepley
1646: Level: advanced
1648: .keywords: TS, set, options, prefix, database
1650: .seealso: TSSetFromOptions()
1652: @*/
1653: PetscErrorCode TSSetOptionsPrefix(TS ts,const char prefix[])
1654: {
1659: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1660: switch(ts->problem_type) {
1661: case TS_NONLINEAR:
1662: if (ts->snes) {
1663: SNESSetOptionsPrefix(ts->snes,prefix);
1664: }
1665: break;
1666: case TS_LINEAR:
1667: if (ts->ksp) {
1668: KSPSetOptionsPrefix(ts->ksp,prefix);
1669: }
1670: break;
1671: }
1672: return(0);
1673: }
1678: /*@C
1679: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1680: TS options in the database.
1682: Collective on TS
1684: Input Parameter:
1685: + ts - The TS context
1686: - prefix - The prefix to prepend to all option names
1688: Notes:
1689: A hyphen (-) must NOT be given at the beginning of the prefix name.
1690: The first character of all runtime options is AUTOMATICALLY the
1691: hyphen.
1693: Contributed by: Matthew Knepley
1695: Level: advanced
1697: .keywords: TS, append, options, prefix, database
1699: .seealso: TSGetOptionsPrefix()
1701: @*/
1702: PetscErrorCode TSAppendOptionsPrefix(TS ts,const char prefix[])
1703: {
1708: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1709: switch(ts->problem_type) {
1710: case TS_NONLINEAR:
1711: if (ts->snes) {
1712: SNESAppendOptionsPrefix(ts->snes,prefix);
1713: }
1714: break;
1715: case TS_LINEAR:
1716: if (ts->ksp) {
1717: KSPAppendOptionsPrefix(ts->ksp,prefix);
1718: }
1719: break;
1720: }
1721: return(0);
1722: }
1726: /*@C
1727: TSGetOptionsPrefix - Sets the prefix used for searching for all
1728: TS options in the database.
1730: Not Collective
1732: Input Parameter:
1733: . ts - The TS context
1735: Output Parameter:
1736: . prefix - A pointer to the prefix string used
1738: Contributed by: Matthew Knepley
1740: Notes: On the fortran side, the user should pass in a string 'prifix' of
1741: sufficient length to hold the prefix.
1743: Level: intermediate
1745: .keywords: TS, get, options, prefix, database
1747: .seealso: TSAppendOptionsPrefix()
1748: @*/
1749: PetscErrorCode TSGetOptionsPrefix(TS ts,const char *prefix[])
1750: {
1756: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1757: return(0);
1758: }
1762: /*@C
1763: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1765: Not Collective, but parallel objects are returned if TS is parallel
1767: Input Parameter:
1768: . ts - The TS context obtained from TSCreate()
1770: Output Parameters:
1771: + J - The Jacobian J of F, where U_t = F(U,t)
1772: . M - The preconditioner matrix, usually the same as J
1773: - ctx - User-defined context for Jacobian evaluation routine
1775: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1777: Level: intermediate
1779: .seealso: TSGetTimeStep(), TSGetMatrices(), TSGetTime(), TSGetTimeStepNumber()
1781: .keywords: TS, timestep, get, matrix, Jacobian
1782: @*/
1783: PetscErrorCode TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1784: {
1786: if (J) *J = ts->Arhs;
1787: if (M) *M = ts->B;
1788: if (ctx) *ctx = ts->jacP;
1789: return(0);
1790: }
1794: /*@C
1795: TSMonitorSolution - Monitors progress of the TS solvers by calling
1796: VecView() for the solution at each timestep
1798: Collective on TS
1800: Input Parameters:
1801: + ts - the TS context
1802: . step - current time-step
1803: . ptime - current time
1804: - dummy - either a viewer or PETSC_NULL
1806: Level: intermediate
1808: .keywords: TS, vector, monitor, view
1810: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
1811: @*/
1812: PetscErrorCode TSMonitorSolution(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1813: {
1815: PetscViewer viewer = (PetscViewer) dummy;
1818: if (!dummy) {
1819: viewer = PETSC_VIEWER_DRAW_(((PetscObject)ts)->comm);
1820: }
1821: VecView(x,viewer);
1822: return(0);
1823: }