3 boost/numeric/odeint/external/gsl/gsl_wrapper.hpp
6 Wrapper for gsl_vector.
9 Copyright 2011-2012 Mario Mulansky
10 Copyright 2011 Karsten Ahnert
12 Distributed under the Boost Software License, Version 1.0.
13 (See accompanying file LICENSE_1_0.txt or
14 copy at http://www.boost.org/LICENSE_1_0.txt)
18 #ifndef BOOST_NUMERIC_ODEINT_EXTERNAL_GSL_GSL_WRAPPER_HPP_INCLUDED
19 #define BOOST_NUMERIC_ODEINT_EXTERNAL_GSL_GSL_WRAPPER_HPP_INCLUDED
23 #include <gsl/gsl_vector.h>
25 #include <boost/type_traits/integral_constant.hpp>
26 #include <boost/range.hpp>
27 #include <boost/iterator/iterator_facade.hpp>
30 #include <boost/numeric/odeint/util/state_wrapper.hpp>
31 #include <boost/numeric/odeint/util/is_resizeable.hpp>
32 #include <boost/numeric/odeint/util/copy.hpp>
34 class const_gsl_vector_iterator;
37 * defines an iterator for gsl_vector
39 class gsl_vector_iterator : public boost::iterator_facade< gsl_vector_iterator , double , boost::random_access_traversal_tag >
43 gsl_vector_iterator( void ): m_p(0) , m_stride( 0 ) { }
44 explicit gsl_vector_iterator( gsl_vector *p ) : m_p( p->data ) , m_stride( p->stride ) { }
45 friend gsl_vector_iterator end_iterator( gsl_vector * );
49 friend class boost::iterator_core_access;
50 friend class const_gsl_vector_iterator;
52 void increment( void ) { m_p += m_stride; }
53 void decrement( void ) { m_p -= m_stride; }
54 void advance( ptrdiff_t n ) { m_p += n*m_stride; }
55 bool equal( const gsl_vector_iterator &other ) const { return this->m_p == other.m_p; }
56 bool equal( const const_gsl_vector_iterator &other ) const;
57 double& dereference( void ) const { return *m_p; }
66 * defines an const iterator for gsl_vector
68 class const_gsl_vector_iterator : public boost::iterator_facade< const_gsl_vector_iterator , const double , boost::random_access_traversal_tag >
72 const_gsl_vector_iterator( void ): m_p(0) , m_stride( 0 ) { }
73 explicit const_gsl_vector_iterator( const gsl_vector *p ) : m_p( p->data ) , m_stride( p->stride ) { }
74 const_gsl_vector_iterator( const gsl_vector_iterator &p ) : m_p( p.m_p ) , m_stride( p.m_stride ) { }
78 friend class boost::iterator_core_access;
79 friend class gsl_vector_iterator;
80 friend const_gsl_vector_iterator end_iterator( const gsl_vector * );
82 void increment( void ) { m_p += m_stride; }
83 void decrement( void ) { m_p -= m_stride; }
84 void advance( ptrdiff_t n ) { m_p += n*m_stride; }
85 bool equal( const const_gsl_vector_iterator &other ) const { return this->m_p == other.m_p; }
86 bool equal( const gsl_vector_iterator &other ) const { return this->m_p == other.m_p; }
87 const double& dereference( void ) const { return *m_p; }
94 bool gsl_vector_iterator::equal( const const_gsl_vector_iterator &other ) const { return this->m_p == other.m_p; }
97 gsl_vector_iterator end_iterator( gsl_vector *x )
99 gsl_vector_iterator iter( x );
100 iter.m_p += iter.m_stride * x->size;
104 const_gsl_vector_iterator end_iterator( const gsl_vector *x )
106 const_gsl_vector_iterator iter( x );
107 iter.m_p += iter.m_stride * x->size;
117 struct range_mutable_iterator< gsl_vector* >
119 typedef gsl_vector_iterator type;
123 struct range_const_iterator< gsl_vector* >
125 typedef const_gsl_vector_iterator type;
133 inline gsl_vector_iterator range_begin( gsl_vector *x )
135 return gsl_vector_iterator( x );
139 inline const_gsl_vector_iterator range_begin( const gsl_vector *x )
141 return const_gsl_vector_iterator( x );
145 inline gsl_vector_iterator range_end( gsl_vector *x )
147 return end_iterator( x );
151 inline const_gsl_vector_iterator range_end( const gsl_vector *x )
153 return end_iterator( x );
168 struct is_resizeable< gsl_vector* >
170 //struct type : public boost::true_type { };
171 typedef boost::true_type type;
172 const static bool value = type::value;
176 struct same_size_impl< gsl_vector* , gsl_vector* >
178 static bool same_size( const gsl_vector* x , const gsl_vector* y )
180 return x->size == y->size;
185 struct resize_impl< gsl_vector* , gsl_vector* >
187 static void resize( gsl_vector* &x , const gsl_vector* y )
189 gsl_vector_free( x );
190 x = gsl_vector_alloc( y->size );
195 struct state_wrapper< gsl_vector* >
197 typedef double value_type;
198 typedef gsl_vector* state_type;
199 typedef state_wrapper< gsl_vector* > state_wrapper_type;
205 m_v = gsl_vector_alloc( 1 );
208 state_wrapper( const state_wrapper_type &x )
210 resize( m_v , x.m_v );
211 gsl_vector_memcpy( m_v , x.m_v );
217 gsl_vector_free( m_v );
229 #endif // BOOST_NUMERIC_ODEINT_EXTERNAL_GSL_GSL_WRAPPER_HPP_INCLUDED