silnikov.cc


//////////////////////////////////////////////////////////////////////////// // Model silnikov.cc SIMLIB // // Silnikov equation: // // dx1/dt = x2 // dx2/dt = x3 // dx3/dt = -a*x3 - x2 + b*x1*(1 - c*x1 - d*x1*x1) // // where a=0.4, b=0.65, c=0, and d=1 are parameters // initial conditions of x1(0) = 0.1234, x2(0)=0.2, x3(0)=0.1 // output for Gnuplot // // SOURCE: Fishwick // // PerPet 1996, 2006 #include "simlib.h" #include <iostream> #include <cmath> using namespace std; const double StepPrn = 0.05; // step of output class Silnikov { public: Parameter a,b,c,d; Integrator x1, x2, x3; Silnikov(double _a, double _b, double _c, double _d) : a(_a), b(_b), c(_c), d(_d), x1(x2, 0.1234), x2(x3, 0.2), x3(-a*x3 - x2 + b*x1*(1 - c*x1 - d*x1*x1), 0.1) {} }; Silnikov e(0.4, 0.65, 0, 1); // output sampling: void Sample() { cout << Time << ' ' << e.x1.Value() << ' ' << e.x2.Value() << '\n'; } Sampler S(Sample, StepPrn); // experiment description: int main() { cout << "# Silnikov equation output \n"; SetStep(1e-8,1e-3); // set step size range Init(0,250); // experiment initialization Run(); // simulation }

Next two plots can be created by Gnuplot script (output data are stored in "silnikov.dat"):

set nokey set style data lines set terminal png large; set output "silnikov.png" set title "Silnikov equation, RKE method, maxstep = 1e-3" plot [][-1.6:1.6] "silnikov.dat" set terminal png large; set output "silnikov-phase.png" set title "Silnikov equation, RKE method, maxstep = 1e-3" plot [-1.6:1.6][-1.6:1.6] "silnikov.dat" using 2:3
[result1]
[result2]

Accumulated errors change the computed results for various integration methods (clearly visible at time > 200; it depends on precision of arithmetic - the result is for 32 bit code on i686 architecture):
[result-methods]


Last modification:

If you have some suggestions to improve this page, please mail to peringer AT fit.vutbr.cz