Chris@16
|
1 /*
|
Chris@16
|
2 [auto_generated]
|
Chris@16
|
3 boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp
|
Chris@16
|
4
|
Chris@16
|
5 [begin_description]
|
Chris@16
|
6 Dense output for Rosenbrock 4.
|
Chris@16
|
7 [end_description]
|
Chris@16
|
8
|
Chris@16
|
9 Copyright 2009-2011 Karsten Ahnert
|
Chris@16
|
10 Copyright 2009-2011 Mario Mulansky
|
Chris@16
|
11
|
Chris@16
|
12 Distributed under the Boost Software License, Version 1.0.
|
Chris@16
|
13 (See accompanying file LICENSE_1_0.txt or
|
Chris@16
|
14 copy at http://www.boost.org/LICENSE_1_0.txt)
|
Chris@16
|
15 */
|
Chris@16
|
16
|
Chris@16
|
17
|
Chris@16
|
18 #ifndef BOOST_NUMERIC_ODEINT_STEPPER_ROSENBROCK4_DENSE_OUTPUT_HPP_INCLUDED
|
Chris@16
|
19 #define BOOST_NUMERIC_ODEINT_STEPPER_ROSENBROCK4_DENSE_OUTPUT_HPP_INCLUDED
|
Chris@16
|
20
|
Chris@16
|
21
|
Chris@16
|
22 #include <utility>
|
Chris@16
|
23
|
Chris@16
|
24 #include <boost/numeric/odeint/util/bind.hpp>
|
Chris@16
|
25
|
Chris@16
|
26 #include <boost/numeric/odeint/stepper/rosenbrock4_controller.hpp>
|
Chris@16
|
27 #include <boost/numeric/odeint/util/is_resizeable.hpp>
|
Chris@16
|
28
|
Chris@16
|
29
|
Chris@16
|
30 namespace boost {
|
Chris@16
|
31 namespace numeric {
|
Chris@16
|
32 namespace odeint {
|
Chris@16
|
33
|
Chris@16
|
34 template< class ControlledStepper >
|
Chris@16
|
35 class rosenbrock4_dense_output
|
Chris@16
|
36 {
|
Chris@16
|
37
|
Chris@16
|
38 public:
|
Chris@16
|
39
|
Chris@16
|
40 typedef ControlledStepper controlled_stepper_type;
|
Chris@16
|
41 typedef typename controlled_stepper_type::stepper_type stepper_type;
|
Chris@16
|
42 typedef typename stepper_type::value_type value_type;
|
Chris@16
|
43 typedef typename stepper_type::state_type state_type;
|
Chris@16
|
44 typedef typename stepper_type::wrapped_state_type wrapped_state_type;
|
Chris@16
|
45 typedef typename stepper_type::time_type time_type;
|
Chris@16
|
46 typedef typename stepper_type::deriv_type deriv_type;
|
Chris@16
|
47 typedef typename stepper_type::wrapped_deriv_type wrapped_deriv_type;
|
Chris@16
|
48 typedef typename stepper_type::resizer_type resizer_type;
|
Chris@16
|
49 typedef dense_output_stepper_tag stepper_category;
|
Chris@16
|
50
|
Chris@16
|
51 typedef rosenbrock4_dense_output< ControlledStepper > dense_output_stepper_type;
|
Chris@16
|
52
|
Chris@16
|
53 rosenbrock4_dense_output( const controlled_stepper_type &stepper = controlled_stepper_type() )
|
Chris@16
|
54 : m_stepper( stepper ) ,
|
Chris@16
|
55 m_x1() , m_x2() ,
|
Chris@16
|
56 m_current_state_x1( true ) ,
|
Chris@16
|
57 m_t() , m_t_old() , m_dt()
|
Chris@16
|
58 {
|
Chris@16
|
59 }
|
Chris@16
|
60
|
Chris@16
|
61
|
Chris@16
|
62
|
Chris@16
|
63 template< class StateType >
|
Chris@16
|
64 void initialize( const StateType &x0 , time_type t0 , time_type dt0 )
|
Chris@16
|
65 {
|
Chris@16
|
66 m_resizer.adjust_size( x0 , detail::bind( &dense_output_stepper_type::template resize_impl< StateType > , detail::ref( *this ) , detail::_1 ) );
|
Chris@16
|
67 get_current_state() = x0;
|
Chris@16
|
68 m_t = t0;
|
Chris@16
|
69 m_dt = dt0;
|
Chris@16
|
70 }
|
Chris@16
|
71
|
Chris@16
|
72 template< class System >
|
Chris@16
|
73 std::pair< time_type , time_type > do_step( System system )
|
Chris@16
|
74 {
|
Chris@16
|
75 const size_t max_count = 1000;
|
Chris@16
|
76
|
Chris@16
|
77 controlled_step_result res = fail;
|
Chris@16
|
78 m_t_old = m_t;
|
Chris@16
|
79 size_t count = 0;
|
Chris@16
|
80 do
|
Chris@16
|
81 {
|
Chris@16
|
82 res = m_stepper.try_step( system , get_current_state() , m_t , get_old_state() , m_dt );
|
Chris@16
|
83 if( count++ == max_count )
|
Chris@16
|
84 throw std::overflow_error( "rosenbrock4 : too much iterations!");
|
Chris@16
|
85 }
|
Chris@16
|
86 while( res == fail );
|
Chris@16
|
87 m_stepper.stepper().prepare_dense_output();
|
Chris@16
|
88 this->toggle_current_state();
|
Chris@16
|
89 return std::make_pair( m_t_old , m_t );
|
Chris@16
|
90 }
|
Chris@16
|
91
|
Chris@16
|
92
|
Chris@16
|
93 /*
|
Chris@16
|
94 * The two overloads are needed in order to solve the forwarding problem.
|
Chris@16
|
95 */
|
Chris@16
|
96 template< class StateOut >
|
Chris@16
|
97 void calc_state( time_type t , StateOut &x )
|
Chris@16
|
98 {
|
Chris@16
|
99 m_stepper.stepper().calc_state( t , x , get_old_state() , m_t_old , get_current_state() , m_t );
|
Chris@16
|
100 }
|
Chris@16
|
101
|
Chris@16
|
102 template< class StateOut >
|
Chris@16
|
103 void calc_state( time_type t , const StateOut &x )
|
Chris@16
|
104 {
|
Chris@16
|
105 m_stepper.stepper().calc_state( t , x , get_old_state() , m_t_old , get_current_state() , m_t );
|
Chris@16
|
106 }
|
Chris@16
|
107
|
Chris@16
|
108
|
Chris@16
|
109 template< class StateType >
|
Chris@16
|
110 void adjust_size( const StateType &x )
|
Chris@16
|
111 {
|
Chris@16
|
112 m_stepper.adjust_size( x );
|
Chris@16
|
113 resize_impl( x );
|
Chris@16
|
114 }
|
Chris@16
|
115
|
Chris@16
|
116
|
Chris@16
|
117
|
Chris@16
|
118
|
Chris@16
|
119 const state_type& current_state( void ) const
|
Chris@16
|
120 {
|
Chris@16
|
121 return get_current_state();
|
Chris@16
|
122 }
|
Chris@16
|
123
|
Chris@16
|
124 time_type current_time( void ) const
|
Chris@16
|
125 {
|
Chris@16
|
126 return m_t;
|
Chris@16
|
127 }
|
Chris@16
|
128
|
Chris@16
|
129 const state_type& previous_state( void ) const
|
Chris@16
|
130 {
|
Chris@16
|
131 return get_old_state();
|
Chris@16
|
132 }
|
Chris@16
|
133
|
Chris@16
|
134 time_type previous_time( void ) const
|
Chris@16
|
135 {
|
Chris@16
|
136 return m_t_old;
|
Chris@16
|
137 }
|
Chris@16
|
138
|
Chris@16
|
139 time_type current_time_step( void ) const
|
Chris@16
|
140 {
|
Chris@16
|
141 return m_dt;
|
Chris@16
|
142 }
|
Chris@16
|
143
|
Chris@16
|
144
|
Chris@16
|
145
|
Chris@16
|
146
|
Chris@16
|
147 private:
|
Chris@16
|
148
|
Chris@16
|
149 state_type& get_current_state( void )
|
Chris@16
|
150 {
|
Chris@16
|
151 return m_current_state_x1 ? m_x1.m_v : m_x2.m_v ;
|
Chris@16
|
152 }
|
Chris@16
|
153
|
Chris@16
|
154 const state_type& get_current_state( void ) const
|
Chris@16
|
155 {
|
Chris@16
|
156 return m_current_state_x1 ? m_x1.m_v : m_x2.m_v ;
|
Chris@16
|
157 }
|
Chris@16
|
158
|
Chris@16
|
159 state_type& get_old_state( void )
|
Chris@16
|
160 {
|
Chris@16
|
161 return m_current_state_x1 ? m_x2.m_v : m_x1.m_v ;
|
Chris@16
|
162 }
|
Chris@16
|
163
|
Chris@16
|
164 const state_type& get_old_state( void ) const
|
Chris@16
|
165 {
|
Chris@16
|
166 return m_current_state_x1 ? m_x2.m_v : m_x1.m_v ;
|
Chris@16
|
167 }
|
Chris@16
|
168
|
Chris@16
|
169 void toggle_current_state( void )
|
Chris@16
|
170 {
|
Chris@16
|
171 m_current_state_x1 = ! m_current_state_x1;
|
Chris@16
|
172 }
|
Chris@16
|
173
|
Chris@16
|
174
|
Chris@16
|
175 template< class StateIn >
|
Chris@16
|
176 bool resize_impl( const StateIn &x )
|
Chris@16
|
177 {
|
Chris@16
|
178 bool resized = false;
|
Chris@16
|
179 resized |= adjust_size_by_resizeability( m_x1 , x , typename is_resizeable<state_type>::type() );
|
Chris@16
|
180 resized |= adjust_size_by_resizeability( m_x2 , x , typename is_resizeable<state_type>::type() );
|
Chris@16
|
181 return resized;
|
Chris@16
|
182 }
|
Chris@16
|
183
|
Chris@16
|
184
|
Chris@16
|
185 controlled_stepper_type m_stepper;
|
Chris@16
|
186 resizer_type m_resizer;
|
Chris@16
|
187 wrapped_state_type m_x1 , m_x2;
|
Chris@16
|
188 bool m_current_state_x1;
|
Chris@16
|
189 time_type m_t , m_t_old , m_dt;
|
Chris@16
|
190 };
|
Chris@16
|
191
|
Chris@16
|
192
|
Chris@16
|
193
|
Chris@16
|
194 } // namespace odeint
|
Chris@16
|
195 } // namespace numeric
|
Chris@16
|
196 } // namespace boost
|
Chris@16
|
197
|
Chris@16
|
198
|
Chris@16
|
199 #endif // BOOST_NUMERIC_ODEINT_STEPPER_ROSENBROCK4_DENSE_OUTPUT_HPP_INCLUDED
|