Actual source code: ts.c

  1: /* $Id: ts.c,v 1.43 2001/09/07 20:12:01 bsmith Exp $ */
 2:  #include src/ts/tsimpl.h

  4: /* Logging support */
  5: int TS_COOKIE;
  6: int TS_Step, TS_PseudoComputeTimeStep, TS_FunctionEval, TS_JacobianEval;

  8: #undef __FUNCT__  
 10: /*
 11:   TSSetTypeFromOptions - Sets the type of ts from user options.

 13:   Collective on TS

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

 18:   Level: intermediate

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

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

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

 49: #undef __FUNCT__  
 51: /*@
 52:    TSSetFromOptions - Sets various TS parameters from user options.

 54:    Collective on TS

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

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

 67:    Level: beginner

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

 71: .seealso: TSGetType
 72: @*/
 73: int TSSetFromOptions(TS ts)
 74: {
 75:   PetscReal  dt;
 76:   PetscTruth opt;
 77:   int        ierr;

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

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

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

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

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

115:   /* Handle subobject options */
116:   switch(ts->problem_type) {
117:     /* Should check for implicit/explicit */
118:   case TS_LINEAR:
119:     if (ts->sles != PETSC_NULL) {
120:       SLESSetFromOptions(ts->sles);
121:     }
122:     break;
123:   case TS_NONLINEAR:
124:     if (ts->snes != PETSC_NULL) {
125:       SNESSetFromOptions(ts->snes);
126:     }
127:     break;
128:   default:
129:     SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid problem type: %d", ts->problem_type);
130:   }

132:   return(0);
133: }

