summaryrefslogtreecommitdiff
path: root/boost/numeric
diff options
context:
space:
mode:
Diffstat (limited to 'boost/numeric')
-rw-r--r--boost/numeric/interval/compare/set.hpp12
-rw-r--r--boost/numeric/interval/detail/bugs.hpp31
-rw-r--r--boost/numeric/interval/detail/msvc_rounding_control.hpp15
-rw-r--r--boost/numeric/interval/limits.hpp2
-rw-r--r--boost/numeric/interval/rounded_arith.hpp8
-rw-r--r--boost/numeric/odeint/algebra/algebra_dispatcher.hpp29
-rw-r--r--boost/numeric/odeint/algebra/array_algebra.hpp309
-rw-r--r--boost/numeric/odeint/algebra/multi_array_algebra.hpp4
-rw-r--r--boost/numeric/odeint/external/eigen/eigen.hpp27
-rw-r--r--boost/numeric/odeint/external/eigen/eigen_algebra.hpp13
-rw-r--r--boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp49
-rw-r--r--boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp55
-rw-r--r--boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp61
-rw-r--r--boost/numeric/odeint/external/thrust/thrust_resize.hpp230
-rw-r--r--boost/numeric/odeint/integrate/integrate.hpp22
-rw-r--r--boost/numeric/odeint/iterator/impl/adaptive_iterator_impl.hpp2
-rw-r--r--boost/numeric/odeint/iterator/impl/const_step_iterator_impl.hpp2
-rw-r--r--boost/numeric/odeint/iterator/integrate/detail/integrate_n_steps.hpp4
-rw-r--r--boost/numeric/odeint/stepper/adams_bashforth.hpp28
-rw-r--r--boost/numeric/odeint/stepper/adams_bashforth_moulton.hpp10
-rw-r--r--boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp35
-rw-r--r--boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp30
-rw-r--r--boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp31
-rw-r--r--boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp2
-rw-r--r--boost/numeric/odeint/stepper/bulirsch_stoer.hpp7
-rw-r--r--boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp5
-rw-r--r--boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp4
-rw-r--r--boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp2
-rw-r--r--boost/numeric/odeint/stepper/euler.hpp4
-rw-r--r--boost/numeric/odeint/stepper/extrapolation_stepper.hpp293
-rw-r--r--boost/numeric/odeint/stepper/implicit_euler.hpp2
-rw-r--r--boost/numeric/odeint/stepper/rosenbrock4_controller.hpp1
-rw-r--r--boost/numeric/odeint/stepper/runge_kutta4_classic.hpp11
-rw-r--r--boost/numeric/odeint/stepper/velocity_verlet.hpp2
-rw-r--r--boost/numeric/odeint/util/bind.hpp3
-rw-r--r--boost/numeric/odeint/util/resize.hpp4
-rw-r--r--boost/numeric/odeint/util/same_instance.hpp2
-rw-r--r--boost/numeric/odeint/util/unwrap_reference.hpp15
-rw-r--r--boost/numeric/ublas/detail/matrix_assign.hpp16
-rw-r--r--boost/numeric/ublas/experimental/sparse_view.hpp1
-rw-r--r--boost/numeric/ublas/functional.hpp3
-rw-r--r--boost/numeric/ublas/matrix_expression.hpp2
-rw-r--r--boost/numeric/ublas/storage.hpp36
-rw-r--r--boost/numeric/ublas/traits.hpp3
-rw-r--r--boost/numeric/ublas/vector.hpp11
45 files changed, 1100 insertions, 338 deletions
diff --git a/boost/numeric/interval/compare/set.hpp b/boost/numeric/interval/compare/set.hpp
index aa4f171..44f4d6e 100644
--- a/boost/numeric/interval/compare/set.hpp
+++ b/boost/numeric/interval/compare/set.hpp
@@ -27,7 +27,7 @@ bool operator<(const interval<T, Policies1>& x, const interval<T, Policies2>& y)
}
template<class T, class Policies> inline
-bool operator<(const interval<T, Policies>& x, const T& y)
+bool operator<(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
@@ -39,7 +39,7 @@ bool operator<=(const interval<T, Policies1>& x, const interval<T, Policies2>& y
}
template<class T, class Policies> inline
-bool operator<=(const interval<T, Policies>& x, const T& y)
+bool operator<=(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
@@ -51,7 +51,7 @@ bool operator>(const interval<T, Policies1>& x, const interval<T, Policies2>& y)
}
template<class T, class Policies> inline
-bool operator>(const interval<T, Policies>& x, const T& y)
+bool operator>(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
@@ -63,7 +63,7 @@ bool operator>=(const interval<T, Policies1>& x, const interval<T, Policies2>& y
}
template<class T, class Policies> inline
-bool operator>=(const interval<T, Policies>& x, const T& y)
+bool operator>=(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
@@ -75,7 +75,7 @@ bool operator==(const interval<T, Policies1>& x, const interval<T, Policies2>& y
}
template<class T, class Policies> inline
-bool operator==(const interval<T, Policies>& x, const T& y)
+bool operator==(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
@@ -87,7 +87,7 @@ bool operator!=(const interval<T, Policies1>& x, const interval<T, Policies2>& y
}
template<class T, class Policies> inline
-bool operator!=(const interval<T, Policies>& x, const T& y)
+bool operator!=(const interval<T, Policies>& , const T& )
{
throw comparison_error();
}
diff --git a/boost/numeric/interval/detail/bugs.hpp b/boost/numeric/interval/detail/bugs.hpp
index cc37988..b9b421a 100644
--- a/boost/numeric/interval/detail/bugs.hpp
+++ b/boost/numeric/interval/detail/bugs.hpp
@@ -45,35 +45,4 @@
# define BOOST_NUMERIC_INTERVAL_using_ahyp(a)
#endif
-#if defined(__GNUC__) && (__GNUC__ <= 2)
-// cf PR c++/1981 for a description of the bug
-#include <algorithm>
-#include <boost/config/no_tr1/cmath.hpp>
-namespace boost {
-namespace numeric {
- using std::min;
- using std::max;
- using std::sqrt;
- using std::exp;
- using std::log;
- using std::cos;
- using std::tan;
- using std::asin;
- using std::acos;
- using std::atan;
- using std::ceil;
- using std::floor;
- using std::sinh;
- using std::cosh;
- using std::tanh;
-# undef BOOST_NUMERIC_INTERVAL_using_max
-# undef BOOST_NUMERIC_INTERVAL_using_math
-# define BOOST_NUMERIC_INTERVAL_using_max(a)
-# define BOOST_NUMERIC_INTERVAL_using_math(a)
-# undef BOOST_NUMERIC_INTERVAL_using_ahyp
-# define BOOST_NUMERIC_INTERVAL_using_ahyp(a)
-} // namespace numeric
-} // namespace boost
-#endif
-
#endif // BOOST_NUMERIC_INTERVAL_DETAIL_BUGS
diff --git a/boost/numeric/interval/detail/msvc_rounding_control.hpp b/boost/numeric/interval/detail/msvc_rounding_control.hpp
index d36e802..ccf1e8b 100644
--- a/boost/numeric/interval/detail/msvc_rounding_control.hpp
+++ b/boost/numeric/interval/detail/msvc_rounding_control.hpp
@@ -88,7 +88,20 @@ struct x86_rounding
static void get_rounding_mode(rounding_mode& mode)
{ mode = msvc2hard(_control87(0, 0)); }
static void set_rounding_mode(const rounding_mode mode)
- { _control87(hard2msvc(mode), _MCW_EM | _MCW_RC | _MCW_PC | _MCW_IC); }
+ {
+ _control87(hard2msvc(mode),
+ _MCW_EM | _MCW_RC
+#if !defined(_M_AMD64) && !defined(_M_ARM)
+ // x64 ignores _MCW_PC and _MCW_IC, and the Debug CRT library actually
+ // asserts when these are passed to _control87.
+ // MSDN says on '_control87' that changing precision (_MCW_PC) or
+ // infinity (_MCW_IC) handling is not supported on the ARM and x64
+ // architectures and that _control87 raises an assertion
+ // and the invalid parameter handler is invoked.
+ | _MCW_PC | _MCW_IC
+#endif
+ );
+ }
static double to_int(const double& x) { return rint(x); }
};
diff --git a/boost/numeric/interval/limits.hpp b/boost/numeric/interval/limits.hpp
index d691cce..973681d 100644
--- a/boost/numeric/interval/limits.hpp
+++ b/boost/numeric/interval/limits.hpp
@@ -11,7 +11,6 @@
#ifndef BOOST_NUMERIC_INTERVAL_LIMITS_HPP
#define BOOST_NUMERIC_INTERVAL_LIMITS_HPP
-#ifndef BOOST_NO_TEMPLATE_PARTIAL_SPECIALIZATION
#include <boost/config.hpp>
#include <boost/limits.hpp>
@@ -46,6 +45,5 @@ private:
} // namespace std
-#endif
#endif // BOOST_NUMERIC_INTERVAL_LIMITS_HPP
diff --git a/boost/numeric/interval/rounded_arith.hpp b/boost/numeric/interval/rounded_arith.hpp
index 8b2d9a7..2520274 100644
--- a/boost/numeric/interval/rounded_arith.hpp
+++ b/boost/numeric/interval/rounded_arith.hpp
@@ -68,8 +68,8 @@ struct rounded_arith_std: Rounding {
{ BOOST_NUMERIC_INTERVAL_using_math(sqrt); BOOST_DN(sqrt(x)); }
T sqrt_up (const T& x)
{ BOOST_NUMERIC_INTERVAL_using_math(sqrt); BOOST_UP(sqrt(x)); }
- T int_down(const T& x) { this->downward(); return to_int(x); }
- T int_up (const T& x) { this->upward(); return to_int(x); }
+ T int_down(const T& x) { this->downward(); return this->to_int(x); }
+ T int_up (const T& x) { this->upward(); return this->to_int(x); }
# undef BOOST_DN
# undef BOOST_NR
# undef BOOST_UP
@@ -105,8 +105,8 @@ struct rounded_arith_opp: Rounding {
{ BOOST_NUMERIC_INTERVAL_using_math(sqrt); BOOST_DN(sqrt(x)); }
T sqrt_up (const T& x)
{ BOOST_NUMERIC_INTERVAL_using_math(sqrt); BOOST_UP(sqrt(x)); }
- T int_down(const T& x) { return -to_int(-x); }
- T int_up (const T& x) { return to_int(x); }
+ T int_down(const T& x) { return -this->to_int(-x); }
+ T int_up (const T& x) { return this->to_int(x); }
# undef BOOST_DN
# undef BOOST_NR
# undef BOOST_UP
diff --git a/boost/numeric/odeint/algebra/algebra_dispatcher.hpp b/boost/numeric/odeint/algebra/algebra_dispatcher.hpp
index df5aaba..88cf159 100644
--- a/boost/numeric/odeint/algebra/algebra_dispatcher.hpp
+++ b/boost/numeric/odeint/algebra/algebra_dispatcher.hpp
@@ -17,8 +17,9 @@
#ifndef BOOST_NUMERIC_ODEINT_ALGEBRA_ALGEBRA_DISPATCHER_HPP_INCLUDED
#define BOOST_NUMERIC_ODEINT_ALGEBRA_ALGEBRA_DISPATCHER_HPP_INCLUDED
-#include <complex>
+#include <boost/numeric/odeint/config.hpp>
+#include <complex>
#include <boost/type_traits/is_floating_point.hpp>
#include <boost/numeric/ublas/vector.hpp>
@@ -38,14 +39,14 @@ namespace odeint {
template< class StateType , class Enabler = void >
struct algebra_dispatcher_sfinae
{
- // range_algebra is the standard algebra^
+ // range_algebra is the standard algebra
typedef range_algebra algebra_type;
};
template< class StateType >
struct algebra_dispatcher : algebra_dispatcher_sfinae< StateType > { };
-//specialize for array
+// specialize for array
template< class T , size_t N >
struct algebra_dispatcher< boost::array< T , N > >
{
@@ -84,4 +85,26 @@ struct algebra_dispatcher< boost::numeric::ublas::matrix< T , L , A > >
}
}
+#ifdef BOOST_NUMERIC_ODEINT_CXX11
+
+// c++11 mode: specialization for std::array if available
+
+#include <array>
+
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+// specialize for std::array
+template< class T , size_t N >
+struct algebra_dispatcher< std::array< T , N > >
+{
+ typedef array_algebra algebra_type;
+};
+
+} } }
+
+#endif
+
+
#endif
diff --git a/boost/numeric/odeint/algebra/array_algebra.hpp b/boost/numeric/odeint/algebra/array_algebra.hpp
index c565d6d..471e866 100644
--- a/boost/numeric/odeint/algebra/array_algebra.hpp
+++ b/boost/numeric/odeint/algebra/array_algebra.hpp
@@ -3,7 +3,11 @@
boost/numeric/odeint/algebra/array_algebra.hpp
[begin_description]
- Algebra for boost::array. Highly specialized for odeint. Const arguments are introduce to work with odeint.
+ Algebra for Arrays. Highly specialized for odeint. Const arguments are
+ introduce to work with odeint.
+ The Array algebra can be used for Array structures with two template
+ parameters:
+ Array<T, N>
[end_description]
Copyright 2011-2013 Mario Mulansky
@@ -29,234 +33,251 @@ namespace odeint {
struct array_algebra
{
- template< typename T , size_t dim , class Op >
- static void for_each1( boost::array< T , dim > &s1 , Op op )
+ //template< typename T , size_t dim , class Op >
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each1( Array< T, dim > &s1, Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] );
}
-
- template< typename T1 , typename T2 , size_t dim , class Op >
- static void for_each2( boost::array< T1 , dim > &s1 ,
- const boost::array< T2 , dim > &s2 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each2( Array< T, dim > &s1, const Array< T, dim > &s2,
+ Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each3( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each3( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] );
}
/* different const signature - required for the scale_sum_swap2 operation */
- template< typename T , size_t dim , class Op >
- static void for_each3( boost::array< T , dim > &s1 ,
- boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each3( Array< T , dim > &s1 ,
+ Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each4( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each4( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each5( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each5( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each6( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each6( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each7( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each7( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each8( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each8( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each9( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each9( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each10( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each10( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each11( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 ,
- const boost::array< T , dim > &s11 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each11( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 ,
+ const Array< T , dim > &s11 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] , s11[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each12( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 ,
- const boost::array< T , dim > &s11 ,
- const boost::array< T , dim > &s12 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each12( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 ,
+ const Array< T , dim > &s11 ,
+ const Array< T , dim > &s12 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] , s11[i] , s12[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each13( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 ,
- const boost::array< T , dim > &s11 ,
- const boost::array< T , dim > &s12 ,
- const boost::array< T , dim > &s13 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each13( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 ,
+ const Array< T , dim > &s11 ,
+ const Array< T , dim > &s12 ,
+ const Array< T , dim > &s13 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] , s11[i] , s12[i] , s13[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each14( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 ,
- const boost::array< T , dim > &s11 ,
- const boost::array< T , dim > &s12 ,
- const boost::array< T , dim > &s13 ,
- const boost::array< T , dim > &s14 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each14( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 ,
+ const Array< T , dim > &s11 ,
+ const Array< T , dim > &s12 ,
+ const Array< T , dim > &s13 ,
+ const Array< T , dim > &s14 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] , s11[i] , s12[i] , s13[i] , s14[i] );
}
- template< typename T , size_t dim , class Op >
- static void for_each15( boost::array< T , dim > &s1 ,
- const boost::array< T , dim > &s2 ,
- const boost::array< T , dim > &s3 ,
- const boost::array< T , dim > &s4 ,
- const boost::array< T , dim > &s5 ,
- const boost::array< T , dim > &s6 ,
- const boost::array< T , dim > &s7 ,
- const boost::array< T , dim > &s8 ,
- const boost::array< T , dim > &s9 ,
- const boost::array< T , dim > &s10 ,
- const boost::array< T , dim > &s11 ,
- const boost::array< T , dim > &s12 ,
- const boost::array< T , dim > &s13 ,
- const boost::array< T , dim > &s14 ,
- const boost::array< T , dim > &s15 , Op op )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim, class Op >
+ static void for_each15( Array< T , dim > &s1 ,
+ const Array< T , dim > &s2 ,
+ const Array< T , dim > &s3 ,
+ const Array< T , dim > &s4 ,
+ const Array< T , dim > &s5 ,
+ const Array< T , dim > &s6 ,
+ const Array< T , dim > &s7 ,
+ const Array< T , dim > &s8 ,
+ const Array< T , dim > &s9 ,
+ const Array< T , dim > &s10 ,
+ const Array< T , dim > &s11 ,
+ const Array< T , dim > &s12 ,
+ const Array< T , dim > &s13 ,
+ const Array< T , dim > &s14 ,
+ const Array< T , dim > &s15 , Op op )
{
for( size_t i=0 ; i<dim ; ++i )
op( s1[i] , s2[i] , s3[i] , s4[i] , s5[i] , s6[i] , s7[i] , s8[i] , s9[i] , s10[i] , s11[i] , s12[i] , s13[i] , s14[i] , s15[i] );
}
- template< typename T , size_t dim >
- static typename norm_result_type< boost::array< T , dim > >::type norm_inf( const boost::array< T , dim > &s )
+ template < template < typename, size_t > class Array, typename T,
+ size_t dim>
+ static typename norm_result_type< Array< T , dim > >::type norm_inf( const Array< T , dim > &s )
{
BOOST_USING_STD_MAX();
using std::abs;
- typedef typename norm_result_type< boost::array< T , dim > >::type result_type;
+ typedef typename norm_result_type< Array< T , dim > >::type result_type;
result_type init = static_cast< result_type >( 0 );
for( size_t i=0 ; i<dim ; ++i )
init = max BOOST_PREVENT_MACRO_SUBSTITUTION ( init , static_cast< result_type >(abs(s[i])) );
diff --git a/boost/numeric/odeint/algebra/multi_array_algebra.hpp b/boost/numeric/odeint/algebra/multi_array_algebra.hpp
index 4d1fc16..0bc476e 100644
--- a/boost/numeric/odeint/algebra/multi_array_algebra.hpp
+++ b/boost/numeric/odeint/algebra/multi_array_algebra.hpp
@@ -113,13 +113,13 @@ struct multi_array_algebra
template< class S1 , class S2 , class S3 , class S4 , class S5 , class S6 ,class S7 , class S8 , class S9 , class S10 , class S11 , class S12 , class S13 , class S14 , class Op >
static void for_each14( S1 &s1 , S2 &s2 , S3 &s3 , S4 &s4 , S5 &s5 , S6 &s6 , S7 &s7 , S8 &s8 , S9 &s9 , S10 &s10 , S11 &s11 , S12 &s12 , S13 &s13 , S14 &s14 , Op op )
{
- detail::for_each14( s1.data() , s1.data() + s1.num_elements() , s2.data() , s3.data() , s4.data() , s5.data() , s6.data() , s7.data() , s8.data() , s9.data() , s10.data() , s11.data() , s12.data() , s14.data() , op );
+ detail::for_each14( s1.data() , s1.data() + s1.num_elements() , s2.data() , s3.data() , s4.data() , s5.data() , s6.data() , s7.data() , s8.data() , s9.data() , s10.data() , s11.data() , s12.data() , s13.data() , s14.data() , op );
}
template< class S1 , class S2 , class S3 , class S4 , class S5 , class S6 ,class S7 , class S8 , class S9 , class S10 , class S11 , class S12 , class S13 , class S14 , class S15 , class Op >
static void for_each15( S1 &s1 , S2 &s2 , S3 &s3 , S4 &s4 , S5 &s5 , S6 &s6 , S7 &s7 , S8 &s8 , S9 &s9 , S10 &s10 , S11 &s11 , S12 &s12 , S13 &s13 , S14 &s14 , S15 &s15 , Op op )
{
- detail::for_each15( s1.data() , s1.data() + s1.num_elements() , s2.data() , s3.data() , s4.data() , s5.data() , s6.data() , s7.data() , s8.data() , s9.data() , s10.data() , s11.data() , s12.data() , s14.data() , s15.data() , op );
+ detail::for_each15( s1.data() , s1.data() + s1.num_elements() , s2.data() , s3.data() , s4.data() , s5.data() , s6.data() , s7.data() , s8.data() , s9.data() , s10.data() , s11.data() , s12.data() , s13.data() , s14.data() , s15.data() , op );
}
template< typename S >
diff --git a/boost/numeric/odeint/external/eigen/eigen.hpp b/boost/numeric/odeint/external/eigen/eigen.hpp
new file mode 100644
index 0000000..620acae
--- /dev/null
+++ b/boost/numeric/odeint/external/eigen/eigen.hpp
@@ -0,0 +1,27 @@
+/*
+ [auto_generated]
+ boost/numeric/odeint/external/eigen/eigen.hpp
+
+ [begin_description]
+ tba.
+ [end_description]
+
+ Copyright 2009-2012 Karsten Ahnert
+ Copyright 2009-2012 Mario Mulansky
+
+ 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)
+*/
+
+
+#ifndef BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_HPP_INCLUDED
+
+
+#include <boost/numeric/odeint/external/eigen/eigen_algebra.hpp>
+#include <boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp>
+#include <boost/numeric/odeint/external/eigen/eigen_resize.hpp>
+
+
+#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_HPP_INCLUDED
diff --git a/boost/numeric/odeint/external/eigen/eigen_algebra.hpp b/boost/numeric/odeint/external/eigen/eigen_algebra.hpp
index 6ca2026..b4ee5c3 100644
--- a/boost/numeric/odeint/external/eigen/eigen_algebra.hpp
+++ b/boost/numeric/odeint/external/eigen/eigen_algebra.hpp
@@ -71,6 +71,19 @@ operator/(const Eigen::MatrixBase<D1> &x1, const Eigen::MatrixBase<D2> &x2) {
return x1.cwiseQuotient(x2);
}
+
+template< typename D >
+inline const
+typename Eigen::CwiseUnaryOp<
+ typename Eigen::internal::scalar_abs_op<
+ typename Eigen::internal::traits< D >::Scalar > ,
+ const D >
+abs( const Eigen::MatrixBase< D > &m ) {
+ return m.cwiseAbs();
+}
+
+
+
} // end Eigen namespace
diff --git a/boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp b/boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp
new file mode 100644
index 0000000..6c8a3a2
--- /dev/null
+++ b/boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp
@@ -0,0 +1,49 @@
+/*
+ [auto_generated]
+ boost/numeric/odeint/external/eigen/eigen_algebra_dispatcher.hpp
+
+ [begin_description]
+ tba.
+ [end_description]
+
+ Copyright 2009-2012 Karsten Ahnert
+ Copyright 2009-2012 Mario Mulansky
+
+ 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)
+*/
+
+
+#ifndef BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_ALGEBRA_DISPATCHER_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_ALGEBRA_DISPATCHER_HPP_INCLUDED
+
+
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+
+template< class Derived >
+struct algebra_dispatcher_sfinae< Derived ,
+ typename boost::enable_if< typename boost::is_base_of< Eigen::MatrixBase< Derived > , Derived >::type >::type >
+{
+ typedef vector_space_algebra algebra_type;
+};
+
+
+template < class Derived >
+struct algebra_dispatcher_sfinae< Derived ,
+ typename boost::enable_if< typename boost::is_base_of< Eigen::ArrayBase< Derived > , Derived >::type >::type >
+{
+ typedef vector_space_algebra algebra_type;
+};
+
+
+
+} // namespace odeint
+} // namespace numeric
+} // namespace boost
+
+
+#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_EIGEN_EIGEN_ALGEBRA_DISPATCHER_HPP_INCLUDED
diff --git a/boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp b/boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp
index 5deba2c..0961346 100644
--- a/boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp
+++ b/boost/numeric/odeint/external/thrust/thrust_algebra_dispatcher.hpp
@@ -24,6 +24,7 @@
#include <boost/numeric/odeint/external/thrust/thrust_algebra.hpp>
#include <boost/numeric/odeint/algebra/algebra_dispatcher.hpp>
+// specializations for the standard thrust containers
namespace boost {
namespace numeric {
@@ -48,5 +49,59 @@ struct algebra_dispatcher< thrust::device_vector< T , A > >
} // namespace boost
+// add support for thrust backend vectors, if available
+
+#include <thrust/version.h>
+
+#if THRUST_VERSION >= 100600
+
+// specialization for thrust cpp vector
+#include <thrust/system/cpp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct algebra_dispatcher< thrust::cpp::vector< T , A > >
+ {
+ typedef thrust_algebra algebra_type;
+ };
+} } }
+
+// specialization for thrust omp vector
+#ifdef _OPENMP
+#include <thrust/system/omp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct algebra_dispatcher< thrust::omp::vector< T , A > >
+ {
+ typedef thrust_algebra algebra_type;
+ };
+} } }
+#endif // _OPENMP
+
+// specialization for thrust tbb vector
+#ifdef TBB_VERSION_MAJOR
+#include <thrust/system/tbb/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct algebra_dispatcher< thrust::tbb::vector< T , A > >
+ {
+ typedef thrust_algebra algebra_type;
+ };
+} } }
+#endif // TBB_VERSION_MAJOR
+
+// specialization for thrust cuda vector
+#ifdef __CUDACC__
+#include <thrust/system/cuda/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct algebra_dispatcher< thrust::cuda::vector< T , A > >
+ {
+ typedef thrust_algebra algebra_type;
+ };
+} } }
+#endif // __CUDACC__
+
+#endif // THRUST_VERSION >= 100600
+
#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_ALGEBRA_DISPATCHER_HPP_DEFINED
diff --git a/boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp b/boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp
index 0105a4e..e9b3a64 100644
--- a/boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp
+++ b/boost/numeric/odeint/external/thrust/thrust_operations_dispatcher.hpp
@@ -6,8 +6,8 @@
operations_dispatcher specialization for thrust
[end_description]
- Copyright 2013 Karsten Ahnert
- Copyright 2013 Mario Mulansky
+ Copyright 2013-2014 Karsten Ahnert
+ Copyright 2013-2014 Mario Mulansky
Distributed under the Boost Software License, Version 1.0.
(See accompanying file LICENSE_1_0.txt or
@@ -24,6 +24,7 @@
#include <boost/numeric/odeint/external/thrust/thrust_operations.hpp>
#include <boost/numeric/odeint/algebra/operations_dispatcher.hpp>
+// support for the standard thrust containers
namespace boost {
namespace numeric {
@@ -47,6 +48,60 @@ struct operations_dispatcher< thrust::device_vector< T , A > >
} // namespace numeric
} // namespace boost
+// add support for thrust backend vectors, if available
-#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_ALGEBRA_DISPATCHER_HPP_DEFINED
+#include <thrust/version.h>
+
+#if THRUST_VERSION >= 100600
+
+// specialization for thrust cpp vector
+#include <thrust/system/cpp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct operations_dispatcher< thrust::cpp::vector< T , A > >
+ {
+ typedef thrust_operations operations_type;
+ };
+} } }
+
+// specialization for thrust omp vector
+#ifdef _OPENMP
+#include <thrust/system/omp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct operations_dispatcher< thrust::omp::vector< T , A > >
+ {
+ typedef thrust_operations operations_type;
+ };
+} } }
+#endif // _OPENMP
+
+// specialization for thrust tbb vector
+#ifdef TBB_VERSION_MAJOR
+#include <thrust/system/tbb/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct operations_dispatcher< thrust::tbb::vector< T , A > >
+ {
+ typedef thrust_operations operations_type;
+ };
+} } }
+#endif // TBB_VERSION_MAJOR
+
+// specialization for thrust cuda vector
+#ifdef __CUDACC__
+#include <thrust/system/cuda/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ template< class T , class A >
+ struct operations_dispatcher< thrust::cuda::vector< T , A > >
+ {
+ typedef thrust_operations operations_type;
+ };
+} } }
+#endif // __CUDACC__
+
+#endif // THRUST_VERSION >= 100600
+
+
+#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_OPERATIONS_DISPATCHER_HPP_DEFINED
diff --git a/boost/numeric/odeint/external/thrust/thrust_resize.hpp b/boost/numeric/odeint/external/thrust/thrust_resize.hpp
index 8f1b181..2f7f7b2 100644
--- a/boost/numeric/odeint/external/thrust/thrust_resize.hpp
+++ b/boost/numeric/odeint/external/thrust/thrust_resize.hpp
@@ -6,7 +6,7 @@
Enable resizing for thrusts device and host_vector.
[end_description]
- Copyright 2010-2012 Mario Mulansky
+ Copyright 2010-2014 Mario Mulansky
Copyright 2010-2011 Karsten Ahnert
Distributed under the Boost Software License, Version 1.0.
@@ -18,102 +18,168 @@
#ifndef BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_RESIZE_HPP_INCLUDED
#define BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_RESIZE_HPP_INCLUDED
+#include <boost/range.hpp>
#include <thrust/device_vector.h>
#include <thrust/host_vector.h>
+#include <thrust/distance.h>
+#include <boost/numeric/odeint/util/resize.hpp>
+#include <boost/numeric/odeint/util/same_size.hpp>
#include <boost/numeric/odeint/util/copy.hpp>
namespace boost {
namespace numeric {
namespace odeint {
-template< class T >
-struct is_resizeable< thrust::device_vector< T > >
-{
- struct type : public boost::true_type { };
- const static bool value = type::value;
-};
-
-template< class T >
-struct same_size_impl< thrust::device_vector< T > , thrust::device_vector< T > >
-{
- static bool same_size( const thrust::device_vector< T > &x , const thrust::device_vector< T > &y )
- {
- return x.size() == y.size();
- }
-};
-
-template< class T >
-struct resize_impl< thrust::device_vector< T > , thrust::device_vector< T > >
-{
- static void resize( thrust::device_vector< T > &x , const thrust::device_vector< T > &y )
- {
- x.resize( y.size() );
- }
-};
-
-
-template< class T >
-struct is_resizeable< thrust::host_vector< T > >
-{
- struct type : public boost::true_type { };
- const static bool value = type::value;
-};
-
-template< class T >
-struct same_size_impl< thrust::host_vector< T > , thrust::host_vector< T > >
-{
- static bool same_size( const thrust::host_vector< T > &x , const thrust::host_vector< T > &y )
- {
- return x.size() == y.size();
- }
-};
-
-template< class T >
-struct resize_impl< thrust::host_vector< T > , thrust::host_vector< T > >
-{
- static void resize( thrust::host_vector< T > &x , const thrust::host_vector< T > &y )
- {
- x.resize( y.size() );
- }
-};
-
-
-
-template< class Container1, class Value >
-struct copy_impl< Container1 , thrust::device_vector< Value > >
-{
- static void copy( const Container1 &from , thrust::device_vector< Value > &to )
- {
- thrust::copy( boost::begin( from ) , boost::end( from ) , boost::begin( to ) );
- }
-};
-
-template< class Value , class Container2 >
-struct copy_impl< thrust::device_vector< Value > , Container2 >
-{
- static void copy( const thrust::device_vector< Value > &from , Container2 &to )
- {
- thrust::copy( boost::begin( from ) , boost::end( from ) , boost::begin( to ) );
- }
-};
-
-template< class Value >
-struct copy_impl< thrust::device_vector< Value > , thrust::device_vector< Value > >
-{
- static void copy( const thrust::device_vector< Value > &from , thrust::device_vector< Value > &to )
- {
- thrust::copy( boost::begin( from ) , boost::end( from ) , boost::begin( to ) );
- }
-};
-
-
+// some macros that define the necessary utilities
+
+#define ODEINT_THRUST_VECTOR_IS_RESIZEABLE( THRUST_VECTOR ) \
+template< class T , class A > \
+struct is_resizeable< THRUST_VECTOR<T,A> > \
+{ \
+ struct type : public boost::true_type { }; \
+ const static bool value = type::value; \
+}; \
+
+#define ODEINT_TRHUST_VECTOR_RESIZE_IMPL( THRUST_VECTOR ) \
+template< class T, class A > \
+struct resize_impl< THRUST_VECTOR<T,A> , THRUST_VECTOR<T,A> > \
+{ \
+ static void resize( THRUST_VECTOR<T,A> &x , \
+ const THRUST_VECTOR<T,A> &y ) \
+ { \
+ x.resize( y.size() ); \
+ } \
+}; \
+template< class T, class A, typename Range > \
+struct resize_impl< THRUST_VECTOR<T,A> , Range > \
+{ \
+ static void resize( THRUST_VECTOR<T,A> &x , \
+ const Range &y ) \
+ { \
+ x.resize( thrust::distance(boost::begin(y), \
+ boost::end(y))); \
+ } \
+}; \
+
+
+#define ODEINT_THRUST_SAME_SIZE_IMPL( THRUST_VECTOR ) \
+template< class T , class A > \
+struct same_size_impl< THRUST_VECTOR<T,A> , THRUST_VECTOR<T,A> > \
+{ \
+ static bool same_size( const THRUST_VECTOR<T,A> &x , \
+ const THRUST_VECTOR<T,A> &y ) \
+ { \
+ return x.size() == y.size(); \
+ } \
+}; \
+template< class T , class A, typename Range > \
+struct same_size_impl< THRUST_VECTOR<T,A> , Range > \
+{ \
+ static bool same_size( const THRUST_VECTOR<T,A> &x , \
+ const Range &y ) \
+ { \
+ return x.size() == thrust::distance(boost::begin(y), \
+ boost::end(y)); \
+ } \
+}; \
+
+
+#define ODEINT_THRUST_COPY_IMPL( THRUST_VECTOR ) \
+template< class Container1 , class T , class A > \
+struct copy_impl< Container1 , THRUST_VECTOR<T,A> > \
+{ \
+ static void copy( const Container1 &from , THRUST_VECTOR<T,A> &to ) \
+ { \
+ thrust::copy( boost::begin( from ) , boost::end( from ) , \
+ boost::begin( to ) ); \
+ } \
+}; \
+ \
+template< class T , class A , class Container2 > \
+struct copy_impl< THRUST_VECTOR<T,A> , Container2 > \
+{ \
+ static void copy( const THRUST_VECTOR<T,A> &from , Container2 &to ) \
+ { \
+ thrust::copy( boost::begin( from ) , boost::end( from ) , \
+ boost::begin( to ) ); \
+ } \
+}; \
+ \
+template< class T , class A > \
+struct copy_impl< THRUST_VECTOR<T,A> , THRUST_VECTOR<T,A> > \
+{ \
+ static void copy( const THRUST_VECTOR<T,A> &from , \
+ THRUST_VECTOR<T,A> &to ) \
+ { \
+ thrust::copy( boost::begin( from ) , boost::end( from ) , \
+ boost::begin( to ) ); \
+ } \
+}; \
+
+// add support for the standard thrust containers
+
+ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::device_vector )
+ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::device_vector )
+ODEINT_THRUST_SAME_SIZE_IMPL( thrust::device_vector )
+ODEINT_THRUST_COPY_IMPL( thrust::device_vector )
+
+
+ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::host_vector )
+ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::host_vector )
+ODEINT_THRUST_SAME_SIZE_IMPL( thrust::host_vector )
+ODEINT_THRUST_COPY_IMPL( thrust::host_vector )
} // odeint
} // numeric
} // boost
+// add support for thrust backend vectors, if available
+
+#include <thrust/version.h>
+
+#if THRUST_VERSION >= 100600
+
+#include <thrust/system/cpp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::cpp::vector )
+ ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::cpp::vector )
+ ODEINT_THRUST_SAME_SIZE_IMPL( thrust::cpp::vector )
+ ODEINT_THRUST_COPY_IMPL( thrust::cpp::vector )
+} } }
+
+#ifdef _OPENMP
+#include <thrust/system/omp/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::omp::vector )
+ ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::omp::vector )
+ ODEINT_THRUST_SAME_SIZE_IMPL( thrust::omp::vector )
+ ODEINT_THRUST_COPY_IMPL( thrust::omp::vector )
+} } }
+#endif // _OPENMP
+
+#ifdef TBB_VERSION_MAJOR
+#include <thrust/system/tbb/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::tbb::vector )
+ ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::tbb::vector )
+ ODEINT_THRUST_SAME_SIZE_IMPL( thrust::tbb::vector )
+ ODEINT_THRUST_COPY_IMPL( thrust::tbb::vector )
+} } }
+#endif // TBB_VERSION_MAJOR
+
+#ifdef __CUDACC__
+#include <thrust/system/cuda/vector.h>
+namespace boost { namespace numeric { namespace odeint {
+ ODEINT_THRUST_VECTOR_IS_RESIZEABLE( thrust::cuda::vector )
+ ODEINT_TRHUST_VECTOR_RESIZE_IMPL( thrust::cuda::vector )
+ ODEINT_THRUST_SAME_SIZE_IMPL( thrust::cuda::vector )
+ ODEINT_THRUST_COPY_IMPL( thrust::cuda::vector )
+} } }
+#endif // __CUDACC__
+
+#endif // THRUST_VERSION >= 100600
#endif // BOOST_NUMERIC_ODEINT_EXTERNAL_THRUST_THRUST_RESIZE_HPP_INCLUDED
diff --git a/boost/numeric/odeint/integrate/integrate.hpp b/boost/numeric/odeint/integrate/integrate.hpp
index 6cb96c0..446656b 100644
--- a/boost/numeric/odeint/integrate/integrate.hpp
+++ b/boost/numeric/odeint/integrate/integrate.hpp
@@ -48,6 +48,15 @@ integrate( System system , State &start_state , Time start_time , Time end_time
return integrate_adaptive( stepper_type() , system , start_state , start_time , end_time , dt , observer );
}
+template< class Value , class System , class State , class Time , class Observer >
+size_t
+integrate( System system , State &start_state , Time start_time , Time end_time , Time dt , Observer observer )
+{
+ typedef controlled_runge_kutta< runge_kutta_dopri5< State , Value , State , Time > > stepper_type;
+ return integrate_adaptive( stepper_type() , system , start_state , start_time , end_time , dt , observer );
+}
+
+
/*
@@ -59,6 +68,13 @@ size_t integrate( System system , State &start_state , Time start_time , Time en
return integrate( system , start_state , start_time , end_time , dt , null_observer() );
}
+template< class Value , class System , class State , class Time >
+size_t integrate( System system , State &start_state , Time start_time , Time end_time , Time dt )
+{
+ return integrate< Value >( system , start_state , start_time , end_time , dt , null_observer() );
+}
+
+
/**
* \fn integrate( System system , State &start_state , Time start_time , Time end_time , Time dt , Observer observer )
@@ -70,6 +86,9 @@ size_t integrate( System system , State &start_state , Time start_time , Time en
* integration with step size control, thus dt changes during the integration.
* This method uses standard error bounds of 1E-6.
* After each step, the observer is called.
+ *
+ * \attention A second version of this function template exists which explicitly
+ * expects the value type as template parameter, i.e. integrate< double >( sys , x , t0 , t1 , dt , obs );
*
* \param system The system function to solve, hence the r.h.s. of the
* ordinary differential equation.
@@ -92,6 +111,9 @@ size_t integrate( System system , State &start_state , Time start_time , Time en
* integration with step size control, thus dt changes during the integration.
* This method uses standard error bounds of 1E-6.
* No observer is called.
+ *
+ * \attention A second version of this function template exists which explicitly
+ * expects the value type as template parameter, i.e. integrate< double >( sys , x , t0 , t1 , dt );
*
* \param system The system function to solve, hence the r.h.s. of the
* ordinary differential equation.
diff --git a/boost/numeric/odeint/iterator/impl/adaptive_iterator_impl.hpp b/boost/numeric/odeint/iterator/impl/adaptive_iterator_impl.hpp
index 643b539..e5b4b3c 100644
--- a/boost/numeric/odeint/iterator/impl/adaptive_iterator_impl.hpp
+++ b/boost/numeric/odeint/iterator/impl/adaptive_iterator_impl.hpp
@@ -205,7 +205,7 @@ namespace odeint {
* \param sys The system function (ODE) to solve.
* \param s The initial state.
*/
- adaptive_iterator_impl( stepper_type stepper , system_type sys , state_type &s )
+ adaptive_iterator_impl( stepper_type stepper , system_type sys , state_type& /* s */ )
: base_type( stepper , sys ) { }
protected:
diff --git a/boost/numeric/odeint/iterator/impl/const_step_iterator_impl.hpp b/boost/numeric/odeint/iterator/impl/const_step_iterator_impl.hpp
index 449acfa..e23474c 100644
--- a/boost/numeric/odeint/iterator/impl/const_step_iterator_impl.hpp
+++ b/boost/numeric/odeint/iterator/impl/const_step_iterator_impl.hpp
@@ -77,7 +77,7 @@ namespace odeint {
* \param sys The system function (ODE) to solve.
* \param s The initial state. const_step_iterator stores a reference of s and changes its value during the iteration.
*/
- const_step_iterator_impl( stepper_type stepper , system_type sys , state_type &s )
+ const_step_iterator_impl( stepper_type stepper , system_type sys , state_type& /* s */ )
: base_type( stepper , sys ) { }
protected:
diff --git a/boost/numeric/odeint/iterator/integrate/detail/integrate_n_steps.hpp b/boost/numeric/odeint/iterator/integrate/detail/integrate_n_steps.hpp
index 1477990..f90f887 100644
--- a/boost/numeric/odeint/iterator/integrate/detail/integrate_n_steps.hpp
+++ b/boost/numeric/odeint/iterator/integrate/detail/integrate_n_steps.hpp
@@ -48,7 +48,7 @@ Time integrate_n_steps(
Observer observer , stepper_tag )
{
// ToDo: is there a better way to extract the final time?
- Time t;
+ Time t = start_time; // Assignment is only here to avoid warnings.
boost::for_each( make_n_step_time_range( stepper , system , start_state ,
start_time , dt , num_of_steps ) ,
obs_caller_time< Observer , Time >( t , observer ) );
@@ -91,7 +91,7 @@ Time integrate_n_steps(
Observer observer , dense_output_stepper_tag )
{
// ToDo: is there a better way to extract the final time?
- Time t;
+ Time t = start_time; // Assignment is only here to avoid warnings.
boost::for_each( make_n_step_time_range( stepper , system , start_state ,
start_time , dt , num_of_steps ) ,
obs_caller_time< Observer , Time >( t , observer ) );
diff --git a/boost/numeric/odeint/stepper/adams_bashforth.hpp b/boost/numeric/odeint/stepper/adams_bashforth.hpp
index 59ed1c9..5ff1e83 100644
--- a/boost/numeric/odeint/stepper/adams_bashforth.hpp
+++ b/boost/numeric/odeint/stepper/adams_bashforth.hpp
@@ -37,6 +37,7 @@
#include <boost/numeric/odeint/stepper/stepper_categories.hpp>
#include <boost/numeric/odeint/stepper/runge_kutta4.hpp>
+#include <boost/numeric/odeint/stepper/extrapolation_stepper.hpp>
#include <boost/numeric/odeint/stepper/base/algebra_stepper_base.hpp>
@@ -44,11 +45,28 @@
#include <boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp>
#include <boost/numeric/odeint/stepper/detail/rotating_buffer.hpp>
+#include <boost/mpl/arithmetic.hpp>
+#include <boost/mpl/min_max.hpp>
+#include <boost/mpl/equal_to.hpp>
+
+namespace mpl = boost::mpl;
+
namespace boost {
namespace numeric {
namespace odeint {
+ using mpl::int_;
+
+ /* if N >= 4, returns the smallest even number > N, otherwise returns 4 */
+ template < int N >
+ struct order_helper
+ : mpl::max< typename mpl::eval_if<
+ mpl::equal_to< mpl::modulus< int_< N >, int_< 2 > >,
+ int_< 0 > >,
+ int_< N >, int_< N + 1 > >::type,
+ int_< 4 > >::type
+ { };
template<
size_t Steps ,
@@ -59,7 +77,9 @@ class Time = Value ,
class Algebra = typename algebra_dispatcher< State >::algebra_type ,
class Operations = typename operations_dispatcher< State >::operations_type ,
class Resizer = initially_resizer ,
-class InitializingStepper = runge_kutta4< State , Value , Deriv , Time , Algebra , Operations, Resizer >
+class InitializingStepper = extrapolation_stepper< order_helper<Steps>::value,
+ State, Value, Deriv, Time,
+ Algebra, Operations, Resizer >
>
class adams_bashforth : public algebra_stepper_base< Algebra , Operations >
{
@@ -181,7 +201,8 @@ public :
{
if( i != 0 ) m_step_storage.rotate();
sys( x , m_step_storage[0].m_v , t );
- stepper.do_step( system , x , m_step_storage[0].m_v , t , dt );
+ stepper.do_step_dxdt_impl( system, x, m_step_storage[0].m_v, t,
+ dt );
t += dt;
}
m_steps_initialized = steps;
@@ -222,7 +243,8 @@ private:
{
if( m_steps_initialized != 0 ) m_step_storage.rotate();
sys( in , m_step_storage[0].m_v , t );
- m_initializing_stepper.do_step( system , in , m_step_storage[0].m_v , t , out , dt );
+ m_initializing_stepper.do_step_dxdt_impl(
+ system, in, m_step_storage[0].m_v, t, out, dt );
++m_steps_initialized;
}
else
diff --git a/boost/numeric/odeint/stepper/adams_bashforth_moulton.hpp b/boost/numeric/odeint/stepper/adams_bashforth_moulton.hpp
index 19b0239..2f7cc4c 100644
--- a/boost/numeric/odeint/stepper/adams_bashforth_moulton.hpp
+++ b/boost/numeric/odeint/stepper/adams_bashforth_moulton.hpp
@@ -49,7 +49,8 @@ class Deriv = State ,
class Time = Value ,
class Algebra = typename algebra_dispatcher< State >::algebra_type ,
class Operations = typename operations_dispatcher< State >::operations_type ,
-class Resizer = initially_resizer
+class Resizer = initially_resizer,
+class InitializingStepper = runge_kutta4< State , Value , Deriv , Time , Algebra , Operations, Resizer >
>
class adams_bashforth_moulton
{
@@ -71,12 +72,13 @@ public :
typedef Operations operations_type;
typedef Resizer resizer_type;
typedef stepper_tag stepper_category;
+ typedef InitializingStepper initializing_stepper_type;
static const size_t steps = Steps;
#ifndef DOXYGEN_SKIP
- typedef adams_bashforth< steps , state_type , value_type , deriv_type , time_type , algebra_type , operations_type , resizer_type > adams_bashforth_type;
+ typedef adams_bashforth< steps , state_type , value_type , deriv_type , time_type , algebra_type , operations_type , resizer_type, initializing_stepper_type > adams_bashforth_type;
typedef adams_moulton< steps , state_type , value_type , deriv_type , time_type , algebra_type , operations_type , resizer_type > adams_moulton_type;
- typedef adams_bashforth_moulton< steps , state_type , value_type , deriv_type , time_type , algebra_type , operations_type , resizer_type > stepper_type;
+ typedef adams_bashforth_moulton< steps , state_type , value_type , deriv_type , time_type , algebra_type , operations_type , resizer_type , initializing_stepper_type> stepper_type;
#endif //DOXYGEN_SKIP
typedef unsigned short order_type;
static const order_type order_value = steps;
@@ -158,7 +160,7 @@ private:
{
m_resizer.adjust_size( x , detail::bind( &stepper_type::template resize_impl< StateInOut > , detail::ref( *this ) , detail::_1 ) );
m_adams_bashforth.do_step( system , x , t , m_x.m_v , dt );
- m_adams_moulton.do_step( system , x , m_x.m_v , t , x , dt , m_adams_bashforth.step_storage() );
+ m_adams_moulton.do_step( system , x , m_x.m_v , t+dt , x , dt , m_adams_bashforth.step_storage() );
}
else
{
diff --git a/boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp b/boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp
index 3ba11a4..08009dc 100644
--- a/boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp
+++ b/boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp
@@ -152,6 +152,22 @@ public:
/*
+ * named Version 2: do_step_dxdt_impl( sys , in , dxdt , t , dt )
+ *
+ * this version is needed when this stepper is used for initializing
+ * multistep stepper like adams-bashforth. Hence we provide an explicitely
+ * named version that is not disabled. Meant for internal use only.
+ */
+ template < class System, class StateInOut, class DerivIn >
+ void do_step_dxdt_impl( System system, StateInOut &x, const DerivIn &dxdt,
+ time_type t, time_type dt )
+ {
+ this->stepper().do_step_impl( system , x , dxdt , t , x , dt );
+ }
+
+
+
+ /*
* Version 3 : do_step( sys , in , t , out , dt )
*
* this version does not solve the forwarding problem, boost.range can not be used
@@ -181,10 +197,21 @@ public:
{
this->stepper().do_step_impl( system , in , dxdt , t , out , dt );
}
-
-
-
-
+
+ /*
+ * named Version 4: do_step_dxdt_impl( sys , in , dxdt , t , out, dt )
+ *
+ * this version is needed when this stepper is used for initializing
+ * multistep stepper like adams-bashforth. Hence we provide an explicitely
+ * named version that is not disabled. Meant for internal use only.
+ */
+ template < class System, class StateIn, class DerivIn, class StateOut >
+ void do_step_dxdt_impl( System system, const StateIn &in,
+ const DerivIn &dxdt, time_type t, StateOut &out,
+ time_type dt )
+ {
+ this->stepper().do_step_impl( system , in , dxdt , t , out , dt );
+ }
/*
* Version 5 :do_step( sys , x , t , dt , xerr )
diff --git a/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp b/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp
index a055c7f..b1d751a 100644
--- a/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp
+++ b/boost/numeric/odeint/stepper/base/explicit_error_stepper_fsal_base.hpp
@@ -151,11 +151,27 @@ public:
/*
+ * named Version 2: do_step_dxdt_impl( sys , in , dxdt , t , dt )
+ *
+ * this version is needed when this stepper is used for initializing
+ * multistep stepper like adams-bashforth. Hence we provide an explicitely
+ * named version that is not disabled. Meant for internal use only.
+ */
+ template< class System , class StateInOut , class DerivInOut >
+ void do_step_dxdt_impl( System system , StateInOut &x , DerivInOut &dxdt , time_type t , time_type dt )
+ {
+ m_first_call = true;
+ this->stepper().do_step_impl( system , x , dxdt , t , x , dxdt , dt );
+ }
+
+ /*
* version 3 : do_step( sys , in , t , out , dt )
*
- * this version does not solve the forwarding problem, boost.range can not be used
+ * this version does not solve the forwarding problem, boost.range can not
+ * be used.
*
- * the disable is needed to avoid ambiguous overloads if state_type = time_type
+ * the disable is needed to avoid ambiguous overloads if
+ * state_type = time_type
*/
template< class System , class StateIn , class StateOut >
typename boost::disable_if< boost::is_same< StateIn , time_type > , void >::type
@@ -174,12 +190,14 @@ public:
*
* this version does not solve the forwarding problem, boost.range can not be used
*/
- template< class System , class StateIn , class DerivIn , class StateOut , class DerivOut >
- void do_step( System system , const StateIn &in , const DerivIn &dxdt_in , time_type t ,
- StateOut &out , DerivOut &dxdt_out , time_type dt )
+ template< class System, class StateIn, class DerivIn, class StateOut,
+ class DerivOut >
+ void do_step( System system, const StateIn &in, const DerivIn &dxdt_in,
+ time_type t, StateOut &out, DerivOut &dxdt_out, time_type dt )
{
m_first_call = true;
- this->stepper().do_step_impl( system , in , dxdt_in , t , out , dxdt_out , dt );
+ this->stepper().do_step_impl( system, in, dxdt_in, t, out, dxdt_out,
+ dt );
}
diff --git a/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp b/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp
index 40aab80..d81c8c7 100644
--- a/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp
+++ b/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp
@@ -139,6 +139,21 @@ public:
/*
+ * named Version 2: do_step_dxdt_impl( sys , in , dxdt , t , dt )
+ *
+ * this version is needed when this stepper is used for initializing
+ * multistep stepper like adams-bashforth. Hence we provide an explicitely
+ * named version that is not disabled. Meant for internal use only.
+ */
+ template < class System, class StateInOut, class DerivIn >
+ void do_step_dxdt_impl( System system, StateInOut &x, const DerivIn &dxdt,
+ time_type t, time_type dt )
+ {
+ this->stepper().do_step_impl( system , x , dxdt , t , x , dt );
+ }
+
+
+ /*
* Version 3 : do_step( sys , in , t , out , dt )
*
* this version does not solve the forwarding problem, boost.range can not be used
@@ -164,6 +179,22 @@ public:
this->stepper().do_step_impl( system , in , dxdt , t , out , dt );
}
+
+ /*
+ * named Version 4: do_step_dxdt_impl( sys , in , dxdt , t , out, dt )
+ *
+ * this version is needed when this stepper is used for initializing
+ * multistep stepper like adams-bashforth. Hence we provide an explicitely
+ * named version. Meant for internal use only.
+ */
+ template < class System, class StateIn, class DerivIn, class StateOut >
+ void do_step_dxdt_impl( System system, const StateIn &in,
+ const DerivIn &dxdt, time_type t, StateOut &out,
+ time_type dt )
+ {
+ this->stepper().do_step_impl( system , in , dxdt , t , out , dt );
+ }
+
template< class StateIn >
void adjust_size( const StateIn &x )
{
diff --git a/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp b/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp
index ea3523e..eb09aef 100644
--- a/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp
+++ b/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp
@@ -182,7 +182,7 @@ private:
// stepper for systems with function for dq/dt = f(p) and dp/dt = -f(q)
template< class System , class StateIn , class StateOut >
- void do_step_impl( System system , const StateIn &in , time_type t , StateOut &out , time_type dt , boost::mpl::true_ )
+ void do_step_impl( System system , const StateIn &in , time_type /* t */ , StateOut &out , time_type dt , boost::mpl::true_ )
{
typedef typename odeint::unwrap_reference< System >::type system_type;
typedef typename odeint::unwrap_reference< typename system_type::first_type >::type coor_deriv_func_type;
diff --git a/boost/numeric/odeint/stepper/bulirsch_stoer.hpp b/boost/numeric/odeint/stepper/bulirsch_stoer.hpp
index b604f86..e71008a 100644
--- a/boost/numeric/odeint/stepper/bulirsch_stoer.hpp
+++ b/boost/numeric/odeint/stepper/bulirsch_stoer.hpp
@@ -194,7 +194,6 @@ public:
static const value_type val1( 1.0 );
- typename odeint::unwrap_reference< System >::type &sys = system;
if( m_resizer.adjust_size( in , detail::bind( &controlled_error_bs_type::template resize_impl< StateIn > , detail::ref( *this ) , detail::_1 ) ) )
{
reset(); // system resized -> reset
@@ -219,12 +218,12 @@ public:
m_midpoint.set_steps( m_interval_sequence[k] );
if( k == 0 )
{
- m_midpoint.do_step( sys , in , dxdt , t , out , dt );
+ m_midpoint.do_step( system , in , dxdt , t , out , dt );
/* the first step, nothing more to do */
}
else
{
- m_midpoint.do_step( sys , in , dxdt , t , m_table[k-1].m_v , dt );
+ m_midpoint.do_step( system , in , dxdt , t , m_table[k-1].m_v , dt );
extrapolate( k , m_table , m_coeff , out );
// get error estimate
m_algebra.for_each3( m_err.m_v , out , m_table[0].m_v ,
@@ -341,7 +340,7 @@ public:
resize_m_dxdt( x );
resize_m_xnew( x );
resize_impl( x );
- m_midpoint.adjust_size();
+ m_midpoint.adjust_size( x );
}
diff --git a/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp b/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp
index c352c96..0036d1d 100644
--- a/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp
+++ b/boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp
@@ -155,8 +155,6 @@ public:
static const value_type val1( 1.0 );
- typename odeint::unwrap_reference< System >::type &sys = system;
-
bool reject( true );
time_vector h_opt( m_k_max+1 );
@@ -261,6 +259,7 @@ public:
{
//calculate dxdt for next step and dense output
+ typename odeint::unwrap_reference< System >::type &sys = system;
sys( out , dxdt_new , t+dt );
//prepare dense output
@@ -367,7 +366,7 @@ public:
void adjust_size( const StateIn &x )
{
resize_impl( x );
- m_midpoint.adjust_size();
+ m_midpoint.adjust_size( x );
}
diff --git a/boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp b/boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp
index 168bc69..cc279ab 100644
--- a/boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp
+++ b/boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp
@@ -117,7 +117,7 @@ struct adams_bashforth_call_algebra< 7 , Algebra , Operations >
{
//BOOST_ASSERT( false ); // not implemented
typedef typename Coefficients::value_type value_type;
- Algebra::for_each9( out , in , steps[0].m_v , steps[1].m_v , steps[2].m_v , steps[3].m_v , steps[4].m_v , steps[5].m_v , steps[6].m_v ,
+ algebra.for_each9( out , in , steps[0].m_v , steps[1].m_v , steps[2].m_v , steps[3].m_v , steps[4].m_v , steps[5].m_v , steps[6].m_v ,
typename Operations::template scale_sum8< value_type , Time , Time , Time , Time , Time , Time >(
1.0 , dt * coef[0] , dt * coef[1] , dt * coef[2] , dt * coef[3] , dt * coef[4] , dt * coef[5] , dt * coef[6] ) );
}
@@ -132,7 +132,7 @@ struct adams_bashforth_call_algebra< 8 , Algebra , Operations >
{
//BOOST_ASSERT( false ); // not implemented
typedef typename Coefficients::value_type value_type;
- Algebra::for_each10( out , in , steps[0].m_v , steps[1].m_v , steps[2].m_v , steps[3].m_v , steps[4].m_v , steps[5].m_v , steps[6].m_v , steps[7].m_v ,
+ algebra.for_each10( out , in , steps[0].m_v , steps[1].m_v , steps[2].m_v , steps[3].m_v , steps[4].m_v , steps[5].m_v , steps[6].m_v , steps[7].m_v ,
typename Operations::template scale_sum9< value_type , Time , Time , Time , Time , Time , Time , Time >(
1.0 , dt * coef[0] , dt * coef[1] , dt * coef[2] , dt * coef[3] , dt * coef[4] , dt * coef[5] , dt * coef[6] , dt * coef[7] ) );
}
diff --git a/boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp b/boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp
index f8a1505..b6f5f2a 100644
--- a/boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp
+++ b/boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp
@@ -32,7 +32,7 @@ template< class Algebra , class Operations >
struct adams_moulton_call_algebra< 1 , Algebra , Operations >
{
template< class StateIn , class StateOut , class DerivIn , class StepStorage , class Coefficients , class Time >
- void operator()( Algebra &algebra , const StateIn &in , StateOut &out , const DerivIn &dxdt , const StepStorage &steps , const Coefficients &coef , Time dt ) const
+ void operator()( Algebra &algebra , const StateIn &in , StateOut &out , const DerivIn &dxdt , const StepStorage& /* steps */ , const Coefficients &coef , Time dt ) const
{
typedef typename Coefficients::value_type value_type;
algebra.for_each3( out , in , dxdt , typename Operations::template scale_sum2< value_type , Time >( 1.0 , dt * coef[0] ) );
diff --git a/boost/numeric/odeint/stepper/euler.hpp b/boost/numeric/odeint/stepper/euler.hpp
index 443f942..1c7c126 100644
--- a/boost/numeric/odeint/stepper/euler.hpp
+++ b/boost/numeric/odeint/stepper/euler.hpp
@@ -76,7 +76,7 @@ public :
{ }
template< class System , class StateIn , class DerivIn , class StateOut >
- void do_step_impl( System system , const StateIn &in , const DerivIn &dxdt , time_type t , StateOut &out , time_type dt )
+ void do_step_impl( System /* system */ , const StateIn &in , const DerivIn &dxdt , time_type /* t */ , StateOut &out , time_type dt )
{
stepper_base_type::m_algebra.for_each3( out , in , dxdt ,
typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , dt ) );
@@ -84,7 +84,7 @@ public :
}
template< class StateOut , class StateIn1 , class StateIn2 >
- void calc_state( StateOut &x , time_type t , const StateIn1 &old_state , time_type t_old , const StateIn2 &current_state , time_type t_new ) const
+ void calc_state( StateOut &x , time_type t , const StateIn1 &old_state , time_type t_old , const StateIn2 & /*current_state*/ , time_type /* t_new */ ) const
{
const time_type delta = t - t_old;
stepper_base_type::m_algebra.for_each3( x , old_state , stepper_base_type::m_dxdt.m_v ,
diff --git a/boost/numeric/odeint/stepper/extrapolation_stepper.hpp b/boost/numeric/odeint/stepper/extrapolation_stepper.hpp
new file mode 100644
index 0000000..3b688e9
--- /dev/null
+++ b/boost/numeric/odeint/stepper/extrapolation_stepper.hpp
@@ -0,0 +1,293 @@
+/*
+ [auto_generated]
+ boost/numeric/odeint/stepper/extrapolation_stepper.hpp
+
+ [begin_description]
+ extrapolation stepper
+ [end_description]
+
+ Copyright 2009-2015 Mario Mulansky
+
+ 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)
+*/
+
+#ifndef BOOST_NUMERIC_ODEINT_STEPPER_EXTRAPOLATION_STEPPER_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_STEPPER_EXTRAPOLATION_STEPPER_HPP_INCLUDED
+
+#include <iostream>
+
+#include <algorithm>
+
+#include <boost/config.hpp> // for min/max guidelines
+#include <boost/static_assert.hpp>
+
+#include <boost/numeric/odeint/util/bind.hpp>
+#include <boost/numeric/odeint/util/unwrap_reference.hpp>
+
+#include <boost/numeric/odeint/stepper/base/explicit_error_stepper_base.hpp>
+#include <boost/numeric/odeint/stepper/modified_midpoint.hpp>
+#include <boost/numeric/odeint/stepper/controlled_step_result.hpp>
+#include <boost/numeric/odeint/algebra/range_algebra.hpp>
+#include <boost/numeric/odeint/algebra/default_operations.hpp>
+#include <boost/numeric/odeint/algebra/algebra_dispatcher.hpp>
+#include <boost/numeric/odeint/algebra/operations_dispatcher.hpp>
+
+#include <boost/numeric/odeint/util/state_wrapper.hpp>
+#include <boost/numeric/odeint/util/is_resizeable.hpp>
+#include <boost/numeric/odeint/util/resizer.hpp>
+#include <boost/numeric/odeint/util/unit_helper.hpp>
+#include <boost/numeric/odeint/util/detail/less_with_sign.hpp>
+
+namespace boost
+{
+namespace numeric
+{
+namespace odeint
+{
+
+template < unsigned short Order, class State, class Value = double,
+ class Deriv = State, class Time = Value,
+ class Algebra = typename algebra_dispatcher< State >::algebra_type,
+ class Operations =
+ typename operations_dispatcher< State >::operations_type,
+ class Resizer = initially_resizer >
+#ifndef DOXYGEN_SKIP
+class extrapolation_stepper
+ : public explicit_error_stepper_base<
+ extrapolation_stepper< Order, State, Value, Deriv, Time, Algebra,
+ Operations, Resizer >,
+ Order, Order, Order - 2, State, Value, Deriv, Time, Algebra,
+ Operations, Resizer >
+#else
+class extrapolation_stepper : public explicit_error_stepper_base
+#endif
+{
+
+ private:
+ // check for Order being odd
+ BOOST_STATIC_ASSERT_MSG(
+ ( ( Order % 2 ) == 0 ) && ( Order > 2 ),
+ "extrapolation_stepper requires even Order larger than 2" );
+
+ public:
+#ifndef DOXYGEN_SKIP
+ typedef explicit_error_stepper_base<
+ extrapolation_stepper< Order, State, Value, Deriv, Time, Algebra,
+ Operations, Resizer >,
+ Order, Order, Order - 2, State, Value, Deriv, Time, Algebra, Operations,
+ Resizer > stepper_base_type;
+#else
+ typedef explicit_error_stepper_base< extrapolation_stepper< ... >, ... >
+ stepper_base_type;
+#endif
+
+ typedef typename stepper_base_type::state_type state_type;
+ typedef typename stepper_base_type::value_type value_type;
+ typedef typename stepper_base_type::deriv_type deriv_type;
+ typedef typename stepper_base_type::time_type time_type;
+ typedef typename stepper_base_type::algebra_type algebra_type;
+ typedef typename stepper_base_type::operations_type operations_type;
+ typedef typename stepper_base_type::resizer_type resizer_type;
+
+#ifndef DOXYGEN_SKIP
+ typedef typename stepper_base_type::stepper_type stepper_type;
+ typedef typename stepper_base_type::wrapped_state_type wrapped_state_type;
+ typedef typename stepper_base_type::wrapped_deriv_type wrapped_deriv_type;
+
+ typedef std::vector< value_type > value_vector;
+ typedef std::vector< value_vector > value_matrix;
+ typedef std::vector< size_t > int_vector;
+ typedef std::vector< wrapped_state_type > state_table_type;
+ typedef modified_midpoint< state_type, value_type, deriv_type, time_type,
+ algebra_type, operations_type,
+ resizer_type > midpoint_stepper_type;
+
+#endif // DOXYGEN_SKIP
+
+ typedef unsigned short order_type;
+ static const order_type order_value = stepper_base_type::order_value;
+ static const order_type stepper_order_value =
+ stepper_base_type::stepper_order_value;
+ static const order_type error_order_value =
+ stepper_base_type::error_order_value;
+
+ const static size_t m_k_max = ( order_value - 2 ) / 2;
+
+ extrapolation_stepper( const algebra_type &algebra = algebra_type() )
+ : stepper_base_type( algebra ), m_interval_sequence( m_k_max + 1 ),
+ m_coeff( m_k_max + 1 ), m_table( m_k_max )
+ {
+ for ( unsigned short i = 0; i < m_k_max + 1; i++ )
+ {
+ m_interval_sequence[i] = 2 * ( i + 1 );
+ m_coeff[i].resize( i );
+ for ( size_t k = 0; k < i; ++k )
+ {
+ const value_type r =
+ static_cast< value_type >( m_interval_sequence[i] ) /
+ static_cast< value_type >( m_interval_sequence[k] );
+ m_coeff[i][k] =
+ static_cast< value_type >( 1 ) /
+ ( r * r - static_cast< value_type >(
+ 1 ) ); // coefficients for extrapolation
+ }
+ }
+ }
+
+ template < class System, class StateIn, class DerivIn, class StateOut,
+ class Err >
+ void do_step_impl( System system, const StateIn &in, const DerivIn &dxdt,
+ time_type t, StateOut &out, time_type dt, Err &xerr )
+ {
+ // std::cout << "dt: " << dt << std::endl;
+ // normal step
+ do_step_impl( system, in, dxdt, t, out, dt );
+
+ static const value_type val1( 1.0 );
+ // additionally, perform the error calculation
+ stepper_base_type::m_algebra.for_each3(
+ xerr, out, m_table[0].m_v,
+ typename operations_type::template scale_sum2<
+ value_type, value_type >( val1, -val1 ) );
+ }
+
+ template < class System, class StateInOut, class DerivIn, class Err >
+ void do_step_impl_io( System system, StateInOut &inout, const DerivIn &dxdt,
+ time_type t, time_type dt, Err &xerr )
+ {
+ // normal step
+ do_step_impl_io( system, inout, dxdt, t, dt );
+
+ static const value_type val1( 1.0 );
+ // additionally, perform the error calculation
+ stepper_base_type::m_algebra.for_each3(
+ xerr, inout, m_table[0].m_v,
+ typename operations_type::template scale_sum2<
+ value_type, value_type >( val1, -val1 ) );
+ }
+
+ template < class System, class StateIn, class DerivIn, class StateOut >
+ void do_step_impl( System system, const StateIn &in, const DerivIn &dxdt,
+ time_type t, StateOut &out, time_type dt )
+ {
+ m_resizer.adjust_size(
+ in, detail::bind( &stepper_type::template resize_impl< StateIn >,
+ detail::ref( *this ), detail::_1 ) );
+ size_t k = 0;
+ m_midpoint.set_steps( m_interval_sequence[k] );
+ m_midpoint.do_step( system, in, dxdt, t, out, dt );
+ for ( k = 1; k <= m_k_max; ++k )
+ {
+ m_midpoint.set_steps( m_interval_sequence[k] );
+ m_midpoint.do_step( system, in, dxdt, t, m_table[k - 1].m_v, dt );
+ extrapolate( k, m_table, m_coeff, out );
+ }
+ }
+
+ template < class System, class StateInOut, class DerivIn >
+ void do_step_impl_io( System system, StateInOut &inout, const DerivIn &dxdt,
+ time_type t, time_type dt )
+ {
+ // special care for inout
+ m_xout_resizer.adjust_size(
+ inout,
+ detail::bind( &stepper_type::template resize_m_xout< StateInOut >,
+ detail::ref( *this ), detail::_1 ) );
+ do_step_impl( system, inout, dxdt, t, m_xout.m_v, dt );
+ boost::numeric::odeint::copy( m_xout.m_v, inout );
+ }
+
+ template < class System, class StateInOut, class DerivIn >
+ void do_step_dxdt_impl( System system, StateInOut &x, const DerivIn &dxdt,
+ time_type t, time_type dt )
+ {
+ do_step_impl_io( system , x , dxdt , t , dt );
+ }
+
+ template < class System, class StateIn, class DerivIn, class StateOut >
+ void do_step_dxdt_impl( System system, const StateIn &in,
+ const DerivIn &dxdt, time_type t, StateOut &out,
+ time_type dt )
+ {
+ do_step_impl( system , in , dxdt , t , out , dt );
+ }
+
+
+ template < class StateIn > void adjust_size( const StateIn &x )
+ {
+ resize_impl( x );
+ m_midpoint.adjust_size( x );
+ }
+
+ private:
+ template < class StateIn > bool resize_impl( const StateIn &x )
+ {
+ bool resized( false );
+ for ( size_t i = 0; i < m_k_max; ++i )
+ resized |= adjust_size_by_resizeability(
+ m_table[i], x, typename is_resizeable< state_type >::type() );
+ return resized;
+ }
+
+ template < class StateIn > bool resize_m_xout( const StateIn &x )
+ {
+ return adjust_size_by_resizeability(
+ m_xout, x, typename is_resizeable< state_type >::type() );
+ }
+
+ template < class StateInOut >
+ void extrapolate( size_t k, state_table_type &table,
+ const value_matrix &coeff, StateInOut &xest )
+ /* polynomial extrapolation, see http://www.nr.com/webnotes/nr3web21.pdf
+ uses the obtained intermediate results to extrapolate to dt->0
+ */
+ {
+ static const value_type val1 = static_cast< value_type >( 1.0 );
+
+ for ( int j = k - 1; j > 0; --j )
+ {
+ stepper_base_type::m_algebra.for_each3(
+ table[j - 1].m_v, table[j].m_v, table[j - 1].m_v,
+ typename operations_type::template scale_sum2<
+ value_type, value_type >( val1 + coeff[k][j],
+ -coeff[k][j] ) );
+ }
+ stepper_base_type::m_algebra.for_each3(
+ xest, table[0].m_v, xest,
+ typename operations_type::template scale_sum2<
+ value_type, value_type >( val1 + coeff[k][0], -coeff[k][0] ) );
+ }
+
+ private:
+ midpoint_stepper_type m_midpoint;
+
+ resizer_type m_resizer;
+ resizer_type m_xout_resizer;
+
+ int_vector m_interval_sequence; // stores the successive interval counts
+ value_matrix m_coeff;
+
+ wrapped_state_type m_xout;
+ state_table_type m_table; // sequence of states for extrapolation
+};
+
+/******** DOXYGEN *******/
+
+/**
+ * \class extrapolation_stepper
+ * \brief Extrapolation stepper with configurable order, and error estimation.
+ *
+ * The extrapolation stepper is a stepper with error estimation and configurable
+ * order. The order is given as template parameter and needs to be an _odd_
+ * number. The stepper is based on several executions of the modified midpoint
+ * method and a Richardson extrapolation. This is essentially the same technique
+ * as for bulirsch_stoer, but without the variable order.
+ *
+ * \note The Order parameter has to be an even number greater 2.
+ */
+}
+}
+}
+#endif
diff --git a/boost/numeric/odeint/stepper/implicit_euler.hpp b/boost/numeric/odeint/stepper/implicit_euler.hpp
index f301988..e1c6416 100644
--- a/boost/numeric/odeint/stepper/implicit_euler.hpp
+++ b/boost/numeric/odeint/stepper/implicit_euler.hpp
@@ -144,7 +144,7 @@ private:
void solve( state_type &x , matrix_type &m )
{
int res = boost::numeric::ublas::lu_factorize( m , m_pm.m_v );
- if( res != 0 ) exit(0);
+ if( res != 0 ) std::exit(0);
boost::numeric::ublas::lu_substitute( m , m_pm.m_v , x );
}
diff --git a/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp b/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp
index 4b3824f..df4e6c4 100644
--- a/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp
+++ b/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp
@@ -66,6 +66,7 @@ public:
{
BOOST_USING_STD_MAX();
using std::abs;
+ using std::sqrt;
const size_t n = x.size();
value_type err = 0.0 , sk = 0.0;
diff --git a/boost/numeric/odeint/stepper/runge_kutta4_classic.hpp b/boost/numeric/odeint/stepper/runge_kutta4_classic.hpp
index 2b8e70b..32bda0b 100644
--- a/boost/numeric/odeint/stepper/runge_kutta4_classic.hpp
+++ b/boost/numeric/odeint/stepper/runge_kutta4_classic.hpp
@@ -121,11 +121,20 @@ public :
// dt * m_dxh = k4
sys( m_x_tmp.m_v , m_dxh.m_v , t + dt );
+
//x += dt/6 * ( m_dxdt + m_dxt + val2*m_dxm )
time_type dt6 = dt / static_cast< value_type >( 6 );
time_type dt3 = dt / static_cast< value_type >( 3 );
stepper_base_type::m_algebra.for_each6( out , in , dxdt , m_dxt.m_v , m_dxm.m_v , m_dxh.m_v ,
- typename operations_type::template scale_sum5< value_type , time_type , time_type , time_type , time_type >( 1.0 , dt6 , dt3 , dt3 , dt6 ) );
+ typename operations_type::template scale_sum5< value_type , time_type , time_type , time_type , time_type >( 1.0 , dt6 , dt3 , dt3 , dt6 ) );
+
+ // x += dt/6 * m_dxdt + dt/3 * m_dxt )
+ // stepper_base_type::m_algebra.for_each4( out , in , dxdt , m_dxt.m_v ,
+ // typename operations_type::template scale_sum3< value_type , time_type , time_type >( 1.0 , dt6 , dt3 ) );
+ // // x += dt/3 * m_dxm + dt/6 * m_dxh )
+ // stepper_base_type::m_algebra.for_each4( out , out , m_dxm.m_v , m_dxh.m_v ,
+ // typename operations_type::template scale_sum3< value_type , time_type , time_type >( 1.0 , dt3 , dt6 ) );
+
}
template< class StateType >
diff --git a/boost/numeric/odeint/stepper/velocity_verlet.hpp b/boost/numeric/odeint/stepper/velocity_verlet.hpp
index 24f3c0e..3a20fc2 100644
--- a/boost/numeric/odeint/stepper/velocity_verlet.hpp
+++ b/boost/numeric/odeint/stepper/velocity_verlet.hpp
@@ -114,7 +114,7 @@ public:
algebra_stepper_base_type::m_algebra.for_each4(
qout , qin , pin , ain ,
- typename operations_type::template scale_sum3< value_type , time_type , time_square_type >( one , one_half * dt , one * dt * dt ) );
+ typename operations_type::template scale_sum3< value_type , time_type , time_square_type >( one , one * dt , one_half * dt * dt ) );
typename odeint::unwrap_reference< System >::type & sys = system;
diff --git a/boost/numeric/odeint/util/bind.hpp b/boost/numeric/odeint/util/bind.hpp
index 966fd03..1201afa 100644
--- a/boost/numeric/odeint/util/bind.hpp
+++ b/boost/numeric/odeint/util/bind.hpp
@@ -41,9 +41,12 @@ using namespace ::std::placeholders;
#else
+// unnamed namespace to avoid multiple declarations (#138)
+namespace {
using ::boost::bind;
boost::arg<1> _1;
boost::arg<2> _2;
+}
// using ::boost::bind;
// using ::_1;
// using ::_2;
diff --git a/boost/numeric/odeint/util/resize.hpp b/boost/numeric/odeint/util/resize.hpp
index db2ce20..75e9eff 100644
--- a/boost/numeric/odeint/util/resize.hpp
+++ b/boost/numeric/odeint/util/resize.hpp
@@ -33,8 +33,8 @@
namespace boost {
namespace numeric {
namespace odeint {
-
-
+
+
template< class StateOut , class StateIn , class Enabler = void >
struct resize_impl_sfinae
{
diff --git a/boost/numeric/odeint/util/same_instance.hpp b/boost/numeric/odeint/util/same_instance.hpp
index dc55db3..a889ee1 100644
--- a/boost/numeric/odeint/util/same_instance.hpp
+++ b/boost/numeric/odeint/util/same_instance.hpp
@@ -25,7 +25,7 @@ namespace odeint {
template< class T1 , class T2 , class Enabler=void >
struct same_instance_impl
{
- static bool same_instance( const T1 &x1 , const T2 &x2 )
+ static bool same_instance( const T1& /* x1 */ , const T2& /* x2 */ )
{
return false;
}
diff --git a/boost/numeric/odeint/util/unwrap_reference.hpp b/boost/numeric/odeint/util/unwrap_reference.hpp
index d622a5c..bc7d423 100644
--- a/boost/numeric/odeint/util/unwrap_reference.hpp
+++ b/boost/numeric/odeint/util/unwrap_reference.hpp
@@ -31,9 +31,9 @@
namespace boost {
#if BOOST_NUMERIC_ODEINT_CXX11
-template<typename T> class reference_wrapper;
+template<typename T> struct reference_wrapper;
-template<typename T> class unwrap_reference;
+template<typename T> struct unwrap_reference;
#endif
namespace numeric {
@@ -43,24 +43,21 @@ namespace odeint {
#if BOOST_NUMERIC_ODEINT_CXX11
template<typename T>
-class unwrap_reference
+struct unwrap_reference
{
-public:
typedef typename std::remove_reference<T>::type type;
};
template<typename T>
-class unwrap_reference< std::reference_wrapper<T> >
+struct unwrap_reference< std::reference_wrapper<T> >
{
-public:
typedef typename std::remove_reference<T>::type type;
};
template<typename T>
-class unwrap_reference< boost::reference_wrapper<T> >
+struct unwrap_reference< boost::reference_wrapper<T> >
{
-public:
- typedef typename boost::unwrap_reference<T>::type type;
+ typedef typename boost::unwrap_reference<T>::type type;
};
#else
diff --git a/boost/numeric/ublas/detail/matrix_assign.hpp b/boost/numeric/ublas/detail/matrix_assign.hpp
index 5dea592..be172dd 100644
--- a/boost/numeric/ublas/detail/matrix_assign.hpp
+++ b/boost/numeric/ublas/detail/matrix_assign.hpp
@@ -1259,20 +1259,20 @@ namespace detail {
void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, row_major_tag) {
typedef F<typename M::iterator2::reference, typename E::reference> functor_type;
// R unnecessary, make_conformant not required
- typedef typename M::size_type size_type;
+ //typedef typename M::size_type size_type; // gcc is complaining that this is not used, although this is not right
typedef typename M::difference_type difference_type;
typename M::iterator1 it1 (m.begin1 ());
typename E::iterator1 it1e (e ().begin1 ());
- difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (e ().end1 () - it1e)));
+ difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), typename M::size_type (e ().end1 () - it1e)));
while (-- size1 >= 0) {
#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
typename M::iterator2 it2 (it1.begin ());
typename E::iterator2 it2e (it1e.begin ());
- difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (it1e.end () - it2e)));
+ difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), typename M::size_type (it1e.end () - it2e)));
#else
typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
- difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (end (it1e, iterator1_tag ()) - it2e)));
+ difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), typename M::size_type (end (it1e, iterator1_tag ()) - it2e)));
#endif
while (-- size2 >= 0)
functor_type::apply (*it2, *it2e), ++ it2, ++ it2e;
@@ -1285,20 +1285,20 @@ namespace detail {
void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, column_major_tag) {
typedef F<typename M::iterator1::reference, typename E::reference> functor_type;
// R unnecessary, make_conformant not required
- typedef typename M::size_type size_type;
+ // typedef typename M::size_type size_type; // gcc is complaining that this is not used, although this is not right
typedef typename M::difference_type difference_type;
typename M::iterator2 it2 (m.begin2 ());
typename E::iterator2 it2e (e ().begin2 ());
- difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (e ().end2 () - it2e)));
+ difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), typename M::size_type (e ().end2 () - it2e)));
while (-- size2 >= 0) {
#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
typename M::iterator1 it1 (it2.begin ());
typename E::iterator1 it1e (it2e.begin ());
- difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (it2e.end () - it1e)));
+ difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), typename M::size_type (it2e.end () - it1e)));
#else
typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
- difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (end (it2e, iterator2_tag ()) - it1e)));
+ difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), typename M::size_type (end (it2e, iterator2_tag ()) - it1e)));
#endif
while (-- size1 >= 0)
functor_type::apply (*it1, *it1e), ++ it1, ++ it1e;
diff --git a/boost/numeric/ublas/experimental/sparse_view.hpp b/boost/numeric/ublas/experimental/sparse_view.hpp
index 52ad509..3a16410 100644
--- a/boost/numeric/ublas/experimental/sparse_view.hpp
+++ b/boost/numeric/ublas/experimental/sparse_view.hpp
@@ -19,6 +19,7 @@
#include <boost/next_prior.hpp>
#include <boost/type_traits/remove_cv.hpp>
+#include <boost/numeric/ublas/storage.hpp>
namespace boost { namespace numeric { namespace ublas {
diff --git a/boost/numeric/ublas/functional.hpp b/boost/numeric/ublas/functional.hpp
index e6f9273..c7c9261 100644
--- a/boost/numeric/ublas/functional.hpp
+++ b/boost/numeric/ublas/functional.hpp
@@ -15,6 +15,8 @@
#include <functional>
+#include <boost/core/ignore_unused.hpp>
+
#include <boost/numeric/ublas/traits.hpp>
#ifdef BOOST_UBLAS_USE_DUFF_DEVICE
#include <boost/numeric/ublas/detail/duff.hpp>
@@ -1161,6 +1163,7 @@ namespace boost { namespace numeric { namespace ublas {
#elif defined(BOOST_UBLAS_HAVE_BINDINGS)
return boost::numeric::bindings::atlas::dot (c1 ().row (i), c2 ().column (j));
#else
+ boost::ignore_unused(j);
return apply (static_cast<const matrix_expression<C1> > (c1), static_cast<const matrix_expression<C2> > (c2, i));
#endif
}
diff --git a/boost/numeric/ublas/matrix_expression.hpp b/boost/numeric/ublas/matrix_expression.hpp
index 3953fad..a363130 100644
--- a/boost/numeric/ublas/matrix_expression.hpp
+++ b/boost/numeric/ublas/matrix_expression.hpp
@@ -3946,7 +3946,9 @@ namespace boost { namespace numeric { namespace ublas {
// (m / t) [i] [j] = m [i] [j] / t
template<class E1, class T2>
BOOST_UBLAS_INLINE
+ typename enable_if< is_convertible<T2, typename E1::value_type>,
typename matrix_binary_scalar2_traits<E1, const T2, scalar_divides<typename E1::value_type, T2> >::result_type
+ >::type
operator / (const matrix_expression<E1> &e1,
const T2 &e2) {
typedef typename matrix_binary_scalar2_traits<E1, const T2, scalar_divides<typename E1::value_type, T2> >::expression_type expression_type;
diff --git a/boost/numeric/ublas/storage.hpp b/boost/numeric/ublas/storage.hpp
index 9dfaaee..8821309 100644
--- a/boost/numeric/ublas/storage.hpp
+++ b/boost/numeric/ublas/storage.hpp
@@ -370,7 +370,7 @@ namespace boost { namespace numeric { namespace ublas {
// Random Access Container
BOOST_UBLAS_INLINE
size_type max_size () const {
- return ALLOC ().max_size();
+ return N;
}
BOOST_UBLAS_INLINE
@@ -537,6 +537,10 @@ namespace boost { namespace numeric { namespace ublas {
BOOST_UBLAS_INLINE
array_adaptor (size_type size, pointer data):
size_ (size), own_ (false), data_ (data) {}
+
+ template <size_t N>
+ BOOST_UBLAS_INLINE array_adaptor (T (&data)[N]):
+ size_ (N), own_ (false), data_ (data) {}
BOOST_UBLAS_INLINE
array_adaptor (const array_adaptor &a):
storage_array<self_type> (),
@@ -602,6 +606,16 @@ namespace boost { namespace numeric { namespace ublas {
resize_internal (size, data, init, true);
}
+ template <size_t N>
+ BOOST_UBLAS_INLINE void resize (T (&data)[N]) {
+ resize_internal (N, data, value_type (), false);
+ }
+
+ template <size_t N>
+ BOOST_UBLAS_INLINE void resize (T (&data)[N], value_type init) {
+ resize_internal (N, data, init, true);
+ }
+
BOOST_UBLAS_INLINE
size_type size () const {
return size_;
@@ -734,7 +748,7 @@ namespace boost { namespace numeric { namespace ublas {
typedef TT *argument_type;
BOOST_UBLAS_INLINE
- result_type operator () (argument_type x) {}
+ result_type operator () (argument_type /* x */) {}
};
public:
@@ -763,6 +777,10 @@ namespace boost { namespace numeric { namespace ublas {
BOOST_UBLAS_INLINE
shallow_array_adaptor (size_type size, pointer data):
size_ (size), own_ (false), data_ (data, leaker<value_type> ()) {}
+ BOOST_UBLAS_INLINE
+ template <size_t N>
+ shallow_array_adaptor (T (&data)[N]):
+ size_ (N), own_ (false), data_ (data, leaker<value_type> ()) {}
BOOST_UBLAS_INLINE
shallow_array_adaptor (const shallow_array_adaptor &a):
@@ -784,6 +802,7 @@ namespace boost { namespace numeric { namespace ublas {
std::fill (data.get () + (std::min) (size, size_), data.get () + size, init);
}
size_ = size;
+ own_ = true;
data_ = data;
}
}
@@ -794,7 +813,8 @@ namespace boost { namespace numeric { namespace ublas {
std::fill (data + (std::min) (size, size_), data + size, init);
}
size_ = size;
- data_ = data;
+ own_ = false;
+ data_.reset(data, leaker<value_type> ());
}
public:
BOOST_UBLAS_INLINE
@@ -813,6 +833,16 @@ namespace boost { namespace numeric { namespace ublas {
void resize (size_type size, pointer data, value_type init) {
resize_internal (size, data, init, true);
}
+ BOOST_UBLAS_INLINE
+ template <size_t N>
+ void resize (T (&data)[N]) {
+ resize_internal (N, data, value_type (), false);
+ }
+ BOOST_UBLAS_INLINE
+ template <size_t N>
+ void resize (T (&data)[N], value_type init) {
+ resize_internal (N, data, init, true);
+ }
BOOST_UBLAS_INLINE
size_type size () const {
diff --git a/boost/numeric/ublas/traits.hpp b/boost/numeric/ublas/traits.hpp
index 86ae9b5..ecd52ca 100644
--- a/boost/numeric/ublas/traits.hpp
+++ b/boost/numeric/ublas/traits.hpp
@@ -20,6 +20,9 @@
#include <boost/numeric/ublas/detail/config.hpp>
#include <boost/numeric/ublas/detail/iterator.hpp>
#include <boost/numeric/ublas/detail/returntype_deduction.hpp>
+#ifdef BOOST_UBLAS_USE_INTERVAL
+#include <boost/numeric/interval.hpp>
+#endif
#include <boost/type_traits.hpp>
#include <complex>
diff --git a/boost/numeric/ublas/vector.hpp b/boost/numeric/ublas/vector.hpp
index c90e429..c1384c1 100644
--- a/boost/numeric/ublas/vector.hpp
+++ b/boost/numeric/ublas/vector.hpp
@@ -27,6 +27,9 @@
#ifdef BOOST_UBLAS_CPP_GE_2011
#include <array>
#include <initializer_list>
+#if defined(BOOST_MSVC) // For std::forward in fixed_vector
+#include <utility>
+#endif
#endif
// Iterators based on ideas of Jeremy Siek
@@ -903,10 +906,18 @@ namespace boost { namespace numeric { namespace ublas {
/// \brief Construct a fixed_vector from a list of values
/// This constructor enables initialization by using any of:
/// fixed_vector<double, 3> v = { 1, 2, 3 } or fixed_vector<double,3> v( {1, 2, 3} ) or fixed_vector<double,3> v( 1, 2, 3 )
+#if defined(BOOST_MSVC)
+ // This may or may not work. Maybe use this for all instead only for MSVC
+ template <typename... U>
+ fixed_vector(U&&... values) :
+ vector_container<self_type> (),
+ data_{{ std::forward<U>(values)... }} {}
+#else
template <typename... Types>
fixed_vector(value_type v0, Types... vrest) :
vector_container<self_type> (),
data_{ { v0, vrest... } } {}
+#endif
// -----------------------
// Random Access Container