Actual source code: ts.c
1: #define PETSCTS_DLL
3: #include include/private/tsimpl.h
5: /* Logging support */
6: PetscCookie 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_monitor_draw - plot information at each timestep
68: Level: beginner
70: .keywords: TS, timestep, set, options, database
72: .seealso: TSGetType
73: @*/
74: PetscErrorCode TSSetFromOptions(TS ts)
75: {
76: PetscReal dt;
77: PetscTruth opt,flg;
78: PetscErrorCode ierr;
79: PetscViewerASCIIMonitor 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","TSMonitorDefault","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
97: if (flg) {
98: PetscViewerASCIIMonitorCreate(ts->comm,monfilename,0,&monviewer);
99: TSMonitorSet(ts,TSMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
100: }
101: PetscOptionsName("-ts_monitor_draw","Monitor timestep size graphically","TSMonitorLG",&opt);
102: if (opt) {
103: TSMonitorSet(ts,TSMonitorLG,PETSC_NULL,PETSC_NULL);
104: }
105: PetscOptionsName("-ts_monitor_solution","Monitor solution graphically","TSMonitorSolution",&opt);
106: if (opt) {
107: TSMonitorSet(ts,TSMonitorSolution,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->Arhs,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->Arhs,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 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 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) {
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;
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: {
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->Arhs,&ts->B,&flg,ts->jacP);
284: PetscStackPop;
285: }
286: MatMult(ts->Arhs,x,y);
287: }
291: return(0);
292: }
296: /*@C
297: TSSetRHSFunction - Sets the routine for evaluating the function,
298: F(t,u), where U_t = F(t,u).
300: Collective on TS
302: Input Parameters:
303: + ts - the TS context obtained from TSCreate()
304: . f - routine for evaluating the right-hand-side function
305: - ctx - [optional] user-defined context for private data for the
306: function evaluation routine (may be PETSC_NULL)
308: Calling sequence of func:
309: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
311: + t - current timestep
312: . u - input vector
313: . F - function vector
314: - ctx - [optional] user-defined function context
316: Important:
317: The user MUST call either this routine or TSSetRHSMatrix().
319: Level: beginner
321: .keywords: TS, timestep, set, right-hand-side, function
323: .seealso: TSSetRHSMatrix()
324: @*/
325: PetscErrorCode TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
326: {
330: if (ts->problem_type == TS_LINEAR) {
331: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
332: }
333: ts->ops->rhsfunction = f;
334: ts->funP = ctx;
335: return(0);
336: }
340: /*@C
341: TSSetMatrices - Sets the functions to compute the matrices Alhs and Arhs,
342: where Alhs(t) U_t = Arhs(t) U.
344: Collective on TS
346: Input Parameters:
347: + ts - the TS context obtained from TSCreate()
348: . Arhs - matrix
349: . frhs - the matrix evaluation routine for Arhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
350: if Arhs is not a function of t.
351: . Alhs - matrix or PETSC_NULL if Alhs is an indentity matrix.
352: . flhs - the matrix evaluation routine for Alhs; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
353: if Alhs is not a function of t.
354: . flag - flag indicating information about the matrix structure of Arhs and Alhs.
355: The available options are
356: SAME_NONZERO_PATTERN - Alhs has the same nonzero structure as Arhs
357: DIFFERENT_NONZERO_PATTERN - Alhs has different nonzero structure as Arhs
358: - ctx - [optional] user-defined context for private data for the
359: matrix evaluation routine (may be PETSC_NULL)
361: Calling sequence of func:
362: $ func(TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
364: + t - current timestep
365: . A - matrix A, where U_t = A(t) U
366: . B - preconditioner matrix, usually the same as A
367: . flag - flag indicating information about the preconditioner matrix
368: structure (same as flag in KSPSetOperators())
369: - ctx - [optional] user-defined context for matrix evaluation routine
371: Notes:
372: The routine func() takes Mat* as the matrix arguments rather than Mat.
373: This allows the matrix evaluation routine to replace Arhs or Alhs with a
374: completely new new matrix structure (not just different matrix elements)
375: when appropriate, for instance, if the nonzero structure is changing
376: throughout the global iterations.
378: Important:
379: The user MUST call either this routine or TSSetRHSFunction().
381: Level: beginner
383: .keywords: TS, timestep, set, matrix
385: .seealso: TSSetRHSFunction()
386: @*/
387: 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)
388: {
394: if (ts->problem_type == TS_NONLINEAR)
395: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
396:
397: ts->ops->rhsmatrix = frhs;
398: ts->ops->lhsmatrix = flhs;
399: ts->jacP = ctx;
400: ts->Arhs = Arhs;
401: ts->Alhs = Alhs;
402: ts->matflg = flag;
403: return(0);
404: }
408: /*@C
409: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
410: Also sets the location to store A.
412: Collective on TS
414: Input Parameters:
415: + ts - the TS context obtained from TSCreate()
416: . A - matrix
417: . B - preconditioner matrix (usually same as A)
418: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
419: if A is not a function of t.
420: - ctx - [optional] user-defined context for private data for the
421: matrix evaluation routine (may be PETSC_NULL)
423: Calling sequence of func:
424: $ func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
426: + t - current timestep
427: . A - matrix A, where U_t = A(t) U
428: . B - preconditioner matrix, usually the same as A
429: . flag - flag indicating information about the preconditioner matrix
430: structure (same as flag in KSPSetOperators())
431: - ctx - [optional] user-defined context for matrix evaluation routine
433: Notes:
434: See KSPSetOperators() for important information about setting the flag
435: output parameter in the routine func(). Be sure to read this information!
437: The routine func() takes Mat * as the matrix arguments rather than Mat.
438: This allows the matrix evaluation routine to replace A and/or B with a
439: completely new new matrix structure (not just different matrix elements)
440: when appropriate, for instance, if the nonzero structure is changing
441: throughout the global iterations.
443: Important:
444: The user MUST call either this routine or TSSetRHSFunction().
446: Level: beginner
448: .keywords: TS, timestep, set, right-hand-side, matrix
450: .seealso: TSSetRHSFunction()
451: @*/
452: PetscErrorCode TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
453: {
460: if (ts->problem_type == TS_NONLINEAR) {
461: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
462: }
464: ts->ops->rhsmatrix = f;
465: ts->jacP = ctx;
466: ts->Arhs = A;
467: ts->B = B;
468: return(0);
469: }
473: /*@C
474: TSSetLHSMatrix - Sets the function to compute the matrix A, where A U_t = F(U).
475: Also sets the location to store A.
477: Collective on TS
479: Input Parameters:
480: + ts - the TS context obtained from TSCreate()
481: . A - matrix
482: . B - ignored
483: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
484: if A is not a function of t.
485: - ctx - [optional] user-defined context for private data for the
486: matrix evaluation routine (may be PETSC_NULL)
488: Calling sequence of func:
489: $ func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
491: + t - current timestep
492: . A - matrix A, where A U_t = F(U)
493: . B - ignored
494: . flag - flag indicating information about the preconditioner matrix
495: structure (same as flag in KSPSetOperators())
496: - ctx - [optional] user-defined context for matrix evaluation routine
498: Notes:
499: See KSPSetOperators() for important information about setting the flag
500: output parameter in the routine func(). Be sure to read this information!
502: The routine func() takes Mat * as the matrix arguments rather than Mat.
503: This allows the matrix evaluation routine to replace A and/or B with a
504: completely new new matrix structure (not just different matrix elements)
505: when appropriate, for instance, if the nonzero structure is changing
506: throughout the global iterations.
508: Notes:
509: Currently, TSSetLHSMatrix() only supports the TS_BEULER type.
511: Level: beginner
513: .keywords: TS, timestep, set, left-hand-side, matrix
515: .seealso: TSSetRHSMatrix()
516: @*/
517: PetscErrorCode TSSetLHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
518: {
519: PetscTruth sametype;
528:
529: if (!ts->type_name) SETERRQ(PETSC_ERR_ARG_NULL,"TS type must be set before calling TSSetLHSMatrix()");
530: PetscTypeCompare((PetscObject)ts,"beuler",&sametype);
531: if (!sametype)
532: SETERRQ1(PETSC_ERR_SUP,"TS type %s not supported for LHSMatrix",ts->type_name);
533: if (ts->problem_type == TS_NONLINEAR) {
534: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems yet");
535: }
537: ts->ops->lhsmatrix = f;
538: ts->jacPlhs = ctx;
539: ts->Alhs = A;
540: return(0);
541: }
545: /*@C
546: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
547: where U_t = F(U,t), as well as the location to store the matrix.
548: Use TSSetRHSMatrix() for linear problems.
550: Collective on TS
552: Input Parameters:
553: + ts - the TS context obtained from TSCreate()
554: . A - Jacobian matrix
555: . B - preconditioner matrix (usually same as A)
556: . f - the Jacobian evaluation routine
557: - ctx - [optional] user-defined context for private data for the
558: Jacobian evaluation routine (may be PETSC_NULL)
560: Calling sequence of func:
561: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
563: + t - current timestep
564: . u - input vector
565: . A - matrix A, where U_t = A(t)u
566: . B - preconditioner matrix, usually the same as A
567: . flag - flag indicating information about the preconditioner matrix
568: structure (same as flag in KSPSetOperators())
569: - ctx - [optional] user-defined context for matrix evaluation routine
571: Notes:
572: See KSPSetOperators() for important information about setting the flag
573: output parameter in the routine func(). Be sure to read this information!
575: The routine func() takes Mat * as the matrix arguments rather than Mat.
576: This allows the matrix evaluation routine to replace A and/or B with a
577: completely new new matrix structure (not just different matrix elements)
578: when appropriate, for instance, if the nonzero structure is changing
579: throughout the global iterations.
581: Level: beginner
582:
583: .keywords: TS, timestep, set, right-hand-side, Jacobian
585: .seealso: TSDefaultComputeJacobianColor(),
586: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()
588: @*/
589: PetscErrorCode TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
590: {
597: if (ts->problem_type != TS_NONLINEAR) {
598: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
599: }
601: ts->ops->rhsjacobian = f;
602: ts->jacP = ctx;
603: ts->Arhs = A;
604: ts->B = B;
605: return(0);
606: }
610: /*@C
611: TSView - Prints the TS data structure.
613: Collective on TS
615: Input Parameters:
616: + ts - the TS context obtained from TSCreate()
617: - viewer - visualization context
619: Options Database Key:
620: . -ts_view - calls TSView() at end of TSStep()
622: Notes:
623: The available visualization contexts include
624: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
625: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
626: output where only the first processor opens
627: the file. All other processors send their
628: data to the first processor to print.
630: The user can open an alternative visualization context with
631: PetscViewerASCIIOpen() - output to a specified file.
633: Level: beginner
635: .keywords: TS, timestep, view
637: .seealso: PetscViewerASCIIOpen()
638: @*/
639: PetscErrorCode TSView(TS ts,PetscViewer viewer)
640: {
642: char *type;
643: PetscTruth iascii,isstring;
647: if (!viewer) {
648: PetscViewerASCIIGetStdout(ts->comm,&viewer);
649: }
653: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
654: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
655: if (iascii) {
656: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
657: TSGetType(ts,(TSType *)&type);
658: if (type) {
659: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
660: } else {
661: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
662: }
663: if (ts->ops->view) {
664: PetscViewerASCIIPushTab(viewer);
665: (*ts->ops->view)(ts,viewer);
666: PetscViewerASCIIPopTab(viewer);
667: }
668: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
669: PetscViewerASCIIPrintf(viewer," maximum time=%G\n",ts->max_time);
670: if (ts->problem_type == TS_NONLINEAR) {
671: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
672: }
673: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
674: } else if (isstring) {
675: TSGetType(ts,(TSType *)&type);
676: PetscViewerStringSPrintf(viewer," %-7.7s",type);
677: }
678: PetscViewerASCIIPushTab(viewer);
679: if (ts->ksp) {KSPView(ts->ksp,viewer);}
680: if (ts->snes) {SNESView(ts->snes,viewer);}
681: PetscViewerASCIIPopTab(viewer);
682: return(0);
683: }
688: /*@C
689: TSSetApplicationContext - Sets an optional user-defined context for
690: the timesteppers.
692: Collective on TS
694: Input Parameters:
695: + ts - the TS context obtained from TSCreate()
696: - usrP - optional user context
698: Level: intermediate
700: .keywords: TS, timestep, set, application, context
702: .seealso: TSGetApplicationContext()
703: @*/
704: PetscErrorCode TSSetApplicationContext(TS ts,void *usrP)
705: {
708: ts->user = usrP;
709: return(0);
710: }
714: /*@C
715: TSGetApplicationContext - Gets the user-defined context for the
716: timestepper.
718: Not Collective
720: Input Parameter:
721: . ts - the TS context obtained from TSCreate()
723: Output Parameter:
724: . usrP - user context
726: Level: intermediate
728: .keywords: TS, timestep, get, application, context
730: .seealso: TSSetApplicationContext()
731: @*/
732: PetscErrorCode TSGetApplicationContext(TS ts,void **usrP)
733: {
736: *usrP = ts->user;
737: return(0);
738: }
742: /*@
743: TSGetTimeStepNumber - Gets the current number of timesteps.
745: Not Collective
747: Input Parameter:
748: . ts - the TS context obtained from TSCreate()
750: Output Parameter:
751: . iter - number steps so far
753: Level: intermediate
755: .keywords: TS, timestep, get, iteration, number
756: @*/
757: PetscErrorCode TSGetTimeStepNumber(TS ts,PetscInt* iter)
758: {
762: *iter = ts->steps;
763: return(0);
764: }
768: /*@
769: TSSetInitialTimeStep - Sets the initial timestep to be used,
770: as well as the initial time.
772: Collective on TS
774: Input Parameters:
775: + ts - the TS context obtained from TSCreate()
776: . initial_time - the initial time
777: - time_step - the size of the timestep
779: Level: intermediate
781: .seealso: TSSetTimeStep(), TSGetTimeStep()
783: .keywords: TS, set, initial, timestep
784: @*/
785: PetscErrorCode TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
786: {
789: ts->time_step = time_step;
790: ts->initial_time_step = time_step;
791: ts->ptime = initial_time;
792: return(0);
793: }
797: /*@
798: TSSetTimeStep - Allows one to reset the timestep at any time,
799: useful for simple pseudo-timestepping codes.
801: Collective on TS
803: Input Parameters:
804: + ts - the TS context obtained from TSCreate()
805: - time_step - the size of the timestep
807: Level: intermediate
809: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
811: .keywords: TS, set, timestep
812: @*/
813: PetscErrorCode TSSetTimeStep(TS ts,PetscReal time_step)
814: {
817: ts->time_step = time_step;
818: return(0);
819: }
823: /*@
824: TSGetTimeStep - Gets the current timestep size.
826: Not Collective
828: Input Parameter:
829: . ts - the TS context obtained from TSCreate()
831: Output Parameter:
832: . dt - the current timestep size
834: Level: intermediate
836: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
838: .keywords: TS, get, timestep
839: @*/
840: PetscErrorCode TSGetTimeStep(TS ts,PetscReal* dt)
841: {
845: *dt = ts->time_step;
846: return(0);
847: }
851: /*@
852: TSGetSolution - Returns the solution at the present timestep. It
853: is valid to call this routine inside the function that you are evaluating
854: in order to move to the new timestep. This vector not changed until
855: the solution at the next timestep has been calculated.
857: Not Collective, but Vec returned is parallel if TS is parallel
859: Input Parameter:
860: . ts - the TS context obtained from TSCreate()
862: Output Parameter:
863: . v - the vector containing the solution
865: Level: intermediate
867: .seealso: TSGetTimeStep()
869: .keywords: TS, timestep, get, solution
870: @*/
871: PetscErrorCode TSGetSolution(TS ts,Vec *v)
872: {
876: *v = ts->vec_sol_always;
877: return(0);
878: }
880: /* ----- Routines to initialize and destroy a timestepper ---- */
883: /*@
884: TSSetProblemType - Sets the type of problem to be solved.
886: Not collective
888: Input Parameters:
889: + ts - The TS
890: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
891: .vb
892: U_t = A U
893: U_t = A(t) U
894: U_t = F(t,U)
895: .ve
897: Level: beginner
899: .keywords: TS, problem type
900: .seealso: TSSetUp(), TSProblemType, TS
901: @*/
902: PetscErrorCode TSSetProblemType(TS ts, TSProblemType type)
903: {
906: ts->problem_type = type;
907: return(0);
908: }
912: /*@C
913: TSGetProblemType - Gets the type of problem to be solved.
915: Not collective
917: Input Parameter:
918: . ts - The TS
920: Output Parameter:
921: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
922: .vb
923: U_t = A U
924: U_t = A(t) U
925: U_t = F(t,U)
926: .ve
928: Level: beginner
930: .keywords: TS, problem type
931: .seealso: TSSetUp(), TSProblemType, TS
932: @*/
933: PetscErrorCode TSGetProblemType(TS ts, TSProblemType *type)
934: {
938: *type = ts->problem_type;
939: return(0);
940: }
944: /*@
945: TSSetUp - Sets up the internal data structures for the later use
946: of a timestepper.
948: Collective on TS
950: Input Parameter:
951: . ts - the TS context obtained from TSCreate()
953: Notes:
954: For basic use of the TS solvers the user need not explicitly call
955: TSSetUp(), since these actions will automatically occur during
956: the call to TSStep(). However, if one wishes to control this
957: phase separately, TSSetUp() should be called after TSCreate()
958: and optional routines of the form TSSetXXX(), but before TSStep().
960: Level: advanced
962: .keywords: TS, timestep, setup
964: .seealso: TSCreate(), TSStep(), TSDestroy()
965: @*/
966: PetscErrorCode TSSetUp(TS ts)
967: {
972: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
973: if (!ts->type_name) {
974: TSSetType(ts,TS_EULER);
975: }
976: (*ts->ops->setup)(ts);
977: ts->setupcalled = 1;
978: return(0);
979: }
983: /*@
984: TSDestroy - Destroys the timestepper context that was created
985: with TSCreate().
987: Collective on TS
989: Input Parameter:
990: . ts - the TS context obtained from TSCreate()
992: Level: beginner
994: .keywords: TS, timestepper, destroy
996: .seealso: TSCreate(), TSSetUp(), TSSolve()
997: @*/
998: PetscErrorCode TSDestroy(TS ts)
999: {
1004: if (--ts->refct > 0) return(0);
1006: /* if memory was published with AMS then destroy it */
1007: PetscObjectDepublish(ts);
1008: if (ts->A) {MatDestroy(ts->A);CHKERRQ(ierr)}
1009: if (ts->ksp) {KSPDestroy(ts->ksp);}
1010: if (ts->snes) {SNESDestroy(ts->snes);}
1011: if (ts->ops->destroy) {(*(ts)->ops->destroy)(ts);}
1012: TSMonitorCancel(ts);
1013: PetscHeaderDestroy(ts);
1014: return(0);
1015: }
1019: /*@
1020: TSGetSNES - Returns the SNES (nonlinear solver) associated with
1021: a TS (timestepper) context. Valid only for nonlinear problems.
1023: Not Collective, but SNES is parallel if TS is parallel
1025: Input Parameter:
1026: . ts - the TS context obtained from TSCreate()
1028: Output Parameter:
1029: . snes - the nonlinear solver context
1031: Notes:
1032: The user can then directly manipulate the SNES context to set various
1033: options, etc. Likewise, the user can then extract and manipulate the
1034: KSP, KSP, and PC contexts as well.
1036: TSGetSNES() does not work for integrators that do not use SNES; in
1037: this case TSGetSNES() returns PETSC_NULL in snes.
1039: Level: beginner
1041: .keywords: timestep, get, SNES
1042: @*/
1043: PetscErrorCode TSGetSNES(TS ts,SNES *snes)
1044: {
1048: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
1049: *snes = ts->snes;
1050: return(0);
1051: }
1055: /*@
1056: TSGetKSP - Returns the KSP (linear solver) associated with
1057: a TS (timestepper) context.
1059: Not Collective, but KSP is parallel if TS is parallel
1061: Input Parameter:
1062: . ts - the TS context obtained from TSCreate()
1064: Output Parameter:
1065: . ksp - the nonlinear solver context
1067: Notes:
1068: The user can then directly manipulate the KSP context to set various
1069: options, etc. Likewise, the user can then extract and manipulate the
1070: KSP and PC contexts as well.
1072: TSGetKSP() does not work for integrators that do not use KSP;
1073: in this case TSGetKSP() returns PETSC_NULL in ksp.
1075: Level: beginner
1077: .keywords: timestep, get, KSP
1078: @*/
1079: PetscErrorCode TSGetKSP(TS ts,KSP *ksp)
1080: {
1084: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1085: *ksp = ts->ksp;
1086: return(0);
1087: }
1089: /* ----------- Routines to set solver parameters ---------- */
1093: /*@
1094: TSGetDuration - Gets the maximum number of timesteps to use and
1095: maximum time for iteration.
1097: Collective on TS
1099: Input Parameters:
1100: + ts - the TS context obtained from TSCreate()
1101: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1102: - maxtime - final time to iterate to, or PETSC_NULL
1104: Level: intermediate
1106: .keywords: TS, timestep, get, maximum, iterations, time
1107: @*/
1108: PetscErrorCode TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1109: {
1112: if (maxsteps) {
1114: *maxsteps = ts->max_steps;
1115: }
1116: if (maxtime ) {
1118: *maxtime = ts->max_time;
1119: }
1120: return(0);
1121: }
1125: /*@
1126: TSSetDuration - Sets the maximum number of timesteps to use and
1127: maximum time for iteration.
1129: Collective on TS
1131: Input Parameters:
1132: + ts - the TS context obtained from TSCreate()
1133: . maxsteps - maximum number of iterations to use
1134: - maxtime - final time to iterate to
1136: Options Database Keys:
1137: . -ts_max_steps <maxsteps> - Sets maxsteps
1138: . -ts_max_time <maxtime> - Sets maxtime
1140: Notes:
1141: The default maximum number of iterations is 5000. Default time is 5.0
1143: Level: intermediate
1145: .keywords: TS, timestep, set, maximum, iterations
1146: @*/
1147: PetscErrorCode TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1148: {
1151: ts->max_steps = maxsteps;
1152: ts->max_time = maxtime;
1153: return(0);
1154: }
1158: /*@
1159: TSSetSolution - Sets the initial solution vector
1160: for use by the TS routines.
1162: Collective on TS and Vec
1164: Input Parameters:
1165: + ts - the TS context obtained from TSCreate()
1166: - x - the solution vector
1168: Level: beginner
1170: .keywords: TS, timestep, set, solution, initial conditions
1171: @*/
1172: PetscErrorCode TSSetSolution(TS ts,Vec x)
1173: {
1177: ts->vec_sol = ts->vec_sol_always = x;
1178: return(0);
1179: }
1183: /*@C
1184: TSSetPreStep - Sets the general-purpose function
1185: called once at the beginning of time stepping.
1187: Collective on TS
1189: Input Parameters:
1190: + ts - The TS context obtained from TSCreate()
1191: - func - The function
1193: Calling sequence of func:
1194: . func (TS ts);
1196: Level: intermediate
1198: .keywords: TS, timestep
1199: @*/
1200: PetscErrorCode TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1201: {
1204: ts->ops->prestep = func;
1205: return(0);
1206: }
1210: /*@
1211: TSDefaultPreStep - The default pre-stepping function which does nothing.
1213: Collective on TS
1215: Input Parameters:
1216: . ts - The TS context obtained from TSCreate()
1218: Level: developer
1220: .keywords: TS, timestep
1221: @*/
1222: PetscErrorCode TSDefaultPreStep(TS ts)
1223: {
1225: return(0);
1226: }
1230: /*@C
1231: TSSetUpdate - Sets the general-purpose update function called
1232: at the beginning of every time step. This function can change
1233: the time step.
1235: Collective on TS
1237: Input Parameters:
1238: + ts - The TS context obtained from TSCreate()
1239: - func - The function
1241: Calling sequence of func:
1242: . func (TS ts, double t, double *dt);
1244: + t - The current time
1245: - dt - The current time step
1247: Level: intermediate
1249: .keywords: TS, update, timestep
1250: @*/
1251: PetscErrorCode TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1252: {
1255: ts->ops->update = func;
1256: return(0);
1257: }
1261: /*@
1262: TSDefaultUpdate - The default update function which does nothing.
1264: Collective on TS
1266: Input Parameters:
1267: + ts - The TS context obtained from TSCreate()
1268: - t - The current time
1270: Output Parameters:
1271: . dt - The current time step
1273: Level: developer
1275: .keywords: TS, update, timestep
1276: @*/
1277: PetscErrorCode TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1278: {
1280: return(0);
1281: }
1285: /*@C
1286: TSSetPostStep - Sets the general-purpose function
1287: called once at the end of time stepping.
1289: Collective on TS
1291: Input Parameters:
1292: + ts - The TS context obtained from TSCreate()
1293: - func - The function
1295: Calling sequence of func:
1296: . func (TS ts);
1298: Level: intermediate
1300: .keywords: TS, timestep
1301: @*/
1302: PetscErrorCode TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1303: {
1306: ts->ops->poststep = func;
1307: return(0);
1308: }
1312: /*@
1313: TSDefaultPostStep - The default post-stepping function which does nothing.
1315: Collective on TS
1317: Input Parameters:
1318: . ts - The TS context obtained from TSCreate()
1320: Level: developer
1322: .keywords: TS, timestep
1323: @*/
1324: PetscErrorCode TSDefaultPostStep(TS ts)
1325: {
1327: return(0);
1328: }
1330: /* ------------ Routines to set performance monitoring options ----------- */
1334: /*@C
1335: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
1336: timestep to display the iteration's progress.
1338: Collective on TS
1340: Input Parameters:
1341: + ts - the TS context obtained from TSCreate()
1342: . func - monitoring routine
1343: . mctx - [optional] user-defined context for private data for the
1344: monitor routine (use PETSC_NULL if no context is desired)
1345: - monitordestroy - [optional] routine that frees monitor context
1346: (may be PETSC_NULL)
1348: Calling sequence of func:
1349: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1351: + ts - the TS context
1352: . steps - iteration number
1353: . time - current time
1354: . x - current iterate
1355: - mctx - [optional] monitoring context
1357: Notes:
1358: This routine adds an additional monitor to the list of monitors that
1359: already has been loaded.
1361: Level: intermediate
1363: .keywords: TS, timestep, set, monitor
1365: .seealso: TSMonitorDefault(), TSMonitorCancel()
1366: @*/
1367: PetscErrorCode TSMonitorSet(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1368: {
1371: if (ts->numbermonitors >= MAXTSMONITORS) {
1372: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1373: }
1374: ts->monitor[ts->numbermonitors] = monitor;
1375: ts->mdestroy[ts->numbermonitors] = mdestroy;
1376: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1377: return(0);
1378: }
1382: /*@C
1383: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
1385: Collective on TS
1387: Input Parameters:
1388: . ts - the TS context obtained from TSCreate()
1390: Notes:
1391: There is no way to remove a single, specific monitor.
1393: Level: intermediate
1395: .keywords: TS, timestep, set, monitor
1397: .seealso: TSMonitorDefault(), TSMonitorSet()
1398: @*/
1399: PetscErrorCode TSMonitorCancel(TS ts)
1400: {
1402: PetscInt i;
1406: for (i=0; i<ts->numbermonitors; i++) {
1407: if (ts->mdestroy[i]) {
1408: (*ts->mdestroy[i])(ts->monitorcontext[i]);
1409: }
1410: }
1411: ts->numbermonitors = 0;
1412: return(0);
1413: }
1417: /*@
1418: TSMonitorDefault - Sets the Default monitor
1419: @*/
1420: PetscErrorCode TSMonitorDefault(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1421: {
1422: PetscErrorCode ierr;
1423: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor)ctx;
1426: if (!ctx) {
1427: PetscViewerASCIIMonitorCreate(ts->comm,"stdout",0,&viewer);
1428: }
1429: PetscViewerASCIIMonitorPrintf(viewer,"timestep %D dt %G time %G\n",step,ts->time_step,ptime);
1430: if (!ctx) {
1431: PetscViewerASCIIMonitorDestroy(viewer);
1432: }
1433: return(0);
1434: }
1438: /*@
1439: TSStep - Steps the requested number of timesteps.
1441: Collective on TS
1443: Input Parameter:
1444: . ts - the TS context obtained from TSCreate()
1446: Output Parameters:
1447: + steps - number of iterations until termination
1448: - ptime - time until termination
1450: Level: beginner
1452: .keywords: TS, timestep, solve
1454: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1455: @*/
1456: PetscErrorCode TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1457: {
1462: if (!ts->setupcalled) {
1463: TSSetUp(ts);
1464: }
1467: (*ts->ops->prestep)(ts);
1468: (*ts->ops->step)(ts, steps, ptime);
1469: (*ts->ops->poststep)(ts);
1472: if (!PetscPreLoadingOn) {
1473: TSViewFromOptions(ts,ts->name);
1474: }
1475: return(0);
1476: }
1480: /*@
1481: TSSolve - Steps the requested number of timesteps.
1483: Collective on TS
1485: Input Parameter:
1486: + ts - the TS context obtained from TSCreate()
1487: - x - the solution vector, or PETSC_NULL if it was set with TSSetSolution()
1489: Level: beginner
1491: .keywords: TS, timestep, solve
1493: .seealso: TSCreate(), TSSetSolution(), TSStep()
1494: @*/
1495: PetscErrorCode TSSolve(TS ts, Vec x)
1496: {
1497: PetscInt steps;
1498: PetscReal ptime;
1502: /* set solution vector if provided */
1503: if (x) { TSSetSolution(ts, x); }
1504: /* reset time step and iteration counters */
1505: ts->steps = 0; ts->linear_its = 0; ts->nonlinear_its = 0;
1506: /* steps the requested number of timesteps. */
1507: TSStep(ts, &steps, &ptime);
1508: return(0);
1509: }
1513: /*
1514: Runs the user provided monitor routines, if they exists.
1515: */
1516: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1517: {
1519: PetscInt i,n = ts->numbermonitors;
1522: for (i=0; i<n; i++) {
1523: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1524: }
1525: return(0);
1526: }
1528: /* ------------------------------------------------------------------------*/
1532: /*@C
1533: TSMonitorLGCreate - Creates a line graph context for use with
1534: TS to monitor convergence of preconditioned residual norms.
1536: Collective on TS
1538: Input Parameters:
1539: + host - the X display to open, or null for the local machine
1540: . label - the title to put in the title bar
1541: . x, y - the screen coordinates of the upper left coordinate of the window
1542: - m, n - the screen width and height in pixels
1544: Output Parameter:
1545: . draw - the drawing context
1547: Options Database Key:
1548: . -ts_monitor_draw - automatically sets line graph monitor
1550: Notes:
1551: Use TSMonitorLGDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1553: Level: intermediate
1555: .keywords: TS, monitor, line graph, residual, seealso
1557: .seealso: TSMonitorLGDestroy(), TSMonitorSet()
1559: @*/
1560: PetscErrorCode TSMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1561: {
1562: PetscDraw win;
1566: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1567: PetscDrawSetType(win,PETSC_DRAW_X);
1568: PetscDrawLGCreate(win,1,draw);
1569: PetscDrawLGIndicateDataPoints(*draw);
1571: PetscLogObjectParent(*draw,win);
1572: return(0);
1573: }
1577: PetscErrorCode TSMonitorLG(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1578: {
1579: PetscDrawLG lg = (PetscDrawLG) monctx;
1580: PetscReal x,y = ptime;
1584: if (!monctx) {
1585: MPI_Comm comm;
1586: PetscViewer viewer;
1588: PetscObjectGetComm((PetscObject)ts,&comm);
1589: viewer = PETSC_VIEWER_DRAW_(comm);
1590: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1591: }
1593: if (!n) {PetscDrawLGReset(lg);}
1594: x = (PetscReal)n;
1595: PetscDrawLGAddPoint(lg,&x,&y);
1596: if (n < 20 || (n % 5)) {
1597: PetscDrawLGDraw(lg);
1598: }
1599: return(0);
1600: }
1604: /*@C
1605: TSMonitorLGDestroy - Destroys a line graph context that was created
1606: with TSMonitorLGCreate().
1608: Collective on PetscDrawLG
1610: Input Parameter:
1611: . draw - the drawing context
1613: Level: intermediate
1615: .keywords: TS, monitor, line graph, destroy
1617: .seealso: TSMonitorLGCreate(), TSMonitorSet(), TSMonitorLG();
1618: @*/
1619: PetscErrorCode TSMonitorLGDestroy(PetscDrawLG drawlg)
1620: {
1621: PetscDraw draw;
1625: PetscDrawLGGetDraw(drawlg,&draw);
1626: PetscDrawDestroy(draw);
1627: PetscDrawLGDestroy(drawlg);
1628: return(0);
1629: }
1633: /*@
1634: TSGetTime - Gets the current time.
1636: Not Collective
1638: Input Parameter:
1639: . ts - the TS context obtained from TSCreate()
1641: Output Parameter:
1642: . t - the current time
1644: Contributed by: Matthew Knepley
1646: Level: beginner
1648: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1650: .keywords: TS, get, time
1651: @*/
1652: PetscErrorCode TSGetTime(TS ts,PetscReal* t)
1653: {
1657: *t = ts->ptime;
1658: return(0);
1659: }
1663: /*@
1664: TSSetTime - Allows one to reset the time.
1666: Collective on TS
1668: Input Parameters:
1669: + ts - the TS context obtained from TSCreate()
1670: - time - the time
1672: Level: intermediate
1674: .seealso: TSGetTime(), TSSetDuration()
1676: .keywords: TS, set, time
1677: @*/
1678: PetscErrorCode TSSetTime(TS ts, PetscReal t)
1679: {
1682: ts->ptime = t;
1683: return(0);
1684: }
1688: /*@C
1689: TSSetOptionsPrefix - Sets the prefix used for searching for all
1690: TS options in the database.
1692: Collective on TS
1694: Input Parameter:
1695: + ts - The TS context
1696: - prefix - The prefix to prepend to all option names
1698: Notes:
1699: A hyphen (-) must NOT be given at the beginning of the prefix name.
1700: The first character of all runtime options is AUTOMATICALLY the
1701: hyphen.
1703: Contributed by: Matthew Knepley
1705: Level: advanced
1707: .keywords: TS, set, options, prefix, database
1709: .seealso: TSSetFromOptions()
1711: @*/
1712: PetscErrorCode TSSetOptionsPrefix(TS ts,const char prefix[])
1713: {
1718: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1719: switch(ts->problem_type) {
1720: case TS_NONLINEAR:
1721: SNESSetOptionsPrefix(ts->snes,prefix);
1722: break;
1723: case TS_LINEAR:
1724: KSPSetOptionsPrefix(ts->ksp,prefix);
1725: break;
1726: }
1727: return(0);
1728: }
1733: /*@C
1734: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1735: TS options in the database.
1737: Collective on TS
1739: Input Parameter:
1740: + ts - The TS context
1741: - prefix - The prefix to prepend to all option names
1743: Notes:
1744: A hyphen (-) must NOT be given at the beginning of the prefix name.
1745: The first character of all runtime options is AUTOMATICALLY the
1746: hyphen.
1748: Contributed by: Matthew Knepley
1750: Level: advanced
1752: .keywords: TS, append, options, prefix, database
1754: .seealso: TSGetOptionsPrefix()
1756: @*/
1757: PetscErrorCode TSAppendOptionsPrefix(TS ts,const char prefix[])
1758: {
1763: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1764: switch(ts->problem_type) {
1765: case TS_NONLINEAR:
1766: SNESAppendOptionsPrefix(ts->snes,prefix);
1767: break;
1768: case TS_LINEAR:
1769: KSPAppendOptionsPrefix(ts->ksp,prefix);
1770: break;
1771: }
1772: return(0);
1773: }
1777: /*@C
1778: TSGetOptionsPrefix - Sets the prefix used for searching for all
1779: TS options in the database.
1781: Not Collective
1783: Input Parameter:
1784: . ts - The TS context
1786: Output Parameter:
1787: . prefix - A pointer to the prefix string used
1789: Contributed by: Matthew Knepley
1791: Notes: On the fortran side, the user should pass in a string 'prifix' of
1792: sufficient length to hold the prefix.
1794: Level: intermediate
1796: .keywords: TS, get, options, prefix, database
1798: .seealso: TSAppendOptionsPrefix()
1799: @*/
1800: PetscErrorCode TSGetOptionsPrefix(TS ts,const char *prefix[])
1801: {
1807: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1808: return(0);
1809: }
1813: /*@C
1814: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1816: Not Collective, but parallel objects are returned if TS is parallel
1818: Input Parameter:
1819: . ts - The TS context obtained from TSCreate()
1821: Output Parameters:
1822: + A - The matrix A, where U_t = A(t) U
1823: . M - The preconditioner matrix, usually the same as A
1824: - ctx - User-defined context for matrix evaluation routine
1826: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1828: Contributed by: Matthew Knepley
1830: Level: intermediate
1832: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1834: .keywords: TS, timestep, get, matrix
1836: @*/
1837: PetscErrorCode TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1838: {
1841: if (A) *A = ts->Arhs;
1842: if (M) *M = ts->B;
1843: if (ctx) *ctx = ts->jacP;
1844: return(0);
1845: }
1849: /*@C
1850: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1852: Not Collective, but parallel objects are returned if TS is parallel
1854: Input Parameter:
1855: . ts - The TS context obtained from TSCreate()
1857: Output Parameters:
1858: + J - The Jacobian J of F, where U_t = F(U,t)
1859: . M - The preconditioner matrix, usually the same as J
1860: - ctx - User-defined context for Jacobian evaluation routine
1862: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1864: Contributed by: Matthew Knepley
1866: Level: intermediate
1868: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1870: .keywords: TS, timestep, get, matrix, Jacobian
1871: @*/
1872: PetscErrorCode TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1873: {
1877: TSGetRHSMatrix(ts,J,M,ctx);
1878: return(0);
1879: }
1883: /*@C
1884: TSMonitorSolution - Monitors progress of the TS solvers by calling
1885: VecView() for the solution at each timestep
1887: Collective on TS
1889: Input Parameters:
1890: + ts - the TS context
1891: . step - current time-step
1892: . ptime - current time
1893: - dummy - either a viewer or PETSC_NULL
1895: Level: intermediate
1897: .keywords: TS, vector, monitor, view
1899: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
1900: @*/
1901: PetscErrorCode TSMonitorSolution(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1902: {
1904: PetscViewer viewer = (PetscViewer) dummy;
1907: if (!dummy) {
1908: viewer = PETSC_VIEWER_DRAW_(ts->comm);
1909: }
1910: VecView(x,viewer);
1911: return(0);
1912: }