Actual source code: ts.c

  1: #include "src/ts/tsimpl.h"        /*I "petscts.h"  I*/

  3: /* Logging support */
  4: PetscCookie TS_COOKIE = 0;
  5: PetscEvent  TS_Step = 0, TS_PseudoComputeTimeStep = 0, TS_FunctionEval = 0, TS_JacobianEval = 0;

  9: /*
 10:   TSSetTypeFromOptions - Sets the type of ts from user options.

 12:   Collective on TS

 14:   Input Parameter:
 15: . ts - The ts

 17:   Level: intermediate

 19: .keywords: TS, set, options, database, type
 20: .seealso: TSSetFromOptions(), TSSetType()
 21: */
 22: static PetscErrorCode TSSetTypeFromOptions(TS ts)
 23: {
 24:   PetscTruth     opt;
 25:   const char     *defaultType;
 26:   char           typeName[256];

 30:   if (ts->type_name != PETSC_NULL) {
 31:     defaultType = ts->type_name;
 32:   } else {
 33:     defaultType = TS_EULER;
 34:   }

 36:   if (!TSRegisterAllCalled) {
 37:     TSRegisterAll(PETSC_NULL);
 38:   }
 39:   PetscOptionsList("-ts_type", "TS method"," TSSetType", TSList, defaultType, typeName, 256, &opt);
 40:   if (opt) {
 41:     TSSetType(ts, typeName);
 42:   } else {
 43:     TSSetType(ts, defaultType);
 44:   }
 45:   return(0);
 46: }

 50: /*@
 51:    TSSetFromOptions - Sets various TS parameters from user options.

 53:    Collective on TS

 55:    Input Parameter:
 56: .  ts - the TS context obtained from TSCreate()

 58:    Options Database Keys:
 59: +  -ts_type <type> - TS_EULER, TS_BEULER, TS_PVODE, TS_PSEUDO, TS_CRANK_NICHOLSON
 60: .  -ts_max_steps maxsteps - maximum number of time-steps to take
 61: .  -ts_max_time time - maximum time to compute to
 62: .  -ts_dt dt - initial time step
 63: .  -ts_monitor - print information at each timestep
 64: -  -ts_xmonitor - plot information at each timestep

 66:    Level: beginner

 68: .keywords: TS, timestep, set, options, database

 70: .seealso: TSGetType
 71: @*/
 72: PetscErrorCode TSSetFromOptions(TS ts)
 73: {
 74:   PetscReal      dt;
 75:   PetscTruth     opt;

 80:   PetscOptionsBegin(ts->comm, ts->prefix, "Time step options", "TS");

 82:   /* Handle generic TS options */
 83:   PetscOptionsInt("-ts_max_steps","Maximum number of time steps","TSSetDuration",ts->max_steps,&ts->max_steps,PETSC_NULL);
 84:   PetscOptionsReal("-ts_max_time","Time to run to","TSSetDuration",ts->max_time,&ts->max_time,PETSC_NULL);
 85:   PetscOptionsReal("-ts_init_time","Initial time","TSSetInitialTime", ts->ptime, &ts->ptime, PETSC_NULL);
 86:   PetscOptionsReal("-ts_dt","Initial time step","TSSetInitialTimeStep",ts->initial_time_step,&dt,&opt);
 87:   if (opt) {
 88:     ts->initial_time_step = ts->time_step = dt;
 89:   }

 91:   /* Monitor options */
 92:     PetscOptionsName("-ts_monitor","Monitor timestep size","TSDefaultMonitor",&opt);
 93:     if (opt) {
 94:       TSSetMonitor(ts,TSDefaultMonitor,PETSC_NULL,PETSC_NULL);
 95:     }
 96:     PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
 97:     if (opt) {
 98:       TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
 99:     }
100:     PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
101:     if (opt) {
102:       TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
103:     }

105:   /* Handle TS type options */
106:   TSSetTypeFromOptions(ts);

108:   /* Handle specific TS options */
109:   if (ts->ops->setfromoptions != PETSC_NULL) {
110:     (*ts->ops->setfromoptions)(ts);
111:   }
112:   PetscOptionsEnd();

114:   /* Handle subobject options */
115:   switch(ts->problem_type) {
116:     /* Should check for implicit/explicit */
117:   case TS_LINEAR:
118:     if (ts->ksp != PETSC_NULL) {
119:       KSPSetOperators(ts->ksp,ts->A,ts->B,DIFFERENT_NONZERO_PATTERN);
120:       KSPSetFromOptions(ts->ksp);
121:     }
122:     break;
123:   case TS_NONLINEAR:
124:     if (ts->snes != PETSC_NULL) {
125:       /* this is a bit of a hack, but it gets the matrix information into SNES earlier
126:          so that SNES and KSP have more information to pick reasonable defaults
127:          before they allow users to set options */
128:       SNESSetJacobian(ts->snes,ts->A,ts->B,0,ts);
129:       SNESSetFromOptions(ts->snes);
130:     }
131:     break;
132:   default:
133:     SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", (int)ts->problem_type);
134:   }

136:   return(0);
137: }

139: #undef  __FUNCT__
141: /*@
142:   TSViewFromOptions - This function visualizes the ts based upon user options.

144:   Collective on TS

146:   Input Parameter:
147: . ts - The ts

149:   Level: intermediate

151: .keywords: TS, view, options, database
152: .seealso: TSSetFromOptions(), TSView()
153: @*/
154: PetscErrorCode TSViewFromOptions(TS ts,const char title[])
155: {
156:   PetscViewer    viewer;
157:   PetscDraw      draw;
158:   PetscTruth     opt;
159:   char           typeName[1024];
160:   char           fileName[PETSC_MAX_PATH_LEN];
161:   size_t         len;

165:   PetscOptionsHasName(ts->prefix, "-ts_view", &opt);
166:   if (opt) {
167:     PetscOptionsGetString(ts->prefix, "-ts_view", typeName, 1024, &opt);
168:     PetscStrlen(typeName, &len);
169:     if (len > 0) {
170:       PetscViewerCreate(ts->comm, &viewer);
171:       PetscViewerSetType(viewer, typeName);
172:       PetscOptionsGetString(ts->prefix, "-ts_view_file", fileName, 1024, &opt);
173:       if (opt) {
174:         PetscViewerSetFilename(viewer, fileName);
175:       } else {
176:         PetscViewerSetFilename(viewer, ts->name);
177:       }
178:       TSView(ts, viewer);
179:       PetscViewerFlush(viewer);
180:       PetscViewerDestroy(viewer);
181:     } else {
182:       TSView(ts, PETSC_NULL);
183:     }
184:   }
185:   PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
186:   if (opt) {
187:     PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
188:     PetscViewerDrawGetDraw(viewer, 0, &draw);
189:     if (title) {
190:       PetscDrawSetTitle(draw, (char *)title);
191:     } else {
192:       PetscObjectName((PetscObject) ts);
193:       PetscDrawSetTitle(draw, ts->name);
194:     }
195:     TSView(ts, viewer);
196:     PetscViewerFlush(viewer);
197:     PetscDrawPause(draw);
198:     PetscViewerDestroy(viewer);
199:   }
200:   return(0);
201: }

