1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
|
/*
* elliptic_functions.cpp
*
* Copyright 2011-2013 Mario Mulansky
* Copyright 2011-2012 Karsten Ahnert
*
* Distributed under the Boost Software License, Version 1.0.
* (See accompanying file LICENSE_1_0.txt or
* copy at http://www.boost.org/LICENSE_1_0.txt)
*/
#include <iostream>
#include <fstream>
#include <cmath>
#include <boost/array.hpp>
#include <boost/numeric/odeint/config.hpp>
#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/stepper/bulirsch_stoer.hpp>
#include <boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp>
using namespace std;
using namespace boost::numeric::odeint;
typedef boost::array< double , 3 > state_type;
/*
* x1' = x2*x3
* x2' = -x1*x3
* x3' = -m*x1*x2
*/
void rhs( const state_type &x , state_type &dxdt , const double t )
{
static const double m = 0.51;
dxdt[0] = x[1]*x[2];
dxdt[1] = -x[0]*x[2];
dxdt[2] = -m*x[0]*x[1];
}
ofstream out;
void write_out( const state_type &x , const double t )
{
out << t << '\t' << x[0] << '\t' << x[1] << '\t' << x[2] << endl;
}
int main()
{
bulirsch_stoer_dense_out< state_type > stepper( 1E-9 , 1E-9 , 1.0 , 0.0 );
state_type x1 = {{ 0.0 , 1.0 , 1.0 }};
double t = 0.0;
double dt = 0.01;
out.open( "elliptic1.dat" );
out.precision(16);
integrate_const( stepper , rhs , x1 , t , 100.0 , dt , write_out );
out.close();
state_type x2 = {{ 0.0 , 1.0 , 1.0 }};
out.open( "elliptic2.dat" );
out.precision(16);
integrate_adaptive( stepper , rhs , x2 , t , 100.0 , dt , write_out );
out.close();
typedef runge_kutta_dopri5< state_type > dopri5_type;
typedef controlled_runge_kutta< dopri5_type > controlled_dopri5_type;
typedef dense_output_runge_kutta< controlled_dopri5_type > dense_output_dopri5_type;
dense_output_dopri5_type dopri5 = make_dense_output( 1E-9 , 1E-9 , dopri5_type() );
//dense_output_dopri5_type dopri5( controlled_dopri5_type( default_error_checker< double >( 1E-9 , 0.0 , 0.0 , 0.0 ) ) );
state_type x3 = {{ 0.0 , 1.0 , 1.0 }};
out.open( "elliptic3.dat" );
out.precision(16);
integrate_adaptive( dopri5 , rhs , x3 , t , 100.0 , dt , write_out );
out.close();
}
|
