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_PVODE, 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;
82: PetscOptionsBegin(ts->comm, 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: PetscOptionsName("-ts_monitor","Monitor timestep size","TSDefaultMonitor",&opt);
95: if (opt) {
96: TSSetMonitor(ts,TSDefaultMonitor,PETSC_NULL,PETSC_NULL);
97: }
98: PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
99: if (opt) {
100: TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
101: }
102: PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
103: if (opt) {
104: TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
105: }
107: /* Handle TS type options */
108: TSSetTypeFromOptions(ts);
110: /* Handle specific TS options */
111: if (ts->ops->setfromoptions) {
112: (*ts->ops->setfromoptions)(ts);
113: }
114: PetscOptionsEnd();
116: /* Handle subobject options */
117: switch(ts->problem_type) {
118: /* Should check for implicit/explicit */
119: case TS_LINEAR:
120: if (ts->ksp) {
121: KSPSetOperators(ts->ksp,ts->A,ts->B,DIFFERENT_NONZERO_PATTERN);
122: KSPSetFromOptions(ts->ksp);
123: }
124: break;
125: case TS_NONLINEAR:
126: if (ts->snes) {
127: /* this is a bit of a hack, but it gets the matrix information into SNES earlier
128: so that SNES and KSP have more information to pick reasonable defaults
129: before they allow users to set options */
130: SNESSetJacobian(ts->snes,ts->A,ts->B,0,ts);
131: SNESSetFromOptions(ts->snes);
132: }
133: break;
134: default:
135: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
136: }
138: return(0);
139: }
141: #undef __FUNCT__
143: /*@
144: TSViewFromOptions - This function visualizes the ts based upon user options.
146: Collective on TS
148: Input Parameter:
149: . ts - The ts
151: Level: intermediate
153: .keywords: TS, view, options, database
154: .seealso: TSSetFromOptions(), TSView()
155: @*/
156: PetscErrorCode PETSCTS_DLLEXPORT TSViewFromOptions(TS ts,const char title[])
157: {
158: PetscViewer viewer;
159: PetscDraw draw;
160: PetscTruth opt;
161: char typeName[1024];
162: char fileName[PETSC_MAX_PATH_LEN];
163: size_t len;
167: PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
168: if (opt) {
169: PetscOptionsGetString(ts->prefix, "-ts_view", typeName, 1024, &opt);
170: PetscStrlen(typeName, &len);
171: if (len > 0) {
172: PetscViewerCreate(ts->comm, &viewer);
173: PetscViewerSetType(viewer, typeName);
174: PetscOptionsGetString(ts->prefix, "-ts_view_file", fileName, 1024, &opt);
175: if (opt) {
176: PetscViewerSetFilename(viewer, fileName);
177: } else {
178: PetscViewerSetFilename(viewer, ts->name);
179: }
180: TSView(ts, viewer);
181: PetscViewerFlush(viewer);
182: PetscViewerDestroy(viewer);
183: } else {
184: TSView(ts, PETSC_NULL);
185: }
186: }
187: PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
188: if (opt) {
189: PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
190: PetscViewerDrawGetDraw(viewer, 0, &draw);
191: if (title) {
192: PetscDrawSetTitle(draw, (char *)title);
193: } else {
194: PetscObjectName((PetscObject) ts);
195: PetscDrawSetTitle(draw, ts->name);
196: }
197: TSView(ts, viewer);
198: PetscViewerFlush(viewer);
199: PetscDrawPause(draw);
200: PetscViewerDestroy(viewer);
201: }
202: return(0);
203: }
207: /*@
208: TSComputeRHSJacobian - Computes the Jacobian matrix that has been
209: set with TSSetRHSJacobian().
211: Collective on TS and Vec
213: Input Parameters:
214: + ts - the SNES context
215: . t - current timestep
216: - x - input vector
218: Output Parameters:
219: + A - Jacobian matrix
220: . B - optional preconditioning matrix
221: - flag - flag indicating matrix structure
223: Notes:
224: Most users should not need to explicitly call this routine, as it
225: is used internally within the nonlinear solvers.
227: See KSPSetOperators() for important information about setting the
228: flag parameter.
230: TSComputeJacobian() is valid only for TS_NONLINEAR
232: Level: developer
234: .keywords: SNES, compute, Jacobian, matrix
236: .seealso: TSSetRHSJacobian(), KSPSetOperators()
237: @*/
238: PetscErrorCode PETSCTS_DLLEXPORT TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
239: {
246: if (ts->problem_type != TS_NONLINEAR) {
247: SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
248: }
249: if (ts->ops->rhsjacobian) {
250: PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
251: *flg = DIFFERENT_NONZERO_PATTERN;
252: PetscStackPush("TS user Jacobian function");
253: (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
254: PetscStackPop;
255: PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
256: /* make sure user returned a correct Jacobian and preconditioner */
259: } else {
260: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
261: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
262: if (*A != *B) {
263: MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
264: MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
265: }
266: }
267: return(0);
268: }
272: /*
273: TSComputeRHSFunction - Evaluates the right-hand-side function.
275: Note: If the user did not provide a function but merely a matrix,
276: this routine applies the matrix.
277: */
278: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
279: {
287: PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
288: if (ts->ops->rhsfunction) {
289: PetscStackPush("TS user right-hand-side function");
290: (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
291: PetscStackPop;
292: } else {
293: if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
294: MatStructure flg;
295: PetscStackPush("TS user right-hand-side matrix function");
296: (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
297: PetscStackPop;
298: }
299: MatMult(ts->A,x,y);
300: }
302: /* apply user-provided boundary conditions (only needed if these are time dependent) */
303: TSComputeRHSBoundaryConditions(ts,t,y);
304: PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);
306: return(0);
307: }
311: /*@C
312: TSSetRHSFunction - Sets the routine for evaluating the function,
313: F(t,u), where U_t = F(t,u).
315: Collective on TS
317: Input Parameters:
318: + ts - the TS context obtained from TSCreate()
319: . f - routine for evaluating the right-hand-side function
320: - ctx - [optional] user-defined context for private data for the
321: function evaluation routine (may be PETSC_NULL)
323: Calling sequence of func:
324: $ func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);
326: + t - current timestep
327: . u - input vector
328: . F - function vector
329: - ctx - [optional] user-defined function context
331: Important:
332: The user MUST call either this routine or TSSetRHSMatrix().
334: Level: beginner
336: .keywords: TS, timestep, set, right-hand-side, function
338: .seealso: TSSetRHSMatrix()
339: @*/
340: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
341: {
345: if (ts->problem_type == TS_LINEAR) {
346: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
347: }
348: ts->ops->rhsfunction = f;
349: ts->funP = ctx;
350: return(0);
351: }
355: /*@C
356: TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
357: Also sets the location to store A.
359: Collective on TS
361: Input Parameters:
362: + ts - the TS context obtained from TSCreate()
363: . A - matrix
364: . B - preconditioner matrix (usually same as A)
365: . f - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
366: if A is not a function of t.
367: - ctx - [optional] user-defined context for private data for the
368: matrix evaluation routine (may be PETSC_NULL)
370: Calling sequence of func:
371: $ func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);
373: + t - current timestep
374: . A - matrix A, where U_t = A(t) U
375: . B - preconditioner matrix, usually the same as A
376: . flag - flag indicating information about the preconditioner matrix
377: structure (same as flag in KSPSetOperators())
378: - ctx - [optional] user-defined context for matrix evaluation routine
380: Notes:
381: See KSPSetOperators() for important information about setting the flag
382: output parameter in the routine func(). Be sure to read this information!
384: The routine func() takes Mat * as the matrix arguments rather than Mat.
385: This allows the matrix evaluation routine to replace A and/or B with a
386: completely new new matrix structure (not just different matrix elements)
387: when appropriate, for instance, if the nonzero structure is changing
388: throughout the global iterations.
390: Important:
391: The user MUST call either this routine or TSSetRHSFunction().
393: Level: beginner
395: .keywords: TS, timestep, set, right-hand-side, matrix
397: .seealso: TSSetRHSFunction()
398: @*/
399: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
400: {
407: if (ts->problem_type == TS_NONLINEAR) {
408: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
409: }
411: ts->ops->rhsmatrix = f;
412: ts->jacP = ctx;
413: ts->A = A;
414: ts->B = B;
416: return(0);
417: }
421: /*@C
422: TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
423: where U_t = F(U,t), as well as the location to store the matrix.
424: Use TSSetRHSMatrix() for linear problems.
426: Collective on TS
428: Input Parameters:
429: + ts - the TS context obtained from TSCreate()
430: . A - Jacobian matrix
431: . B - preconditioner matrix (usually same as A)
432: . f - the Jacobian evaluation routine
433: - ctx - [optional] user-defined context for private data for the
434: Jacobian evaluation routine (may be PETSC_NULL)
436: Calling sequence of func:
437: $ func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);
439: + t - current timestep
440: . u - input vector
441: . A - matrix A, where U_t = A(t)u
442: . B - preconditioner matrix, usually the same as A
443: . flag - flag indicating information about the preconditioner matrix
444: structure (same as flag in KSPSetOperators())
445: - ctx - [optional] user-defined context for matrix evaluation routine
447: Notes:
448: See KSPSetOperators() for important information about setting the flag
449: output parameter in the routine func(). Be sure to read this information!
451: The routine func() takes Mat * as the matrix arguments rather than Mat.
452: This allows the matrix evaluation routine to replace A and/or B with a
453: completely new new matrix structure (not just different matrix elements)
454: when appropriate, for instance, if the nonzero structure is changing
455: throughout the global iterations.
457: Level: beginner
458:
459: .keywords: TS, timestep, set, right-hand-side, Jacobian
461: .seealso: TSDefaultComputeJacobianColor(),
462: SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()
464: @*/
465: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
466: {
473: if (ts->problem_type != TS_NONLINEAR) {
474: SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
475: }
477: ts->ops->rhsjacobian = f;
478: ts->jacP = ctx;
479: ts->A = A;
480: ts->B = B;
481: return(0);
482: }
486: /*
487: TSComputeRHSBoundaryConditions - Evaluates the boundary condition function.
489: Note: If the user did not provide a function but merely a matrix,
490: this routine applies the matrix.
491: */
492: PetscErrorCode TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
493: {
501: if (ts->ops->rhsbc) {
502: PetscStackPush("TS user boundary condition function");
503: (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
504: PetscStackPop;
505: return(0);
506: }
508: return(0);
509: }
513: /*@C
514: TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
515: boundary conditions for the function F.
517: Collective on TS
519: Input Parameters:
520: + ts - the TS context obtained from TSCreate()
521: . f - routine for evaluating the boundary condition function
522: - ctx - [optional] user-defined context for private data for the
523: function evaluation routine (may be PETSC_NULL)
525: Calling sequence of func:
526: $ func (TS ts,PetscReal t,Vec F,void *ctx);
528: + t - current timestep
529: . F - function vector
530: - ctx - [optional] user-defined function context
532: Level: intermediate
534: .keywords: TS, timestep, set, boundary conditions, function
535: @*/
536: PetscErrorCode PETSCTS_DLLEXPORT TSSetRHSBoundaryConditions(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
537: {
541: if (ts->problem_type != TS_LINEAR) {
542: SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
543: }
544: ts->ops->rhsbc = f;
545: ts->bcP = ctx;
546: return(0);
547: }
551: /*@C
552: TSView - Prints the TS data structure.
554: Collective on TS
556: Input Parameters:
557: + ts - the TS context obtained from TSCreate()
558: - viewer - visualization context
560: Options Database Key:
561: . -ts_view - calls TSView() at end of TSStep()
563: Notes:
564: The available visualization contexts include
565: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
566: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
567: output where only the first processor opens
568: the file. All other processors send their
569: data to the first processor to print.
571: The user can open an alternative visualization context with
572: PetscViewerASCIIOpen() - output to a specified file.
574: Level: beginner
576: .keywords: TS, timestep, view
578: .seealso: PetscViewerASCIIOpen()
579: @*/
580: PetscErrorCode PETSCTS_DLLEXPORT TSView(TS ts,PetscViewer viewer)
581: {
583: char *type;
584: PetscTruth iascii,isstring;
588: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
592: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
593: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
594: if (iascii) {
595: PetscViewerASCIIPrintf(viewer,"TS Object:\n");
596: TSGetType(ts,(TSType *)&type);
597: if (type) {
598: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
599: } else {
600: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
601: }
602: if (ts->ops->view) {
603: PetscViewerASCIIPushTab(viewer);
604: (*ts->ops->view)(ts,viewer);
605: PetscViewerASCIIPopTab(viewer);
606: }
607: PetscViewerASCIIPrintf(viewer," maximum steps=%D\n",ts->max_steps);
608: PetscViewerASCIIPrintf(viewer," maximum time=%g\n",ts->max_time);
609: if (ts->problem_type == TS_NONLINEAR) {
610: PetscViewerASCIIPrintf(viewer," total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
611: }
612: PetscViewerASCIIPrintf(viewer," total number of linear solver iterations=%D\n",ts->linear_its);
613: } else if (isstring) {
614: TSGetType(ts,(TSType *)&type);
615: PetscViewerStringSPrintf(viewer," %-7.7s",type);
616: }
617: PetscViewerASCIIPushTab(viewer);
618: if (ts->ksp) {KSPView(ts->ksp,viewer);}
619: if (ts->snes) {SNESView(ts->snes,viewer);}
620: PetscViewerASCIIPopTab(viewer);
621: return(0);
622: }
627: /*@C
628: TSSetApplicationContext - Sets an optional user-defined context for
629: the timesteppers.
631: Collective on TS
633: Input Parameters:
634: + ts - the TS context obtained from TSCreate()
635: - usrP - optional user context
637: Level: intermediate
639: .keywords: TS, timestep, set, application, context
641: .seealso: TSGetApplicationContext()
642: @*/
643: PetscErrorCode PETSCTS_DLLEXPORT TSSetApplicationContext(TS ts,void *usrP)
644: {
647: ts->user = usrP;
648: return(0);
649: }
653: /*@C
654: TSGetApplicationContext - Gets the user-defined context for the
655: timestepper.
657: Not Collective
659: Input Parameter:
660: . ts - the TS context obtained from TSCreate()
662: Output Parameter:
663: . usrP - user context
665: Level: intermediate
667: .keywords: TS, timestep, get, application, context
669: .seealso: TSSetApplicationContext()
670: @*/
671: PetscErrorCode PETSCTS_DLLEXPORT TSGetApplicationContext(TS ts,void **usrP)
672: {
675: *usrP = ts->user;
676: return(0);
677: }
681: /*@
682: TSGetTimeStepNumber - Gets the current number of timesteps.
684: Not Collective
686: Input Parameter:
687: . ts - the TS context obtained from TSCreate()
689: Output Parameter:
690: . iter - number steps so far
692: Level: intermediate
694: .keywords: TS, timestep, get, iteration, number
695: @*/
696: PetscErrorCode PETSCTS_DLLEXPORT TSGetTimeStepNumber(TS ts,PetscInt* iter)
697: {
701: *iter = ts->steps;
702: return(0);
703: }
707: /*@
708: TSSetInitialTimeStep - Sets the initial timestep to be used,
709: as well as the initial time.
711: Collective on TS
713: Input Parameters:
714: + ts - the TS context obtained from TSCreate()
715: . initial_time - the initial time
716: - time_step - the size of the timestep
718: Level: intermediate
720: .seealso: TSSetTimeStep(), TSGetTimeStep()
722: .keywords: TS, set, initial, timestep
723: @*/
724: PetscErrorCode PETSCTS_DLLEXPORT TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
725: {
728: ts->time_step = time_step;
729: ts->initial_time_step = time_step;
730: ts->ptime = initial_time;
731: return(0);
732: }
736: /*@
737: TSSetTimeStep - Allows one to reset the timestep at any time,
738: useful for simple pseudo-timestepping codes.
740: Collective on TS
742: Input Parameters:
743: + ts - the TS context obtained from TSCreate()
744: - time_step - the size of the timestep
746: Level: intermediate
748: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
750: .keywords: TS, set, timestep
751: @*/
752: PetscErrorCode PETSCTS_DLLEXPORT TSSetTimeStep(TS ts,PetscReal time_step)
753: {
756: ts->time_step = time_step;
757: return(0);
758: }
762: /*@
763: TSGetTimeStep - Gets the current timestep size.
765: Not Collective
767: Input Parameter:
768: . ts - the TS context obtained from TSCreate()
770: Output Parameter:
771: . dt - the current timestep size
773: Level: intermediate
775: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
777: .keywords: TS, get, timestep
778: @*/
779: PetscErrorCode PETSCTS_DLLEXPORT TSGetTimeStep(TS ts,PetscReal* dt)
780: {
784: *dt = ts->time_step;
785: return(0);
786: }
790: /*@C
791: TSGetSolution - Returns the solution at the present timestep. It
792: is valid to call this routine inside the function that you are evaluating
793: in order to move to the new timestep. This vector not changed until
794: the solution at the next timestep has been calculated.
796: Not Collective, but Vec returned is parallel if TS is parallel
798: Input Parameter:
799: . ts - the TS context obtained from TSCreate()
801: Output Parameter:
802: . v - the vector containing the solution
804: Level: intermediate
806: .seealso: TSGetTimeStep()
808: .keywords: TS, timestep, get, solution
809: @*/
810: PetscErrorCode PETSCTS_DLLEXPORT TSGetSolution(TS ts,Vec *v)
811: {
815: *v = ts->vec_sol_always;
816: return(0);
817: }
819: /* ----- Routines to initialize and destroy a timestepper ---- */
822: /*@C
823: TSSetProblemType - Sets the type of problem to be solved.
825: Not collective
827: Input Parameters:
828: + ts - The TS
829: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
830: .vb
831: U_t = A U
832: U_t = A(t) U
833: U_t = F(t,U)
834: .ve
836: Level: beginner
838: .keywords: TS, problem type
839: .seealso: TSSetUp(), TSProblemType, TS
840: @*/
841: PetscErrorCode PETSCTS_DLLEXPORT TSSetProblemType(TS ts, TSProblemType type)
842: {
845: ts->problem_type = type;
846: return(0);
847: }
851: /*@C
852: TSGetProblemType - Gets the type of problem to be solved.
854: Not collective
856: Input Parameter:
857: . ts - The TS
859: Output Parameter:
860: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
861: .vb
862: U_t = A U
863: U_t = A(t) U
864: U_t = F(t,U)
865: .ve
867: Level: beginner
869: .keywords: TS, problem type
870: .seealso: TSSetUp(), TSProblemType, TS
871: @*/
872: PetscErrorCode PETSCTS_DLLEXPORT TSGetProblemType(TS ts, TSProblemType *type)
873: {
877: *type = ts->problem_type;
878: return(0);
879: }
883: /*@
884: TSSetUp - Sets up the internal data structures for the later use
885: of a timestepper.
887: Collective on TS
889: Input Parameter:
890: . ts - the TS context obtained from TSCreate()
892: Notes:
893: For basic use of the TS solvers the user need not explicitly call
894: TSSetUp(), since these actions will automatically occur during
895: the call to TSStep(). However, if one wishes to control this
896: phase separately, TSSetUp() should be called after TSCreate()
897: and optional routines of the form TSSetXXX(), but before TSStep().
899: Level: advanced
901: .keywords: TS, timestep, setup
903: .seealso: TSCreate(), TSStep(), TSDestroy()
904: @*/
905: PetscErrorCode PETSCTS_DLLEXPORT TSSetUp(TS ts)
906: {
911: if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
912: if (!ts->type_name) {
913: TSSetType(ts,TS_EULER);
914: }
915: (*ts->ops->setup)(ts);
916: ts->setupcalled = 1;
917: return(0);
918: }
922: /*@C
923: TSDestroy - Destroys the timestepper context that was created
924: with TSCreate().
926: Collective on TS
928: Input Parameter:
929: . ts - the TS context obtained from TSCreate()
931: Level: beginner
933: .keywords: TS, timestepper, destroy
935: .seealso: TSCreate(), TSSetUp(), TSSolve()
936: @*/
937: PetscErrorCode PETSCTS_DLLEXPORT TSDestroy(TS ts)
938: {
943: if (--ts->refct > 0) return(0);
945: /* if memory was published with AMS then destroy it */
946: PetscObjectDepublish(ts);
948: if (ts->ksp) {KSPDestroy(ts->ksp);}
949: if (ts->snes) {SNESDestroy(ts->snes);}
950: (*(ts)->ops->destroy)(ts);
951: TSClearMonitor(ts);
952: PetscHeaderDestroy(ts);
953: return(0);
954: }
958: /*@C
959: TSGetSNES - Returns the SNES (nonlinear solver) associated with
960: a TS (timestepper) context. Valid only for nonlinear problems.
962: Not Collective, but SNES is parallel if TS is parallel
964: Input Parameter:
965: . ts - the TS context obtained from TSCreate()
967: Output Parameter:
968: . snes - the nonlinear solver context
970: Notes:
971: The user can then directly manipulate the SNES context to set various
972: options, etc. Likewise, the user can then extract and manipulate the
973: KSP, KSP, and PC contexts as well.
975: TSGetSNES() does not work for integrators that do not use SNES; in
976: this case TSGetSNES() returns PETSC_NULL in snes.
978: Level: beginner
980: .keywords: timestep, get, SNES
981: @*/
982: PetscErrorCode PETSCTS_DLLEXPORT TSGetSNES(TS ts,SNES *snes)
983: {
987: if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
988: *snes = ts->snes;
989: return(0);
990: }
994: /*@C
995: TSGetKSP - Returns the KSP (linear solver) associated with
996: a TS (timestepper) context.
998: Not Collective, but KSP is parallel if TS is parallel
1000: Input Parameter:
1001: . ts - the TS context obtained from TSCreate()
1003: Output Parameter:
1004: . ksp - the nonlinear solver context
1006: Notes:
1007: The user can then directly manipulate the KSP context to set various
1008: options, etc. Likewise, the user can then extract and manipulate the
1009: KSP and PC contexts as well.
1011: TSGetKSP() does not work for integrators that do not use KSP;
1012: in this case TSGetKSP() returns PETSC_NULL in ksp.
1014: Level: beginner
1016: .keywords: timestep, get, KSP
1017: @*/
1018: PetscErrorCode PETSCTS_DLLEXPORT TSGetKSP(TS ts,KSP *ksp)
1019: {
1023: if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1024: *ksp = ts->ksp;
1025: return(0);
1026: }
1028: /* ----------- Routines to set solver parameters ---------- */
1032: /*@
1033: TSGetDuration - Gets the maximum number of timesteps to use and
1034: maximum time for iteration.
1036: Collective on TS
1038: Input Parameters:
1039: + ts - the TS context obtained from TSCreate()
1040: . maxsteps - maximum number of iterations to use, or PETSC_NULL
1041: - maxtime - final time to iterate to, or PETSC_NULL
1043: Level: intermediate
1045: .keywords: TS, timestep, get, maximum, iterations, time
1046: @*/
1047: PetscErrorCode PETSCTS_DLLEXPORT TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1048: {
1051: if (maxsteps) {
1053: *maxsteps = ts->max_steps;
1054: }
1055: if (maxtime ) {
1057: *maxtime = ts->max_time;
1058: }
1059: return(0);
1060: }
1064: /*@
1065: TSSetDuration - Sets the maximum number of timesteps to use and
1066: maximum time for iteration.
1068: Collective on TS
1070: Input Parameters:
1071: + ts - the TS context obtained from TSCreate()
1072: . maxsteps - maximum number of iterations to use
1073: - maxtime - final time to iterate to
1075: Options Database Keys:
1076: . -ts_max_steps <maxsteps> - Sets maxsteps
1077: . -ts_max_time <maxtime> - Sets maxtime
1079: Notes:
1080: The default maximum number of iterations is 5000. Default time is 5.0
1082: Level: intermediate
1084: .keywords: TS, timestep, set, maximum, iterations
1085: @*/
1086: PetscErrorCode PETSCTS_DLLEXPORT TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1087: {
1090: ts->max_steps = maxsteps;
1091: ts->max_time = maxtime;
1092: return(0);
1093: }
1097: /*@
1098: TSSetSolution - Sets the initial solution vector
1099: for use by the TS routines.
1101: Collective on TS and Vec
1103: Input Parameters:
1104: + ts - the TS context obtained from TSCreate()
1105: - x - the solution vector
1107: Level: beginner
1109: .keywords: TS, timestep, set, solution, initial conditions
1110: @*/
1111: PetscErrorCode PETSCTS_DLLEXPORT TSSetSolution(TS ts,Vec x)
1112: {
1116: ts->vec_sol = ts->vec_sol_always = x;
1117: return(0);
1118: }
1122: /*@
1123: TSDefaultRhsBC - The default boundary condition function which does nothing.
1125: Collective on TS
1127: Input Parameters:
1128: + ts - The TS context obtained from TSCreate()
1129: . rhs - The Rhs
1130: - ctx - The user-context
1132: Level: developer
1134: .keywords: TS, Rhs, boundary conditions
1135: @*/
1136: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultRhsBC(TS ts, Vec rhs, void *ctx)
1137: {
1139: return(0);
1140: }
1144: /*@C
1145: TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1146: to the system matrix and preconditioner of each system.
1148: Collective on TS
1150: Input Parameters:
1151: + ts - The TS context obtained from TSCreate()
1152: - func - The function
1154: Calling sequence of func:
1155: . func (TS ts, Mat A, Mat B, void *ctx);
1157: + A - The current system matrix
1158: . B - The current preconditioner
1159: - ctx - The user-context
1161: Level: intermediate
1163: .keywords: TS, System matrix, boundary conditions
1164: @*/
1165: PetscErrorCode PETSCTS_DLLEXPORT TSSetSystemMatrixBC(TS ts, PetscErrorCode (*func)(TS, Mat, Mat, void *))
1166: {
1169: ts->ops->applymatrixbc = func;
1170: return(0);
1171: }
1175: /*@
1176: TSDefaultSystemMatrixBC - The default boundary condition function which
1177: does nothing.
1179: Collective on TS
1181: Input Parameters:
1182: + ts - The TS context obtained from TSCreate()
1183: . A - The system matrix
1184: . B - The preconditioner
1185: - ctx - The user-context
1187: Level: developer
1189: .keywords: TS, System matrix, boundary conditions
1190: @*/
1191: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1192: {
1194: return(0);
1195: }
1199: /*@C
1200: TSSetSolutionBC - Sets the function which applies boundary conditions
1201: to the solution of each system. This is necessary in nonlinear systems
1202: which time dependent boundary conditions.
1204: Collective on TS
1206: Input Parameters:
1207: + ts - The TS context obtained from TSCreate()
1208: - func - The function
1210: Calling sequence of func:
1211: . func (TS ts, Vec rsol, void *ctx);
1213: + sol - The current solution vector
1214: - ctx - The user-context
1216: Level: intermediate
1218: .keywords: TS, solution, boundary conditions
1219: @*/
1220: PetscErrorCode PETSCTS_DLLEXPORT TSSetSolutionBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1221: {
1224: ts->ops->applysolbc = func;
1225: return(0);
1226: }
1230: /*@
1231: TSDefaultSolutionBC - The default boundary condition function which
1232: does nothing.
1234: Collective on TS
1236: Input Parameters:
1237: + ts - The TS context obtained from TSCreate()
1238: . sol - The solution
1239: - ctx - The user-context
1241: Level: developer
1243: .keywords: TS, solution, boundary conditions
1244: @*/
1245: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1246: {
1248: return(0);
1249: }
1253: /*@C
1254: TSSetPreStep - Sets the general-purpose function
1255: called once at the beginning of time stepping.
1257: Collective on TS
1259: Input Parameters:
1260: + ts - The TS context obtained from TSCreate()
1261: - func - The function
1263: Calling sequence of func:
1264: . func (TS ts);
1266: Level: intermediate
1268: .keywords: TS, timestep
1269: @*/
1270: PetscErrorCode PETSCTS_DLLEXPORT TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1271: {
1274: ts->ops->prestep = func;
1275: return(0);
1276: }
1280: /*@
1281: TSDefaultPreStep - The default pre-stepping function which does nothing.
1283: Collective on TS
1285: Input Parameters:
1286: . ts - The TS context obtained from TSCreate()
1288: Level: developer
1290: .keywords: TS, timestep
1291: @*/
1292: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultPreStep(TS ts)
1293: {
1295: return(0);
1296: }
1300: /*@C
1301: TSSetUpdate - Sets the general-purpose update function called
1302: at the beginning of every time step. This function can change
1303: the time step.
1305: Collective on TS
1307: Input Parameters:
1308: + ts - The TS context obtained from TSCreate()
1309: - func - The function
1311: Calling sequence of func:
1312: . func (TS ts, double t, double *dt);
1314: + t - The current time
1315: - dt - The current time step
1317: Level: intermediate
1319: .keywords: TS, update, timestep
1320: @*/
1321: PetscErrorCode PETSCTS_DLLEXPORT TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1322: {
1325: ts->ops->update = func;
1326: return(0);
1327: }
1331: /*@
1332: TSDefaultUpdate - The default update function which does nothing.
1334: Collective on TS
1336: Input Parameters:
1337: + ts - The TS context obtained from TSCreate()
1338: - t - The current time
1340: Output Parameters:
1341: . dt - The current time step
1343: Level: developer
1345: .keywords: TS, update, timestep
1346: @*/
1347: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1348: {
1350: return(0);
1351: }
1355: /*@C
1356: TSSetPostStep - Sets the general-purpose function
1357: called once at the end of time stepping.
1359: Collective on TS
1361: Input Parameters:
1362: + ts - The TS context obtained from TSCreate()
1363: - func - The function
1365: Calling sequence of func:
1366: . func (TS ts);
1368: Level: intermediate
1370: .keywords: TS, timestep
1371: @*/
1372: PetscErrorCode PETSCTS_DLLEXPORT TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1373: {
1376: ts->ops->poststep = func;
1377: return(0);
1378: }
1382: /*@
1383: TSDefaultPostStep - The default post-stepping function which does nothing.
1385: Collective on TS
1387: Input Parameters:
1388: . ts - The TS context obtained from TSCreate()
1390: Level: developer
1392: .keywords: TS, timestep
1393: @*/
1394: PetscErrorCode PETSCTS_DLLEXPORT TSDefaultPostStep(TS ts)
1395: {
1397: return(0);
1398: }
1400: /* ------------ Routines to set performance monitoring options ----------- */
1404: /*@C
1405: TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1406: timestep to display the iteration's progress.
1408: Collective on TS
1410: Input Parameters:
1411: + ts - the TS context obtained from TSCreate()
1412: . func - monitoring routine
1413: . mctx - [optional] user-defined context for private data for the
1414: monitor routine (use PETSC_NULL if no context is desired)
1415: - monitordestroy - [optional] routine that frees monitor context
1416: (may be PETSC_NULL)
1418: Calling sequence of func:
1419: $ int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)
1421: + ts - the TS context
1422: . steps - iteration number
1423: . time - current time
1424: . x - current iterate
1425: - mctx - [optional] monitoring context
1427: Notes:
1428: This routine adds an additional monitor to the list of monitors that
1429: already has been loaded.
1431: Level: intermediate
1433: .keywords: TS, timestep, set, monitor
1435: .seealso: TSDefaultMonitor(), TSClearMonitor()
1436: @*/
1437: PetscErrorCode PETSCTS_DLLEXPORT TSSetMonitor(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1438: {
1441: if (ts->numbermonitors >= MAXTSMONITORS) {
1442: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1443: }
1444: ts->monitor[ts->numbermonitors] = monitor;
1445: ts->mdestroy[ts->numbermonitors] = mdestroy;
1446: ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
1447: return(0);
1448: }
1452: /*@C
1453: TSClearMonitor - Clears all the monitors that have been set on a time-step object.
1455: Collective on TS
1457: Input Parameters:
1458: . ts - the TS context obtained from TSCreate()
1460: Notes:
1461: There is no way to remove a single, specific monitor.
1463: Level: intermediate
1465: .keywords: TS, timestep, set, monitor
1467: .seealso: TSDefaultMonitor(), TSSetMonitor()
1468: @*/
1469: PetscErrorCode PETSCTS_DLLEXPORT TSClearMonitor(TS ts)
1470: {
1472: PetscInt i;
1476: for (i=0; i<ts->numbermonitors; i++) {
1477: if (ts->mdestroy[i]) {
1478: (*ts->mdestroy[i])(ts->monitorcontext[i]);
1479: }
1480: }
1481: ts->numbermonitors = 0;
1482: return(0);
1483: }
1487: PetscErrorCode TSDefaultMonitor(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1488: {
1492: PetscPrintf(ts->comm,"timestep %D dt %g time %g\n",step,ts->time_step,ptime);
1493: return(0);
1494: }
1498: /*@
1499: TSStep - Steps the requested number of timesteps.
1501: Collective on TS
1503: Input Parameter:
1504: . ts - the TS context obtained from TSCreate()
1506: Output Parameters:
1507: + steps - number of iterations until termination
1508: - ptime - time until termination
1510: Level: beginner
1512: .keywords: TS, timestep, solve
1514: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1515: @*/
1516: PetscErrorCode PETSCTS_DLLEXPORT TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1517: {
1522: if (!ts->setupcalled) {
1523: TSSetUp(ts);
1524: }
1526: PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1527: (*ts->ops->prestep)(ts);
1528: (*ts->ops->step)(ts, steps, ptime);
1529: (*ts->ops->poststep)(ts);
1530: PetscLogEventEnd(TS_Step, ts, 0, 0, 0);
1532: if (!PetscPreLoadingOn) {
1533: TSViewFromOptions(ts,ts->name);
1534: }
1535: return(0);
1536: }
1540: /*
1541: Runs the user provided monitor routines, if they exists.
1542: */
1543: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1544: {
1546: PetscInt i,n = ts->numbermonitors;
1549: for (i=0; i<n; i++) {
1550: (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1551: }
1552: return(0);
1553: }
1555: /* ------------------------------------------------------------------------*/
1559: /*@C
1560: TSLGMonitorCreate - Creates a line graph context for use with
1561: TS to monitor convergence of preconditioned residual norms.
1563: Collective on TS
1565: Input Parameters:
1566: + host - the X display to open, or null for the local machine
1567: . label - the title to put in the title bar
1568: . x, y - the screen coordinates of the upper left coordinate of the window
1569: - m, n - the screen width and height in pixels
1571: Output Parameter:
1572: . draw - the drawing context
1574: Options Database Key:
1575: . -ts_xmonitor - automatically sets line graph monitor
1577: Notes:
1578: Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().
1580: Level: intermediate
1582: .keywords: TS, monitor, line graph, residual, seealso
1584: .seealso: TSLGMonitorDestroy(), TSSetMonitor()
1586: @*/
1587: PetscErrorCode PETSCTS_DLLEXPORT TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1588: {
1589: PetscDraw win;
1593: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1594: PetscDrawSetType(win,PETSC_DRAW_X);
1595: PetscDrawLGCreate(win,1,draw);
1596: PetscDrawLGIndicateDataPoints(*draw);
1598: PetscLogObjectParent(*draw,win);
1599: return(0);
1600: }
1604: PetscErrorCode TSLGMonitor(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1605: {
1606: PetscDrawLG lg = (PetscDrawLG) monctx;
1607: PetscReal x,y = ptime;
1611: if (!monctx) {
1612: MPI_Comm comm;
1613: PetscViewer viewer;
1615: PetscObjectGetComm((PetscObject)ts,&comm);
1616: viewer = PETSC_VIEWER_DRAW_(comm);
1617: PetscViewerDrawGetDrawLG(viewer,0,&lg);
1618: }
1620: if (!n) {PetscDrawLGReset(lg);}
1621: x = (PetscReal)n;
1622: PetscDrawLGAddPoint(lg,&x,&y);
1623: if (n < 20 || (n % 5)) {
1624: PetscDrawLGDraw(lg);
1625: }
1626: return(0);
1627: }
1631: /*@C
1632: TSLGMonitorDestroy - Destroys a line graph context that was created
1633: with TSLGMonitorCreate().
1635: Collective on PetscDrawLG
1637: Input Parameter:
1638: . draw - the drawing context
1640: Level: intermediate
1642: .keywords: TS, monitor, line graph, destroy
1644: .seealso: TSLGMonitorCreate(), TSSetMonitor(), TSLGMonitor();
1645: @*/
1646: PetscErrorCode PETSCTS_DLLEXPORT TSLGMonitorDestroy(PetscDrawLG drawlg)
1647: {
1648: PetscDraw draw;
1652: PetscDrawLGGetDraw(drawlg,&draw);
1653: PetscDrawDestroy(draw);
1654: PetscDrawLGDestroy(drawlg);
1655: return(0);
1656: }
1660: /*@
1661: TSGetTime - Gets the current time.
1663: Not Collective
1665: Input Parameter:
1666: . ts - the TS context obtained from TSCreate()
1668: Output Parameter:
1669: . t - the current time
1671: Contributed by: Matthew Knepley
1673: Level: beginner
1675: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()
1677: .keywords: TS, get, time
1678: @*/
1679: PetscErrorCode PETSCTS_DLLEXPORT TSGetTime(TS ts,PetscReal* t)
1680: {
1684: *t = ts->ptime;
1685: return(0);
1686: }
1690: /*@C
1691: TSSetOptionsPrefix - Sets the prefix used for searching for all
1692: TS options in the database.
1694: Collective on TS
1696: Input Parameter:
1697: + ts - The TS context
1698: - prefix - The prefix to prepend to all option names
1700: Notes:
1701: A hyphen (-) must NOT be given at the beginning of the prefix name.
1702: The first character of all runtime options is AUTOMATICALLY the
1703: hyphen.
1705: Contributed by: Matthew Knepley
1707: Level: advanced
1709: .keywords: TS, set, options, prefix, database
1711: .seealso: TSSetFromOptions()
1713: @*/
1714: PetscErrorCode PETSCTS_DLLEXPORT TSSetOptionsPrefix(TS ts,const char prefix[])
1715: {
1720: PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1721: switch(ts->problem_type) {
1722: case TS_NONLINEAR:
1723: SNESSetOptionsPrefix(ts->snes,prefix);
1724: break;
1725: case TS_LINEAR:
1726: KSPSetOptionsPrefix(ts->ksp,prefix);
1727: break;
1728: }
1729: return(0);
1730: }
1735: /*@C
1736: TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1737: TS options in the database.
1739: Collective on TS
1741: Input Parameter:
1742: + ts - The TS context
1743: - prefix - The prefix to prepend to all option names
1745: Notes:
1746: A hyphen (-) must NOT be given at the beginning of the prefix name.
1747: The first character of all runtime options is AUTOMATICALLY the
1748: hyphen.
1750: Contributed by: Matthew Knepley
1752: Level: advanced
1754: .keywords: TS, append, options, prefix, database
1756: .seealso: TSGetOptionsPrefix()
1758: @*/
1759: PetscErrorCode PETSCTS_DLLEXPORT TSAppendOptionsPrefix(TS ts,const char prefix[])
1760: {
1765: PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1766: switch(ts->problem_type) {
1767: case TS_NONLINEAR:
1768: SNESAppendOptionsPrefix(ts->snes,prefix);
1769: break;
1770: case TS_LINEAR:
1771: KSPAppendOptionsPrefix(ts->ksp,prefix);
1772: break;
1773: }
1774: return(0);
1775: }
1779: /*@C
1780: TSGetOptionsPrefix - Sets the prefix used for searching for all
1781: TS options in the database.
1783: Not Collective
1785: Input Parameter:
1786: . ts - The TS context
1788: Output Parameter:
1789: . prefix - A pointer to the prefix string used
1791: Contributed by: Matthew Knepley
1793: Notes: On the fortran side, the user should pass in a string 'prifix' of
1794: sufficient length to hold the prefix.
1796: Level: intermediate
1798: .keywords: TS, get, options, prefix, database
1800: .seealso: TSAppendOptionsPrefix()
1801: @*/
1802: PetscErrorCode PETSCTS_DLLEXPORT TSGetOptionsPrefix(TS ts,const char *prefix[])
1803: {
1809: PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1810: return(0);
1811: }
1815: /*@C
1816: TSGetRHSMatrix - Returns the matrix A at the present timestep.
1818: Not Collective, but parallel objects are returned if TS is parallel
1820: Input Parameter:
1821: . ts - The TS context obtained from TSCreate()
1823: Output Parameters:
1824: + A - The matrix A, where U_t = A(t) U
1825: . M - The preconditioner matrix, usually the same as A
1826: - ctx - User-defined context for matrix evaluation routine
1828: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1830: Contributed by: Matthew Knepley
1832: Level: intermediate
1834: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()
1836: .keywords: TS, timestep, get, matrix
1838: @*/
1839: PetscErrorCode PETSCTS_DLLEXPORT TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1840: {
1843: if (A) *A = ts->A;
1844: if (M) *M = ts->B;
1845: if (ctx) *ctx = ts->jacP;
1846: return(0);
1847: }
1851: /*@C
1852: TSGetRHSJacobian - Returns the Jacobian J at the present timestep.
1854: Not Collective, but parallel objects are returned if TS is parallel
1856: Input Parameter:
1857: . ts - The TS context obtained from TSCreate()
1859: Output Parameters:
1860: + J - The Jacobian J of F, where U_t = F(U,t)
1861: . M - The preconditioner matrix, usually the same as J
1862: - ctx - User-defined context for Jacobian evaluation routine
1864: Notes: You can pass in PETSC_NULL for any return argument you do not need.
1866: Contributed by: Matthew Knepley
1868: Level: intermediate
1870: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()
1872: .keywords: TS, timestep, get, matrix, Jacobian
1873: @*/
1874: PetscErrorCode PETSCTS_DLLEXPORT TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1875: {
1879: TSGetRHSMatrix(ts,J,M,ctx);
1880: return(0);
1881: }
1885: /*@C
1886: TSVecViewMonitor - Monitors progress of the TS solvers by calling
1887: VecView() for the solution at each timestep
1889: Collective on TS
1891: Input Parameters:
1892: + ts - the TS context
1893: . step - current time-step
1894: . ptime - current time
1895: - dummy - either a viewer or PETSC_NULL
1897: Level: intermediate
1899: .keywords: TS, vector, monitor, view
1901: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1902: @*/
1903: PetscErrorCode PETSCTS_DLLEXPORT TSVecViewMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1904: {
1906: PetscViewer viewer = (PetscViewer) dummy;
1909: if (!viewer) {
1910: MPI_Comm comm;
1911: PetscObjectGetComm((PetscObject)ts,&comm);
1912: viewer = PETSC_VIEWER_DRAW_(comm);
1913: }
1914: VecView(x,viewer);
1915: return(0);
1916: }