Actual source code: euler.c

petsc-3.4.5 2014-06-29
  1: /*
  2:        Code for Timestepping with explicit Euler.
  3: */
  4: #include <petsc-private/tsimpl.h>                /*I   "petscts.h"   I*/

  6: typedef struct {
  7:   Vec update;     /* work vector where new solution is formed  */
  8: } TS_Euler;

 12: static PetscErrorCode TSStep_Euler(TS ts)
 13: {
 14:   TS_Euler       *euler = (TS_Euler*)ts->data;
 15:   Vec            sol    = ts->vec_sol,update = euler->update;

 19:   TSPreStep(ts);
 20:   TSPreStage(ts,ts->ptime);
 21:   TSComputeRHSFunction(ts,ts->ptime,sol,update);
 22:   VecAXPY(sol,ts->time_step,update);
 23:   ts->ptime += ts->time_step;
 24:   ts->steps++;
 25:   return(0);
 26: }
 27: /*------------------------------------------------------------*/

 31: static PetscErrorCode TSSetUp_Euler(TS ts)
 32: {
 33:   TS_Euler       *euler = (TS_Euler*)ts->data;

 37:   VecDuplicate(ts->vec_sol,&euler->update);
 38:   return(0);
 39: }

 43: static PetscErrorCode TSReset_Euler(TS ts)
 44: {
 45:   TS_Euler       *euler = (TS_Euler*)ts->data;

 49:   VecDestroy(&euler->update);
 50:   return(0);
 51: }

 55: static PetscErrorCode TSDestroy_Euler(TS ts)
 56: {

 60:   TSReset_Euler(ts);
 61:   PetscFree(ts->data);
 62:   return(0);
 63: }
 64: /*------------------------------------------------------------*/

 68: static PetscErrorCode TSSetFromOptions_Euler(TS ts)
 69: {
 71:   return(0);
 72: }

 76: static PetscErrorCode TSView_Euler(TS ts,PetscViewer viewer)
 77: {
 79:   return(0);
 80: }

 84: static PetscErrorCode TSInterpolate_Euler(TS ts,PetscReal t,Vec X)
 85: {
 86:   PetscReal      alpha = (ts->ptime - t)/ts->time_step;

 90:   VecAXPBY(ts->vec_sol,1.0-alpha,alpha,X);
 91:   return(0);
 92: }

 96: PetscErrorCode TSComputeLinearStability_Euler(TS ts,PetscReal xr,PetscReal xi,PetscReal *yr,PetscReal *yi)
 97: {
 99:   *yr = 1.0 + xr;
100:   *yi = xi;
101:   return(0);
102: }
103: /* ------------------------------------------------------------ */

105: /*MC
106:       TSEULER - ODE solver using the explicit forward Euler method

108:   Level: beginner

110: .seealso:  TSCreate(), TS, TSSetType(), TSBEULER

112: M*/
115: PETSC_EXTERN PetscErrorCode TSCreate_Euler(TS ts)
116: {
117:   TS_Euler       *euler;

121:   ts->ops->setup           = TSSetUp_Euler;
122:   ts->ops->step            = TSStep_Euler;
123:   ts->ops->reset           = TSReset_Euler;
124:   ts->ops->destroy         = TSDestroy_Euler;
125:   ts->ops->setfromoptions  = TSSetFromOptions_Euler;
126:   ts->ops->view            = TSView_Euler;
127:   ts->ops->interpolate     = TSInterpolate_Euler;
128:   ts->ops->linearstability = TSComputeLinearStability_Euler;

130:   PetscNewLog(ts,TS_Euler,&euler);
131:   ts->data = (void*)euler;
132:   return(0);
133: }