1 /* 2 Code for Timestepping with explicit Euler. 3 */ 4 #include "src/ts/tsimpl.h" /*I "petscts.h" I*/ 5 6 typedef struct { 7 Vec update; /* work vector where F(t[i],u[i]) is stored */ 8 } TS_Euler; 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "TSSetUp_Euler" 12 static int TSSetUp_Euler(TS ts) 13 { 14 TS_Euler *euler = (TS_Euler*)ts->data; 15 int ierr; 16 17 PetscFunctionBegin; 18 ierr = VecDuplicate(ts->vec_sol,&euler->update);CHKERRQ(ierr); 19 PetscFunctionReturn(0); 20 } 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "TSStep_Euler" 24 static int TSStep_Euler(TS ts,int *steps,PetscReal *ptime) 25 { 26 TS_Euler *euler = (TS_Euler*)ts->data; 27 Vec sol = ts->vec_sol,update = euler->update; 28 int ierr,i,max_steps = ts->max_steps; 29 PetscScalar dt = ts->time_step; 30 31 PetscFunctionBegin; 32 *steps = -ts->steps; 33 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 34 35 for (i=0; i<max_steps; i++) { 36 ts->ptime += ts->time_step; 37 ierr = TSComputeRHSFunction(ts,ts->ptime,sol,update);CHKERRQ(ierr); 38 ierr = VecAXPY(&dt,update,sol);CHKERRQ(ierr); 39 ts->steps++; 40 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 41 if (ts->ptime > ts->max_time) break; 42 } 43 44 *steps += ts->steps; 45 *ptime = ts->ptime; 46 PetscFunctionReturn(0); 47 } 48 /*------------------------------------------------------------*/ 49 #undef __FUNCT__ 50 #define __FUNCT__ "TSDestroy_Euler" 51 static int TSDestroy_Euler(TS ts) 52 { 53 TS_Euler *euler = (TS_Euler*)ts->data; 54 int ierr; 55 56 PetscFunctionBegin; 57 if (euler->update) {ierr = VecDestroy(euler->update);CHKERRQ(ierr);} 58 ierr = PetscFree(euler);CHKERRQ(ierr); 59 PetscFunctionReturn(0); 60 } 61 /*------------------------------------------------------------*/ 62 63 #undef __FUNCT__ 64 #define __FUNCT__ "TSSetFromOptions_Euler" 65 static int TSSetFromOptions_Euler(TS ts) 66 { 67 PetscFunctionBegin; 68 PetscFunctionReturn(0); 69 } 70 71 #undef __FUNCT__ 72 #define __FUNCT__ "TSView_Euler" 73 static int TSView_Euler(TS ts,PetscViewer viewer) 74 { 75 PetscFunctionBegin; 76 PetscFunctionReturn(0); 77 } 78 79 /* ------------------------------------------------------------ */ 80 81 /*MC 82 TS_EULER - ODE solver using the explicit forward Euler method 83 84 Level: beginner 85 86 .seealso: TSCreate(), TS, TSSetType(), TS_BEULER 87 88 M*/ 89 EXTERN_C_BEGIN 90 #undef __FUNCT__ 91 #define __FUNCT__ "TSCreate_Euler" 92 int TSCreate_Euler(TS ts) 93 { 94 TS_Euler *euler; 95 int ierr; 96 97 PetscFunctionBegin; 98 ts->ops->setup = TSSetUp_Euler; 99 ts->ops->step = TSStep_Euler; 100 ts->ops->destroy = TSDestroy_Euler; 101 ts->ops->setfromoptions = TSSetFromOptions_Euler; 102 ts->ops->view = TSView_Euler; 103 104 ierr = PetscNew(TS_Euler,&euler);CHKERRQ(ierr); 105 PetscLogObjectMemory(ts,sizeof(TS_Euler)); 106 ts->data = (void*)euler; 107 108 PetscFunctionReturn(0); 109 } 110 EXTERN_C_END 111 112 113 114 115