Actual source code: euler.c
petsc-3.5.4 2015-05-23
1: /*
2: Code for Timestepping with explicit Euler.
3: */
4: #include <petsc-private/tsimpl.h> /*I "petscts.h" I*/
6: typedef struct {
7: Vec update; /* work vector where new solution is formed */
8: } TS_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;
19: TSPreStep(ts);
20: TSPreStage(ts,ts->ptime);
21: TSComputeRHSFunction(ts,ts->ptime,sol,update);
22: VecAXPY(sol,ts->time_step,update);
23: TSPostStage(ts,ts->ptime,0,&sol);
24: ts->ptime += ts->time_step;
25: ts->steps++;
26: return(0);
27: }
28: /*------------------------------------------------------------*/
32: static PetscErrorCode TSSetUp_Euler(TS ts)
33: {
34: TS_Euler *euler = (TS_Euler*)ts->data;
38: VecDuplicate(ts->vec_sol,&euler->update);
39: return(0);
40: }
44: static PetscErrorCode TSReset_Euler(TS ts)
45: {
46: TS_Euler *euler = (TS_Euler*)ts->data;
50: VecDestroy(&euler->update);
51: return(0);
52: }
56: static PetscErrorCode TSDestroy_Euler(TS ts)
57: {
61: TSReset_Euler(ts);
62: PetscFree(ts->data);
63: return(0);
64: }
65: /*------------------------------------------------------------*/
69: static PetscErrorCode TSSetFromOptions_Euler(TS ts)
70: {
72: return(0);
73: }
77: static PetscErrorCode TSView_Euler(TS ts,PetscViewer viewer)
78: {
80: return(0);
81: }
85: static PetscErrorCode TSInterpolate_Euler(TS ts,PetscReal t,Vec X)
86: {
87: PetscReal alpha = (ts->ptime - t)/ts->time_step;
91: VecAXPBY(ts->vec_sol,1.0-alpha,alpha,X);
92: return(0);
93: }
97: PetscErrorCode TSComputeLinearStability_Euler(TS ts,PetscReal xr,PetscReal xi,PetscReal *yr,PetscReal *yi)
98: {
100: *yr = 1.0 + xr;
101: *yi = xi;
102: return(0);
103: }
104: /* ------------------------------------------------------------ */
106: /*MC
107: TSEULER - ODE solver using the explicit forward Euler method
109: Level: beginner
111: .seealso: TSCreate(), TS, TSSetType(), TSBEULER
113: M*/
116: PETSC_EXTERN PetscErrorCode TSCreate_Euler(TS ts)
117: {
118: TS_Euler *euler;
122: ts->ops->setup = TSSetUp_Euler;
123: ts->ops->step = TSStep_Euler;
124: ts->ops->reset = TSReset_Euler;
125: ts->ops->destroy = TSDestroy_Euler;
126: ts->ops->setfromoptions = TSSetFromOptions_Euler;
127: ts->ops->view = TSView_Euler;
128: ts->ops->interpolate = TSInterpolate_Euler;
129: ts->ops->linearstability = TSComputeLinearStability_Euler;
131: PetscNewLog(ts,&euler);
132: ts->data = (void*)euler;
133: return(0);
134: }