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 PetscErrorCode TSSetUp_Euler(TS ts) 13 { 14 TS_Euler *euler = (TS_Euler*)ts->data; 15 PetscErrorCode 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 PetscErrorCode TSStep_Euler(TS ts,PetscInt *steps,PetscReal *ptime) 25 { 26 TS_Euler *euler = (TS_Euler*)ts->data; 27 Vec sol = ts->vec_sol,update = euler->update; 28 PetscErrorCode ierr; 29 PetscInt i,max_steps = ts->max_steps; 30 PetscScalar dt = ts->time_step; 31 32 PetscFunctionBegin; 33 *steps = -ts->steps; 34 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 35 36 for (i=0; i<max_steps; i++) { 37 ts->ptime += ts->time_step; 38 ierr = TSComputeRHSFunction(ts,ts->ptime,sol,update);CHKERRQ(ierr); 39 ierr = VecAXPY(&dt,update,sol);CHKERRQ(ierr); 40 ts->steps++; 41 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 42 if (ts->ptime > ts->max_time) break; 43 } 44 45 *steps += ts->steps; 46 *ptime = ts->ptime; 47 PetscFunctionReturn(0); 48 } 49 /*------------------------------------------------------------*/ 50 #undef __FUNCT__ 51 #define __FUNCT__ "TSDestroy_Euler" 52 static PetscErrorCode TSDestroy_Euler(TS ts) 53 { 54 TS_Euler *euler = (TS_Euler*)ts->data; 55 PetscErrorCode ierr; 56 57 PetscFunctionBegin; 58 if (euler->update) {ierr = VecDestroy(euler->update);CHKERRQ(ierr);} 59 ierr = PetscFree(euler);CHKERRQ(ierr); 60 PetscFunctionReturn(0); 61 } 62 /*------------------------------------------------------------*/ 63 64 #undef __FUNCT__ 65 #define __FUNCT__ "TSSetFromOptions_Euler" 66 static PetscErrorCode TSSetFromOptions_Euler(TS ts) 67 { 68 PetscFunctionBegin; 69 PetscFunctionReturn(0); 70 } 71 72 #undef __FUNCT__ 73 #define __FUNCT__ "TSView_Euler" 74 static PetscErrorCode TSView_Euler(TS ts,PetscViewer viewer) 75 { 76 PetscFunctionBegin; 77 PetscFunctionReturn(0); 78 } 79 80 /* ------------------------------------------------------------ */ 81 82 /*MC 83 TS_EULER - ODE solver using the explicit forward Euler method 84 85 Level: beginner 86 87 .seealso: TSCreate(), TS, TSSetType(), TS_BEULER 88 89 M*/ 90 EXTERN_C_BEGIN 91 #undef __FUNCT__ 92 #define __FUNCT__ "TSCreate_Euler" 93 PetscErrorCode TSCreate_Euler(TS ts) 94 { 95 TS_Euler *euler; 96 PetscErrorCode ierr; 97 98 PetscFunctionBegin; 99 ts->ops->setup = TSSetUp_Euler; 100 ts->ops->step = TSStep_Euler; 101 ts->ops->destroy = TSDestroy_Euler; 102 ts->ops->setfromoptions = TSSetFromOptions_Euler; 103 ts->ops->view = TSView_Euler; 104 105 ierr = PetscNew(TS_Euler,&euler);CHKERRQ(ierr); 106 PetscLogObjectMemory(ts,sizeof(TS_Euler)); 107 ts->data = (void*)euler; 108 109 PetscFunctionReturn(0); 110 } 111 EXTERN_C_END 112 113 114 115 116