Actual source code: ts.c

  1: #define PETSCTS_DLL

 3:  #include src/ts/tsimpl.h

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

 11: /*
 12:   TSSetTypeFromOptions - Sets the type of ts from user options.

 14:   Collective on TS

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

 19:   Level: intermediate

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

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

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

 52: /*@
 53:    TSSetFromOptions - Sets various TS parameters from user options.

 55:    Collective on TS

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

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

 68:    Level: beginner

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

 72: .seealso: TSGetType
 73: @*/
 74: PetscErrorCode  TSSetFromOptions(TS ts)
 75: {
 76:   PetscReal      dt;
 77:   PetscTruth     opt,flg;
 79:   PetscViewer    monviewer;
 80:   char           monfilename[PETSC_MAX_PATH_LEN];

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

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

 95:     /* Monitor options */
 96:     PetscOptionsString("-ts_monitor","Monitor timestep size","TSDefaultMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
 97:     if (flg) {
 98:       PetscViewerASCIIOpen(ts->comm,monfilename,&monviewer);
 99:       TSSetMonitor(ts,TSDefaultMonitor,monviewer,(PetscErrorCode (*)(void*))PetscViewerDestroy);
100:     }
101:     PetscOptionsName("-ts_xmonitor","Monitor timestep size graphically","TSLGMonitor",&opt);
102:     if (opt) {
103:       TSSetMonitor(ts,TSLGMonitor,PETSC_NULL,PETSC_NULL);
104:     }
105:     PetscOptionsName("-ts_vecmonitor","Monitor solution graphically","TSVecViewMonitor",&opt);
106:     if (opt) {
107:       TSSetMonitor(ts,TSVecViewMonitor,PETSC_NULL,PETSC_NULL);
108:     }

110:     /* Handle TS type options */
111:     TSSetTypeFromOptions(ts);

113:     /* Handle specific TS options */
114:     if (ts->ops->setfromoptions) {
115:       (*ts->ops->setfromoptions)(ts);
116:     }
117:   PetscOptionsEnd();

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

141:   return(0);
142: }

144: #undef  __FUNCT__
146: /*@
147:   TSViewFromOptions - This function visualizes the ts based upon user options.

149:   Collective on TS

151:   Input Parameter:
152: . ts - The ts

154:   Level: intermediate

156: .keywords: TS, view, options, database
157: .seealso: TSSetFromOptions(), TSView()
158: @*/
159: PetscErrorCode  TSViewFromOptions(TS ts,const char title[])
160: {
161:   PetscViewer    viewer;
162:   PetscDraw      draw;
163:   PetscTruth     opt;
164:   char           fileName[PETSC_MAX_PATH_LEN];

168:   PetscOptionsGetString(ts->prefix, "-ts_view", fileName, PETSC_MAX_PATH_LEN, &opt);
169:   if (opt && !PetscPreLoadingOn) {
170:     PetscViewerASCIIOpen(ts->comm,fileName,&viewer);
171:     TSView(ts, viewer);
172:     PetscViewerDestroy(viewer);
173:   }
174:   PetscOptionsHasName(ts->prefix, "-ts_view_draw", &opt);
175:   if (opt) {
176:     PetscViewerDrawOpen(ts->comm, 0, 0, 0, 0, 300, 300, &viewer);
177:     PetscViewerDrawGetDraw(viewer, 0, &draw);
178:     if (title) {
179:       PetscDrawSetTitle(draw, (char *)title);
180:     } else {
181:       PetscObjectName((PetscObject) ts);
182:       PetscDrawSetTitle(draw, ts->name);
183:     }
184:     TSView(ts, viewer);
185:     PetscViewerFlush(viewer);
186:     PetscDrawPause(draw);
187:     PetscViewerDestroy(viewer);
188:   }
189:   return(0);
190: }

194: /*@
195:    TSComputeRHSJacobian - Computes the Jacobian matrix that has been
196:       set with TSSetRHSJacobian().

198:    Collective on TS and Vec

200:    Input Parameters:
201: +  ts - the SNES context
202: .  t - current timestep
203: -  x - input vector

205:    Output Parameters:
206: +  A - Jacobian matrix
207: .  B - optional preconditioning matrix
208: -  flag - flag indicating matrix structure

210:    Notes: 
211:    Most users should not need to explicitly call this routine, as it 
212:    is used internally within the nonlinear solvers. 

214:    See KSPSetOperators() for important information about setting the
215:    flag parameter.

217:    TSComputeJacobian() is valid only for TS_NONLINEAR

219:    Level: developer

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

223: .seealso:  TSSetRHSJacobian(), KSPSetOperators()
224: @*/
225: PetscErrorCode  TSComputeRHSJacobian(TS ts,PetscReal t,Vec X,Mat *A,Mat *B,MatStructure *flg)
226: {

233:   if (ts->problem_type != TS_NONLINEAR) {
234:     SETERRQ(PETSC_ERR_ARG_WRONG,"For TS_NONLINEAR only");
235:   }
236:   if (ts->ops->rhsjacobian) {
238:     *flg = DIFFERENT_NONZERO_PATTERN;
239:     PetscStackPush("TS user Jacobian function");
240:     (*ts->ops->rhsjacobian)(ts,t,X,A,B,flg,ts->jacP);
241:     PetscStackPop;
243:     /* make sure user returned a correct Jacobian and preconditioner */
246:   } else {
247:     MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
248:     MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
249:     if (*A != *B) {
250:       MatAssemblyBegin(*B,MAT_FINAL_ASSEMBLY);
251:       MatAssemblyEnd(*B,MAT_FINAL_ASSEMBLY);
252:     }
253:   }
254:   return(0);
255: }

259: /*
260:    TSComputeRHSFunction - Evaluates the right-hand-side function. 

262:    Note: If the user did not provide a function but merely a matrix,
263:    this routine applies the matrix.
264: */
265: PetscErrorCode TSComputeRHSFunction(TS ts,PetscReal t,Vec x,Vec y)
266: {


275:   if (ts->ops->rhsfunction) {
276:     PetscStackPush("TS user right-hand-side function");
277:     (*ts->ops->rhsfunction)(ts,t,x,y,ts->funP);
278:     PetscStackPop;
279:   } else {
280:     if (ts->ops->rhsmatrix) { /* assemble matrix for this timestep */
281:       MatStructure flg;
282:       PetscStackPush("TS user right-hand-side matrix function");
283:       (*ts->ops->rhsmatrix)(ts,t,&ts->A,&ts->B,&flg,ts->jacP);
284:       PetscStackPop;
285:     }
286:     MatMult(ts->A,x,y);
287:   }


291:   return(0);
292: }

296: /*@C
297:     TSSetRHSFunction - Sets the routine for evaluating the function,
298:     F(t,u), where U_t = F(t,u).

300:     Collective on TS

302:     Input Parameters:
303: +   ts - the TS context obtained from TSCreate()
304: .   f - routine for evaluating the right-hand-side function
305: -   ctx - [optional] user-defined context for private data for the 
306:           function evaluation routine (may be PETSC_NULL)

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

311: +   t - current timestep
312: .   u - input vector
313: .   F - function vector
314: -   ctx - [optional] user-defined function context 

316:     Important: 
317:     The user MUST call either this routine or TSSetRHSMatrix().

319:     Level: beginner

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

323: .seealso: TSSetRHSMatrix()
324: @*/
325: PetscErrorCode  TSSetRHSFunction(TS ts,PetscErrorCode (*f)(TS,PetscReal,Vec,Vec,void*),void *ctx)
326: {

330:   if (ts->problem_type == TS_LINEAR) {
331:     SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot set function for linear problem");
332:   }
333:   ts->ops->rhsfunction = f;
334:   ts->funP             = ctx;
335:   return(0);
336: }

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

344:    Collective on TS

346:    Input Parameters:
347: +  ts  - the TS context obtained from TSCreate()
348: .  A   - matrix
349: .  B   - preconditioner matrix (usually same as A)
350: .  f   - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
351:          if A is not a function of t.
352: -  ctx - [optional] user-defined context for private data for the 
353:           matrix evaluation routine (may be PETSC_NULL)

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

358: +  t - current timestep
359: .  A - matrix A, where U_t = A(t) U
360: .  B - preconditioner matrix, usually the same as A
361: .  flag - flag indicating information about the preconditioner matrix
362:           structure (same as flag in KSPSetOperators())
363: -  ctx - [optional] user-defined context for matrix evaluation routine

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

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

375:    Important: 
376:    The user MUST call either this routine or TSSetRHSFunction().

378:    Level: beginner

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

382: .seealso: TSSetRHSFunction()
383: @*/
384: PetscErrorCode  TSSetRHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
385: {
392:   if (ts->problem_type == TS_NONLINEAR) {
393:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems; use TSSetRHSJacobian()");
394:   }

396:   ts->ops->rhsmatrix = f;
397:   ts->jacP           = ctx;
398:   ts->A              = A;
399:   ts->B              = B;

401:   return(0);
402: }

406: /*@C
407:    TSSetLHSMatrix - Sets the function to compute the matrix A, where A U_t = F(U).
408:    Also sets the location to store A.

410:    Collective on TS

412:    Input Parameters:
413: +  ts  - the TS context obtained from TSCreate()
414: .  A   - matrix
415: .  B   - preconditioner matrix (usually same as A)
416: .  f   - the matrix evaluation routine; use PETSC_NULL (PETSC_NULL_FUNCTION in fortran)
417:          if A is not a function of t.
418: -  ctx - [optional] user-defined context for private data for the 
419:           matrix evaluation routine (may be PETSC_NULL)

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

424: +  t - current timestep
425: .  A - matrix A, where A U_t = F(U)
426: .  B - preconditioner matrix, usually the same as A
427: .  flag - flag indicating information about the preconditioner matrix
428:           structure (same as flag in KSPSetOperators())
429: -  ctx - [optional] user-defined context for matrix evaluation routine

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

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

441:    Notes: 
442:    Currently, TSSetLHSMatrix() only supports the TS_BEULER type.

444:    Level: beginner

446: .keywords: TS, timestep, set, left-hand-side, matrix

448: .seealso: TSSetRHSMatrix()
449: @*/
450: PetscErrorCode  TSSetLHSMatrix(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void *ctx)
451: {
452:   PetscTruth     sametype;

461: 
462:   if (!ts->type_name) SETERRQ(PETSC_ERR_ARG_NULL,"TS type must be set before calling TSSetLHSMatrix()");
463:   PetscTypeCompare((PetscObject)ts,"beuler",&sametype);
464:   if (!sametype)
465:     SETERRQ1(PETSC_ERR_SUP,"TS type %s not supported for LHSMatrix",ts->type_name);
466:   if (ts->problem_type == TS_NONLINEAR) {
467:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for nonlinear problems yet");
468:   }

470:   ts->ops->lhsmatrix = f;
471:   ts->jacPlhs        = ctx;
472:   ts->Alhs           = A;
473:   ts->Blhs           = B;
474:   return(0);
475: }

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

484:    Collective on TS

486:    Input Parameters:
487: +  ts  - the TS context obtained from TSCreate()
488: .  A   - Jacobian matrix
489: .  B   - preconditioner matrix (usually same as A)
490: .  f   - the Jacobian evaluation routine
491: -  ctx - [optional] user-defined context for private data for the 
492:          Jacobian evaluation routine (may be PETSC_NULL)

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

497: +  t - current timestep
498: .  u - input vector
499: .  A - matrix A, where U_t = A(t)u
500: .  B - preconditioner matrix, usually the same as A
501: .  flag - flag indicating information about the preconditioner matrix
502:           structure (same as flag in KSPSetOperators())
503: -  ctx - [optional] user-defined context for matrix evaluation routine

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

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

515:    Level: beginner
516:    
517: .keywords: TS, timestep, set, right-hand-side, Jacobian

519: .seealso: TSDefaultComputeJacobianColor(),
520:           SNESDefaultComputeJacobianColor(), TSSetRHSFunction(), TSSetRHSMatrix()

522: @*/
523: PetscErrorCode  TSSetRHSJacobian(TS ts,Mat A,Mat B,PetscErrorCode (*f)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
524: {
531:   if (ts->problem_type != TS_NONLINEAR) {
532:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not for linear problems; use TSSetRHSMatrix()");
533:   }

535:   ts->ops->rhsjacobian = f;
536:   ts->jacP             = ctx;
537:   ts->A                = A;
538:   ts->B                = B;
539:   return(0);
540: }

544: /*@C
545:     TSView - Prints the TS data structure.

547:     Collective on TS

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

553:     Options Database Key:
554: .   -ts_view - calls TSView() at end of TSStep()

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

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

567:     Level: beginner

569: .keywords: TS, timestep, view

571: .seealso: PetscViewerASCIIOpen()
572: @*/
573: PetscErrorCode  TSView(TS ts,PetscViewer viewer)
574: {
576:   char           *type;
577:   PetscTruth     iascii,isstring;

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

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


620: /*@C
621:    TSSetApplicationContext - Sets an optional user-defined context for 
622:    the timesteppers.

624:    Collective on TS

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

630:    Level: intermediate

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

634: .seealso: TSGetApplicationContext()
635: @*/
636: PetscErrorCode  TSSetApplicationContext(TS ts,void *usrP)
637: {
640:   ts->user = usrP;
641:   return(0);
642: }

646: /*@C
647:     TSGetApplicationContext - Gets the user-defined context for the 
648:     timestepper.

650:     Not Collective

652:     Input Parameter:
653: .   ts - the TS context obtained from TSCreate()

655:     Output Parameter:
656: .   usrP - user context

658:     Level: intermediate

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

662: .seealso: TSSetApplicationContext()
663: @*/
664: PetscErrorCode  TSGetApplicationContext(TS ts,void **usrP)
665: {
668:   *usrP = ts->user;
669:   return(0);
670: }

674: /*@
675:    TSGetTimeStepNumber - Gets the current number of timesteps.

677:    Not Collective

679:    Input Parameter:
680: .  ts - the TS context obtained from TSCreate()

682:    Output Parameter:
683: .  iter - number steps so far

685:    Level: intermediate

687: .keywords: TS, timestep, get, iteration, number
688: @*/
689: PetscErrorCode  TSGetTimeStepNumber(TS ts,PetscInt* iter)
690: {
694:   *iter = ts->steps;
695:   return(0);
696: }

700: /*@
701:    TSSetInitialTimeStep - Sets the initial timestep to be used, 
702:    as well as the initial time.

704:    Collective on TS

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

711:    Level: intermediate

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

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

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

733:    Collective on TS

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

739:    Level: intermediate

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

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

755: /*@
756:    TSGetTimeStep - Gets the current timestep size.

758:    Not Collective

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

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

766:    Level: intermediate

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

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

783: /*@
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: PetscErrorCode  TSGetSolution(TS ts,Vec *v)
804: {
808:   *v = ts->vec_sol_always;
809:   return(0);
810: }

812: /* ----- Routines to initialize and destroy a timestepper ---- */
815: /*@
816:   TSSetProblemType - Sets the type of problem to be solved.

818:   Not collective

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

829:    Level: beginner

831: .keywords: TS, problem type
832: .seealso: TSSetUp(), TSProblemType, TS
833: @*/
834: PetscErrorCode  TSSetProblemType(TS ts, TSProblemType type)
835: {
838:   ts->problem_type = type;
839:   return(0);
840: }

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

847:   Not collective

849:   Input Parameter:
850: . ts   - The TS

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

860:    Level: beginner

862: .keywords: TS, problem type
863: .seealso: TSSetUp(), TSProblemType, TS
864: @*/
865: PetscErrorCode  TSGetProblemType(TS ts, TSProblemType *type)
866: {
870:   *type = ts->problem_type;
871:   return(0);
872: }

876: /*@
877:    TSSetUp - Sets up the internal data structures for the later use
878:    of a timestepper.  

880:    Collective on TS

882:    Input Parameter:
883: .  ts - the TS context obtained from TSCreate()

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

892:    Level: advanced

894: .keywords: TS, timestep, setup

896: .seealso: TSCreate(), TSStep(), TSDestroy()
897: @*/
898: PetscErrorCode  TSSetUp(TS ts)
899: {

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

915: /*@
916:    TSDestroy - Destroys the timestepper context that was created
917:    with TSCreate().

919:    Collective on TS

921:    Input Parameter:
922: .  ts - the TS context obtained from TSCreate()

924:    Level: beginner

926: .keywords: TS, timestepper, destroy

928: .seealso: TSCreate(), TSSetUp(), TSSolve()
929: @*/
930: PetscErrorCode  TSDestroy(TS ts)
931: {

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

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

941:   if (ts->ksp) {KSPDestroy(ts->ksp);}
942:   if (ts->snes) {SNESDestroy(ts->snes);}
943:   if (ts->ops->destroy) {(*(ts)->ops->destroy)(ts);}
944:   TSClearMonitor(ts);
945:   PetscHeaderDestroy(ts);
946:   return(0);
947: }

951: /*@
952:    TSGetSNES - Returns the SNES (nonlinear solver) associated with 
953:    a TS (timestepper) context. Valid only for nonlinear problems.

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

957:    Input Parameter:
958: .  ts - the TS context obtained from TSCreate()

960:    Output Parameter:
961: .  snes - the nonlinear solver context

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

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

971:    Level: beginner

973: .keywords: timestep, get, SNES
974: @*/
975: PetscErrorCode  TSGetSNES(TS ts,SNES *snes)
976: {
980:   if (ts->problem_type == TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Nonlinear only; use TSGetKSP()");
981:   *snes = ts->snes;
982:   return(0);
983: }

987: /*@
988:    TSGetKSP - Returns the KSP (linear solver) associated with 
989:    a TS (timestepper) context.

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

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

996:    Output Parameter:
997: .  ksp - the nonlinear solver context

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

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

1007:    Level: beginner

1009: .keywords: timestep, get, KSP
1010: @*/
1011: PetscErrorCode  TSGetKSP(TS ts,KSP *ksp)
1012: {
1016:   if (ts->problem_type != TS_LINEAR) SETERRQ(PETSC_ERR_ARG_WRONG,"Linear only; use TSGetSNES()");
1017:   *ksp = ts->ksp;
1018:   return(0);
1019: }

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

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: PetscErrorCode  TSGetDuration(TS ts, PetscInt *maxsteps, PetscReal *maxtime)
1041: {
1044:   if (maxsteps) {
1046:     *maxsteps = ts->max_steps;
1047:   }
1048:   if (maxtime ) {
1050:     *maxtime  = ts->max_time;
1051:   }
1052:   return(0);
1053: }

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: PetscErrorCode  TSSetDuration(TS ts,PetscInt maxsteps,PetscReal maxtime)
1080: {
1083:   ts->max_steps = maxsteps;
1084:   ts->max_time  = maxtime;
1085:   return(0);
1086: }

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: PetscErrorCode  TSSetSolution(TS ts,Vec x)
1105: {
1109:   ts->vec_sol        = ts->vec_sol_always = x;
1110:   return(0);
1111: }

1115: /*@C
1116:   TSSetPreStep - Sets the general-purpose function
1117:   called once at the beginning of time stepping.

1119:   Collective on TS

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

1125:   Calling sequence of func:
1126: . func (TS ts);

1128:   Level: intermediate

1130: .keywords: TS, timestep
1131: @*/
1132: PetscErrorCode  TSSetPreStep(TS ts, PetscErrorCode (*func)(TS))
1133: {
1136:   ts->ops->prestep = func;
1137:   return(0);
1138: }

1142: /*@
1143:   TSDefaultPreStep - The default pre-stepping function which does nothing.

1145:   Collective on TS

1147:   Input Parameters:
1148: . ts  - The TS context obtained from TSCreate()

1150:   Level: developer

1152: .keywords: TS, timestep
1153: @*/
1154: PetscErrorCode  TSDefaultPreStep(TS ts)
1155: {
1157:   return(0);
1158: }

1162: /*@C
1163:   TSSetUpdate - Sets the general-purpose update function called
1164:   at the beginning of every time step. This function can change
1165:   the time step.

1167:   Collective on TS

1169:   Input Parameters:
1170: + ts   - The TS context obtained from TSCreate()
1171: - func - The function

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

1176: + t   - The current time
1177: - dt  - The current time step

1179:   Level: intermediate

1181: .keywords: TS, update, timestep
1182: @*/
1183: PetscErrorCode  TSSetUpdate(TS ts, PetscErrorCode (*func)(TS, PetscReal, PetscReal *))
1184: {
1187:   ts->ops->update = func;
1188:   return(0);
1189: }

1193: /*@
1194:   TSDefaultUpdate - The default update function which does nothing.

1196:   Collective on TS

1198:   Input Parameters:
1199: + ts  - The TS context obtained from TSCreate()
1200: - t   - The current time

1202:   Output Parameters:
1203: . dt  - The current time step

1205:   Level: developer

1207: .keywords: TS, update, timestep
1208: @*/
1209: PetscErrorCode  TSDefaultUpdate(TS ts, PetscReal t, PetscReal *dt)
1210: {
1212:   return(0);
1213: }

1217: /*@C
1218:   TSSetPostStep - Sets the general-purpose function
1219:   called once at the end of time stepping.

1221:   Collective on TS

1223:   Input Parameters:
1224: + ts   - The TS context obtained from TSCreate()
1225: - func - The function

1227:   Calling sequence of func:
1228: . func (TS ts);

1230:   Level: intermediate

1232: .keywords: TS, timestep
1233: @*/
1234: PetscErrorCode  TSSetPostStep(TS ts, PetscErrorCode (*func)(TS))
1235: {
1238:   ts->ops->poststep = func;
1239:   return(0);
1240: }

1244: /*@
1245:   TSDefaultPostStep - The default post-stepping function which does nothing.

1247:   Collective on TS

1249:   Input Parameters:
1250: . ts  - The TS context obtained from TSCreate()

1252:   Level: developer

1254: .keywords: TS, timestep
1255: @*/
1256: PetscErrorCode  TSDefaultPostStep(TS ts)
1257: {
1259:   return(0);
1260: }

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

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

1270:    Collective on TS

1272:    Input Parameters:
1273: +  ts - the TS context obtained from TSCreate()
1274: .  func - monitoring routine
1275: .  mctx - [optional] user-defined context for private data for the 
1276:              monitor routine (use PETSC_NULL if no context is desired)
1277: -  monitordestroy - [optional] routine that frees monitor context
1278:           (may be PETSC_NULL)

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

1283: +    ts - the TS context
1284: .    steps - iteration number
1285: .    time - current time
1286: .    x - current iterate
1287: -    mctx - [optional] monitoring context

1289:    Notes:
1290:    This routine adds an additional monitor to the list of monitors that 
1291:    already has been loaded.

1293:    Level: intermediate

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

1297: .seealso: TSDefaultMonitor(), TSClearMonitor()
1298: @*/
1299: PetscErrorCode  TSSetMonitor(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void*))
1300: {
1303:   if (ts->numbermonitors >= MAXTSMONITORS) {
1304:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
1305:   }
1306:   ts->monitor[ts->numbermonitors]           = monitor;
1307:   ts->mdestroy[ts->numbermonitors]          = mdestroy;
1308:   ts->monitorcontext[ts->numbermonitors++]  = (void*)mctx;
1309:   return(0);
1310: }

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

1317:    Collective on TS

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

1322:    Notes:
1323:    There is no way to remove a single, specific monitor.

1325:    Level: intermediate

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

1329: .seealso: TSDefaultMonitor(), TSSetMonitor()
1330: @*/
1331: PetscErrorCode  TSClearMonitor(TS ts)
1332: {
1334:   PetscInt       i;

1338:   for (i=0; i<ts->numbermonitors; i++) {
1339:     if (ts->mdestroy[i]) {
1340:       (*ts->mdestroy[i])(ts->monitorcontext[i]);
1341:     }
1342:   }
1343:   ts->numbermonitors = 0;
1344:   return(0);
1345: }

1349: /*@
1350:    TSDefaultMonitor - Sets the Default monitor
1351: @*/
1352: PetscErrorCode TSDefaultMonitor(TS ts,PetscInt step,PetscReal ptime,Vec v,void *ctx)
1353: {
1355:   PetscViewer    viewer = (PetscViewer)ctx;

1358:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ts->comm);
1359:   PetscViewerASCIIPrintf(viewer,"timestep %D dt %G time %G\n",step,ts->time_step,ptime);
1360:   return(0);
1361: }

1365: /*@
1366:    TSStep - Steps the requested number of timesteps.

1368:    Collective on TS

1370:    Input Parameter:
1371: .  ts - the TS context obtained from TSCreate()

1373:    Output Parameters:
1374: +  steps - number of iterations until termination
1375: -  ptime - time until termination

1377:    Level: beginner

1379: .keywords: TS, timestep, solve

1381: .seealso: TSCreate(), TSSetUp(), TSDestroy()
1382: @*/
1383: PetscErrorCode  TSStep(TS ts,PetscInt *steps,PetscReal *ptime)
1384: {

1389:   if (!ts->setupcalled) {
1390:     TSSetUp(ts);
1391:   }

1394:   (*ts->ops->prestep)(ts);
1395:   (*ts->ops->step)(ts, steps, ptime);
1396:   (*ts->ops->poststep)(ts);

1399:   if (!PetscPreLoadingOn) {
1400:     TSViewFromOptions(ts,ts->name);
1401:   }
1402:   return(0);
1403: }

1407: /*
1408:      Runs the user provided monitor routines, if they exists.
1409: */
1410: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x)
1411: {
1413:   PetscInt i,n = ts->numbermonitors;

1416:   for (i=0; i<n; i++) {
1417:     (*ts->monitor[i])(ts,step,ptime,x,ts->monitorcontext[i]);
1418:   }
1419:   return(0);
1420: }

1422: /* ------------------------------------------------------------------------*/

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

1430:    Collective on TS

1432:    Input Parameters:
1433: +  host - the X display to open, or null for the local machine
1434: .  label - the title to put in the title bar
1435: .  x, y - the screen coordinates of the upper left coordinate of the window
1436: -  m, n - the screen width and height in pixels

1438:    Output Parameter:
1439: .  draw - the drawing context

1441:    Options Database Key:
1442: .  -ts_xmonitor - automatically sets line graph monitor

1444:    Notes: 
1445:    Use TSLGMonitorDestroy() to destroy this line graph, not PetscDrawLGDestroy().

1447:    Level: intermediate

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

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

1453: @*/
1454: PetscErrorCode  TSLGMonitorCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
1455: {
1456:   PetscDraw      win;

1460:   PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
1461:   PetscDrawSetType(win,PETSC_DRAW_X);
1462:   PetscDrawLGCreate(win,1,draw);
1463:   PetscDrawLGIndicateDataPoints(*draw);

1465:   PetscLogObjectParent(*draw,win);
1466:   return(0);
1467: }

1471: PetscErrorCode TSLGMonitor(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1472: {
1473:   PetscDrawLG    lg = (PetscDrawLG) monctx;
1474:   PetscReal      x,y = ptime;

1478:   if (!monctx) {
1479:     MPI_Comm    comm;
1480:     PetscViewer viewer;

1482:     PetscObjectGetComm((PetscObject)ts,&comm);
1483:     viewer = PETSC_VIEWER_DRAW_(comm);
1484:     PetscViewerDrawGetDrawLG(viewer,0,&lg);
1485:   }

1487:   if (!n) {PetscDrawLGReset(lg);}
1488:   x = (PetscReal)n;
1489:   PetscDrawLGAddPoint(lg,&x,&y);
1490:   if (n < 20 || (n % 5)) {
1491:     PetscDrawLGDraw(lg);
1492:   }
1493:   return(0);
1494: }

1498: /*@C
1499:    TSLGMonitorDestroy - Destroys a line graph context that was created 
1500:    with TSLGMonitorCreate().

1502:    Collective on PetscDrawLG

1504:    Input Parameter:
1505: .  draw - the drawing context

1507:    Level: intermediate

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

1511: .seealso: TSLGMonitorCreate(),  TSSetMonitor(), TSLGMonitor();
1512: @*/
1513: PetscErrorCode  TSLGMonitorDestroy(PetscDrawLG drawlg)
1514: {
1515:   PetscDraw      draw;

1519:   PetscDrawLGGetDraw(drawlg,&draw);
1520:   PetscDrawDestroy(draw);
1521:   PetscDrawLGDestroy(drawlg);
1522:   return(0);
1523: }

1527: /*@
1528:    TSGetTime - Gets the current time.

1530:    Not Collective

1532:    Input Parameter:
1533: .  ts - the TS context obtained from TSCreate()

1535:    Output Parameter:
1536: .  t  - the current time

1538:    Contributed by: Matthew Knepley

1540:    Level: beginner

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

1544: .keywords: TS, get, time
1545: @*/
1546: PetscErrorCode  TSGetTime(TS ts,PetscReal* t)
1547: {
1551:   *t = ts->ptime;
1552:   return(0);
1553: }

1557: /*@C
1558:    TSSetOptionsPrefix - Sets the prefix used for searching for all
1559:    TS options in the database.

1561:    Collective on TS

1563:    Input Parameter:
1564: +  ts     - The TS context
1565: -  prefix - The prefix to prepend to all option names

1567:    Notes:
1568:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1569:    The first character of all runtime options is AUTOMATICALLY the
1570:    hyphen.

1572:    Contributed by: Matthew Knepley

1574:    Level: advanced

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

1578: .seealso: TSSetFromOptions()

1580: @*/
1581: PetscErrorCode  TSSetOptionsPrefix(TS ts,const char prefix[])
1582: {

1587:   PetscObjectSetOptionsPrefix((PetscObject)ts,prefix);
1588:   switch(ts->problem_type) {
1589:     case TS_NONLINEAR:
1590:       SNESSetOptionsPrefix(ts->snes,prefix);
1591:       break;
1592:     case TS_LINEAR:
1593:       KSPSetOptionsPrefix(ts->ksp,prefix);
1594:       break;
1595:   }
1596:   return(0);
1597: }


1602: /*@C
1603:    TSAppendOptionsPrefix - Appends to the prefix used for searching for all
1604:    TS options in the database.

1606:    Collective on TS

1608:    Input Parameter:
1609: +  ts     - The TS context
1610: -  prefix - The prefix to prepend to all option names

1612:    Notes:
1613:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1614:    The first character of all runtime options is AUTOMATICALLY the
1615:    hyphen.

1617:    Contributed by: Matthew Knepley

1619:    Level: advanced

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

1623: .seealso: TSGetOptionsPrefix()

1625: @*/
1626: PetscErrorCode  TSAppendOptionsPrefix(TS ts,const char prefix[])
1627: {

1632:   PetscObjectAppendOptionsPrefix((PetscObject)ts,prefix);
1633:   switch(ts->problem_type) {
1634:     case TS_NONLINEAR:
1635:       SNESAppendOptionsPrefix(ts->snes,prefix);
1636:       break;
1637:     case TS_LINEAR:
1638:       KSPAppendOptionsPrefix(ts->ksp,prefix);
1639:       break;
1640:   }
1641:   return(0);
1642: }

1646: /*@C
1647:    TSGetOptionsPrefix - Sets the prefix used for searching for all
1648:    TS options in the database.

1650:    Not Collective

1652:    Input Parameter:
1653: .  ts - The TS context

1655:    Output Parameter:
1656: .  prefix - A pointer to the prefix string used

1658:    Contributed by: Matthew Knepley

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

1663:    Level: intermediate

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

1667: .seealso: TSAppendOptionsPrefix()
1668: @*/
1669: PetscErrorCode  TSGetOptionsPrefix(TS ts,const char *prefix[])
1670: {

1676:   PetscObjectGetOptionsPrefix((PetscObject)ts,prefix);
1677:   return(0);
1678: }

1682: /*@C
1683:    TSGetRHSMatrix - Returns the matrix A at the present timestep.

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

1687:    Input Parameter:
1688: .  ts  - The TS context obtained from TSCreate()

1690:    Output Parameters:
1691: +  A   - The matrix A, where U_t = A(t) U
1692: .  M   - The preconditioner matrix, usually the same as A
1693: -  ctx - User-defined context for matrix evaluation routine

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

1697:    Contributed by: Matthew Knepley

1699:    Level: intermediate

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

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

1705: @*/
1706: PetscErrorCode  TSGetRHSMatrix(TS ts,Mat *A,Mat *M,void **ctx)
1707: {
1710:   if (A)   *A = ts->A;
1711:   if (M)   *M = ts->B;
1712:   if (ctx) *ctx = ts->jacP;
1713:   return(0);
1714: }

1718: /*@C
1719:    TSGetRHSJacobian - Returns the Jacobian J at the present timestep.

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

1723:    Input Parameter:
1724: .  ts  - The TS context obtained from TSCreate()

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

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

1733:    Contributed by: Matthew Knepley

1735:    Level: intermediate

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

1739: .keywords: TS, timestep, get, matrix, Jacobian
1740: @*/
1741: PetscErrorCode  TSGetRHSJacobian(TS ts,Mat *J,Mat *M,void **ctx)
1742: {

1746:   TSGetRHSMatrix(ts,J,M,ctx);
1747:   return(0);
1748: }

1752: /*@C
1753:    TSVecViewMonitor - Monitors progress of the TS solvers by calling 
1754:    VecView() for the solution at each timestep

1756:    Collective on TS

1758:    Input Parameters:
1759: +  ts - the TS context
1760: .  step - current time-step
1761: .  ptime - current time
1762: -  dummy - either a viewer or PETSC_NULL

1764:    Level: intermediate

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

1768: .seealso: TSSetMonitor(), TSDefaultMonitor(), VecView()
1769: @*/
1770: PetscErrorCode  TSVecViewMonitor(TS ts,PetscInt step,PetscReal ptime,Vec x,void *dummy)
1771: {
1773:   PetscViewer    viewer = (PetscViewer) dummy;

1776:   if (!viewer) {
1777:     MPI_Comm comm;
1778:     PetscObjectGetComm((PetscObject)ts,&comm);
1779:     viewer = PETSC_VIEWER_DRAW_(comm);
1780:   }
1781:   VecView(x,viewer);
1782:   return(0);
1783: }