205: /*@
206:    TSComputeRHSJacobian - Computes the Jacobian matrix that has been
207:       set with TSSetRHSJacobian().

209:    Collective on TS and Vec

211:    Input Parameters:
212: +  ts - the SNES context
213: .  t - current timestep
214: -  x - input vector

216:    Output Parameters:
217: +  A - Jacobian matrix
218: .  B - optional preconditioning matrix
219: -  flag - flag indicating matrix structure

221:    Notes: 
222:    Most users should not need to explicitly call this routine, as it 
223:    is used internally within the nonlinear solvers. 

225:    See KSPSetOperators() for important information about setting the
226:    flag parameter.

228:    TSComputeJacobian() is valid only for TS_NONLINEAR

230:    Level: developer

232: .keywords: SNES, compute, Jacobian, matrix

234: .seealso:  TSSetRHSJacobian(), KSPSetOperators()
235: @*/
236: PetscErrorCode TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
237: {

244:   if (ts->problem_type != TS_NONLINEAR) {
245:     SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
246:   }
247:   if (ts->ops->rhsjacobian) {
248:     PetscLogEventBegin(TS_JacobianEval,ts,X,*A,*B);
249:     *flg = DIFFERENT_NONZERO_PATTERN;
250:     PetscStackPush("TS user Jacobian function");
251:     (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
252:     PetscStackPop;
253:     PetscLogEventEnd(TS_JacobianEval,ts,X,*A,*B);
254:     /* make sure user returned a correct Jacobian and preconditioner */
257:   } else {
258:     MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
259:     MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
260:     if (*A != *B) {
261:       MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
262:       MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
263:     }
264:   }
265:   return(0);
266: }

270: /*
271:    TSComputeRHSFunction - Evaluates the right-hand-side function. 

273:    Note: If the user did not provide a function but merely a matrix,
274:    this routine applies the matrix.
275: */
276: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
277: {


285:   PetscLogEventBegin(TS_FunctionEval,ts,x,y,0);
286:   if (ts->ops->rhsfunction) {
287:     PetscStackPush("TS user right-hand-side function");
288:     (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
289:     PetscStackPop;
290:   } else {
291:     if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
292:       MatStructure flg;
293:       PetscStackPush("TS user right-hand-side matrix function");
294:       (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
295:       PetscStackPop;
296:     }
297:     MatMult(ts->A,x,y);
298:   }

300:   /* apply user-provided boundary conditions (only needed if these are time dependent) */
301:   TSComputeRHSBoundaryConditions(ts,t,y);
302:   PetscLogEventEnd(TS_FunctionEval,ts,x,y,0);

304:   return(0);
305: }

309: /*@C
310:     TSSetRHSFunction - Sets the routine for evaluating the function,
311:     F(t,u), where U_t = F(t,u).

313:     Collective on TS

315:     Input Parameters:
316: +   ts - the TS context obtained from TSCreate()
317: .   f - routine for evaluating the right-hand-side function
318: -   ctx - [optional] user-defined context for private data for the 
319:           function evaluation routine (may be PETSC_NULL)

321:     Calling sequence of func:
322: $     func (TS ts,PetscReal t,Vec u,Vec F,void *ctx);

324: +   t - current timestep
325: .   u - input vector
326: .   F - function vector
327: -   ctx - [optional] user-defined function context 

329:     Important: 
330:     The user MUST call either this routine or TSSetRHSMatrix().

332:     Level: beginner

334: .keywords: TS, timestep, set, right-hand-side, function

336: .seealso: TSSetRHSMatrix()
337: @*/
338: PetscErrorCode TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
339: {

343:   if (ts->problem_type == TS_LINEAR) {
344:     SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
345:   }
346:   ts->ops->rhsfunction = f;
347:   ts->funP             = ctx;
348:   return(0);
349: }

353: /*@C
354:    TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
355:    Also sets the location to store A.

357:    Collective on TS

359:    Input Parameters:
360: +  ts  - the TS context obtained from TSCreate()
361: .  A   - matrix
362: .  B   - preconditioner matrix (usually same as A)
363: .  f   - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
364:          if A is not a function of t.
365: -  ctx - [optional] user-defined context for private data for the 
366:           matrix evaluation routine (may be PETSC_NULL)

368:    Calling sequence of func:
369: $     func (TS ts,PetscReal t,Mat *A,Mat *B,PetscInt *flag,void *ctx);

371: +  t - current timestep
372: .  A - matrix A, where U_t = A(t) U
373: .  B - preconditioner matrix, usually the same as A
374: .  flag - flag indicating information about the preconditioner matrix
375:           structure (same as flag in KSPSetOperators())
376: -  ctx - [optional] user-defined context for matrix evaluation routine

378:    Notes: 
379:    See KSPSetOperators() for important information about setting the flag
380:    output parameter in the routine func().  Be sure to read this information!

382:    The routine func() takes Mat * as the matrix arguments rather than Mat.  
383:    This allows the matrix evaluation routine to replace A and/or B with a 
384:    completely new new matrix structure (not just different matrix elements)
385:    when appropriate, for instance, if the nonzero structure is changing
386:    throughout the global iterations.

388:    Important: 
389:    The user MUST call either this routine or TSSetRHSFunction().

391:    Level: beginner

393: .keywords: TS, timestep, set, right-hand-side, matrix

395: .seealso: TSSetRHSFunction()
396: @*/
397: PetscErrorCode TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
398: {
405:   if (ts->problem_type == TS_NONLINEAR) {
406:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
407:   }

409:   ts->ops->rhsmatrix = f;
410:   ts->jacP           = ctx;
411:   ts->A              = A;
412:   ts->B              = B;

414:   return(0);
415: }

419: /*@C
420:    TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
421:    where U_t = F(U,t), as well as the location to store the matrix.
422:    Use TSSetRHSMatrix() for linear problems.

424:    Collective on TS

426:    Input Parameters:
427: +  ts  - the TS context obtained from TSCreate()
428: .  A   - Jacobian matrix
429: .  B   - preconditioner matrix (usually same as A)
430: .  f   - the Jacobian evaluation routine
431: -  ctx - [optional] user-defined context for private data for the 
432:          Jacobian evaluation routine (may be PETSC_NULL)

434:    Calling sequence of func:
435: $     func (TS ts,PetscReal t,Vec u,Mat *A,Mat *B,MatStructure *flag,void *ctx);

437: +  t - current timestep
438: .  u - input vector
439: .  A - matrix A, where U_t = A(t)u
440: .  B - preconditioner matrix, usually the same as A
441: .  flag - flag indicating information about the preconditioner matrix
442:           structure (same as flag in KSPSetOperators())
443: -  ctx - [optional] user-defined context for matrix evaluation routine

445:    Notes: 
446:    See KSPSetOperators() for important information about setting the flag
447:    output parameter in the routine func().  Be sure to read this information!

449:    The routine func() takes Mat * as the matrix arguments rather than Mat.  
450:    This allows the matrix evaluation routine to replace A and/or B with a 
451:    completely new new matrix structure (not just different matrix elements)
452:    when appropriate, for instance, if the nonzero structure is changing
453:    throughout the global iterations.

455:    Level: beginner
456:    
457: .keywords: TS, timestep, set, right-hand-side, Jacobian

459: .seealso: TSDefaultComputeJacobianColor(),
460:           SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()

462: @*/
463: PetscErrorCode TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
464: {
471:   if (ts->problem_type != TS_NONLINEAR) {
472:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
473:   }

475:   ts->ops->rhsjacobian = f;
476:   ts->jacP             = ctx;
477:   ts->A                = A;
478:   ts->B                = B;
479:   return(0);
480: }

484: /*
485:    TSComputeRHSBoundaryConditions - Evaluates the boundary condition function. 

487:    Note: If the user did not provide a function but merely a matrix,
488:    this routine applies the matrix.
489: */
490: PetscErrorCode TSComputeRHSBoundaryConditions(TS ts,PetscReal t,Vec x)
491: {


499:   if (ts->ops->rhsbc) {
500:     PetscStackPush("TS user boundary condition function");
501:     (*ts->ops->rhsbc)(ts,t,x,ts->bcP);
502:     PetscStackPop;
503:     return(0);
504:   }

506:   return(0);
507: }

511: /*@C
512:     TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
513:     boundary conditions for the function F.

515:     Collective on TS

517:     Input Parameters:
518: +   ts - the TS context obtained from TSCreate()
519: .   f - routine for evaluating the boundary condition function
520: -   ctx - [optional] user-defined context for private data for the 
521:           function evaluation routine (may be PETSC_NULL)

523:     Calling sequence of func:
524: $     func (TS ts,PetscReal t,Vec F,void *ctx);

526: +   t - current timestep
527: .   F - function vector
528: -   ctx - [optional] user-defined function context 

530:     Level: intermediate

532: .keywords: TS, timestep, set, boundary conditions, function
533: @*/
534: PetscErrorCode TSSetRHSBoundaryConditions(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,void*),void *ctx)
535: {

539:   if (ts->problem_type != TS_LINEAR) {
540:     SETERRQ(PETSC_ERR_ARG_WRONG,"For linear problems only");
541:   }
542:   ts->ops->rhsbc = f;
543:   ts->bcP        = ctx;
544:   return(0);
545: }

549: /*@C
550:     TSView - Prints the TS data structure.

552:     Collective on TS

554:     Input Parameters:
555: +   ts - the TS context obtained from TSCreate()
556: -   viewer - visualization context

558:     Options Database Key:
559: .   -ts_view - calls TSView() at end of TSStep()

561:     Notes:
562:     The available visualization contexts include
563: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
564: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
565:          output where only the first processor opens
566:          the file.  All other processors send their 
567:          data to the first processor to print. 

569:     The user can open an alternative visualization context with
570:     PetscViewerASCIIOpen() - output to a specified file.

572:     Level: beginner

574: .keywords: TS, timestep, view

576: .seealso: PetscViewerASCIIOpen()
577: @*/
578: PetscErrorCode TSView(TS ts,PetscViewer viewer)
579: {
581:   char           *type;
582:   PetscTruth     iascii,isstring;

586:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);

590:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
591:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);
592:   if (iascii) {
593:     PetscViewerASCIIPrintf(viewer,"TS Object:\n");
594:     TSGetType(ts,(TSType *)&type);
595:     if (type) {
596:       PetscViewerASCIIPrintf(viewer,"  type: %s\n",type);
597:     } else {
598:       PetscViewerASCIIPrintf(viewer,"  type: not yet set\n");
599:     }
600:     if (ts->ops->view) {
601:       PetscViewerASCIIPushTab(viewer);
602:       (*ts->ops->view)(ts,viewer);
603:       PetscViewerASCIIPopTab(viewer);
604:     }
605:     PetscViewerASCIIPrintf(viewer,"  maximum steps=%D\n",ts->max_steps);
606:     PetscViewerASCIIPrintf(viewer,"  maximum time=%g\n",ts->max_time);
607:     if (ts->problem_type == TS_NONLINEAR) {
608:       PetscViewerASCIIPrintf(viewer,"  total number of nonlinear solver iterations=%D\n",ts->nonlinear_its);
609:     }
610:     PetscViewerASCIIPrintf(viewer,"  total number of linear solver iterations=%D\n",ts->linear_its);
611:   } else if (isstring) {
612:     TSGetType(ts,(TSType *)&type);
613:     PetscViewerStringSPrintf(viewer," %-7.7s",type);
614:   }
615:   PetscViewerASCIIPushTab(viewer);
616:   if (ts->ksp) {KSPView(ts->ksp,viewer);}
617:   if (ts->snes) {SNESView(ts->snes,viewer);}
618:   PetscViewerASCIIPopTab(viewer);
619:   return(0);
620: }


