Chris@16: ///////////////////////////////////////////////////////////////////////////////
Chris@16: // sum_kahan.hpp
Chris@16: //
Chris@16: //  Copyright 2010 Gaetano Mendola, 2011 Simon West. Distributed under the Boost
Chris@16: //  Software License, Version 1.0. (See accompanying file
Chris@16: //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16: 
Chris@16: #ifndef BOOST_ACCUMULATORS_STATISTICS_SUM_KAHAN_HPP_EAN_26_07_2010
Chris@16: #define BOOST_ACCUMULATORS_STATISTICS_SUM_KAHAN_HPP_EAN_26_07_2010
Chris@16: 
Chris@16: #include <boost/accumulators/framework/accumulator_base.hpp>
Chris@16: #include <boost/accumulators/framework/parameters/sample.hpp>
Chris@16: #include <boost/accumulators/statistics_fwd.hpp>
Chris@16: #include <boost/accumulators/statistics/sum.hpp>
Chris@16: #include <boost/accumulators/statistics/weighted_sum_kahan.hpp>
Chris@16: #include <boost/numeric/conversion/cast.hpp>
Chris@16: 
Chris@16: namespace boost { namespace accumulators
Chris@16: {
Chris@16: 
Chris@16: namespace impl
Chris@16: {
Chris@16: 
Chris@16: #if _MSC_VER > 1400
Chris@16: # pragma float_control(push)
Chris@16: # pragma float_control(precise, on)
Chris@16: #endif
Chris@16: 
Chris@16: template<typename Sample, typename Tag>
Chris@16: struct sum_kahan_impl
Chris@16:   : accumulator_base
Chris@16: {
Chris@16:     typedef Sample result_type;
Chris@16: 
Chris@16:     ////////////////////////////////////////////////////////////////////////////
Chris@16:     // sum_kahan_impl
Chris@16:     /**
Chris@16:         @brief Kahan summation algorithm
Chris@16: 
Chris@16:         The Kahan summation algorithm reduces the numerical error obtained with standard
Chris@16:         sequential sum.
Chris@16: 
Chris@16:     */
Chris@16:     template<typename Args>
Chris@16:     sum_kahan_impl(Args const & args)
Chris@16:       : sum(args[parameter::keyword<Tag>::get() | Sample()]),
Chris@16:         compensation(boost::numeric_cast<Sample>(0.0))
Chris@16:     {
Chris@16:     }
Chris@16: 
Chris@16:     template<typename Args>
Chris@16:     void 
Chris@16: #if BOOST_ACCUMULATORS_GCC_VERSION > 40305
Chris@16:     __attribute__((__optimize__("no-associative-math")))
Chris@16: #endif
Chris@16:     operator ()(Args const & args)
Chris@16:     {
Chris@16:         const Sample myTmp1 = args[parameter::keyword<Tag>::get()] - this->compensation;
Chris@16:         const Sample myTmp2 = this->sum + myTmp1;
Chris@16:         this->compensation = (myTmp2 - this->sum) - myTmp1;
Chris@16:         this->sum = myTmp2;
Chris@16:     }
Chris@16: 
Chris@16:     result_type result(dont_care) const
Chris@16:     {
Chris@16:       return this->sum;
Chris@16:     }
Chris@16: 
Chris@16: private:
Chris@16:     Sample sum;
Chris@16:     Sample compensation;
Chris@16: };
Chris@16: 
Chris@16: #if _MSC_VER > 1400
Chris@16: # pragma float_control(pop)
Chris@16: #endif
Chris@16: 
Chris@16: } // namespace impl
Chris@16: 
Chris@16: ///////////////////////////////////////////////////////////////////////////////
Chris@16: // tag::sum_kahan
Chris@16: // tag::sum_of_weights_kahan
Chris@16: // tag::sum_of_variates_kahan
Chris@16: //
Chris@16: namespace tag
Chris@16: {
Chris@16: 
Chris@16:     struct sum_kahan
Chris@16:       : depends_on<>
Chris@16:     {
Chris@16:         /// INTERNAL ONLY
Chris@16:         ///
Chris@16:         typedef impl::sum_kahan_impl< mpl::_1, tag::sample > impl;
Chris@16:     };
Chris@16: 
Chris@16:     struct sum_of_weights_kahan
Chris@16:       : depends_on<>
Chris@16:     {
Chris@16:         typedef mpl::true_ is_weight_accumulator;
Chris@16:         /// INTERNAL ONLY
Chris@16:         ///
Chris@16:         typedef accumulators::impl::sum_kahan_impl<mpl::_2, tag::weight> impl;
Chris@16:     };
Chris@16: 
Chris@16:     template<typename VariateType, typename VariateTag>
Chris@16:     struct sum_of_variates_kahan
Chris@16:       : depends_on<>
Chris@16:     {
Chris@16:         /// INTERNAL ONLY
Chris@16:         ///
Chris@16:         typedef mpl::always<accumulators::impl::sum_kahan_impl<VariateType, VariateTag> > impl;
Chris@16:     };
Chris@16: 
Chris@16: } // namespace tag
Chris@16: 
Chris@16: ///////////////////////////////////////////////////////////////////////////////
Chris@16: // extract::sum_kahan
Chris@16: // extract::sum_of_weights_kahan
Chris@16: // extract::sum_of_variates_kahan
Chris@16: //
Chris@16: namespace extract
Chris@16: {
Chris@16:     extractor<tag::sum_kahan> const sum_kahan = {};
Chris@16:     extractor<tag::sum_of_weights_kahan> const sum_of_weights_kahan = {};
Chris@16:     extractor<tag::abstract_sum_of_variates> const sum_of_variates_kahan = {};
Chris@16: 
Chris@16:     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_kahan)
Chris@16:     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_of_weights_kahan)
Chris@16:     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_of_variates_kahan)
Chris@16: } // namespace extract
Chris@16: 
Chris@16: using extract::sum_kahan;
Chris@16: using extract::sum_of_weights_kahan;
Chris@16: using extract::sum_of_variates_kahan;
Chris@16: 
Chris@16: // sum(kahan) -> sum_kahan
Chris@16: template<>
Chris@16: struct as_feature<tag::sum(kahan)>
Chris@16: {
Chris@16:     typedef tag::sum_kahan type;
Chris@16: };
Chris@16: 
Chris@16: // sum_of_weights(kahan) -> sum_of_weights_kahan
Chris@16: template<>
Chris@16: struct as_feature<tag::sum_of_weights(kahan)>
Chris@16: {
Chris@16:     typedef tag::sum_of_weights_kahan type;
Chris@16: };
Chris@16: 
Chris@16: // So that sum_kahan can be automatically substituted with
Chris@16: // weighted_sum_kahan when the weight parameter is non-void.
Chris@16: template<>
Chris@16: struct as_weighted_feature<tag::sum_kahan>
Chris@16: {
Chris@16:     typedef tag::weighted_sum_kahan type;
Chris@16: };
Chris@16: 
Chris@16: template<>
Chris@16: struct feature_of<tag::weighted_sum_kahan>
Chris@16:   : feature_of<tag::sum>
Chris@16: {};
Chris@16: 
Chris@16: // for the purposes of feature-based dependency resolution,
Chris@16: // sum_kahan provides the same feature as sum
Chris@16: template<>
Chris@16: struct feature_of<tag::sum_kahan>
Chris@16:   : feature_of<tag::sum>
Chris@16: {
Chris@16: };
Chris@16: 
Chris@16: // for the purposes of feature-based dependency resolution,
Chris@16: // sum_of_weights_kahan provides the same feature as sum_of_weights
Chris@16: template<>
Chris@16: struct feature_of<tag::sum_of_weights_kahan>
Chris@16:   : feature_of<tag::sum_of_weights>
Chris@16: {
Chris@16: };
Chris@16: 
Chris@16: template<typename VariateType, typename VariateTag>
Chris@16: struct feature_of<tag::sum_of_variates_kahan<VariateType, VariateTag> >
Chris@16:   : feature_of<tag::abstract_sum_of_variates>
Chris@16: {
Chris@16: };
Chris@16: 
Chris@16: }} // namespace boost::accumulators
Chris@16: 
Chris@16: #endif
Chris@16: