Actual source code: euler.c

      1: #define PETSCTS_DLL

      3: /*
      4:        Code for Timestepping with explicit Euler.
      5: */
     6:  #include src/ts/tsimpl.h

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

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

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

     26: static PetscErrorCode TSStep_Euler(TS ts,PetscInt *steps,PetscReal *ptime)
     27: {
     28:   TS_Euler       *euler = (TS_Euler*)ts->data;
     29:   Vec            sol = ts->vec_sol,update = euler->update;
     31:   PetscInt       i,max_steps = ts->max_steps;
     32:   PetscScalar    dt = ts->time_step;
     33:
     35:   *steps = -ts->steps;
     36:   TSMonitor(ts,ts->steps,ts->ptime,sol);

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

     47:   *steps += ts->steps;
     48:   *ptime  = ts->ptime;
     49:   return(0);
     50: }
     51: /*------------------------------------------------------------*/
     54: static PetscErrorCode TSDestroy_Euler(TS ts)
     55: {
     56:   TS_Euler       *euler = (TS_Euler*)ts->data;

     60:   if (euler->update) {VecDestroy(euler->update);}
     61:   PetscFree(euler);
     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: }

     82: /* ------------------------------------------------------------ */

     84: /*MC
     85:       TS_EULER - ODE solver using the explicit forward Euler method

     87:   Level: beginner

     89: .seealso:  TSCreate(), TS, TSSetType(), TS_BEULER

     91: M*/
     95: PetscErrorCode PETSCTS_DLLEXPORT TSCreate_Euler(TS ts)
     96: {
     97:   TS_Euler       *euler;

    101:   ts->ops->setup           = TSSetUp_Euler;
    102:   ts->ops->step            = TSStep_Euler;
    103:   ts->ops->destroy         = TSDestroy_Euler;
    104:   ts->ops->setfromoptions  = TSSetFromOptions_Euler;
    105:   ts->ops->view            = TSView_Euler;

    107:   PetscNew(TS_Euler,&euler);
    108:   ts->data = (void*)euler;

    110:   return(0);
    111: }