625: /*@C
626:    TSSetApplicationContext - Sets an optional user-defined context for 
627:    the timesteppers.

629:    Collective on TS

631:    Input Parameters:
632: +  ts - the TS context obtained from TSCreate()
633: -  usrP - optional user context

635:    Level: intermediate

637: .keywords: TS, timestep, set, application, context

639: .seealso: TSGetApplicationContext()
640: @*/
641: PetscErrorCode TSSetApplicationContext(TS ts,void *usrP)
642: {
645:   ts->user = usrP;
646:   return(0);
647: }

651: /*@C
652:     TSGetApplicationContext - Gets the user-defined context for the 
653:     timestepper.

655:     Not Collective

657:     Input Parameter:
658: .   ts - the TS context obtained from TSCreate()

660:     Output Parameter:
661: .   usrP - user context

663:     Level: intermediate

665: .keywords: TS, timestep, get, application, context

667: .seealso: TSSetApplicationContext()
668: @*/
669: PetscErrorCode TSGetApplicationContext(TS ts,void **usrP)
670: {
673:   *usrP = ts->user;
674:   return(0);
675: }

679: /*@
680:    TSGetTimeStepNumber - Gets the current number of timesteps.

682:    Not Collective

684:    Input Parameter:
685: .  ts - the TS context obtained from TSCreate()

687:    Output Parameter:
688: .  iter - number steps so far

690:    Level: intermediate

692: .keywords: TS, timestep, get, iteration, number
693: @*/
694: PetscErrorCode TSGetTimeStepNumber(TS ts,PetscInt* iter)
695: {
699:   *iter = ts->steps;
700:   return(0);
701: }

