}
template<class F>
- Real integrate(const F f, Real* error, Real* L1, const char* function, Real left_min_complement, Real right_min_complement, Real tolerance, std::size_t* levels) const;
+ decltype(std::declval<F>()(std::declval<Real>(), std::declval<Real>())) integrate(const F f, Real* error, Real* L1, const char* function, Real left_min_complement, Real right_min_complement, Real tolerance, std::size_t* levels) const;
private:
const std::vector<Real>& get_abscissa_row(std::size_t n)const
template<class Real, class Policy>
template<class F>
-Real tanh_sinh_detail<Real, Policy>::integrate(const F f, Real* error, Real* L1, const char* function, Real left_min_complement, Real right_min_complement, Real tolerance, std::size_t* levels) const
+decltype(std::declval<F>()(std::declval<Real>(), std::declval<Real>())) tanh_sinh_detail<Real, Policy>::integrate(const F f, Real* error, Real* L1, const char* function, Real left_min_complement, Real right_min_complement, Real tolerance, std::size_t* levels) const
{
using std::abs;
using std::fabs;
//
BOOST_ASSERT(m_abscissas[0][max_left_position] < 0);
BOOST_ASSERT(m_abscissas[0][max_right_position] < 0);
-
+ //
+ // The type of the result:
+ typedef decltype(std::declval<F>()(std::declval<Real>(), std::declval<Real>())) result_type;
Real h = m_t_max / m_inital_row_length;
- Real I0 = half_pi<Real>()*f(0, 1);
+ result_type I0 = half_pi<Real>()*f(0, 1);
Real L1_I0 = abs(I0);
for(size_t i = 1; i < m_abscissas[0].size(); ++i)
{
}
else
xc = x - 1;
- Real yp, ym;
+ result_type yp, ym;
yp = i <= max_right_position ? f(x, -xc) : 0;
ym = i <= max_left_position ? f(-x, xc) : 0;
I0 += (yp + ym)*w;
// L1_I0 and L1_I1 are the absolute integral values.
//
size_t k = 1;
- Real I1 = I0;
+ result_type I1 = I0;
Real L1_I1 = L1_I0;
Real err = 0;
//
I1 = half<Real>()*I0;
L1_I1 = half<Real>()*L1_I0;
h *= half<Real>();
- Real sum = 0;
+ result_type sum = 0;
Real absum = 0;
auto const& abscissa_row = this->get_abscissa_row(k);
auto const& weight_row = this->get_weight_row(k);
xc = x - 1;
}
- Real yp = j > max_right_index ? 0 : f(x, -xc);
- Real ym = j > max_left_index ? 0 : f(-x, xc);
- Real term = (yp + ym)*w;
+ result_type yp = j > max_right_index ? 0 : f(x, -xc);
+ result_type ym = j > max_left_index ? 0 : f(-x, xc);
+ result_type term = (yp + ym)*w;
sum += term;
// A question arises as to how accurately we actually need to estimate the L1 integral.
// parameters. We could keep hunting until we find something, but that would handicap
// integrals which really are zero.... so a compromise then!
//
- if (err <= tolerance*L1_I1)
+ if (err <= abs(tolerance*L1_I1))
{
break;
}
m_first_complements = {
1, 0, 1, 1, 3, 5, 11, 22,
};
-
+#ifndef BOOST_MATH_NO_ATOMIC_INT
m_committed_refinements = static_cast<boost::math::detail::atomic_unsigned_integer_type>(m_abscissas.size() - 1);
+#else
+ m_committed_refinements = m_abscissas.size() - 1;
+#endif
if (m_max_refinements >= m_abscissas.size())
{
m_first_complements = {
1, 0, 1, 1, 3, 5, 11, 22,
};
-
+#ifndef BOOST_MATH_NO_ATOMIC_INT
m_committed_refinements = static_cast<boost::math::detail::atomic_unsigned_integer_type>(m_abscissas.size() - 1);
+#else
+ m_committed_refinements = m_abscissas.size() - 1;
+#endif
if (m_max_refinements >= m_abscissas.size())
{
m_first_complements = {
1, 0, 1, 1, 3, 5, 11, 22,
};
-
+#ifndef BOOST_MATH_NO_ATOMIC_INT
m_committed_refinements = static_cast<boost::math::detail::atomic_unsigned_integer_type>(m_abscissas.size() - 1);
+#else
+ m_committed_refinements = m_abscissas.size() - 1;
+#endif
if (m_max_refinements >= m_abscissas.size())
{
m_first_complements = {
1, 0, 1, 1, 3, 5, 11, 22,
};
-
+#ifndef BOOST_MATH_NO_ATOMIC_INT
m_committed_refinements = static_cast<boost::math::detail::atomic_unsigned_integer_type>(m_abscissas.size() - 1);
+#else
+ m_committed_refinements = m_abscissas.size() - 1;
+#endif
if (m_max_refinements >= m_abscissas.size())
{