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