705: /*@
706:    TSSetInitialTimeStep - Sets the initial timestep to be used, 
707:    as well as the initial time.

709:    Collective on TS

711:    Input Parameters:
712: +  ts - the TS context obtained from TSCreate()
713: .  initial_time - the initial time
714: -  time_step - the size of the timestep

716:    Level: intermediate

718: .seealso: TSSetTimeStep(), TSGetTimeStep()

720: .keywords: TS, set, initial, timestep
721: @*/
722: PetscErrorCode TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
723: {
726:   ts->time_step         = time_step;
727:   ts->initial_time_step = time_step;
728:   ts->ptime             = initial_time;
729:   return(0);
730: }

734: /*@
735:    TSSetTimeStep - Allows one to reset the timestep at any time,
736:    useful for simple pseudo-timestepping codes.

738:    Collective on TS

740:    Input Parameters:
741: +  ts - the TS context obtained from TSCreate()
742: -  time_step - the size of the timestep

744:    Level: intermediate

746: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()

748: .keywords: TS, set, timestep
749: @*/
750: PetscErrorCode TSSetTimeStep(TS ts,PetscReal time_step)
751: {
754:   ts->time_step = time_step;
755:   return(0);
756: }

760: /*@
761:    TSGetTimeStep - Gets the current timestep size.

763:    Not Collective

765:    Input Parameter:
766: .  ts - the TS context obtained from TSCreate()

768:    Output Parameter:
769: .  dt - the current timestep size

771:    Level: intermediate

773: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()

775: .keywords: TS, get, timestep
776: @*/
777: PetscErrorCode TSGetTimeStep(TS ts,PetscReal* dt)
778: {
782:   *dt = ts->time_step;
783:   return(0);
784: }

788: /*@C
789:    TSGetSolution - Returns the solution at the present timestep. It
790:    is valid to call this routine inside the function that you are evaluating
791:    in order to move to the new timestep. This vector not changed until
792:    the solution at the next timestep has been calculated.

794:    Not Collective, but Vec returned is parallel if TS is parallel

796:    Input Parameter:
797: .  ts - the TS context obtained from TSCreate()

799:    Output Parameter:
800: .  v - the vector containing the solution

802:    Level: intermediate

804: .seealso: TSGetTimeStep()

806: .keywords: TS, timestep, get, solution
807: @*/
808: PetscErrorCode TSGetSolution(TS ts,Vec *v)
809: {
813:   *v = ts->vec_sol_always;
814:   return(0);
815: }

817: /* ----- Routines to initialize and destroy a timestepper ---- */
820: /*@C
821:   TSSetProblemType - Sets the type of problem to be solved.

823:   Not collective

825:   Input Parameters:
826: + ts   - The TS
827: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
828: .vb
829:          U_t = A U    
830:          U_t = A(t) U 
831:          U_t = F(t,U) 
832: .ve

834:    Level: beginner

836: .keywords: TS, problem type
837: .seealso: TSSetUp(), TSProblemType, TS
838: @*/
839: PetscErrorCode TSSetProblemType(TS ts, TSProblemType type)
840: {
843:   ts->problem_type = type;
844:   return(0);
845: }

849: /*@C
850:   TSGetProblemType - Gets the type of problem to be solved.

852:   Not collective

854:   Input Parameter:
855: . ts   - The TS

857:   Output Parameter:
858: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
859: .vb
860:          U_t = A U    
861:          U_t = A(t) U 
862:          U_t = F(t,U) 
863: .ve

865:    Level: beginner

867: .keywords: TS, problem type
868: .seealso: TSSetUp(), TSProblemType, TS
869: @*/
870: PetscErrorCode TSGetProblemType(TS ts, TSProblemType *type)
871: {
875:   *type = ts->problem_type;
876:   return(0);
877: }

881: /*@
882:    TSSetUp - Sets up the internal data structures for the later use
883:    of a timestepper.  

885:    Collective on TS

887:    Input Parameter:
888: .  ts - the TS context obtained from TSCreate()

890:    Notes:
891:    For basic use of the TS solvers the user need not explicitly call
892:    TSSetUp(), since these actions will automatically occur during
893:    the call to TSStep().  However, if one wishes to control this
894:    phase separately, TSSetUp() should be called after TSCreate()
895:    and optional routines of the form TSSetXXX(), but before TSStep().  

897:    Level: advanced

899: .keywords: TS, timestep, setup

901: .seealso: TSCreate(), TSStep(), TSDestroy()
902: @*/
903: PetscErrorCode TSSetUp(TS ts)
904: {

909:   if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
910:   if (!ts->type_name) {
911:     TSSetType(ts,TS_EULER);
912:   }
913:   (*ts->ops->setup)(ts);
914:   ts->setupcalled = 1;
915:   return(0);
916: }

