// Copyright (C) 2006 Fokko Beekhof // Email contact: Fokko.Beekhof@cui.unige.ch // The OMPTL library is free software; you can redistribute it and/or // modify it under the terms of the GNU Lesser General Public // License as published by the Free Software Foundation; either // version 2.1 of the License, or (at your option) any later version. // This library is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU // Lesser General Public License for more details. // You should have received a copy of the GNU Lesser General Public // License along with this library; if not, write to the Free Software // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA #ifndef OMPTL_NUMERIC_H #define OMPTL_NUMERIC_H 1 #include <utility> #include <functional> #include <iterator> #include "omptl_algorithm" #include "omptl_tools.h" namespace omptl { template <class IteratorTag> struct _Accumulate { template <class InputIterator, class T, class BinaryFunction> static T accumulate(InputIterator first, InputIterator last, T init, T par_init, BinaryFunction binary_op, const unsigned P) { assert(P > 0u); if (_linear_serial_is_faster(first, last, P)) return ::std::accumulate(first, last, init, binary_op); assert(2*(int)P <= std::distance(first, last)); ::std::vector< ::std::pair<InputIterator, InputIterator> > partitions(P); ::omptl::_partition_range(first, last, partitions, P); ::std::vector<T> results(P); #pragma omp parallel for for (int t = 0; t < int(P); ++t) results[t] = ::std::accumulate( partitions[t].first, partitions[t].second, par_init, binary_op); return ::std::accumulate(results.begin(), results.end(), init, binary_op); } }; template <> struct _Accumulate< ::std::input_iterator_tag > { template <class InputIterator, class T, class BinaryFunction> static T accumulate(InputIterator first, InputIterator last, T init, T par_init, BinaryFunction binary_op, const unsigned P) { return ::std::accumulate(first, last, init, binary_op); } }; template <class InputIterator, class BinaryFunction> struct _AccumulateOp { template <class T> static T accumulate(InputIterator first, InputIterator last, T init, BinaryFunction binary_op, const unsigned P) { assert(P > 0u); return ::std::accumulate(first, last, init, binary_op); } }; template <class InputIterator> struct _AccumulateOp<InputIterator, ::std::plus<typename ::std::iterator_traits<InputIterator>::value_type> > { typedef ::std::plus<typename ::std::iterator_traits<InputIterator>::value_type> BinaryFunction; template <class T> static T accumulate(InputIterator first, InputIterator last, T init, BinaryFunction binary_op, const unsigned P) { assert(P > 0u); return ::omptl::_Accumulate< typename ::std::iterator_traits<InputIterator>::iterator_category>:: accumulate(first, last, init, T(0), binary_op, P); } }; template <class InputIterator> struct _AccumulateOp<InputIterator, ::std::multiplies<typename ::std::iterator_traits<InputIterator>::value_type> > { typedef ::std::multiplies<typename ::std::iterator_traits<InputIterator>::value_type> BinaryFunction; template <class T> static T accumulate(InputIterator first, InputIterator last, T init, BinaryFunction binary_op, const unsigned P) { assert(P > 0u); return ::omptl::_Accumulate<typename ::std::iterator_traits<InputIterator>::iterator_category>:: accumulate(first, last, init, T(1), binary_op, P); } }; template <class InputIterator, class T, class BinaryFunction> T accumulate(InputIterator first, InputIterator last, T init, BinaryFunction binary_op, const unsigned P) { assert(P > 0u); return ::omptl::_AccumulateOp<InputIterator, BinaryFunction>:: accumulate(first, last, init, binary_op, P); } template <class InputIterator, class T> T accumulate(InputIterator first, InputIterator last, T init, const unsigned P) { assert(P > 0u); typedef typename std::iterator_traits<InputIterator>::value_type VT; return ::omptl::accumulate(first, last, init, std::plus<VT>(), P); } template <class InputIterator, class OutputIterator, class BinaryFunction> OutputIterator adjacent_difference(InputIterator first, InputIterator last, OutputIterator result, BinaryFunction binary_op, const unsigned P) { return ::std::adjacent_difference(first, last, result, binary_op); } template <class InputIterator, class OutputIterator> OutputIterator adjacent_difference(InputIterator first, InputIterator last, OutputIterator result, const unsigned P) { // ::std::minus<typename ::std::iterator_traits<InputIterator>::value_type>() return ::std::adjacent_difference(first, last, result); } template <class Iterator1Tag, class Iterator2Tag> struct _InnerProduct { template <class InputIterator1, class InputIterator2, class T, class BinaryFunction1, class BinaryFunction2> static T inner_product( InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, BinaryFunction1 binary_op1, BinaryFunction2 binary_op2, const unsigned P) { return ::std::inner_product(first1, last1, first2, init, binary_op1, binary_op2); } template <class Iterator1, class Iterator2, class T, class BinaryFunction2> static T inner_product(Iterator1 first1, Iterator1 last1, Iterator2 first2, T init, ::std::plus<T> binary_op1, BinaryFunction2 binary_op2, const unsigned P) { assert(P > 0u); if (_linear_serial_is_faster(first1, last1, P)) return ::std::inner_product(first1, last1, first2, init, binary_op1, binary_op2); assert(2*(int)P <= std::distance(first1, last1)); ::std::vector< ::std::pair<Iterator1, Iterator1> > partitions1(P); ::omptl::_partition_range(first1, last1, partitions1, P); ::std::vector<Iterator2> partitions2(P); ::omptl::_copy_partitions(partitions1, first2, partitions2, P); #pragma omp parallel for reduction(+:init) for (int t = 0; t < int(P); ++t) init += ::std::inner_product( partitions1[t].first, partitions1[t].second, partitions2[t], T(0.0), binary_op1, binary_op2); return init; } }; template <class Iterator2Tag> struct _InnerProduct< ::std::input_iterator_tag, Iterator2Tag> { template <class InputIterator1, class InputIterator2, class T, class BinaryFunction1, class BinaryFunction2> static T inner_product(InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, BinaryFunction1 binary_op1, BinaryFunction2 binary_op2, const unsigned P) { return ::std::inner_product(first1, last1, first2, init, binary_op1, binary_op2); } }; template <class Iterator1Tag> struct _InnerProduct<Iterator1Tag, ::std::input_iterator_tag> { template <class InputIterator1, class InputIterator2, class T, class BinaryFunction1, class BinaryFunction2> static T inner_product(InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, BinaryFunction1 binary_op1, BinaryFunction2 binary_op2, const unsigned P) { return ::std::inner_product(first1, last1, first2, init, binary_op1, binary_op2); } }; template <> struct _InnerProduct< ::std::input_iterator_tag, ::std::input_iterator_tag > { template <class InputIterator1, class InputIterator2, class T, class BinaryFunction1, class BinaryFunction2> static T inner_product(InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, BinaryFunction1 binary_op1, BinaryFunction2 binary_op2, const unsigned P) { return ::std::inner_product(first1, last1, first2, init, binary_op1, binary_op2); } }; template <class InputIterator1, class InputIterator2, class T, class BinaryFunction1, class BinaryFunction2> T inner_product(InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, BinaryFunction1 binary_op1, BinaryFunction2 binary_op2, const unsigned P) { return _InnerProduct< typename ::std::iterator_traits<InputIterator1>::iterator_category, typename ::std::iterator_traits<InputIterator2>::iterator_category> ::inner_product(first1, last1, first2, init, binary_op1, binary_op2, P); } template <class InputIterator1, class InputIterator2, class T> T inner_product(InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, T init, const unsigned P) { return ::omptl::inner_product(first1, last1, first2, init, ::std::plus<T>(), ::std::multiplies<T>(), P); } template <class InputIterator, class OutputIterator, class BinaryOperation> OutputIterator partial_sum(InputIterator first, InputIterator last, OutputIterator result, BinaryOperation binary_op, const unsigned P) { return ::std::partial_sum(first, last, result, binary_op); } template <class InputIterator, class OutputIterator> OutputIterator partial_sum(InputIterator first, InputIterator last, OutputIterator result, const unsigned P) { // ::std::plus<typename ::std::iterator_traits<InputIterator>::value_type>(), return ::std::partial_sum(first, last, result); } } // namespace omptl #endif /* OMPTL_NUMERIC_H */