135: #undef  __FUNCT__
137: /*@
138:   TSViewFromOptions - This function visualizes the ts based upon user options.

140:   Collective on TS

142:   Input Parameter:
143: . ts - The ts

145:   Level: intermediate

147: .keywords: TS, view, options, database
148: .seealso: TSSetFromOptions(), TSView()
149: @*/
150: int TSViewFromOptions(TS ts, char *title)
151: {
152:   PetscViewer viewer;
153:   PetscDraw   draw;
154:   PetscTruth  opt;
155:   char       *titleStr;
156:   char        typeName[1024];
157:   char        fileName[PETSC_MAX_PATH_LEN];
158:   int         len;
159:   int         ierr;

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

201: #undef __FUNCT__  
203: /*@
204:    TSComputeRHSJacobian - Computes the Jacobian matrix that has been
205:       set with TSSetRHSJacobian().

207:    Collective on TS and Vec

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

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

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

223:    See SLESSetOperators() for important information about setting the
224:    flag parameter.

226:    TSComputeJacobian() is valid only for TS_NONLINEAR

228:    Level: developer

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

232: .seealso:  TSSetRHSJacobian(), SLESSetOperators()
233: @*/
234: int TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
235: {

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

266: #undef __FUNCT__  
268: /*
269:    TSComputeRHSFunction - Evaluates the right-hand-side function. 

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


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

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

302:   return(0);
303: }

305: #undef __FUNCT__  
307: /*@C
308:     TSSetRHSFunction - Sets the routine for evaluating the function,
309:     F(t,u), where U_t = F(t,u).

311:     Collective on TS

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

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

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

327:     Important: 
328:     The user MUST call either this routine or TSSetRHSMatrix().

330:     Level: beginner

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

334: .seealso: TSSetRHSMatrix()
335: @*/
336: int TSSetRHSFunction(TS ts,int (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
337: {

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

349: #undef __FUNCT__  
351: /*@C
352:    TSSetRHSMatrix - Sets the function to compute the matrix A, where U_t = A(t) U.
353:    Also sets the location to store A.

355:    Collective on TS

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

366:    Calling sequence of func:
367: $     func (TS ts,PetscReal t,Mat *A,Mat *B,int *flag,void *ctx);

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

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

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

386:    Important: 
387:    The user MUST call either this routine or TSSetRHSFunction().

389:    Level: beginner

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

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

407:   ts->ops->rhsmatrix = f;
408:   ts->jacP           = ctx;
409:   ts->A              = A;
410:   ts->B              = B;

412:   return(0);
413: }

415: #undef __FUNCT__  
417: /*@C
418:    TSSetRHSJacobian - Sets the function to compute the Jacobian of F,
419:    where U_t = F(U,t), as well as the location to store the matrix.

421:    Collective on TS

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

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

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

442:    Notes: 
443:    See SLESSetOperators() for important information about setting the flag
444:    output parameter in the routine func().  Be sure to read this information!

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

452:    Level: beginner
453:    
454: .keywords: TS, timestep, set, right-hand-side, Jacobian

456: .seealso: TSDefaultComputeJacobianColor(),
457:           SNESDefaultComputeJacobianColor()

459: @*/
460: int TSSetRHSJacobian(TS ts,Mat A,Mat B,int (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
461: {
468:   if (ts->problem_type != TS_NONLINEAR) {
469:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
470:   }

472:   ts->ops->rhsjacobian = f;
473:   ts->jacP             = ctx;
474:   ts->A                = A;
475:   ts->B                = B;
476:   return(0);
477: }

479: #undef __FUNCT__  
481: /*
482:    TSComputeRHSBoundaryConditions - Evaluates the boundary condition function. 

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


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

503:   return(0);
504: }

506: #undef __FUNCT__  
508: /*@C
509:     TSSetRHSBoundaryConditions - Sets the routine for evaluating the function,
510:     boundary conditions for the function F.

512:     Collective on TS

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

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

523: +   t - current timestep
524: .   F - function vector
525: -   ctx - [optional] user-defined function context 

527:     Level: intermediate

529: .keywords: TS, timestep, set, boundary conditions, function
530: @*/
531: int TSSetRHSBoundaryConditions(TS ts,int (*f)(TS,PetscReal,Vec,void*),void *ctx)
532: {

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

544: #undef __FUNCT__  
546: /*@C
547:     TSView - Prints the TS data structure.

549:     Collective on TS

551:     Input Parameters:
552: +   ts - the TS context obtained from TSCreate()
553: -   viewer - visualization context

555:     Options Database Key:
556: .   -ts_view - calls TSView() at end of TSStep()

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

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

569:     Level: beginner

571: .keywords: TS, timestep, view

573: .seealso: PetscViewerASCIIOpen()
574: @*/
575: int TSView(TS ts,PetscViewer viewer)
576: {
577:   int        ierr;
578:   char       *type;
579:   PetscTruth isascii,isstring;

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

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


620: #undef __FUNCT__  
622: /*@C
623:    TSSetApplicationContext - Sets an optional user-defined context for 
624:    the timesteppers.

626:    Collective on TS

628:    Input Parameters:
629: +  ts - the TS context obtained from TSCreate()
630: -  usrP - optional user context

632:    Level: intermediate

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

636: .seealso: TSGetApplicationContext()
637: @*/
638: int TSSetApplicationContext(TS ts,void *usrP)
639: {
642:   ts->user = usrP;
643:   return(0);
644: }

646: #undef __FUNCT__  
648: /*@C
649:     TSGetApplicationContext - Gets the user-defined context for the 
650:     timestepper.

652:     Not Collective

654:     Input Parameter:
655: .   ts - the TS context obtained from TSCreate()

657:     Output Parameter:
658: .   usrP - user context

660:     Level: intermediate

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

664: .seealso: TSSetApplicationContext()
665: @*/
666: int TSGetApplicationContext(TS ts,void **usrP)
667: {
670:   *usrP = ts->user;
671:   return(0);
672: }

674: #undef __FUNCT__  
676: /*@
677:    TSGetTimeStepNumber - Gets the current number of timesteps.

679:    Not Collective

681:    Input Parameter:
682: .  ts - the TS context obtained from TSCreate()

684:    Output Parameter:
685: .  iter - number steps so far

687:    Level: intermediate

689: .keywords: TS, timestep, get, iteration, number
690: @*/
691: int TSGetTimeStepNumber(TS ts,int* iter)
692: {
695:   *iter = ts->steps;
696:   return(0);
697: }

699: #undef __FUNCT__  
701: /*@
702:    TSSetInitialTimeStep - Sets the initial timestep to be used, 
703:    as well as the initial time.

705:    Collective on TS

707:    Input Parameters:
708: +  ts - the TS context obtained from TSCreate()
709: .  initial_time - the initial time
710: -  time_step - the size of the timestep

712:    Level: intermediate

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

716: .keywords: TS, set, initial, timestep
717: @*/
718: int TSSetInitialTimeStep(TS ts,PetscReal initial_time,PetscReal time_step)
719: {
722:   ts->time_step         = time_step;
723:   ts->initial_time_step = time_step;
724:   ts->ptime             = initial_time;
725:   return(0);
726: }

728: #undef __FUNCT__  
730: /*@
731:    TSSetTimeStep - Allows one to reset the timestep at any time,
732:    useful for simple pseudo-timestepping codes.

734:    Collective on TS

736:    Input Parameters:
737: +  ts - the TS context obtained from TSCreate()
738: -  time_step - the size of the timestep

740:    Level: intermediate

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

744: .keywords: TS, set, timestep
745: @*/
746: int TSSetTimeStep(TS ts,PetscReal time_step)
747: {
750:   ts->time_step = time_step;
751:   return(0);
752: }

754: #undef __FUNCT__  
756: /*@
757:    TSGetTimeStep - Gets the current timestep size.

759:    Not Collective

761:    Input Parameter:
762: .  ts - the TS context obtained from TSCreate()

764:    Output Parameter:
765: .  dt - the current timestep size

767:    Level: intermediate

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

771: .keywords: TS, get, timestep
772: @*/
773: int TSGetTimeStep(TS ts,PetscReal* dt)
774: {
777:   *dt = ts->time_step;
778:   return(0);
779: }

781: #undef __FUNCT__  
783: /*@C
784:    TSGetSolution - Returns the solution at the present timestep. It
785:    is valid to call this routine inside the function that you are evaluating
786:    in order to move to the new timestep. This vector not changed until
787:    the solution at the next timestep has been calculated.

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

791:    Input Parameter:
792: .  ts - the TS context obtained from TSCreate()

794:    Output Parameter:
795: .  v - the vector containing the solution

797:    Level: intermediate

799: .seealso: TSGetTimeStep()

801: .keywords: TS, timestep, get, solution
802: @*/
803: int TSGetSolution(TS ts,Vec *v)
804: {
807:   *v = ts->vec_sol_always;
808:   return(0);
809: }

811: /* ----- Routines to initialize and destroy a timestepper ---- */
812: #undef __FUNCT__  
814: /*@C
815:   TSSetProblemType - Sets the type of problem to be solved.

817:   Not collective

819:   Input Parameters:
820: + ts   - The TS
821: - type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
822: .vb
823:          U_t = A U    
824:          U_t = A(t) U 
825:          U_t = F(t,U) 
826: .ve

828:    Level: beginner

830: .keywords: TS, problem type
831: .seealso: TSSetUp(), TSProblemType, TS
832: @*/
833: int TSSetProblemType(TS ts, TSProblemType type) {
836:   ts->problem_type = type;
837:   return(0);
838: }

840: #undef __FUNCT__  
842: /*@C
843:   TSGetProblemType - Gets the type of problem to be solved.

845:   Not collective

847:   Input Parameter:
848: . ts   - The TS

850:   Output Parameter:
851: . type - One of TS_LINEAR, TS_NONLINEAR where these types refer to problems of the forms
852: .vb
853:          U_t = A U    
854:          U_t = A(t) U 
855:          U_t = F(t,U) 
856: .ve

858:    Level: beginner

860: .keywords: TS, problem type
861: .seealso: TSSetUp(), TSProblemType, TS
862: @*/
863: int TSGetProblemType(TS ts, TSProblemType *type) {
867:   *type = ts->problem_type;
868:   return(0);
869: }

871: #undef __FUNCT__  
873: /*@
874:    TSSetUp - Sets up the internal data structures for the later use
875:    of a timestepper.  

877:    Collective on TS

879:    Input Parameter:
880: .  ts - the TS context obtained from TSCreate()

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

889:    Level: advanced

891: .keywords: TS, timestep, setup

893: .seealso: TSCreate(), TSStep(), TSDestroy()
894: @*/
895: int TSSetUp(TS ts)
896: {

901:   if (!ts->vec_sol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call TSSetSolution() first");
902:   if (!ts->type_name) {
903:     TSSetType(ts,TS_EULER);
904:   }
905:   (*ts->ops->setup)(ts);
906:   ts->setupcalled = 1;
907:   return(0);
908: }

910: #undef __FUNCT__  
912: /*@C
913:    TSDestroy - Destroys the timestepper context that was created
914:    with TSCreate().

916:    Collective on TS

918:    Input Parameter:
919: .  ts - the TS context obtained from TSCreate()

921:    Level: beginner

923: .keywords: TS, timestepper, destroy

925: .seealso: TSCreate(), TSSetUp(), TSSolve()
926: @*/
927: int TSDestroy(TS ts)
928: {
929:   int ierr,i;

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

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

938:   if (ts->sles) {SLESDestroy(ts->sles);}
939:   if (ts->snes) {SNESDestroy(ts->snes);}
940:   (*(ts)->ops->destroy)(ts);
941:   for (i=0; i<ts->numbermonitors; i++) {
942:     if (ts->mdestroy[i]) {
943:       (*ts->mdestroy[i])(ts->monitorcontext[i]);
944:     }
945:   }
946:   PetscLogObjectDestroy((PetscObject)ts);
947:   PetscHeaderDestroy((PetscObject)ts);
948:   return(0);
949: }

951: #undef __FUNCT__  
953: /*@C
954:    TSGetSNES - Returns the SNES (nonlinear solver) associated with 
955:    a TS (timestepper) context. Valid only for nonlinear problems.

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

959:    Input Parameter:
960: .  ts - the TS context obtained from TSCreate()

962:    Output Parameter:
963: .  snes - the nonlinear solver context

965:    Notes:
966:    The user can then directly manipulate the SNES context to set various
967:    options, etc.  Likewise, the user can then extract and manipulate the 
968:    SLES, KSP, and PC contexts as well.

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

973:    Level: beginner

975: .keywords: timestep, get, SNES
976: @*/
977: int TSGetSNES(TS ts,SNES *snes)
978: {
981:   if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetSLES()");
982:   *snes = ts->snes;
983:   return(0);
984: }

986: #undef __FUNCT__  
988: /*@C
989:    TSGetSLES - Returns the SLES (linear solver) associated with 
990:    a TS (timestepper) context.

992:    Not Collective, but SLES is parallel if TS is parallel

994:    Input Parameter:
995: .  ts - the TS context obtained from TSCreate()

997:    Output Parameter:
998: .  sles - the nonlinear solver context

1000:    Notes:
1001:    The user can then directly manipulate the SLES context to set various
1002:    options, etc.  Likewise, the user can then extract and manipulate the 
1003:    KSP and PC contexts as well.

1005:    TSGetSLES() does not work for integrators that do not use SLES;
1006:    in this case TSGetSLES() returns PETSC_NULL in sles.

1008:    Level: beginner

1010: .keywords: timestep, get, SLES
1011: @*/
1012: int TSGetSLES(TS ts,SLES *sles)
1013: {
1016:   if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1017:   *sles = ts->sles;
1018:   return(0);
1019: }

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

1023: #undef __FUNCT__  
1025: /*@
1026:    TSGetDuration - Gets the maximum number of timesteps to use and 
1027:    maximum time for iteration.

1029:    Collective on TS

1031:    Input Parameters:
1032: +  ts       - the TS context obtained from TSCreate()
1033: .  maxsteps - maximum number of iterations to use, or PETSC_NULL
1034: -  maxtime  - final time to iterate to, or PETSC_NULL

1036:    Level: intermediate

1038: .keywords: TS, timestep, get, maximum, iterations, time
1039: @*/
1040: int TSGetDuration(TS ts, int *maxsteps, PetscReal *maxtime)
1041: {
1044:   if (maxsteps != PETSC_NULL) {
1046:     *maxsteps = ts->max_steps;
1047:   }
1048:   if (maxtime  != PETSC_NULL) {
1050:     *maxtime  = ts->max_time;
1051:   }
1052:   return(0);
1053: }

1055: #undef __FUNCT__  
1057: /*@
1058:    TSSetDuration - Sets the maximum number of timesteps to use and 
1059:    maximum time for iteration.

1061:    Collective on TS

1063:    Input Parameters:
1064: +  ts - the TS context obtained from TSCreate()
1065: .  maxsteps - maximum number of iterations to use
1066: -  maxtime - final time to iterate to

1068:    Options Database Keys:
1069: .  -ts_max_steps <maxsteps> - Sets maxsteps
1070: .  -ts_max_time <maxtime> - Sets maxtime

1072:    Notes:
1073:    The default maximum number of iterations is 5000. Default time is 5.0

1075:    Level: intermediate

1077: .keywords: TS, timestep, set, maximum, iterations
1078: @*/
1079: int TSSetDuration(TS ts,int maxsteps,PetscReal maxtime)
1080: {
1083:   ts->max_steps = maxsteps;
1084:   ts->max_time  = maxtime;
1085:   return(0);
1086: }

1088: #undef __FUNCT__  
1090: /*@
1091:    TSSetSolution - Sets the initial solution vector
1092:    for use by the TS routines.

1094:    Collective on TS and Vec

1096:    Input Parameters:
1097: +  ts - the TS context obtained from TSCreate()
1098: -  x - the solution vector

1100:    Level: beginner

1102: .keywords: TS, timestep, set, solution, initial conditions
1103: @*/
1104: int TSSetSolution(TS ts,Vec x)
1105: {
1108:   ts->vec_sol        = ts->vec_sol_always = x;
1109:   return(0);
1110: }

1112: #undef __FUNCT__  
1114: /*@
1115:   TSSetRhsBC - Sets the function which applies boundary conditions
1116:   to the Rhs of each system.

1118:   Collective on TS

1120:   Input Parameters:
1121: + ts   - The TS context obtained from TSCreate()
1122: - func - The function

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

1127: + rhs - The current rhs vector
1128: - ctx - The user-context

1130:   Level: intermediate

1132: .keywords: TS, Rhs, boundary conditions
1133: @*/
1134: int TSSetRhsBC(TS ts, int (*func)(TS, Vec, void *))
1135: {
1138:   ts->ops->applyrhsbc = func;
1139:   return(0);
1140: }

1142: #undef __FUNCT__  
1144: /*@
1145:   TSDefaultRhsBC - The default boundary condition function which does nothing.

1147:   Collective on TS

1149:   Input Parameters:
1150: + ts  - The TS context obtained from TSCreate()
1151: . rhs - The Rhs
1152: - ctx - The user-context

1154:   Level: developer

1156: .keywords: TS, Rhs, boundary conditions
1157: @*/
1158: int TSDefaultRhsBC(TS ts,  Vec rhs, void *ctx)
1159: {
1161:   return(0);
1162: }

1164: #undef __FUNCT__  
1166: /*@
1167:   TSSetSystemMatrixBC - Sets the function which applies boundary conditions
1168:   to the system matrix and preconditioner of each system.

1170:   Collective on TS

1172:   Input Parameters:
1173: + ts   - The TS context obtained from TSCreate()
1174: - func - The function

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

1179: + A   - The current system matrix
1180: . B   - The current preconditioner
1181: - ctx - The user-context

1183:   Level: intermediate

1185: .keywords: TS, System matrix, boundary conditions
1186: @*/
1187: int TSSetSystemMatrixBC(TS ts, int (*func)(TS, Mat, Mat, void *))
1188: {
1191:   ts->ops->applymatrixbc = func;
1192:   return(0);
1193: }

1195: #undef __FUNCT__  
1197: /*@
1198:   TSDefaultSystemMatrixBC - The default boundary condition function which
1199:   does nothing.

1201:   Collective on TS

1203:   Input Parameters:
1204: + ts  - The TS context obtained from TSCreate()
1205: . A   - The system matrix
1206: . B   - The preconditioner
1207: - ctx - The user-context

1209:   Level: developer

1211: .keywords: TS, System matrix, boundary conditions
1212: @*/
1213: int TSDefaultSystemMatrixBC(TS ts, Mat A, Mat B, void *ctx)
1214: {
1216:   return(0);
1217: }

1219: #undef __FUNCT__  
1221: /*@
1222:   TSSetSolutionBC - Sets the function which applies boundary conditions
1223:   to the solution of each system. This is necessary in nonlinear systems
1224:   which time dependent boundary conditions.

1226:   Collective on TS

1228:   Input Parameters:
1229: + ts   - The TS context obtained from TSCreate()
1230: - func - The function

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

1235: + sol - The current solution vector
1236: - ctx - The user-context

1238:   Level: intermediate

1240: .keywords: TS, solution, boundary conditions
1241: @*/
1242: int TSSetSolutionBC(TS ts, int (*func)(TS, Vec, void *))
1243: {
1246:   ts->ops->applysolbc = func;
1247:   return(0);
1248: }

1250: #undef __FUNCT__  
1252: /*@
1253:   TSDefaultSolutionBC - The default boundary condition function which
1254:   does nothing.

1256:   Collective on TS

1258:   Input Parameters:
1259: + ts  - The TS context obtained from TSCreate()
1260: . sol - The solution
1261: - ctx - The user-context

1263:   Level: developer

1265: .keywords: TS, solution, boundary conditions
1266: @*/
1267: int TSDefaultSolutionBC(TS ts, Vec sol, void *ctx)
1268: {
1270:   return(0);
1271: }

1273: #undef __FUNCT__  
1275: /*@
1276:   TSSetPreStep - Sets the general-purpose function
1277:   called once at the beginning of time stepping.

1279:   Collective on TS

1281:   Input Parameters:
1282: + ts   - The TS context obtained from TSCreate()
1283: - func - The function

1285:   Calling sequence of func:
1286: . func (TS ts);

1288:   Level: intermediate

1290: .keywords: TS, timestep
1291: @*/
1292: int TSSetPreStep(TS ts, int (*func)(TS))
1293: {
1296:   ts->ops->prestep = func;
1297:   return(0);
1298: }

1300: #undef __FUNCT__  
1302: /*@
1303:   TSDefaultPreStep - The default pre-stepping function which does nothing.

1305:   Collective on TS

1307:   Input Parameters:
1308: . ts  - The TS context obtained from TSCreate()

1310:   Level: developer

1312: .keywords: TS, timestep
1313: @*/
1314: int TSDefaultPreStep(TS ts)
1315: {
1317:   return(0);
1318: }

1320: #undef __FUNCT__  
1322: /*@
1323:   TSSetUpdate - Sets the general-purpose update function called
1324:   at the beginning of every time step. This function can change
1325:   the time step.

1327:   Collective on TS

1329:   Input Parameters:
1330: + ts   - The TS context obtained from TSCreate()
1331: - func - The function

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

1336: + t   - The current time
1337: - dt  - The current time step

1339:   Level: intermediate

1341: .keywords: TS, update, timestep
1342: @*/
1343: int TSSetUpdate(TS ts, int (*func)(TS, PetscReal, PetscReal *))
1344: {
1347:   ts->ops->update = func;
1348:   return(0);
1349: }

1351: #undef __FUNCT__  
1353: /*@
1354:   TSDefaultUpdate - The default update function which does nothing.

1356:   Collective on TS

1358:   Input Parameters:
1359: + ts  - The TS context obtained from TSCreate()
1360: - t   - The current time

1362:   Output Parameters:
1363: . dt  - The current time step

1365:   Level: developer

1367: .keywords: TS, update, timestep
1368: @*/
1369: int TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1370: {
1372:   return(0);
1373: }

1375: #undef __FUNCT__  
1377: /*@
1378:   TSSetPostStep - Sets the general-purpose function
1379:   called once at the end of time stepping.

1381:   Collective on TS

1383:   Input Parameters:
1384: + ts   - The TS context obtained from TSCreate()
1385: - func - The function

1387:   Calling sequence of func:
1388: . func (TS ts);

1390:   Level: intermediate

1392: .keywords: TS, timestep
1393: @*/
1394: int TSSetPostStep(TS ts, int (*func)(TS))
1395: {
1398:   ts->ops->poststep = func;
1399:   return(0);
1400: }

1402: #undef __FUNCT__  
1404: /*@
1405:   TSDefaultPostStep - The default post-stepping function which does nothing.

1407:   Collective on TS

1409:   Input Parameters:
1410: . ts  - The TS context obtained from TSCreate()

1412:   Level: developer

1414: .keywords: TS, timestep
1415: @*/
1416: int TSDefaultPostStep(TS ts)
1417: {
1419:   return(0);
1420: }

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

1424: #undef __FUNCT__  
1426: /*@C
1427:    TSSetMonitor - Sets an ADDITIONAL function that is to be used at every
1428:    timestep to display the iteration's  progress.   

1430:    Collective on TS

1432:    Input Parameters:
1433: +  ts - the TS context obtained from TSCreate()
1434: .  func - monitoring routine
1435: .  mctx - [optional] user-defined context for private data for the 
1436:              monitor routine (use PETSC_NULL if no context is desired)
1437: -  monitordestroy - [optional] routine that frees monitor context
1438:           (may be PETSC_NULL)

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

1443: +    ts - the TS context
1444: .    steps - iteration number
1445: .    time - current time
1446: .    x - current iterate
1447: -    mctx - [optional] monitoring context

1449:    Notes:
1450:    This routine adds an additional monitor to the list of monitors that 
1451:    already has been loaded.

1453:    Level: intermediate

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

1457: .seealso: TSDefaultMonitor(), TSClearMonitor()
1458: @*/
1459: int TSSetMonitor(TS ts,int (*monitor)(TS,int,PetscReal,Vec,void*),void *mctx,int (*mdestroy)(void*))
1460: {
1463:   if (ts->numbermonitors >= MAXTSMONITORS) {
1464:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1465:   }
1466:   ts->monitor[ts->numbermonitors]           = monitor;
1467:   ts->mdestroy[ts->numbermonitors]          = mdestroy;
1468:   ts->monitorcontext[ts->numbermonitors++]  = (void*)mctx;
1469:   return(0);
1470: }

1472: #undef __FUNCT__  
1474: /*@C
1475:    TSClearMonitor - Clears all the monitors that have been set on a time-step object.   

1477:    Collective on TS

1479:    Input Parameters:
1480: .  ts - the TS context obtained from TSCreate()

1482:    Notes:
1483:    There is no way to remove a single, specific monitor.

1485:    Level: intermediate

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

1489: .seealso: TSDefaultMonitor(), TSSetMonitor()
1490: @*/
1491: int TSClearMonitor(TS ts)
1492: {
1495:   ts->numbermonitors = 0;
1496:   return(0);
1497: }

1499: #undef __FUNCT__  
1501: int TSDefaultMonitor(TS ts,int step,PetscReal ptime,Vec v,void *ctx)
1502: {

1506:   PetscPrintf(ts->comm,"timestep %d dt %g time %gn",step,ts->time_step,ptime);
1507:   return(0);
1508: }

1510: #undef __FUNCT__  
1512: /*@
1513:    TSStep - Steps the requested number of timesteps.

1515:    Collective on TS

1517:    Input Parameter:
1518: .  ts - the TS context obtained from TSCreate()

1520:    Output Parameters:
1521: +  steps - number of iterations until termination
1522: -  ptime - time until termination

1524:    Level: beginner

1526: .keywords: TS, timestep, solve

1528: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1529: @*/
1530: int TSStep(TS ts,int *steps,PetscReal *ptime)
1531: {
1532:   int        ierr;

1536:   if (!ts->setupcalled) {
1537:     TSSetUp(ts);
1538:   }

1540:   PetscLogEventBegin(TS_Step, ts, 0, 0, 0);
1541:   (*ts->ops->prestep)(ts);
1542:   (*ts->ops->step)(ts, steps, ptime);
1543:   (*ts->ops->poststep)(ts);
1544:   PetscLogEventEnd(TS_Step, ts, 0, 0, 0);

1546:   if (!PetscPreLoadingOn) {
1547:     TSViewFromOptions(ts,ts->name);
1548:   }
1549:   return(0);
1550: }

1552: #undef __FUNCT__  
1554: /*
1555:      Runs the user provided monitor routines, if they exists.
1556: */
1557: int TSMonitor(TS ts,int step,PetscReal ptime,Vec x)
1558: {
1559:   int i,ierr,n = ts->numbermonitors;

1562:   for (i=0; i<n; i++) {
1563:     (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1564:   }
1565:   return(0);
1566: }

1568: /* ------------------------------------------------------------------------*/

1570: #undef __FUNCT__  
1572: /*@C
1573:    TSLGMonitorCreate - Creates a line graph context for use with 
1574:    TS to monitor convergence of preconditioned residual norms.

1576:    Collective on TS

1578:    Input Parameters:
1579: +  host - the X display to open, or null for the local machine
1580: .  label - the title to put in the title bar
1581: .  x, y - the screen coordinates of the upper left coordinate of the window
1582: -  m, n - the screen width and height in pixels

1584:    Output Parameter:
1585: .  draw - the drawing context

1587:    Options Database Key:
1588: .  -ts_xmonitor - automatically sets line graph monitor

1590:    Notes: 
1591:    Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().

1593:    Level: intermediate

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

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

1599: @*/
1600: int TSLGMonitorCreate(char *host,char *label,int x,int y,int m,int n,PetscDrawLG *draw)
1601: {
1602:   PetscDraw win;
1603:   int       ierr;

1606:   PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1607:   PetscDrawSetType(win,PETSC_DRAW_X);
1608:   PetscDrawLGCreate(win,1,draw);
1609:   PetscDrawLGIndicateDataPoints(*draw);

1611:   PetscLogObjectParent(*draw,win);
1612:   return(0);
1613: }

1615: #undef __FUNCT__  
1617: int TSLGMonitor(TS ts,int n,PetscReal ptime,Vec v,void *monctx)
1618: {
1619:   PetscDrawLG lg = (PetscDrawLG) monctx;
1620:   PetscReal      x,y = ptime;
1621:   int         ierr;

1624:   if (!monctx) {
1625:     MPI_Comm    comm;
1626:     PetscViewer viewer;

1628:     ierr   = PetscObjectGetComm((PetscObject)ts,&comm);
1629:     viewer = PETSC_VIEWER_DRAW_(comm);
1630:     ierr   = PetscViewerDrawGetDrawLG(viewer,0,&lg);
1631:   }

1633:   if (!n) {PetscDrawLGReset(lg);}
1634:   x = (PetscReal)n;
1635:   PetscDrawLGAddPoint(lg,&x,&y);
1636:   if (n < 20 || (n % 5)) {
1637:     PetscDrawLGDraw(lg);
1638:   }
1639:   return(0);
1640: }

1642: #undef __FUNCT__  
1644: /*@C
1645:    TSLGMonitorDestroy - Destroys a line graph context that was created 
1646:    with TSLGMonitorCreate().

1648:    Collective on PetscDrawLG

1650:    Input Parameter:
1651: .  draw - the drawing context

1653:    Level: intermediate

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

1657: .seealso: TSLGMonitorCreate(),  TSSetMonitor(), TSLGMonitor();
1658: @*/
1659: int TSLGMonitorDestroy(PetscDrawLG drawlg)
1660: {
1661:   PetscDraw draw;
1662:   int       ierr;

1665:   PetscDrawLGGetDraw(drawlg,&draw);
1666:   PetscDrawDestroy(draw);
1667:   PetscDrawLGDestroy(drawlg);
1668:   return(0);
1669: }

1671: #undef __FUNCT__  
1673: /*@
1674:    TSGetTime - Gets the current time.

1676:    Not Collective

1678:    Input Parameter:
1679: .  ts - the TS context obtained from TSCreate()

1681:    Output Parameter:
1682: .  t  - the current time

1684:    Contributed by: Matthew Knepley

1686:    Level: beginner

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

1690: .keywords: TS, get, time
1691: @*/
1692: int TSGetTime(TS ts,PetscReal* t)
1693: {
1696:   *t = ts->ptime;
1697:   return(0);
1698: }

1700: #undef __FUNCT__
1702: /*@C
1703:    TSSetOptionsPrefix - Sets the prefix used for searching for all
1704:    TS options in the database.

1706:    Collective on TS

1708:    Input Parameter:
1709: +  ts     - The TS context
1710: -  prefix - The prefix to prepend to all option names

1712:    Notes:
1713:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1714:    The first character of all runtime options is AUTOMATICALLY the
1715:    hyphen.

1717:    Contributed by: Matthew Knepley

1719:    Level: advanced

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

1723: .seealso: TSSetFromOptions()

1725: @*/
1726: int TSSetOptionsPrefix(TS ts,char *prefix)
1727: {

1732:   PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1733:   switch(ts->problem_type) {
1734:     case TS_NONLINEAR:
1735:       SNESSetOptionsPrefix(ts->snes,prefix);
1736:       break;
1737:     case TS_LINEAR:
1738:       SLESSetOptionsPrefix(ts->sles,prefix);
1739:       break;
1740:   }
1741:   return(0);
1742: }


1745: #undef __FUNCT__
1747: /*@C
1748:    TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1749:    TS options in the database.

1751:    Collective on TS

1753:    Input Parameter:
1754: +  ts     - The TS context
1755: -  prefix - The prefix to prepend to all option names

1757:    Notes:
1758:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1759:    The first character of all runtime options is AUTOMATICALLY the
1760:    hyphen.

1762:    Contributed by: Matthew Knepley

1764:    Level: advanced

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

1768: .seealso: TSGetOptionsPrefix()

1770: @*/
1771: int TSAppendOptionsPrefix(TS ts,char *prefix)
1772: {

1777:   PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1778:   switch(ts->problem_type) {
1779:     case TS_NONLINEAR:
1780:       SNESAppendOptionsPrefix(ts->snes,prefix);
1781:       break;
1782:     case TS_LINEAR:
1783:       SLESAppendOptionsPrefix(ts->sles,prefix);
1784:       break;
1785:   }
1786:   return(0);
1787: }

1789: #undef __FUNCT__
1791: /*@C
1792:    TSGetOptionsPrefix - Sets the prefix used for searching for all
1793:    TS options in the database.

1795:    Not Collective

1797:    Input Parameter:
1798: .  ts - The TS context

1800:    Output Parameter:
1801: .  prefix - A pointer to the prefix string used

1803:    Contributed by: Matthew Knepley

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

1808:    Level: intermediate

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

1812: .seealso: TSAppendOptionsPrefix()
1813: @*/
1814: int TSGetOptionsPrefix(TS ts,char **prefix)
1815: {

1820:   PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1821:   return(0);
1822: }

1824: #undef __FUNCT__
1826: /*@C
1827:    TSGetRHSMatrix - Returns the matrix A at the present timestep.

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

1831:    Input Parameter:
1832: .  ts  - The TS context obtained from TSCreate()

1834:    Output Parameters:
1835: +  A   - The matrix A, where U_t = A(t) U
1836: .  M   - The preconditioner matrix, usually the same as A
1837: -  ctx - User-defined context for matrix evaluation routine

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

1841:    Contributed by: Matthew Knepley

1843:    Level: intermediate

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

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

1849: @*/
1850: int TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1851: {
1854:   if (A)   *A = ts->A;
1855:   if (M)   *M = ts->B;
1856:   if (ctx) *ctx = ts->jacP;
1857:   return(0);
1858: }

1860: #undef __FUNCT__
1862: /*@C
1863:    TSGetRHSJacobian - Returns the Jacobian J at the present timestep.

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

1867:    Input Parameter:
1868: .  ts  - The TS context obtained from TSCreate()

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

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

1877:    Contributed by: Matthew Knepley

1879:    Level: intermediate

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

1883: .keywords: TS, timestep, get, matrix, Jacobian
1884: @*/
1885: int TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1886: {

1890:   TSGetRHSMatrix(ts,J,M,ctx);
1891:   return(0);
1892: }

1894: #undef __FUNCT__  
1896: /*@C
1897:    TSVecViewMonitor - Monitors progress of the TS solvers by calling 
1898:    VecView() for the solution at each timestep

1900:    Collective on TS

1902:    Input Parameters:
1903: +  ts - the TS context
1904: .  step - current time-step
1905: .  ptime - current time
1906: -  dummy - either a viewer or PETSC_NULL

1908:    Level: intermediate

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

1912: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1913: @*/
1914: int TSVecViewMonitor(TS ts,int step,PetscReal ptime,Vec x,void *dummy)
1915: {
1916:   int         ierr;
1917:   PetscViewer viewer = (PetscViewer) dummy;

1920:   if (!viewer) {
1921:     MPI_Comm comm;
1922:     ierr   = PetscObjectGetComm((PetscObject)ts,&comm);
1923:     viewer = PETSC_VIEWER_DRAW_(comm);
1924:   }
1925:   VecView(x,viewer);
1926:   return(0);
1927: }