920: /*@C
921:    TSDestroy - Destroys the timestepper context that was created
922:    with TSCreate().

924:    Collective on TS

926:    Input Parameter:
927: .  ts - the TS context obtained from TSCreate()

929:    Level: beginner

931: .keywords: TS, timestepper, destroy

933: .seealso: TSCreate(), TSSetUp(), TSSolve()
934: @*/
935: PetscErrorCode TSDestroy(TS ts)
936: {
938:   PetscInt       i;

942:   if (--ts->refct > 0) return(0);

944:   /* if memory was published with AMS then destroy it */
945:   PetscObjectDepublish(ts);

947:   if (ts->ksp) {KSPDestroy(ts->ksp);}
948:   if (ts->snes) {SNESDestroy(ts->snes);}
949:   (*(ts)->ops->destroy)(ts);
950:   for (i=0; i<ts->numbermonitors; i++) {
951:     if (ts->mdestroy[i]) {
952:       (*ts->mdestroy[i])(ts->monitorcontext[i]);
953:     }
954:   }
955:   PetscLogObjectDestroy((PetscObject)ts);
956:   PetscHeaderDestroy((PetscObject)ts);
957:   return(0);
958: }

962: /*@C
963:    TSGetSNES - Returns the SNES (nonlinear solver) associated with 
964:    a TS (timestepper) context. Valid only for nonlinear problems.

966:    Not Collective, but SNES is parallel if TS is parallel

968:    Input Parameter:
969: .  ts - the TS context obtained from TSCreate()

971:    Output Parameter:
972: .  snes - the nonlinear solver context

974:    Notes:
975:    The user can then directly manipulate the SNES context to set various
976:    options, etc.  Likewise, the user can then extract and manipulate the 
977:    KSP, KSP, and PC contexts as well.

979:    TSGetSNES() does not work for integrators that do not use SNES; in
980:    this case TSGetSNES() returns PETSC_NULL in snes.

982:    Level: beginner

984: .keywords: timestep, get, SNES
985: @*/
986: PetscErrorCode TSGetSNES(TS ts,SNES *snes)
987: {
991:   if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
992:   *snes = ts->snes;
993:   return(0);
994: }

998: /*@C
999:    TSGetKSP - Returns the KSP (linear solver) associated with 
1000:    a TS (timestepper) context.

1002:    Not Collective, but KSP is parallel if TS is parallel

1004:    Input Parameter:
1005: .  ts - the TS context obtained from TSCreate()

1007:    Output Parameter:
1008: .  ksp - the nonlinear solver context

1010:    Notes:
1011:    The user can then directly manipulate the KSP context to set various
1012:    options, etc.  Likewise, the user can then extract and manipulate the 
1013:    KSP and PC contexts as well.

1015:    TSGetKSP() does not work for integrators that do not use KSP;
1016:    in this case TSGetKSP() returns PETSC_NULL in ksp.

1018:    Level: beginner

1020: .keywords: timestep, get, KSP
1021: @*/
1022: PetscErrorCode TSGetKSP(TS ts,KSP *ksp)
1023: {
1027:   if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1028:   *ksp = ts->ksp;
1029:   return(0);
1030: }

1032: /* ----------- Routines to set solver parameters ---------- */

1036: /*@
1037:    TSGetDuration - Gets the maximum number of timesteps to use and 
1038:    maximum time for iteration.

1040:    Collective on TS

1042:    Input Parameters:
1043: +  ts       - the TS context obtained from TSCreate()
1044: .  maxsteps - maximum number of iterations to use, or PETSC_NULL
1045: -  maxtime  - final time to iterate to, or PETSC_NULL

1047:    Level: intermediate

1049: .keywords: TS, timestep, get, maximum, iterations, time
1050: @*/
1051: PetscErrorCode TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1052: {
1055:   if (maxsteps != PETSC_NULL) {
1057:     *maxsteps = ts->max_steps;
1058:   }
1059:   if (maxtime  != PETSC_NULL) {
1061:     *maxtime  = ts->max_time;
1062:   }
1063:   return(0);
1064: }

1068: /*@
1069:    TSSetDuration - Sets the maximum number of timesteps to use and 
1070:    maximum time for iteration.

1072:    Collective on TS

1074:    Input Parameters:
1075: +  ts - the TS context obtained from TSCreate()
1076: .  maxsteps - maximum number of iterations to use
1077: -  maxtime - final time to iterate to

1079:    Options Database Keys:
1080: .  -ts_max_steps <maxsteps> - Sets maxsteps
1081: .  -ts_max_time <maxtime> - Sets maxtime

1083:    Notes:
1084:    The default maximum number of iterations is 5000. Default time is 5.0

1086:    Level: intermediate

1088: .keywords: TS, timestep, set, maximum, iterations
1089: @*/
1090: PetscErrorCode TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1091: {
1094:   ts->max_steps = maxsteps;
1095:   ts->max_time  = maxtime;
1096:   return(0);
1097: }

1101: /*@
1102:    TSSetSolution - Sets the initial solution vector
1103:    for use by the TS routines.

1105:    Collective on TS and Vec

1107:    Input Parameters:
1108: +  ts - the TS context obtained from TSCreate()
1109: -  x - the solution vector

1111:    Level: beginner

1113: .keywords: TS, timestep, set, solution, initial conditions
1114: @*/
1115: PetscErrorCode TSSetSolution(TS ts,Vec x)
1116: {
1120:   ts->vec_sol        = ts->vec_sol_always = x;
1121:   return(0);
1122: }

1126: /*@C
1127:   TSSetRhsBC - Sets the function which applies boundary conditions
1128:   to the Rhs of each system.

1130:   Collective on TS

1132:   Input Parameters:
1133: + ts   - The TS context obtained from TSCreate()
1134: - func - The function

1136:   Calling sequence of func:
1137: . func (TS ts, Vec rhs, void *ctx);

1139: + rhs - The current rhs vector
1140: - ctx - The user-context

1142:   Level: intermediate

1144: .keywords: TS, Rhs, boundary conditions
1145: @*/
1146: PetscErrorCode TSSetRhsBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1147: {
1150:   ts->ops->applyrhsbc = func;
1151:   return(0);
1152: }

