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;

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

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

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

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

 45:   *steps += ts->steps;
 46:   *ptime  = ts->ptime;
 47:   return(0);
 48: }
 49: /*------------------------------------------------------------*/
 52: static int TSDestroy_Euler(TS ts)
 53: {
 54:   TS_Euler *euler = (TS_Euler*)ts->data;
 55:   int      ierr;

 58:   if (euler->update) {VecDestroy(euler->update);}
 59:   PetscFree(euler);
 60:   return(0);
 61: }
 62: /*------------------------------------------------------------*/

 66: static int TSSetFromOptions_Euler(TS ts)
 67: {
 69:   return(0);
 70: }

 74: static int TSView_Euler(TS ts,PetscViewer viewer)
 75: {
 77:   return(0);
 78: }

 80: /* ------------------------------------------------------------ */
 81: EXTERN_C_BEGIN
 84: int TSCreate_Euler(TS ts)
 85: {
 86:   TS_Euler *euler;
 87:   int      ierr;

 90:   ts->ops->setup           = TSSetUp_Euler;
 91:   ts->ops->step            = TSStep_Euler;
 92:   ts->ops->destroy         = TSDestroy_Euler;
 93:   ts->ops->setfromoptions  = TSSetFromOptions_Euler;
 94:   ts->ops->view            = TSView_Euler;

 96:   PetscNew(TS_Euler,&euler);
 97:   PetscLogObjectMemory(ts,sizeof(TS_Euler));
 98:   ts->data = (void*)euler;

100:   return(0);
101: }
102: EXTERN_C_END