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 new solution is formed */ 9 } TS_Euler; 10 11 #undef __FUNCT__ 12 #define __FUNCT__ "TSStep_Euler" 13 static PetscErrorCode TSStep_Euler(TS ts,PetscInt *steps,PetscReal *ptime) 14 { 15 TS_Euler *euler = (TS_Euler*)ts->data; 16 Vec sol = ts->vec_sol,update = euler->update; 17 PetscInt i,max_steps = ts->max_steps; 18 PetscErrorCode ierr; 19 20 PetscFunctionBegin; 21 *steps = -ts->steps; 22 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 23 24 for (i=0; i<max_steps; i++) { 25 PetscReal dt = ts->time_step; 26 27 ierr = TSPreStep(ts);CHKERRQ(ierr); 28 ts->ptime += dt; 29 ierr = TSComputeRHSFunction(ts,ts->ptime,sol,update);CHKERRQ(ierr); 30 ierr = VecAXPY(sol,dt,update);CHKERRQ(ierr); 31 ts->steps++; 32 ierr = TSPostStep(ts);CHKERRQ(ierr); 33 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol);CHKERRQ(ierr); 34 if (ts->ptime > ts->max_time) break; 35 } 36 37 *steps += ts->steps; 38 *ptime = ts->ptime; 39 PetscFunctionReturn(0); 40 } 41 /*------------------------------------------------------------*/ 42 43 #undef __FUNCT__ 44 #define __FUNCT__ "TSSetUp_Euler" 45 static PetscErrorCode TSSetUp_Euler(TS ts) 46 { 47 TS_Euler *euler = (TS_Euler*)ts->data; 48 PetscErrorCode ierr; 49 50 PetscFunctionBegin; 51 ierr = VecDuplicate(ts->vec_sol,&euler->update);CHKERRQ(ierr); 52 PetscFunctionReturn(0); 53 } 54 55 #undef __FUNCT__ 56 #define __FUNCT__ "TSReset_Euler" 57 static PetscErrorCode TSReset_Euler(TS ts) 58 { 59 TS_Euler *euler = (TS_Euler*)ts->data; 60 PetscErrorCode ierr; 61 62 PetscFunctionBegin; 63 if (euler->update) {ierr = VecDestroy(euler->update);CHKERRQ(ierr);} 64 PetscFunctionReturn(0); 65 } 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "TSDestroy_Euler" 69 static PetscErrorCode TSDestroy_Euler(TS ts) 70 { 71 PetscErrorCode ierr; 72 73 PetscFunctionBegin; 74 ierr = TSReset_Euler(ts);CHKERRQ(ierr); 75 ierr = PetscFree(ts->data);CHKERRQ(ierr); 76 PetscFunctionReturn(0); 77 } 78 /*------------------------------------------------------------*/ 79 80 #undef __FUNCT__ 81 #define __FUNCT__ "TSSetFromOptions_Euler" 82 static PetscErrorCode TSSetFromOptions_Euler(TS ts) 83 { 84 PetscFunctionBegin; 85 PetscFunctionReturn(0); 86 } 87 88 #undef __FUNCT__ 89 #define __FUNCT__ "TSView_Euler" 90 static PetscErrorCode TSView_Euler(TS ts,PetscViewer viewer) 91 { 92 PetscFunctionBegin; 93 PetscFunctionReturn(0); 94 } 95 96 /* ------------------------------------------------------------ */ 97 98 /*MC 99 TSEULER - ODE solver using the explicit forward Euler method 100 101 Level: beginner 102 103 .seealso: TSCreate(), TS, TSSetType(), TSBEULER 104 105 M*/ 106 EXTERN_C_BEGIN 107 #undef __FUNCT__ 108 #define __FUNCT__ "TSCreate_Euler" 109 PetscErrorCode TSCreate_Euler(TS ts) 110 { 111 TS_Euler *euler; 112 PetscErrorCode ierr; 113 114 PetscFunctionBegin; 115 ts->ops->setup = TSSetUp_Euler; 116 ts->ops->step = TSStep_Euler; 117 ts->ops->reset = TSReset_Euler; 118 ts->ops->destroy = TSDestroy_Euler; 119 ts->ops->setfromoptions = TSSetFromOptions_Euler; 120 ts->ops->view = TSView_Euler; 121 122 ierr = PetscNewLog(ts,TS_Euler,&euler);CHKERRQ(ierr); 123 ts->data = (void*)euler; 124 125 PetscFunctionReturn(0); 126 } 127 EXTERN_C_END 128 129 130 131 132