1156: /*@
1157:   TSDefaultRhsBC - The default boundary condition function which does nothing.

1159:   Collective on TS

1161:   Input Parameters:
1162: + ts  - The TS context obtained from TSCreate()
1163: . rhs - The Rhs
1164: - ctx - The user-context

1166:   Level: developer

1168: .keywords: TS, Rhs, boundary conditions
1169: @*/
1170: PetscErrorCode TSDefaultRhsBC(TS ts,  Vec rhs, void *ctx)
1171: {
1173:   return(0);
1174: }

1178: /*@C
1179:   TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1180:   to the system matrix and preconditioner of each system.

1182:   Collective on TS

1184:   Input Parameters:
1185: + ts   - The TS context obtained from TSCreate()
1186: - func - The function

1188:   Calling sequence of func:
1189: . func (TS ts, Mat A, Mat B, void *ctx);

1191: + A   - The current system matrix
1192: . B   - The current preconditioner
1193: - ctx - The user-context

1195:   Level: intermediate

1197: .keywords: TS, System matrix, boundary conditions
1198: @*/
1199: PetscErrorCode TSSetSystemMatrixBC(TS ts, PetscErrorCode (*func)(TS, Mat, Mat, void *))
1200: {
1203:   ts->ops->applymatrixbc = func;
1204:   return(0);
1205: }

1209: /*@
1210:   TSDefaultSystemMatrixBC - The default boundary condition function which
1211:   does nothing.

1213:   Collective on TS

1215:   Input Parameters:
1216: + ts  - The TS context obtained from TSCreate()
1217: . A   - The system matrix
1218: . B   - The preconditioner
1219: - ctx - The user-context

1221:   Level: developer

1223: .keywords: TS, System matrix, boundary conditions
1224: @*/
1225: PetscErrorCode TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1226: {
1228:   return(0);
1229: }

1233: /*@C
1234:   TSSetSolutionBC - Sets the function which applies boundary conditions
1235:   to the solution of each system. This is necessary in nonlinear systems
1236:   which time dependent boundary conditions.

1238:   Collective on TS

1240:   Input Parameters:
1241: + ts   - The TS context obtained from TSCreate()
1242: - func - The function

1244:   Calling sequence of func:
1245: . func (TS ts, Vec rsol, void *ctx);

1247: + sol - The current solution vector
1248: - ctx - The user-context

1250:   Level: intermediate

1252: .keywords: TS, solution, boundary conditions
1253: @*/
1254: PetscErrorCode TSSetSolutionBC(TS ts, PetscErrorCode (*func)(TS, Vec, void *))
1255: {
1258:   ts->ops->applysolbc = func;
1259:   return(0);
1260: }

1264: /*@
1265:   TSDefaultSolutionBC - The default boundary condition function which
1266:   does nothing.

1268:   Collective on TS

1270:   Input Parameters:
1271: + ts  - The TS context obtained from TSCreate()
1272: . sol - The solution
1273: - ctx - The user-context

1275:   Level: developer

1277: .keywords: TS, solution, boundary conditions
1278: @*/
1279: PetscErrorCode TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1280: {
1282:   return(0);
1283: }

1287: /*@C
1288:   TSSetPreStep - Sets the general-purpose function
1289:   called once at the beginning of time stepping.

1291:   Collective on TS

1293:   Input Parameters:
1294: + ts   - The TS context obtained from TSCreate()
1295: - func - The function

1297:   Calling sequence of func:
1298: . func (TS ts);

1300:   Level: intermediate

1302: .keywords: TS, timestep
1303: @*/
1304: PetscErrorCode TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1305: {
1308:   ts->ops->prestep = func;
1309:   return(0);
1310: }

1314: /*@
1315:   TSDefaultPreStep - The default pre-stepping function which does nothing.

1317:   Collective on TS

1319:   Input Parameters:
1320: . ts  - The TS context obtained from TSCreate()

1322:   Level: developer

1324: .keywords: TS, timestep
1325: @*/
1326: PetscErrorCode TSDefaultPreStep(TS ts)
1327: {
1329:   return(0);
1330: }

1334: /*@C
1335:   TSSetUpdate - Sets the general-purpose update function called
1336:   at the beginning of every time step. This function can change
1337:   the time step.

1339:   Collective on TS

1341:   Input Parameters:
1342: + ts   - The TS context obtained from TSCreate()
1343: - func - The function

1345:   Calling sequence of func:
1346: . func (TS ts, double t, double *dt);

1348: + t   - The current time
1349: - dt  - The current time step

1351:   Level: intermediate

1353: .keywords: TS, update, timestep
1354: @*/
1355: PetscErrorCode TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1356: {
1359:   ts->ops->update = func;
1360:   return(0);
1361: }

1365: /*@
1366:   TSDefaultUpdate - The default update function which does nothing.

1368:   Collective on TS

1370:   Input Parameters:
1371: + ts  - The TS context obtained from TSCreate()
1372: - t   - The current time

1374:   Output Parameters:
1375: . dt  - The current time step

1377:   Level: developer

1379: .keywords: TS, update, timestep
1380: @*/
1381: PetscErrorCode TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1382: {
1384:   return(0);
1385: }

1389: /*@C
1390:   TSSetPostStep - Sets the general-purpose function
1391:   called once at the end of time stepping.

1393:   Collective on TS

1395:   Input Parameters:
1396: + ts   - The TS context obtained from TSCreate()
1397: - func - The function

1399:   Calling sequence of func:
1400: . func (TS ts);

1402:   Level: intermediate

1404: .keywords: TS, timestep
1405: @*/
1406: PetscErrorCode TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1407: {
1410:   ts->ops->poststep = func;
1411:   return(0);
1412: }

1416: /*@
1417:   TSDefaultPostStep - The default post-stepping function which does nothing.

1419:   Collective on TS

1421:   Input Parameters:
1422: . ts  - The TS context obtained from TSCreate()

1424:   Level: developer

1426: .keywords: TS, timestep
1427: @*/
1428: PetscErrorCode TSDefaultPostStep(TS ts)
1429: {
1431:   return(0);
1432: }

1434: /* ------------ Routines to set performance monitoring options ----------- */

