Actual source code: euler.c

  1: /*$Id: euler.c,v 1.30 2001/08/07 03:04:22 balay Exp $*/
  2: /*
  3:        Code for Timestepping with explicit Euler.
  4: */
 5:  #include src/ts/tsimpl.h

  7: typedef struct {
  8:   Vec update;     /* work vector where F(t[i],u[i]) is stored */
  9: } TS_Euler;

 11: static int TSSetUp_Euler(TS ts)
 12: {
 13:   TS_Euler *euler = (TS_Euler*)ts->data;
 14:   int      ierr;

 17:   VecDuplicate(ts->vec_sol,&euler->update);
 18:   return(0);
 19: }

 21: static int TSStep_Euler(TS ts,int *steps,PetscReal *ptime)
 22: {
 23:   TS_Euler *euler = (TS_Euler*)ts->data;
 24:   Vec      sol = ts->vec_sol,update = euler->update;
 25:   int      ierr,i,max_steps = ts->max_steps;
 26:   PetscScalar   dt = ts->time_step;
 27: 
 29:   *steps = -ts->steps;
 30:   TSMonitor(ts,ts->steps,ts->ptime,sol);

 32:   for (i=0; i<max_steps; i++) {
 33:     ts->ptime += ts->time_step;
 34:     TSComputeRHSFunction(ts,ts->ptime,sol,update);
 35:     VecAXPY(&dt,update,sol);
 36:     ts->steps++;
 37:     TSMonitor(ts,ts->steps,ts->ptime,sol);
 38:     if (ts->ptime > ts->max_time) break;
 39:   }

 41:   *steps += ts->steps;
 42:   *ptime  = ts->ptime;
 43:   return(0);
 44: }
 45: /*------------------------------------------------------------*/
 46: static int TSDestroy_Euler(TS ts)
 47: {
 48:   TS_Euler *euler = (TS_Euler*)ts->data;
 49:   int      ierr;

 52:   if (euler->update) {VecDestroy(euler->update);}
 53:   PetscFree(euler);
 54:   return(0);
 55: }
 56: /*------------------------------------------------------------*/

 58: static int TSSetFromOptions_Euler(TS ts)
 59: {
 61:   return(0);
 62: }

 64: static int TSView_Euler(TS ts,PetscViewer viewer)
 65: {
 67:   return(0);
 68: }

 70: /* ------------------------------------------------------------ */
 71: EXTERN_C_BEGIN
 72: int TSCreate_Euler(TS ts)
 73: {
 74:   TS_Euler *euler;
 75:   int      ierr;

 78:   ts->ops->setup           = TSSetUp_Euler;
 79:   ts->ops->step            = TSStep_Euler;
 80:   ts->ops->destroy         = TSDestroy_Euler;
 81:   ts->ops->setfromoptions  = TSSetFromOptions_Euler;
 82:   ts->ops->view            = TSView_Euler;

 84:   PetscNew(TS_Euler,&euler);
 85:   PetscLogObjectMemory(ts,sizeof(TS_Euler));
 86:   ts->data = (void*)euler;

 88:   return(0);
 89: }
 90: EXTERN_C_END