1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: euler.c,v 1.11 1997/10/19 03:28:22 bsmith Exp bsmith $"; 3 #endif 4 /* 5 Code for Timestepping with explicit Euler. 6 */ 7 #include <math.h> 8 #include "src/ts/tsimpl.h" /*I "ts.h" I*/ 9 #include "pinclude/pviewer.h" 10 11 12 typedef struct { 13 Vec update; /* work vector where F(t[i],u[i]) is stored */ 14 } TS_Euler; 15 16 #undef __FUNC__ 17 #define __FUNC__ "TSSetUp_Euler" 18 static int TSSetUp_Euler(TS ts) 19 { 20 TS_Euler *euler = (TS_Euler*) ts->data; 21 int ierr; 22 23 PetscFunctionBegin; 24 ierr = VecDuplicate(ts->vec_sol,&euler->update); CHKERRQ(ierr); 25 PetscFunctionReturn(0); 26 } 27 28 #undef __FUNC__ 29 #define __FUNC__ "TSStep_Euler" 30 static int TSStep_Euler(TS ts,int *steps,double *time) 31 { 32 TS_Euler *euler = (TS_Euler*) ts->data; 33 Vec sol = ts->vec_sol,update = euler->update; 34 int ierr,i,max_steps = ts->max_steps; 35 Scalar dt = ts->time_step; 36 37 PetscFunctionBegin; 38 *steps = -ts->steps; 39 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol); CHKERRQ(ierr); 40 41 for ( i=0; i<max_steps; i++ ) { 42 ts->ptime += ts->time_step; 43 ierr = TSComputeRHSFunction(ts,ts->ptime,sol,update); CHKERRQ(ierr); 44 ierr = VecAXPY(&dt,update,sol); CHKERRQ(ierr); 45 ts->steps++; 46 ierr = TSMonitor(ts,ts->steps,ts->ptime,sol); CHKERRQ(ierr); 47 if (ts->ptime > ts->max_time) break; 48 } 49 50 *steps += ts->steps; 51 *time = ts->ptime; 52 PetscFunctionReturn(0); 53 } 54 /*------------------------------------------------------------*/ 55 #undef __FUNC__ 56 #define __FUNC__ "TSDestroy_Euler" 57 static int TSDestroy_Euler(PetscObject obj ) 58 { 59 TS ts = (TS) obj; 60 TS_Euler *euler = (TS_Euler*) ts->data; 61 62 PetscFunctionBegin; 63 VecDestroy(euler->update); 64 PetscFree(euler); 65 PetscFunctionReturn(0); 66 } 67 /*------------------------------------------------------------*/ 68 69 #undef __FUNC__ 70 #define __FUNC__ "TSSetFromOptions_Euler" 71 static int TSSetFromOptions_Euler(TS ts) 72 { 73 PetscFunctionBegin; 74 PetscFunctionReturn(0); 75 } 76 77 #undef __FUNC__ 78 #define __FUNC__ "TSPrintHelp_Euler" 79 static int TSPrintHelp_Euler(TS ts,char *p) 80 { 81 PetscFunctionBegin; 82 PetscFunctionReturn(0); 83 } 84 85 #undef __FUNC__ 86 #define __FUNC__ "TSView_Euler" 87 static int TSView_Euler(PetscObject obj,Viewer viewer) 88 { 89 PetscFunctionBegin; 90 PetscFunctionReturn(0); 91 } 92 93 /* ------------------------------------------------------------ */ 94 #undef __FUNC__ 95 #define __FUNC__ "TSCreate_Euler" 96 int TSCreate_Euler(TS ts ) 97 { 98 TS_Euler *euler; 99 100 PetscFunctionBegin; 101 ts->setup = TSSetUp_Euler; 102 ts->step = TSStep_Euler; 103 ts->destroy = TSDestroy_Euler; 104 ts->printhelp = TSPrintHelp_Euler; 105 ts->setfromoptions = TSSetFromOptions_Euler; 106 ts->view = TSView_Euler; 107 108 euler = PetscNew(TS_Euler); CHKPTRQ(euler); 109 PLogObjectMemory(ts,sizeof(TS_Euler)); 110 ts->data = (void *) euler; 111 112 PetscFunctionReturn(0); 113 } 114 115 116 117 118 119