1438: /*@C
1439:    TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1440:    timestep to display the iteration's  progress.   

1442:    Collective on TS

1444:    Input Parameters:
1445: +  ts - the TS context obtained from TSCreate()
1446: .  func - monitoring routine
1447: .  mctx - [optional] user-defined context for private data for the 
1448:              monitor routine (use PETSC_NULL if no context is desired)
1449: -  monitordestroy - [optional] routine that frees monitor context
1450:           (may be PETSC_NULL)

1452:    Calling sequence of func:
1453: $    int func(TS ts,PetscInt steps,PetscReal time,Vec x,void *mctx)

1455: +    ts - the TS context
1456: .    steps - iteration number
1457: .    time - current time
1458: .    x - current iterate
1459: -    mctx - [optional] monitoring context

1461:    Notes:
1462:    This routine adds an additional monitor to the list of monitors that 
1463:    already has been loaded.

1465:    Level: intermediate

1467: .keywords: TS, timestep, set, monitor

1469: .seealso: TSDefaultMonitor(), TSClearMonitor()
1470: @*/
1471: PetscErrorCode TSSetMonitor(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1472: {
1475:   if (ts->numbermonitors >= MAXTSMONITORS) {
1476:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1477:   }
1478:   ts->monitor[ts->numbermonitors]           = monitor;
1479:   ts->mdestroy[ts->numbermonitors]          = mdestroy;
1480:   ts->monitorcontext[ts->numbermonitors++]  = (void*)mctx;
1481:   return(0);
1482: }

1486: /*@C
1487:    TSClearMonitor - Clears all the monitors that have been set on a time-step object.   

1489:    Collective on TS

1491:    Input Parameters:
1492: .  ts - the TS context obtained from TSCreate()

1494:    Notes:
1495:    There is no way to remove a single, specific monitor.

1497:    Level: intermediate

1499: .keywords: TS, timestep, set, monitor

1501: .seealso: TSDefaultMonitor(), TSSetMonitor()
1502: @*/
1503: PetscErrorCode TSClearMonitor(TS ts)
1504: {
1507:   ts->numbermonitors = 0;
1508:   return(0);
1509: }

1513: PetscErrorCode TSDefaultMonitor(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1514: {

1518:   PetscPrintf(ts->comm,"timestep %D dt %g time %g\n",step,ts->time_step,ptime);
1519:   return(0);
1520: }

1524: /*@
1525:    TSStep - Steps the requested number of timesteps.

1527:    Collective on TS

1529:    Input Parameter:
1530: .  ts - the TS context obtained from TSCreate()

1532:    Output Parameters:
1533: +  steps - number of iterations until termination
1534: -  ptime - time until termination

1536:    Level: beginner

1538: .keywords: TS, timestep, solve

1540: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1541: @*/
1542: PetscErrorCode TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1543: {

1548:   if (!ts->setupcalled) {
1549:     TSSetUp(ts);
1550:   }

1552:   PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1553:   (*ts->ops->prestep)(ts);
1554:   (*ts->ops->step)(ts, steps, ptime);
1555:   (*ts->ops->poststep)(ts);
1556:   PetscLogEventEnd(TS_Step, ts, 0, 0, 0);

1558:   if (!PetscPreLoadingOn) {
1559:     TSViewFromOptions(ts,ts->name);
1560:   }
1561:   return(0);
1562: }

1566: /*
1567:      Runs the user provided monitor routines, if they exists.
1568: */
1569: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1570: {
1572:   PetscInt i,n = ts->numbermonitors;

1575:   for (i=0; i<n; i++) {
1576:     (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1577:   }
1578:   return(0);
1579: }

1581: /* ------------------------------------------------------------------------*/

1585: /*@C
1586:    TSLGMonitorCreate - Creates a line graph context for use with 
1587:    TS to monitor convergence of preconditioned residual norms.

1589:    Collective on TS

1591:    Input Parameters:
1592: +  host - the X display to open, or null for the local machine
1593: .  label - the title to put in the title bar
1594: .  x, y - the screen coordinates of the upper left coordinate of the window
1595: -  m, n - the screen width and height in pixels

1597:    Output Parameter:
1598: .  draw - the drawing context

1600:    Options Database Key:
1601: .  -ts_xmonitor - automatically sets line graph monitor

1603:    Notes: 
1604:    Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().

1606:    Level: intermediate

1608: .keywords: TS, monitor, line graph, residual, seealso

1610: .seealso: TSLGMonitorDestroy(), TSSetMonitor()

1612: @*/
1613: PetscErrorCode TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1614: {
1615:   PetscDraw      win;

1619:   PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1620:   PetscDrawSetType(win,PETSC_DRAW_X);
1621:   PetscDrawLGCreate(win,1,draw);
1622:   PetscDrawLGIndicateDataPoints(*draw);

1624:   PetscLogObjectParent(*draw,win);
1625:   return(0);
1626: }

1630: PetscErrorCode TSLGMonitor(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1631: {
1632:   PetscDrawLG    lg = (PetscDrawLG) monctx;
1633:   PetscReal      x,y = ptime;

1637:   if (!monctx) {
1638:     MPI_Comm    comm;
1639:     PetscViewer viewer;

1641:     PetscObjectGetComm((PetscObject)ts,&comm);
1642:     viewer = PETSC_VIEWER_DRAW_(comm);
1643:     PetscViewerDrawGetDrawLG(viewer,0,&lg);
1644:   }

1646:   if (!n) {PetscDrawLGReset(lg);}
1647:   x = (PetscReal)n;
1648:   PetscDrawLGAddPoint(lg,&x,&y);
1649:   if (n < 20 || (n % 5)) {
1650:     PetscDrawLGDraw(lg);
1651:   }
1652:   return(0);
1653: }

1657: /*@C
1658:    TSLGMonitorDestroy - Destroys a line graph context that was created 
1659:    with TSLGMonitorCreate().

1661:    Collective on PetscDrawLG

1663:    Input Parameter:
1664: .  draw - the drawing context

1666:    Level: intermediate

1668: .keywords: TS, monitor, line graph, destroy

1670: .seealso: TSLGMonitorCreate(),  TSSetMonitor(), TSLGMonitor();
1671: @*/
1672: PetscErrorCode TSLGMonitorDestroy(PetscDrawLG drawlg)
1673: {
1674:   PetscDraw      draw;

1678:   PetscDrawLGGetDraw(drawlg,&draw);
1679:   PetscDrawDestroy(draw);
1680:   PetscDrawLGDestroy(drawlg);
1681:   return(0);
1682: }

1686: /*@
1687:    TSGetTime - Gets the current time.

1689:    Not Collective

1691:    Input Parameter:
1692: .  ts - the TS context obtained from TSCreate()

1694:    Output Parameter:
1695: .  t  - the current time

1697:    Contributed by: Matthew Knepley

1699:    Level: beginner

1701: .seealso: TSSetInitialTimeStep(), TSGetTimeStep()

1703: .keywords: TS, get, time
1704: @*/
1705: PetscErrorCode TSGetTime(TS ts,PetscReal* t)
1706: {
1710:   *t = ts->ptime;
1711:   return(0);
1712: }

1716: /*@C
1717:    TSSetOptionsPrefix - Sets the prefix used for searching for all
1718:    TS options in the database.

1720:    Collective on TS

1722:    Input Parameter:
1723: +  ts     - The TS context
1724: -  prefix - The prefix to prepend to all option names

1726:    Notes:
1727:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1728:    The first character of all runtime options is AUTOMATICALLY the
1729:    hyphen.

1731:    Contributed by: Matthew Knepley

1733:    Level: advanced

1735: .keywords: TS, set, options, prefix, database

1737: .seealso: TSSetFromOptions()

1739: @*/
1740: PetscErrorCode TSSetOptionsPrefix(TS ts,const char prefix[])
1741: {

1746:   PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1747:   switch(ts->problem_type) {
1748:     case TS_NONLINEAR:
1749:       SNESSetOptionsPrefix(ts->snes,prefix);
1750:       break;
1751:     case TS_LINEAR:
1752:       KSPSetOptionsPrefix(ts->ksp,prefix);
1753:       break;
1754:   }
1755:   return(0);
1756: }


1761: /*@C
1762:    TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1763:    TS options in the database.

1765:    Collective on TS

1767:    Input Parameter:
1768: +  ts     - The TS context
1769: -  prefix - The prefix to prepend to all option names

1771:    Notes:
1772:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1773:    The first character of all runtime options is AUTOMATICALLY the
1774:    hyphen.

1776:    Contributed by: Matthew Knepley

1778:    Level: advanced

1780: .keywords: TS, append, options, prefix, database

1782: .seealso: TSGetOptionsPrefix()

1784: @*/
1785: PetscErrorCode TSAppendOptionsPrefix(TS ts,const char prefix[])
1786: {

1791:   PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1792:   switch(ts->problem_type) {
1793:     case TS_NONLINEAR:
1794:       SNESAppendOptionsPrefix(ts->snes,prefix);
1795:       break;
1796:     case TS_LINEAR:
1797:       KSPAppendOptionsPrefix(ts->ksp,prefix);
1798:       break;
1799:   }
1800:   return(0);
1801: }

1805: /*@C
1806:    TSGetOptionsPrefix - Sets the prefix used for searching for all
1807:    TS options in the database.

1809:    Not Collective

1811:    Input Parameter:
1812: .  ts - The TS context

1814:    Output Parameter:
1815: .  prefix - A pointer to the prefix string used

1817:    Contributed by: Matthew Knepley

1819:    Notes: On the fortran side, the user should pass in a string 'prifix' of
1820:    sufficient length to hold the prefix.

1822:    Level: intermediate

1824: .keywords: TS, get, options, prefix, database

1826: .seealso: TSAppendOptionsPrefix()
1827: @*/
1828: PetscErrorCode TSGetOptionsPrefix(TS ts,char *prefix[])
1829: {

1835:   PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1836:   return(0);
1837: }

1841: /*@C
1842:    TSGetRHSMatrix - Returns the matrix A at the present timestep.

1844:    Not Collective, but parallel objects are returned if TS is parallel

1846:    Input Parameter:
1847: .  ts  - The TS context obtained from TSCreate()

1849:    Output Parameters:
1850: +  A   - The matrix A, where U_t = A(t) U
1851: .  M   - The preconditioner matrix, usually the same as A
1852: -  ctx - User-defined context for matrix evaluation routine

1854:    Notes: You can pass in PETSC_NULL for any return argument you do not need.

1856:    Contributed by: Matthew Knepley

1858:    Level: intermediate

1860: .seealso: TSGetTimeStep(), TSGetTime(), TSGetTimeStepNumber(), TSGetRHSJacobian()

1862: .keywords: TS, timestep, get, matrix

1864: @*/
1865: PetscErrorCode TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1866: {
1869:   if (A)   *A = ts->A;
1870:   if (M)   *M = ts->B;
1871:   if (ctx) *ctx = ts->jacP;
1872:   return(0);
1873: }

1877: /*@C
1878:    TSGetRHSJacobian - Returns the Jacobian J at the present timestep.

1880:    Not Collective, but parallel objects are returned if TS is parallel

1882:    Input Parameter:
1883: .  ts  - The TS context obtained from TSCreate()

1885:    Output Parameters:
1886: +  J   - The Jacobian J of F, where U_t = F(U,t)
1887: .  M   - The preconditioner matrix, usually the same as J
1888: - ctx - User-defined context for Jacobian evaluation routine

1890:    Notes: You can pass in PETSC_NULL for any return argument you do not need.

1892:    Contributed by: Matthew Knepley

1894:    Level: intermediate

1896: .seealso: TSGetTimeStep(), TSGetRHSMatrix(), TSGetTime(), TSGetTimeStepNumber()

1898: .keywords: TS, timestep, get, matrix, Jacobian
1899: @*/
1900: PetscErrorCode TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1901: {

1905:   TSGetRHSMatrix(ts,J,M,ctx);
1906:   return(0);
1907: }

1911: /*@C
1912:    TSVecViewMonitor - Monitors progress of the TS solvers by calling 
1913:    VecView() for the solution at each timestep

1915:    Collective on TS

1917:    Input Parameters:
1918: +  ts - the TS context
1919: .  step - current time-step
1920: .  ptime - current time
1921: -  dummy - either a viewer or PETSC_NULL

1923:    Level: intermediate

1925: .keywords: TS,  vector, monitor, view

1927: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1928: @*/
1929: PetscErrorCode TSVecViewMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1930: {
1932:   PetscViewer    viewer = (PetscViewer) dummy;

1935:   if (!viewer) {
1936:     MPI_Comm comm;
1937:     PetscObjectGetComm((PetscObject)ts,&comm);
1938:     viewer = PETSC_VIEWER_DRAW_(comm);
1939:   }
1940:   VecView(x,viewer);
1941:   return(0);
1942: }