1 /* 2 Code for Timestepping with explicit Euler. 3 */ 4 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 5 6 typedef struct { 7 Vec update; /* work vector where new solution is formed */ 8 } TS_Euler; 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "TSStep_Euler" 12 static PetscErrorCode TSStep_Euler(TS ts) 13 { 14 TS_Euler *euler = (TS_Euler*)ts->data; 15 Vec sol = ts->vec_sol,update = euler->update; 16 PetscErrorCode ierr; 17 18 PetscFunctionBegin; 19 ierr = TSPreStep(ts);CHKERRQ(ierr); 20 ierr = TSPreStage(ts,ts->ptime);CHKERRQ(ierr); 21 ierr = TSComputeRHSFunction(ts,ts->ptime,sol,update);CHKERRQ(ierr); 22 ierr = VecAXPY(sol,ts->time_step,update);CHKERRQ(ierr); 23 ierr = TSPostStage(ts,ts->ptime,0,&sol);CHKERRQ(ierr); 24 ts->ptime += ts->time_step; 25 ts->steps++; 26 PetscFunctionReturn(0); 27 } 28 /*------------------------------------------------------------*/ 29 30 #undef __FUNCT__ 31 #define __FUNCT__ "TSSetUp_Euler" 32 static PetscErrorCode TSSetUp_Euler(TS ts) 33 { 34 TS_Euler *euler = (TS_Euler*)ts->data; 35 PetscErrorCode ierr; 36 37 PetscFunctionBegin; 38 ierr = VecDuplicate(ts->vec_sol,&euler->update);CHKERRQ(ierr); 39 PetscFunctionReturn(0); 40 } 41 42 #undef __FUNCT__ 43 #define __FUNCT__ "TSReset_Euler" 44 static PetscErrorCode TSReset_Euler(TS ts) 45 { 46 TS_Euler *euler = (TS_Euler*)ts->data; 47 PetscErrorCode ierr; 48 49 PetscFunctionBegin; 50 ierr = VecDestroy(&euler->update);CHKERRQ(ierr); 51 PetscFunctionReturn(0); 52 } 53 54 #undef __FUNCT__ 55 #define __FUNCT__ "TSDestroy_Euler" 56 static PetscErrorCode TSDestroy_Euler(TS ts) 57 { 58 PetscErrorCode ierr; 59 60 PetscFunctionBegin; 61 ierr = TSReset_Euler(ts);CHKERRQ(ierr); 62 ierr = PetscFree(ts->data);CHKERRQ(ierr); 63 PetscFunctionReturn(0); 64 } 65 /*------------------------------------------------------------*/ 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "TSSetFromOptions_Euler" 69 static PetscErrorCode TSSetFromOptions_Euler(PetscOptions *PetscOptionsObject,TS ts) 70 { 71 PetscFunctionBegin; 72 PetscFunctionReturn(0); 73 } 74 75 #undef __FUNCT__ 76 #define __FUNCT__ "TSView_Euler" 77 static PetscErrorCode TSView_Euler(TS ts,PetscViewer viewer) 78 { 79 PetscFunctionBegin; 80 PetscFunctionReturn(0); 81 } 82 83 #undef __FUNCT__ 84 #define __FUNCT__ "TSInterpolate_Euler" 85 static PetscErrorCode TSInterpolate_Euler(TS ts,PetscReal t,Vec X) 86 { 87 TS_Euler *euler = (TS_Euler*)ts->data; 88 Vec update = euler->update; 89 PetscReal alpha = (ts->ptime - t)/ts->time_step; 90 PetscErrorCode ierr; 91 92 PetscFunctionBegin; 93 ierr = VecWAXPY(X,-ts->time_step,update,ts->vec_sol);CHKERRQ(ierr); 94 ierr = VecAXPBY(X,1.0-alpha,alpha,ts->vec_sol);CHKERRQ(ierr); 95 PetscFunctionReturn(0); 96 } 97 98 #undef __FUNCT__ 99 #define __FUNCT__ "TSComputeLinearStability_Euler" 100 PetscErrorCode TSComputeLinearStability_Euler(TS ts,PetscReal xr,PetscReal xi,PetscReal *yr,PetscReal *yi) 101 { 102 PetscFunctionBegin; 103 *yr = 1.0 + xr; 104 *yi = xi; 105 PetscFunctionReturn(0); 106 } 107 /* ------------------------------------------------------------ */ 108 109 /*MC 110 TSEULER - ODE solver using the explicit forward Euler method 111 112 Level: beginner 113 114 .seealso: TSCreate(), TS, TSSetType(), TSBEULER 115 116 M*/ 117 #undef __FUNCT__ 118 #define __FUNCT__ "TSCreate_Euler" 119 PETSC_EXTERN PetscErrorCode TSCreate_Euler(TS ts) 120 { 121 TS_Euler *euler; 122 PetscErrorCode ierr; 123 124 PetscFunctionBegin; 125 ts->ops->setup = TSSetUp_Euler; 126 ts->ops->step = TSStep_Euler; 127 ts->ops->reset = TSReset_Euler; 128 ts->ops->destroy = TSDestroy_Euler; 129 ts->ops->setfromoptions = TSSetFromOptions_Euler; 130 ts->ops->view = TSView_Euler; 131 ts->ops->interpolate = TSInterpolate_Euler; 132 ts->ops->linearstability = TSComputeLinearStability_Euler; 133 134 ierr = PetscNewLog(ts,&euler);CHKERRQ(ierr); 135 ts->data = (void*)euler; 136 PetscFunctionReturn(0); 137 } 138