accumulators/statistics/sum_kahan.hpp

0001 ///////////////////////////////////////////////////////////////////////////////
0002 // sum_kahan.hpp
0003 //
0004 //  Copyright 2010 Gaetano Mendola, 2011 Simon West. Distributed under the Boost
0005 //  Software License, Version 1.0. (See accompanying file
0006 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
0007
0008 #ifndef BOOST_ACCUMULATORS_STATISTICS_SUM_KAHAN_HPP_EAN_26_07_2010
0009 #define BOOST_ACCUMULATORS_STATISTICS_SUM_KAHAN_HPP_EAN_26_07_2010
0010
0011 #include <boost/accumulators/framework/accumulator_base.hpp>
0012 #include <boost/accumulators/framework/parameters/sample.hpp>
0013 #include <boost/accumulators/statistics_fwd.hpp>
0014 #include <boost/accumulators/statistics/sum.hpp>
0015 #include <boost/accumulators/statistics/weighted_sum_kahan.hpp>
0016 #include <boost/numeric/conversion/cast.hpp>
0017
0018 namespace boost { namespace accumulators
0019 {
0020
0021 namespace impl
0022 {
0023
0024 #if _MSC_VER > 1400
0025 # pragma float_control(push)
0026 # pragma float_control(precise, on)
0027 #endif
0028
0029 template<typename Sample, typename Tag>
0030 struct sum_kahan_impl
0031   : accumulator_base
0032 {
0033     typedef Sample result_type;
0034
0035     ////////////////////////////////////////////////////////////////////////////
0036     // sum_kahan_impl
0037     /**
0038         @brief Kahan summation algorithm
0039
0040         The Kahan summation algorithm reduces the numerical error obtained with standard
0041         sequential sum.
0042
0043     */
0044     template<typename Args>
0045     sum_kahan_impl(Args const & args)
0046       : sum(args[parameter::keyword<Tag>::get() | Sample()]),
0047         compensation(boost::numeric_cast<Sample>(0.0))
0048     {
0049     }
0050
0051     template<typename Args>
0052     void
0053 #if BOOST_ACCUMULATORS_GCC_VERSION > 40305
0054     __attribute__((__optimize__("no-associative-math")))
0055 #endif
0056     operator ()(Args const & args)
0057     {
0058         const Sample myTmp1 = args[parameter::keyword<Tag>::get()] - this->compensation;
0059         const Sample myTmp2 = this->sum + myTmp1;
0060         this->compensation = (myTmp2 - this->sum) - myTmp1;
0061         this->sum = myTmp2;
0062     }
0063
0064     result_type result(dont_care) const
0065     {
0066       return this->sum;
0067     }
0068
0069     // make this accumulator serializeable
0070     template<class Archive>
0071     void serialize(Archive & ar, const unsigned int file_version)
0072     {
0073         ar & sum;
0074         ar & compensation;
0075     }
0076
0077 private:
0078     Sample sum;
0079     Sample compensation;
0080 };
0081
0082 #if _MSC_VER > 1400
0083 # pragma float_control(pop)
0084 #endif
0085
0086 } // namespace impl
0087
0088 ///////////////////////////////////////////////////////////////////////////////
0089 // tag::sum_kahan
0090 // tag::sum_of_weights_kahan
0091 // tag::sum_of_variates_kahan
0092 //
0093 namespace tag
0094 {
0095
0096     struct sum_kahan
0097       : depends_on<>
0098     {
0099         /// INTERNAL ONLY
0100         ///
0101         typedef impl::sum_kahan_impl< mpl::_1, tag::sample > impl;
0102     };
0103
0104     struct sum_of_weights_kahan
0105       : depends_on<>
0106     {
0107         typedef mpl::true_ is_weight_accumulator;
0108         /// INTERNAL ONLY
0109         ///
0110         typedef accumulators::impl::sum_kahan_impl<mpl::_2, tag::weight> impl;
0111     };
0112
0113     template<typename VariateType, typename VariateTag>
0114     struct sum_of_variates_kahan
0115       : depends_on<>
0116     {
0117         /// INTERNAL ONLY
0118         ///
0119         typedef mpl::always<accumulators::impl::sum_kahan_impl<VariateType, VariateTag> > impl;
0120     };
0121
0122 } // namespace tag
0123
0124 ///////////////////////////////////////////////////////////////////////////////
0125 // extract::sum_kahan
0126 // extract::sum_of_weights_kahan
0127 // extract::sum_of_variates_kahan
0128 //
0129 namespace extract
0130 {
0131     extractor<tag::sum_kahan> const sum_kahan = {};
0132     extractor<tag::sum_of_weights_kahan> const sum_of_weights_kahan = {};
0133     extractor<tag::abstract_sum_of_variates> const sum_of_variates_kahan = {};
0134
0135     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_kahan)
0136     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_of_weights_kahan)
0137     BOOST_ACCUMULATORS_IGNORE_GLOBAL(sum_of_variates_kahan)
0138 } // namespace extract
0139
0140 using extract::sum_kahan;
0141 using extract::sum_of_weights_kahan;
0142 using extract::sum_of_variates_kahan;
0143
0144 // sum(kahan) -> sum_kahan
0145 template<>
0146 struct as_feature<tag::sum(kahan)>
0147 {
0148     typedef tag::sum_kahan type;
0149 };
0150
0151 // sum_of_weights(kahan) -> sum_of_weights_kahan
0152 template<>
0153 struct as_feature<tag::sum_of_weights(kahan)>
0154 {
0155     typedef tag::sum_of_weights_kahan type;
0156 };
0157
0158 // So that sum_kahan can be automatically substituted with
0159 // weighted_sum_kahan when the weight parameter is non-void.
0160 template<>
0161 struct as_weighted_feature<tag::sum_kahan>
0162 {
0163     typedef tag::weighted_sum_kahan type;
0164 };
0165
0166 template<>
0167 struct feature_of<tag::weighted_sum_kahan>
0168   : feature_of<tag::sum>
0169 {};
0170
0171 // for the purposes of feature-based dependency resolution,
0172 // sum_kahan provides the same feature as sum
0173 template<>
0174 struct feature_of<tag::sum_kahan>
0175   : feature_of<tag::sum>
0176 {
0177 };
0178
0179 // for the purposes of feature-based dependency resolution,
0180 // sum_of_weights_kahan provides the same feature as sum_of_weights
0181 template<>
0182 struct feature_of<tag::sum_of_weights_kahan>
0183   : feature_of<tag::sum_of_weights>
0184 {
0185 };
0186
0187 template<typename VariateType, typename VariateTag>
0188 struct feature_of<tag::sum_of_variates_kahan<VariateType, VariateTag> >
0189   : feature_of<tag::abstract_sum_of_variates>
0190 {
0191 };
0192
0193 }} // namespace boost::accumulators
0194
0195 #endif
0196