48 const double b = 0.5/a;
52 add(x, (1. - b)*dt, dxdt, x1);
57 add(x1, b*dt, dxdt, x);
83 add(3./4, x, 1./4, y, y);
89 add(1./3, x, 2./3, y, x);
147 for (
int i = 0; i < s; i++)
163 for (
int l = 0, i = 1; i < s; i++)
165 add(x, a[l++]*dt, k[0], y);
166 for (
int j = 1; j < i; j++)
167 y.
Add(a[l++]*dt, k[j]);
172 for (
int i = 0; i < s; i++)
173 x.
Add(b[i]*dt, k[i]);
182 const double RK6Solver::a[] = {
184 .1923996296296296296296296296296296296296e-1,
185 .7669337037037037037037037037037037037037e-1,
189 1.318683415233148260919747276431735612861,
191 -5.042058063628562225427761634715637693344,
192 4.220674648395413964508014358283902080483,
193 -41.87259166432751461803757780644346812905,
195 159.4325621631374917700365669070346830453,
196 -122.1192135650100309202516203389242140663,
197 5.531743066200053768252631238332999150076,
198 -54.43015693531650433250642051294142461271,
200 207.0672513650184644273657173866509835987,
201 -158.6108137845899991828742424365058599469,
202 6.991816585950242321992597280791793907096,
203 -.1859723106220323397765171799549294623692e-1,
204 -54.66374178728197680241215648050386959351,
206 207.9528062553893734515824816699834244238,
207 -159.2889574744995071508959805871426654216,
208 7.018743740796944434698170760964252490817,
209 -.1833878590504572306472782005141738268361e-1,
210 -.5119484997882099077875432497245168395840e-3
212 const double RK6Solver::b[] = {
213 .3438957868357036009278820124728322386520e-1,
216 .2582624555633503404659558098586120858767,
217 .4209371189673537150642551514069801967032,
218 4.405396469669310170148836816197095664891,
219 -176.4831190242986576151740942499002125029,
220 172.3641334014150730294022582711902413315
222 const double RK6Solver::c[] = {
224 .9593333333333333333333333333333333333333e-1,
232 const double RK8Solver::a[] = {
239 .3613975628004575124052940721184028345129,
241 -1.341524066700492771819987788202715834917,
242 1.370126503900035259414693716084313000404,
243 .490472027972027972027972027972027972028e-1,
246 .2350972042214404739862988335493427143122,
247 .180855592981356728810903963653454488485,
248 .6169289044289044289044289044289044289044e-1,
251 .1123656831464027662262557035130015442303,
252 -.3885046071451366767049048108111244567456e-1,
253 .1979188712522045855379188712522045855379e-1,
254 -1.767630240222326875735597119572145586714,
258 -6.061889377376669100821361459659331999758,
259 5.650823198222763138561298030600840174201,
260 65.62169641937623283799566054863063741227,
261 -1.180945066554970799825116282628297957882,
264 -41.50473441114320841606641502701994225874,
265 -4.434438319103725011225169229846100211776,
266 4.260408188586133024812193710744693240761,
267 43.75364022446171584987676829438379303004,
268 .787142548991231068744647504422630755086e-2,
269 -1.281405999441488405459510291182054246266,
272 -45.04713996013986630220754257136007322267,
273 -4.731362069449576477311464265491282810943,
274 4.514967016593807841185851584597240996214,
275 47.44909557172985134869022392235929015114,
276 .1059228297111661135687393955516542875228e-1,
277 -.5746842263844616254432318478286296232021e-2,
278 -1.724470134262485191756709817484481861731,
281 -60.92349008483054016518434619253765246063,
282 -5.95151837622239245520283276706185486829,
283 5.556523730698456235979791650843592496839,
284 63.98301198033305336837536378635995939281,
285 .1464202825041496159275921391759452676003e-1,
286 .6460408772358203603621865144977650714892e-1,
287 -.7930323169008878984024452548693373291447e-1,
288 -3.301622667747079016353994789790983625569,
291 -118.011272359752508566692330395789886851,
292 -10.14142238845611248642783916034510897595,
293 9.139311332232057923544012273556827000619,
294 123.3759428284042683684847180986501894364,
295 4.623244378874580474839807625067630924792,
296 -3.383277738068201923652550971536811240814,
297 4.527592100324618189451265339351129035325,
298 -5.828495485811622963193088019162985703755
300 const double RK8Solver::b[] = {
301 .4427989419007951074716746668098518862111e-1,
306 .3541049391724448744815552028733568354121,
307 .2479692154956437828667629415370663023884,
308 -15.69420203883808405099207034271191213468,
309 25.08406496555856261343930031237186278518,
310 -31.73836778626027646833156112007297739997,
311 22.93828327398878395231483560344797018313,
312 -.2361324633071542145259900641263517600737
314 const double RK8Solver::c[] = {
322 .901802041735856958259707940678372149956,
362 gamma = (3. - sqrt(3.))/6.;
363 else if (gamma_opt == 2)
364 gamma = (2. - sqrt(2.))/2.;
365 else if (gamma_opt == 3)
366 gamma = (2. + sqrt(2.))/2.;
368 gamma = (3. + sqrt(3.))/6.;
414 const double a = 1./sqrt(3.)*cos(M_PI/18.) + 0.5;
415 const double b = 1./(6.*(2.*a-1.)*(2.*a-1.));
419 add(x, (0.5-a)*dt,
k,
y);
420 add(x, (2.*a)*dt,
k,
z);
425 z.
Add((1.-4.*a)*dt,
k);
426 x.
Add((1.-2.*b)*dt,
k);
449 const double a = 0.435866521508458999416019;
450 const double b = 1.20849664917601007033648;
451 const double c = 0.717933260754229499708010;
455 add(x, (c-a)*dt,
k,
y);
460 x.
Add((1.-a-b)*dt,
k);
virtual void Init(TimeDependentOperator &_f)
virtual void Init(TimeDependentOperator &_f)
virtual void Step(Vector &x, double &t, double &dt)
void SetSize(int s)
Resizes the vector if the new size is different.
Base abstract class for time dependent operators: (x,t) -> f(x,t)
int Width() const
Get the width (size of input) of the Operator. Synonym with NumCols.
virtual void SetTime(const double _t)
virtual void Step(Vector &x, double &t, double &dt)
virtual void Mult(const Vector &x, Vector &y) const =0
Operator application.
virtual void Init(TimeDependentOperator &_f)
virtual void Step(Vector &x, double &t, double &dt)
void add(const Vector &v1, const Vector &v2, Vector &v)
virtual void Step(Vector &x, double &t, double &dt)
virtual void Init(TimeDependentOperator &_f)
virtual void Init(TimeDependentOperator &_f)
virtual ~ExplicitRKSolver()
virtual void Step(Vector &x, double &t, double &dt)
virtual void Init(TimeDependentOperator &_f)
virtual void ImplicitSolve(const double dt, const Vector &x, Vector &k)
virtual void Step(Vector &x, double &t, double &dt)
ExplicitRKSolver(int _s, const double *_a, const double *_b, const double *_c)
virtual void Step(Vector &x, double &t, double &dt)
virtual void Init(TimeDependentOperator &_f)
virtual void Init(TimeDependentOperator &_f)
virtual void Init(TimeDependentOperator &_f)
Vector & Add(const double a, const Vector &Va)
(*this) += a * Va
virtual void Step(Vector &x, double &t, double &dt)
virtual void Init(TimeDependentOperator &_f)
virtual void Step(Vector &x, double &t, double &dt)
TimeDependentOperator * f
virtual void Step(Vector &x, double &t, double &dt)
SDIRK23Solver(int gamma_opt=1)
virtual void Init(TimeDependentOperator &_f)