Actual source code: ts.c
1: #define PETSCTS_DLL
3: #include src/ts/tsimpl.h
5: /* Logging support */
6: PetscCookie PETSCTS_DLLEXPORT TS_COOKIE = 0;
7: PetscEvent TS_Step = 0, TS_PseudoComputeTimeStep = 0, TS_FunctionEval = 0, TS_JacobianEval = 0;
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 (ts->type_name) {
33: defaultType = ts->type_name;
34: } else {
35: defaultType = TS_EULER;
36: }
38: if (!TSRegisterAllCalled) {
39: TSRegisterAll(PETSC_NULL);
40: }
41: PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
42: if (opt) {
43: TSSetType(ts, typeName);
44: } else {
45: TSSetType(ts, defaultType);
46: }
47: return(0);
48: }
52: /*@
53: TSSetFromOptions - Sets various TS parameters from user options.
55: Collective on TS
57: Input Parameter:
58: . ts - the TS context obtained from TSCreate()
60: Options Database Keys:
61: + -ts_type <type> - TS_EULER, TS_BEULER, TS_SUNDIALS, TS_PSEUDO, TS_CRANK_NICHOLSON
62: . -ts_max_steps maxsteps - maximum number of time-steps to take
63: . -ts_max_time time - maximum time to compute to
64: . -ts_dt dt - initial time step
65: . -ts_monitor - print information at each timestep
66: - -ts_xmonitor - plot information at each timestep
68: Level: beginner
70: .keywords: TS, timestep, set, options, database
72: .seealso: TSGetType
73: @*/
74: PetscErrorCode PETSCTS_DLLEXPORT TSSetFromOptions(TS ts)
75: {
76: PetscReal dt;
77: PetscTruth opt,flg;
79: PetscViewer monviewer;
80: char monfilename[PETSC_MAX_PATH_LEN];
84: PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");
86: /* Handle generic TS options */
87: PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
88: PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
89: PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
90: PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
91: if (opt) {
92: ts->initial_time_step = ts->time_step = dt;
93: }
95: /* Monitor options */
96: PetscOptionsString("-ts_monitor","Monitor timestep size","TSDefaultMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
97: if (flg) {
98: PetscViewerASCIIOpen(ts->comm,monfilename,&monviewer);
99: TSSetMonitor(ts,TSDefaultMonitor,monviewer,(PetscErrorCode (*)(void*))PetscViewerDestroy);
100: }
101: PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
102: if (opt) {
103: TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
104: }
105: PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
106: if (opt) {
107: TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
108: }
110: /* Handle TS type options */
111: TSSetTypeFromOptions(ts);
113: /* Handle specific TS options */
114: if (ts->ops->setfromoptions) {
115: (*ts->ops->setfromoptions)(ts);
116: }
117: PetscOptionsEnd();
119: /* Handle subobject options */
120: switch(ts->problem_type) {
121: /* Should check for implicit/explicit */
122: case TS_LINEAR:
123: if (ts->ksp) {
124: KSPSetOperators(ts->ksp,ts->A,ts->B,DIFFERENT_NONZERO_PATTERN);
125: KSPSetFromOptions(ts->ksp);
126: }
127: break;
128: case TS_NONLINEAR:
129: if (ts->snes) {
130: /* this is a bit of a hack, but it gets the matrix information into SNES earlier
131: so that SNES and KSP have more information to pick reasonable defaults
132: before they allow users to set options */
133: SNESSetJacobian(ts->snes,ts->A,ts->B,0,ts);
134: SNESSetFromOptions(ts->snes);
135: }
136: break;
137: default:
138: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
139: }
141: return(0);
142: }
144: #undef __FUNCT__
146: /*@
147: TSViewFromOptions - This function visualizes the ts based upon user options.
149: Collective on TS
151: Input Parameter:
152: . ts - The ts
154: Level: intermediate
156: .keywords: TS, view, options, database
157: .seealso: TSSetFromOptions(), TSView()
158: @*/
159: PetscErrorCode PETSCTS_DLLEXPORT TSViewFromOptions(TS ts,const char title[])
160: {
161: PetscViewer viewer;
162: PetscDraw draw;
163: PetscTruth opt;
164: char fileName[PETSC_MAX_PATH_LEN];
168: PetscOptionsGetString(ts->prefix, "-ts_view", fileName, PETSC_MAX_PATH_LEN, &opt);
169: if (opt && !PetscPreLoadingOn) {
170: PetscViewerASCIIOpen(ts->comm,fileName,&viewer);
171: TSView(ts, viewer);
172: PetscViewerDestroy(viewer);
173: }
174: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
175: if (opt) {
176: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
177: PetscViewerDrawGetDraw(viewer, 0, &draw);
178: if (title) {
179: PetscDrawSetTitle(draw, (char *)title);
180: } else {
181: PetscObjectName((PetscObject) ts);
182: PetscDrawSetTitle(draw, ts->name);
183: }
184: TSView(ts, viewer);
185: PetscViewerFlush(viewer);
186: PetscDrawPause(draw);
187: PetscViewerDestroy(viewer);
188: }
189: return(0);
190: }
194: /*@
195: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
196: set with TSSetRHSJacobian().
198: Collective on TS and Vec
200: Input Parameters:
201: + ts - the SNES context
202: . t - current timestep
203: - x - input vector
205: Output Parameters:
206: + A - Jacobian matrix
207: . B - optional preconditioning matrix
208: - flag - flag indicating matrix structure
210: Notes:
211: Most users should not need to explicitly call this routine, as it
212: is used internally within the nonlinear solvers.
214: See KSPSetOperators() for important information about setting the
215: flag parameter.
217: TSComputeJacobian() is valid only for TS_NONLINEAR
219: Level: developer
221: .keywords: SNES, compute, Jacobian, matrix
223: .seealso: TSSetRHSJacobian(), KSPSetOperators()
224: @*/
225: PetscErrorCode PETSCTS_DLLEXPORT TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
226: {
233: if (ts->problem_type != TS_NONLINEAR) {
234: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
235: }
236: if (ts->ops->rhsjacobian) {
237: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
238: *flg = DIFFERENT_NONZERO_PATTERN;
239: PetscStackPush("TS user Jacobian function");
240: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
241: PetscStackPop;
242: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
243: /* make sure user returned a correct Jacobian and preconditioner */
246: } else {
247: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
248: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
249: if (*A != *B) {
250: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
251: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
252: }
253: }
254: return(0);
255: }
259: /*
260: TSComputeRHSFunction - Evaluates the right-hand-side function.
262: Note: If the user did not provide a function but merely a matrix,
263: this routine applies the matrix.
264: */
265: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
266: {
274: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
275: if (ts->ops->rhsfunction) {
276: PetscStackPush("TS user right-hand-side function");
277: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
278: PetscStackPop;
279: } else {
280: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
281: MatStructure flg;
282: PetscStackPush("TS user right-hand-side matrix function");
283: (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
284: PetscStackPop;
285: }
286: MatMult(ts->A,x,y);
287: }
289: /* apply user-provided boundary conditions (only needed if these are time dependent) */
290: TSComputeRHSBoundaryConditions(ts,t,y);
291: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
293: return(0);
294: }
298: /*@C
299: TSSetRHSFunction - Sets the routine for evaluating the function,
300: F(t,u), where U_t = F(t,u).
302: Collective on TS
304: Input Parameters:
305: + ts - the TS context obtained from TSCreate()
306: . f - routine for evaluating the right-hand-side function
307: - ctx - [optional] user-defined context for private data for the
308: function evaluation routine (may be PETSC_NULL)
310: Calling sequence of func:
311: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
313: + t - current timestep
314: . u - input vector
315: . F - function vector
316: - ctx - [optional] user-defined function context
318: Important:
319: The user MUST call either this routine or TSSetRHSMatrix().
321: Level: beginner
323: .keywords: TS, timestep, set, right-hand-side, function
325: .seealso: TSSetRHSMatrix()
326: @*/
327: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
328: {
332: if (ts->problem_type == TS_LINEAR) {
333: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
334: }
335: ts->ops->rhsfunction = f;
336: ts->funP = ctx;
337: return(0);
338: }
342: /*@C
343: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
344: Also sets the location to store A.
346: Collective on TS
348: Input Parameters:
349: + ts - the TS context obtained from TSCreate()
350: . A - matrix
351: . B - preconditioner matrix (usually same as A)
352: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
353: if A is not a function of t.
354: - ctx - [optional] user-defined context for private data for the
355: matrix evaluation routine (may be PETSC_NULL)
357: Calling sequence of func:
358: $ func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
360: + t - current timestep
361: . A - matrix A, where U_t = A(t) U
362: . B - preconditioner matrix, usually the same as A
363: . flag - flag indicating information about the preconditioner matrix
364: structure (same as flag in KSPSetOperators())
365: - ctx - [optional] user-defined context for matrix evaluation routine
367: Notes:
368: See KSPSetOperators() for important information about setting the flag
369: output parameter in the routine func(). Be sure to read this information!
371: The routine func() takes Mat * as the matrix arguments rather than Mat.
372: This allows the matrix evaluation routine to replace A and/or B with a
373: completely new new matrix structure (not just different matrix elements)
374: when appropriate, for instance, if the nonzero structure is changing
375: throughout the global iterations.
377: Important:
378: The user MUST call either this routine or TSSetRHSFunction().
380: Level: beginner
382: .keywords: TS, timestep, set, right-hand-side, matrix
384: .seealso: TSSetRHSFunction()
385: @*/
386: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
387: {
394: if (ts->problem_type == TS_NONLINEAR) {
395: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
396: }
398: ts->ops->rhsmatrix = f;
399: ts->jacP = ctx;
400: ts->A = A;
401: ts->B = B;
403: return(0);
404: }
408: /*@C
409: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
410: where U_t = F(U,t), as well as the location to store the matrix.
411: Use TSSetRHSMatrix() for linear problems.
413: Collective on TS
415: Input Parameters:
416: + ts - the TS context obtained from TSCreate()
417: . A - Jacobian matrix
418: . B - preconditioner matrix (usually same as A)
419: . f - the Jacobian evaluation routine
420: - ctx - [optional] user-defined context for private data for the
421: Jacobian evaluation routine (may be PETSC_NULL)
423: Calling sequence of func:
424: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
426: + t - current timestep
427: . u - input vector
428: . A - matrix A, where U_t = A(t)u
429: . B - preconditioner matrix, usually the same as A
430: . flag - flag indicating information about the preconditioner matrix
431: structure (same as flag in KSPSetOperators())
432: - ctx - [optional] user-defined context for matrix evaluation routine
434: Notes:
435: See KSPSetOperators() for important information about setting the flag
436: output parameter in the routine func(). Be sure to read this information!
438: The routine func() takes Mat * as the matrix arguments rather than Mat.
439: This allows the matrix evaluation routine to replace A and/or B with a
440: completely new new matrix structure (not just different matrix elements)
441: when appropriate, for instance, if the nonzero structure is changing
442: throughout the global iterations.
444: Level: beginner
445:
446: .keywords: TS, timestep, set, right-hand-side, Jacobian
448: .seealso: TSDefaultComputeJacobianColor(),
449: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()
451: @*/
452: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
453: {
460: if (ts->problem_type != TS_NONLINEAR) {
461: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
462: }
464: ts->ops->rhsjacobian = f;
465: ts->jacP = ctx;
466: ts->A = A;
467: ts->B = B;
468: return(0);
469: }
473: /*
474: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
476: Note: If the user did not provide a function but merely a matrix,
477: this routine applies the matrix.
478: */
479: PetscErrorCode TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
480: {
488: if (ts->ops->rhsbc) {
489: PetscStackPush("TS user boundary condition function");
490: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
491: PetscStackPop;
492: return(0);
493: }
495: return(0);
496: }
500: /*@C
501: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
502: boundary conditions for the function F.
504: Collective on TS
506: Input Parameters:
507: + ts - the TS context obtained from TSCreate()
508: . f - routine for evaluating the boundary condition function
509: - ctx - [optional] user-defined context for private data for the
510: function evaluation routine (may be PETSC_NULL)
512: Calling sequence of func:
513: $ func (TS ts,PetscReal t,Vec F,void *ctx);
515: + t - current timestep
516: . F - function vector
517: - ctx - [optional] user-defined function context
519: Level: intermediate
521: .keywords: TS, timestep, set, boundary conditions, function
522: @*/
523: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSBoundaryConditions(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
524: {
528: if (ts->problem_type != TS_LINEAR) {
529: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
530: }
531: ts->ops->rhsbc = f;
532: ts->bcP = ctx;
533: return(0);
534: }
538: /*@C
539: TSView - Prints the TS data structure.
541: Collective on TS
543: Input Parameters:
544: + ts - the TS context obtained from TSCreate()
545: - viewer - visualization context
547: Options Database Key:
548: . -ts_view - calls TSView() at end of TSStep()
550: Notes:
551: The available visualization contexts include
552: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
553: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
554: output where only the first processor opens
555: the file. All other processors send their
556: data to the first processor to print.
558: The user can open an alternative visualization context with
559: PetscViewerASCIIOpen() - output to a specified file.
561: Level: beginner
563: .keywords: TS, timestep, view
565: .seealso: PetscViewerASCIIOpen()
566: @*/
567: PetscErrorCode PETSCTS_DLLEXPORT TSView(TS ts,PetscViewer viewer)
568: {
570: char *type;
571: PetscTruth iascii,isstring;
575: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
579: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
580: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
581: if (iascii) {
582: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
583: TSGetType(ts,(TSType *)&type);
584: if (type) {
585: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
586: } else {
587: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
588: }
589: if (ts->ops->view) {
590: PetscViewerASCIIPushTab(viewer);
591: (*ts->ops->view)(ts,viewer);
592: PetscViewerASCIIPopTab(viewer);
593: }
594: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
595: PetscViewerASCIIPrintf(viewer," maximum time=%G\n",ts->max_time);
596: if (ts->problem_type == TS_NONLINEAR) {
597: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
598: }
599: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
600: } else if (isstring) {
601: TSGetType(ts,(TSType *)&type);
602: PetscViewerStringSPrintf(viewer," %-7.7s",type);
603: }
604: PetscViewerASCIIPushTab(viewer);
605: if (ts->ksp) {KSPView(ts->ksp,viewer);}
606: if (ts->snes) {SNESView(ts->snes,viewer);}
607: PetscViewerASCIIPopTab(viewer);
608: return(0);
609: }
614: /*@C
615: TSSetApplicationContext - Sets an optional user-defined context for
616: the timesteppers.
618: Collective on TS
620: Input Parameters:
621: + ts - the TS context obtained from TSCreate()
622: - usrP - optional user context
624: Level: intermediate
626: .keywords: TS, timestep, set, application, context
628: .seealso: TSGetApplicationContext()
629: @*/
630: PetscErrorCode PETSCTS_DLLEXPORT TSSetApplicationContext(TS ts,void *usrP)
631: {
634: ts->user = usrP;
635: return(0);
636: }
640: /*@C
641: TSGetApplicationContext - Gets the user-defined context for the
642: timestepper.
644: Not Collective
646: Input Parameter:
647: . ts - the TS context obtained from TSCreate()
649: Output Parameter:
650: . usrP - user context
652: Level: intermediate
654: .keywords: TS, timestep, get, application, context
656: .seealso: TSSetApplicationContext()
657: @*/
658: PetscErrorCode PETSCTS_DLLEXPORT TSGetApplicationContext(TS ts,void **usrP)
659: {
662: *usrP = ts->user;
663: return(0);
664: }
668: /*@
669: TSGetTimeStepNumber - Gets the current number of timesteps.
671: Not Collective
673: Input Parameter:
674: . ts - the TS context obtained from TSCreate()
676: Output Parameter:
677: . iter - number steps so far
679: Level: intermediate
681: .keywords: TS, timestep, get, iteration, number
682: @*/
683: PetscErrorCode PETSCTS_DLLEXPORT TSGetTimeStepNumber(TS ts,PetscInt* iter)
684: {
688: *iter = ts->steps;
689: return(0);
690: }
694: /*@
695: TSSetInitialTimeStep - Sets the initial timestep to be used,
696: as well as the initial time.
698: Collective on TS
700: Input Parameters:
701: + ts - the TS context obtained from TSCreate()
702: . initial_time - the initial time
703: - time_step - the size of the timestep
705: Level: intermediate
707: .seealso: TSSetTimeStep(), TSGetTimeStep()
709: .keywords: TS, set, initial, timestep
710: @*/
711: PetscErrorCode PETSCTS_DLLEXPORT TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
712: {
715: ts->time_step = time_step;
716: ts->initial_time_step = time_step;
717: ts->ptime = initial_time;
718: return(0);
719: }
723: /*@
724: TSSetTimeStep - Allows one to reset the timestep at any time,
725: useful for simple pseudo-timestepping codes.
727: Collective on TS
729: Input Parameters:
730: + ts - the TS context obtained from TSCreate()
731: - time_step - the size of the timestep
733: Level: intermediate
735: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
737: .keywords: TS, set, timestep
738: @*/
739: PetscErrorCode PETSCTS_DLLEXPORT TSSetTimeStep(TS ts,PetscReal time_step)
740: {
743: ts->time_step = time_step;
744: return(0);
745: }
749: /*@
750: TSGetTimeStep - Gets the current timestep size.
752: Not Collective
754: Input Parameter:
755: . ts - the TS context obtained from TSCreate()
757: Output Parameter:
758: . dt - the current timestep size
760: Level: intermediate
762: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
764: .keywords: TS, get, timestep
765: @*/
766: PetscErrorCode PETSCTS_DLLEXPORT TSGetTimeStep(TS ts,PetscReal* dt)
767: {
771: *dt = ts->time_step;
772: return(0);
773: }
777: /*@
778: TSGetSolution - Returns the solution at the present timestep. It
779: is valid to call this routine inside the function that you are evaluating
780: in order to move to the new timestep. This vector not changed until
781: the solution at the next timestep has been calculated.
783: Not Collective, but Vec returned is parallel if TS is parallel
785: Input Parameter:
786: . ts - the TS context obtained from TSCreate()
788: Output Parameter:
789: . v - the vector containing the solution
791: Level: intermediate
793: .seealso: TSGetTimeStep()
795: .keywords: TS, timestep, get, solution
796: @*/
797: PetscErrorCode PETSCTS_DLLEXPORT TSGetSolution(TS ts,Vec *v)
798: {
802: *v = ts->vec_sol_always;
803: return(0);
804: }
806: /* ----- Routines to initialize and destroy a timestepper ---- */
809: /*@
810: TSSetProblemType - Sets the type of problem to be solved.
812: Not collective
814: Input Parameters:
815: + ts - The TS
816: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
817: .vb
818: U_t = A U
819: U_t = A(t) U
820: U_t = F(t,U)
821: .ve
823: Level: beginner
825: .keywords: TS, problem type
826: .seealso: TSSetUp(), TSProblemType, TS
827: @*/
828: PetscErrorCode PETSCTS_DLLEXPORT TSSetProblemType(TS ts, TSProblemType type)
829: {
832: ts->problem_type = type;
833: return(0);
834: }
838: /*@C
839: TSGetProblemType - Gets the type of problem to be solved.
841: Not collective
843: Input Parameter:
844: . ts - The TS
846: Output Parameter:
847: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
848: .vb
849: U_t = A U
850: U_t = A(t) U
851: U_t = F(t,U)
852: .ve
854: Level: beginner
856: .keywords: TS, problem type
857: .seealso: TSSetUp(), TSProblemType, TS
858: @*/
859: PetscErrorCode PETSCTS_DLLEXPORT TSGetProblemType(TS ts, TSProblemType *type)
860: {
864: *type = ts->problem_type;
865: return(0);
866: }
870: /*@
871: TSSetUp - Sets up the internal data structures for the later use
872: of a timestepper.
874: Collective on TS
876: Input Parameter:
877: . ts - the TS context obtained from TSCreate()
879: Notes:
880: For basic use of the TS solvers the user need not explicitly call
881: TSSetUp(), since these actions will automatically occur during
882: the call to TSStep(). However, if one wishes to control this
883: phase separately, TSSetUp() should be called after TSCreate()
884: and optional routines of the form TSSetXXX(), but before TSStep().
886: Level: advanced
888: .keywords: TS, timestep, setup
890: .seealso: TSCreate(), TSStep(), TSDestroy()
891: @*/
892: PetscErrorCode PETSCTS_DLLEXPORT TSSetUp(TS ts)
893: {
898: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
899: if (!ts->type_name) {
900: TSSetType(ts,TS_EULER);
901: }
902: (*ts->ops->setup)(ts);
903: ts->setupcalled = 1;
904: return(0);
905: }
909: /*@
910: TSDestroy - Destroys the timestepper context that was created
911: with TSCreate().
913: Collective on TS
915: Input Parameter:
916: . ts - the TS context obtained from TSCreate()
918: Level: beginner
920: .keywords: TS, timestepper, destroy
922: .seealso: TSCreate(), TSSetUp(), TSSolve()
923: @*/
924: PetscErrorCode PETSCTS_DLLEXPORT TSDestroy(TS ts)
925: {
930: if (--ts->refct > 0) return(0);
932: /* if memory was published with AMS then destroy it */
933: PetscObjectDepublish(ts);
935: if (ts->ksp) {KSPDestroy(ts->ksp);}
936: if (ts->snes) {SNESDestroy(ts->snes);}
937: (*(ts)->ops->destroy)(ts);
938: TSClearMonitor(ts);
939: PetscHeaderDestroy(ts);
940: return(0);
941: }
945: /*@
946: TSGetSNES - Returns the SNES (nonlinear solver) associated with
947: a TS (timestepper) context. Valid only for nonlinear problems.
949: Not Collective, but SNES is parallel if TS is parallel
951: Input Parameter:
952: . ts - the TS context obtained from TSCreate()
954: Output Parameter:
955: . snes - the nonlinear solver context
957: Notes:
958: The user can then directly manipulate the SNES context to set various
959: options, etc. Likewise, the user can then extract and manipulate the
960: KSP, KSP, and PC contexts as well.
962: TSGetSNES() does not work for integrators that do not use SNES; in
963: this case TSGetSNES() returns PETSC_NULL in snes.
965: Level: beginner
967: .keywords: timestep, get, SNES
968: @*/
969: PetscErrorCode PETSCTS_DLLEXPORT TSGetSNES(TS ts,SNES *snes)
970: {
974: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
975: *snes = ts->snes;
976: return(0);
977: }
981: /*@
982: TSGetKSP - Returns the KSP (linear solver) associated with
983: a TS (timestepper) context.
985: Not Collective, but KSP is parallel if TS is parallel
987: Input Parameter:
988: . ts - the TS context obtained from TSCreate()
990: Output Parameter:
991: . ksp - the nonlinear solver context
993: Notes:
994: The user can then directly manipulate the KSP context to set various
995: options, etc. Likewise, the user can then extract and manipulate the
996: KSP and PC contexts as well.
998: TSGetKSP() does not work for integrators that do not use KSP;
999: in this case TSGetKSP() returns PETSC_NULL in ksp.
1001: Level: beginner
1003: .keywords: timestep, get, KSP
1004: @*/
1005: PetscErrorCode PETSCTS_DLLEXPORT TSGetKSP(TS ts,KSP *ksp)
1006: {
1010: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1011: *ksp = ts->ksp;
1012: return(0);
1013: }
1015: /* ----------- Routines to set solver parameters ---------- */
1019: /*@
1020: TSGetDuration - Gets the maximum number of timesteps to use and
1021: maximum time for iteration.
1023: Collective on TS
1025: Input Parameters:
1026: + ts - the TS context obtained from TSCreate()
1027: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1028: - maxtime - final time to iterate to, or PETSC_NULL
1030: Level: intermediate
1032: .keywords: TS, timestep, get, maximum, iterations, time
1033: @*/
1034: PetscErrorCode PETSCTS_DLLEXPORT TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1035: {
1038: if (maxsteps) {
1040: *maxsteps = ts->max_steps;
1041: }
1042: if (maxtime ) {
1044: *maxtime = ts->max_time;
1045: }
1046: return(0);
1047: }
1051: /*@
1052: TSSetDuration - Sets the maximum number of timesteps to use and
1053: maximum time for iteration.
1055: Collective on TS
1057: Input Parameters:
1058: + ts - the TS context obtained from TSCreate()
1059: . maxsteps - maximum number of iterations to use
1060: - maxtime - final time to iterate to
1062: Options Database Keys:
1063: . -ts_max_steps <maxsteps> - Sets maxsteps
1064: . -ts_max_time <maxtime> - Sets maxtime
1066: Notes:
1067: The default maximum number of iterations is 5000. Default time is 5.0
1069: Level: intermediate
1071: .keywords: TS, timestep, set, maximum, iterations
1072: @*/
1073: PetscErrorCode PETSCTS_DLLEXPORT TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1074: {
1077: ts->max_steps = maxsteps;
1078: ts->max_time = maxtime;
1079: return(0);
1080: }
1084: /*@
1085: TSSetSolution - Sets the initial solution vector
1086: for use by the TS routines.
1088: Collective on TS and Vec
1090: Input Parameters:
1091: + ts - the TS context obtained from TSCreate()
1092: - x - the solution vector
1094: Level: beginner
1096: .keywords: TS, timestep, set, solution, initial conditions
1097: @*/
1098: PetscErrorCode PETSCTS_DLLEXPORT TSSetSolution(TS ts,Vec x)
1099: {
1103: ts->vec_sol = ts->vec_sol_always = x;
1104: return(0);
1105: }
1109: /*@
1110: TSDefaultRhsBC - The default boundary condition function which does nothing.
1112: Collective on TS
1114: Input Parameters:
1115: + ts - The TS context obtained from TSCreate()
1116: . rhs - The Rhs
1117: - ctx - The user-context
1119: Level: developer
1121: .keywords: TS, Rhs, boundary conditions
1122: @*/
1123: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1124: {
1126: return(0);
1127: }
1131: /*@C
1132: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1133: to the system matrix and preconditioner of each system.
1135: Collective on TS
1137: Input Parameters:
1138: + ts - The TS context obtained from TSCreate()
1139: - func - The function
1141: Calling sequence of func:
1142: . func (TS ts, Mat A, Mat B, void *ctx);
1144: + A - The current system matrix
1145: . B - The current preconditioner
1146: - ctx - The user-context
1148: Level: intermediate
1150: .keywords: TS, System matrix, boundary conditions
1151: @*/
1152: PetscErrorCode PETSCTS_DLLEXPORT TSSetSystemMatrixBC(TS ts, PetscErrorCode (*func)(TS, Mat, Mat, void *))
1153: {
1156: ts->ops->applymatrixbc = func;
1157: return(0);
1158: }
1162: /*@
1163: TSDefaultSystemMatrixBC - The default boundary condition function which
1164: does nothing.
1166: Collective on TS
1168: Input Parameters:
1169: + ts - The TS context obtained from TSCreate()
1170: . A - The system matrix
1171: . B - The preconditioner
1172: - ctx - The user-context
1174: Level: developer
1176: .keywords: TS, System matrix, boundary conditions
1177: @*/
1178: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1179: {
1181: return(0);
1182: }
1186: /*@C
1187: TSSetSolutionBC - Sets the function which applies boundary conditions
1188: to the solution of each system. This is necessary in nonlinear systems
1189: which time dependent boundary conditions.
1191: Collective on TS
1193: Input Parameters:
1194: + ts - The TS context obtained from TSCreate()
1195: - func - The function
1197: Calling sequence of func:
1198: . func (TS ts, Vec rsol, void *ctx);
1200: + sol - The current solution vector
1201: - ctx - The user-context
1203: Level: intermediate
1205: .keywords: TS, solution, boundary conditions
1206: @*/
1207: PetscErrorCode PETSCTS_DLLEXPORT TSSetSolutionBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1208: {
1211: ts->ops->applysolbc = func;
1212: return(0);
1213: }
1217: /*@
1218: TSDefaultSolutionBC - The default boundary condition function which
1219: does nothing.
1221: Collective on TS
1223: Input Parameters:
1224: + ts - The TS context obtained from TSCreate()
1225: . sol - The solution
1226: - ctx - The user-context
1228: Level: developer
1230: .keywords: TS, solution, boundary conditions
1231: @*/
1232: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1233: {
1235: return(0);
1236: }
1240: /*@C
1241: TSSetPreStep - Sets the general-purpose function
1242: called once at the beginning of time stepping.
1244: Collective on TS
1246: Input Parameters:
1247: + ts - The TS context obtained from TSCreate()
1248: - func - The function
1250: Calling sequence of func:
1251: . func (TS ts);
1253: Level: intermediate
1255: .keywords: TS, timestep
1256: @*/
1257: PetscErrorCode PETSCTS_DLLEXPORT TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1258: {
1261: ts->ops->prestep = func;
1262: return(0);
1263: }
1267: /*@
1268: TSDefaultPreStep - The default pre-stepping function which does nothing.
1270: Collective on TS
1272: Input Parameters:
1273: . ts - The TS context obtained from TSCreate()
1275: Level: developer
1277: .keywords: TS, timestep
1278: @*/
1279: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultPreStep(TS ts)
1280: {
1282: return(0);
1283: }
1287: /*@C
1288: TSSetUpdate - Sets the general-purpose update function called
1289: at the beginning of every time step. This function can change
1290: the time step.
1292: Collective on TS
1294: Input Parameters:
1295: + ts - The TS context obtained from TSCreate()
1296: - func - The function
1298: Calling sequence of func:
1299: . func (TS ts, double t, double *dt);
1301: + t - The current time
1302: - dt - The current time step
1304: Level: intermediate
1306: .keywords: TS, update, timestep
1307: @*/
1308: PetscErrorCode PETSCTS_DLLEXPORT TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1309: {
1312: ts->ops->update = func;
1313: return(0);
1314: }
1318: /*@
1319: TSDefaultUpdate - The default update function which does nothing.
1321: Collective on TS
1323: Input Parameters:
1324: + ts - The TS context obtained from TSCreate()
1325: - t - The current time
1327: Output Parameters:
1328: . dt - The current time step
1330: Level: developer
1332: .keywords: TS, update, timestep
1333: @*/
1334: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1335: {
1337: return(0);
1338: }
1342: /*@C
1343: TSSetPostStep - Sets the general-purpose function
1344: called once at the end of time stepping.
1346: Collective on TS
1348: Input Parameters:
1349: + ts - The TS context obtained from TSCreate()
1350: - func - The function
1352: Calling sequence of func:
1353: . func (TS ts);
1355: Level: intermediate
1357: .keywords: TS, timestep
1358: @*/
1359: PetscErrorCode PETSCTS_DLLEXPORT TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1360: {
1363: ts->ops->poststep = func;
1364: return(0);
1365: }
1369: /*@
1370: TSDefaultPostStep - The default post-stepping function which does nothing.
1372: Collective on TS
1374: Input Parameters:
1375: . ts - The TS context obtained from TSCreate()
1377: Level: developer
1379: .keywords: TS, timestep
1380: @*/
1381: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultPostStep(TS ts)
1382: {
1384: return(0);
1385: }
1387: /* ------------ Routines to set performance monitoring options ----------- */
1391: /*@C
1392: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1393: timestep to display the iteration's progress.
1395: Collective on TS
1397: Input Parameters:
1398: + ts - the TS context obtained from TSCreate()
1399: . func - monitoring routine
1400: . mctx - [optional] user-defined context for private data for the
1401: monitor routine (use PETSC_NULL if no context is desired)
1402: - monitordestroy - [optional] routine that frees monitor context
1403: (may be PETSC_NULL)
1405: Calling sequence of func:
1406: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1408: + ts - the TS context
1409: . steps - iteration number
1410: . time - current time
1411: . x - current iterate
1412: - mctx - [optional] monitoring context
1414: Notes:
1415: This routine adds an additional monitor to the list of monitors that
1416: already has been loaded.
1418: Level: intermediate
1420: .keywords: TS, timestep, set, monitor
1422: .seealso: TSDefaultMonitor(), TSClearMonitor()
1423: @*/
1424: PetscErrorCode PETSCTS_DLLEXPORT TSSetMonitor(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1425: {
1428: if (ts->numbermonitors >= MAXTSMONITORS) {
1429: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1430: }
1431: ts->monitor[ts->numbermonitors] = monitor;
1432: ts->mdestroy[ts->numbermonitors] = mdestroy;
1433: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1434: return(0);
1435: }
1439: /*@C
1440: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1442: Collective on TS
1444: Input Parameters:
1445: . ts - the TS context obtained from TSCreate()
1447: Notes:
1448: There is no way to remove a single, specific monitor.
1450: Level: intermediate
1452: .keywords: TS, timestep, set, monitor
1454: .seealso: TSDefaultMonitor(), TSSetMonitor()
1455: @*/
1456: PetscErrorCode PETSCTS_DLLEXPORT TSClearMonitor(TS ts)
1457: {
1459: PetscInt i;
1463: for (i=0; i<ts->numbermonitors; i++) {
1464: if (ts->mdestroy[i]) {
1465: (*ts->mdestroy[i])(ts->monitorcontext[i]);
1466: }
1467: }
1468: ts->numbermonitors = 0;
1469: return(0);
1470: }
1474: /*@
1475: TSDefaultMonitor - Sets the Default monitor
1476: @*/
1477: PetscErrorCode TSDefaultMonitor(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1478: {
1480: PetscViewer viewer = (PetscViewer)ctx;
1483: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
1484: PetscViewerASCIIPrintf(viewer,"timestep %D dt %G time %G\n",step,ts->time_step,ptime);
1485: return(0);
1486: }
1490: /*@
1491: TSStep - Steps the requested number of timesteps.
1493: Collective on TS
1495: Input Parameter:
1496: . ts - the TS context obtained from TSCreate()
1498: Output Parameters:
1499: + steps - number of iterations until termination
1500: - ptime - time until termination
1502: Level: beginner
1504: .keywords: TS, timestep, solve
1506: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1507: @*/
1508: PetscErrorCode PETSCTS_DLLEXPORT TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1509: {
1514: if (!ts->setupcalled) {
1515: TSSetUp(ts);
1516: }
1518: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1519: (*ts->ops->prestep)(ts);
1520: (*ts->ops->step)(ts, steps, ptime);
1521: (*ts->ops->poststep)(ts);
1522: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1524: if (!PetscPreLoadingOn) {
1525: TSViewFromOptions(ts,ts->name);
1526: }
1527: return(0);
1528: }
1532: /*
1533: Runs the user provided monitor routines, if they exists.
1534: */
1535: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1536: {
1538: PetscInt i,n = ts->numbermonitors;
1541: for (i=0; i<n; i++) {
1542: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1543: }
1544: return(0);
1545: }
1547: /* ------------------------------------------------------------------------*/
1551: /*@C
1552: TSLGMonitorCreate - Creates a line graph context for use with
1553: TS to monitor convergence of preconditioned residual norms.
1555: Collective on TS
1557: Input Parameters:
1558: + host - the X display to open, or null for the local machine
1559: . label - the title to put in the title bar
1560: . x, y - the screen coordinates of the upper left coordinate of the window
1561: - m, n - the screen width and height in pixels
1563: Output Parameter:
1564: . draw - the drawing context
1566: Options Database Key:
1567: . -ts_xmonitor - automatically sets line graph monitor
1569: Notes:
1570: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1572: Level: intermediate
1574: .keywords: TS, monitor, line graph, residual, seealso
1576: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1578: @*/
1579: PetscErrorCode PETSCTS_DLLEXPORT TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1580: {
1581: PetscDraw win;
1585: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1586: PetscDrawSetType(win,PETSC_DRAW_X);
1587: PetscDrawLGCreate(win,1,draw);
1588: PetscDrawLGIndicateDataPoints(*draw);
1590: PetscLogObjectParent(*draw,win);
1591: return(0);
1592: }
1596: PetscErrorCode TSLGMonitor(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1597: {
1598: PetscDrawLG lg = (PetscDrawLG) monctx;
1599: PetscReal x,y = ptime;
1603: if (!monctx) {
1604: MPI_Comm comm;
1605: PetscViewer viewer;
1607: PetscObjectGetComm((PetscObject)ts,&comm);
1608: viewer = PETSC_VIEWER_DRAW_(comm);
1609: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1610: }
1612: if (!n) {PetscDrawLGReset(lg);}
1613: x = (PetscReal)n;
1614: PetscDrawLGAddPoint(lg,&x,&y);
1615: if (n < 20 || (n % 5)) {
1616: PetscDrawLGDraw(lg);
1617: }
1618: return(0);
1619: }
1623: /*@C
1624: TSLGMonitorDestroy - Destroys a line graph context that was created
1625: with TSLGMonitorCreate().
1627: Collective on PetscDrawLG
1629: Input Parameter:
1630: . draw - the drawing context
1632: Level: intermediate
1634: .keywords: TS, monitor, line graph, destroy
1636: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1637: @*/
1638: PetscErrorCode PETSCTS_DLLEXPORT TSLGMonitorDestroy(PetscDrawLG drawlg)
1639: {
1640: PetscDraw draw;
1644: PetscDrawLGGetDraw(drawlg,&draw);
1645: PetscDrawDestroy(draw);
1646: PetscDrawLGDestroy(drawlg);
1647: return(0);
1648: }
1652: /*@
1653: TSGetTime - Gets the current time.
1655: Not Collective
1657: Input Parameter:
1658: . ts - the TS context obtained from TSCreate()
1660: Output Parameter:
1661: . t - the current time
1663: Contributed by: Matthew Knepley
1665: Level: beginner
1667: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1669: .keywords: TS, get, time
1670: @*/
1671: PetscErrorCode PETSCTS_DLLEXPORT TSGetTime(TS ts,PetscReal* t)
1672: {
1676: *t = ts->ptime;
1677: return(0);
1678: }
1682: /*@C
1683: TSSetOptionsPrefix - Sets the prefix used for searching for all
1684: TS options in the database.
1686: Collective on TS
1688: Input Parameter:
1689: + ts - The TS context
1690: - prefix - The prefix to prepend to all option names
1692: Notes:
1693: A hyphen (-) must NOT be given at the beginning of the prefix name.
1694: The first character of all runtime options is AUTOMATICALLY the
1695: hyphen.
1697: Contributed by: Matthew Knepley
1699: Level: advanced
1701: .keywords: TS, set, options, prefix, database
1703: .seealso: TSSetFromOptions()
1705: @*/
1706: PetscErrorCode PETSCTS_DLLEXPORT TSSetOptionsPrefix(TS ts,const char prefix[])
1707: {
1712: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1713: switch(ts->problem_type) {
1714: case TS_NONLINEAR:
1715: SNESSetOptionsPrefix(ts->snes,prefix);
1716: break;
1717: case TS_LINEAR:
1718: KSPSetOptionsPrefix(ts->ksp,prefix);
1719: break;
1720: }
1721: return(0);
1722: }
1727: /*@C
1728: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1729: TS options in the database.
1731: Collective on TS
1733: Input Parameter:
1734: + ts - The TS context
1735: - prefix - The prefix to prepend to all option names
1737: Notes:
1738: A hyphen (-) must NOT be given at the beginning of the prefix name.
1739: The first character of all runtime options is AUTOMATICALLY the
1740: hyphen.
1742: Contributed by: Matthew Knepley
1744: Level: advanced
1746: .keywords: TS, append, options, prefix, database
1748: .seealso: TSGetOptionsPrefix()
1750: @*/
1751: PetscErrorCode PETSCTS_DLLEXPORT TSAppendOptionsPrefix(TS ts,const char prefix[])
1752: {
1757: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1758: switch(ts->problem_type) {
1759: case TS_NONLINEAR:
1760: SNESAppendOptionsPrefix(ts->snes,prefix);
1761: break;
1762: case TS_LINEAR:
1763: KSPAppendOptionsPrefix(ts->ksp,prefix);
1764: break;
1765: }
1766: return(0);
1767: }
1771: /*@C
1772: TSGetOptionsPrefix - Sets the prefix used for searching for all
1773: TS options in the database.
1775: Not Collective
1777: Input Parameter:
1778: . ts - The TS context
1780: Output Parameter:
1781: . prefix - A pointer to the prefix string used
1783: Contributed by: Matthew Knepley
1785: Notes: On the fortran side, the user should pass in a string 'prifix' of
1786: sufficient length to hold the prefix.
1788: Level: intermediate
1790: .keywords: TS, get, options, prefix, database
1792: .seealso: TSAppendOptionsPrefix()
1793: @*/
1794: PetscErrorCode PETSCTS_DLLEXPORT TSGetOptionsPrefix(TS ts,const char *prefix[])
1795: {
1801: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1802: return(0);
1803: }
1807: /*@C
1808: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1810: Not Collective, but parallel objects are returned if TS is parallel
1812: Input Parameter:
1813: . ts - The TS context obtained from TSCreate()
1815: Output Parameters:
1816: + A - The matrix A, where U_t = A(t) U
1817: . M - The preconditioner matrix, usually the same as A
1818: - ctx - User-defined context for matrix evaluation routine
1820: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1822: Contributed by: Matthew Knepley
1824: Level: intermediate
1826: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1828: .keywords: TS, timestep, get, matrix
1830: @*/
1831: PetscErrorCode PETSCTS_DLLEXPORT TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1832: {
1835: if (A) *A = ts->A;
1836: if (M) *M = ts->B;
1837: if (ctx) *ctx = ts->jacP;
1838: return(0);
1839: }
1843: /*@C
1844: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1846: Not Collective, but parallel objects are returned if TS is parallel
1848: Input Parameter:
1849: . ts - The TS context obtained from TSCreate()
1851: Output Parameters:
1852: + J - The Jacobian J of F, where U_t = F(U,t)
1853: . M - The preconditioner matrix, usually the same as J
1854: - ctx - User-defined context for Jacobian evaluation routine
1856: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1858: Contributed by: Matthew Knepley
1860: Level: intermediate
1862: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1864: .keywords: TS, timestep, get, matrix, Jacobian
1865: @*/
1866: PetscErrorCode PETSCTS_DLLEXPORT TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1867: {
1871: TSGetRHSMatrix(ts,J,M,ctx);
1872: return(0);
1873: }
1877: /*@C
1878: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1879: VecView() for the solution at each timestep
1881: Collective on TS
1883: Input Parameters:
1884: + ts - the TS context
1885: . step - current time-step
1886: . ptime - current time
1887: - dummy - either a viewer or PETSC_NULL
1889: Level: intermediate
1891: .keywords: TS, vector, monitor, view
1893: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1894: @*/
1895: PetscErrorCode PETSCTS_DLLEXPORT TSVecViewMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1896: {
1898: PetscViewer viewer = (PetscViewer) dummy;
1901: if (!viewer) {
1902: MPI_Comm comm;
1903: PetscObjectGetComm((PetscObject)ts,&comm);
1904: viewer = PETSC_VIEWER_DRAW_(comm);
1905: }
1906: VecView(x,viewer);
1907: return(0);
1908: }