Chris@16: // Copyright (C) 2005-2006 Douglas Gregor <doug.gregor@gmail.com>.
Chris@16: // Copyright (C) 2004 The Trustees of Indiana University
Chris@16: 
Chris@16: // Use, modification and distribution is subject to the Boost Software
Chris@16: // License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
Chris@16: // http://www.boost.org/LICENSE_1_0.txt)
Chris@16: 
Chris@16: //   Authors: Douglas Gregor
Chris@16: //            Andrew Lumsdaine
Chris@16: 
Chris@16: // Message Passing Interface 1.1 -- Section 4.9.1. Scan
Chris@16: #ifndef BOOST_MPI_SCAN_HPP
Chris@16: #define BOOST_MPI_SCAN_HPP
Chris@16: 
Chris@16: #include <boost/mpi/exception.hpp>
Chris@16: #include <boost/mpi/datatype.hpp>
Chris@16: 
Chris@16: // For (de-)serializing sends and receives
Chris@16: #include <boost/mpi/packed_oarchive.hpp>
Chris@16: #include <boost/mpi/packed_iarchive.hpp>
Chris@16: 
Chris@16: // For packed_[io]archive sends and receives
Chris@16: #include <boost/mpi/detail/point_to_point.hpp>
Chris@16: 
Chris@16: #include <boost/mpi/communicator.hpp>
Chris@16: #include <boost/mpi/environment.hpp>
Chris@16: #include <boost/mpi/detail/computation_tree.hpp>
Chris@16: #include <boost/mpi/operations.hpp>
Chris@16: #include <algorithm>
Chris@16: #include <exception>
Chris@16: #include <boost/assert.hpp>
Chris@16: 
Chris@16: namespace boost { namespace mpi {
Chris@16: 
Chris@16: 
Chris@16: /************************************************************************
Chris@16:  * Implementation details                                               *
Chris@16:  ************************************************************************/
Chris@16: namespace detail {
Chris@16:   /**********************************************************************
Chris@16:    * Simple prefix reduction with MPI_Scan                              *
Chris@16:    **********************************************************************/
Chris@16: 
Chris@16:   // We are performing prefix reduction for a type that has an
Chris@16:   // associated MPI datatype and operation, so we'll use MPI_Scan
Chris@16:   // directly.
Chris@16:   template<typename T, typename Op>
Chris@16:   void
Chris@16:   scan_impl(const communicator& comm, const T* in_values, int n, T* out_values,
Chris@16:             Op op, mpl::true_ /*is_mpi_op*/, mpl::true_ /*is_mpi_datatype*/)
Chris@16:   {
Chris@16:     BOOST_MPI_CHECK_RESULT(MPI_Scan,
Chris@16:                            (const_cast<T*>(in_values), out_values, n,
Chris@16:                             boost::mpi::get_mpi_datatype<T>(*in_values),
Chris@16:                             (is_mpi_op<Op, T>::op()), comm));
Chris@16:   }
Chris@16: 
Chris@16:   /**********************************************************************
Chris@16:    * User-defined prefix reduction with MPI_Scan                        *
Chris@16:    **********************************************************************/
Chris@16: 
Chris@16:   // We are performing prefix reduction for a type that has an
Chris@16:   // associated MPI datatype but with a custom operation. We'll use
Chris@16:   // MPI_Scan directly, but we'll need to create an MPI_Op manually.
Chris@16:   template<typename T, typename Op>
Chris@16:   void
Chris@16:   scan_impl(const communicator& comm, const T* in_values, int n, T* out_values,
Chris@16:             Op op, mpl::false_ /*is_mpi_op*/, mpl::true_ /*is_mpi_datatype*/)
Chris@16:   {
Chris@16:     user_op<Op, T> mpi_op(op);
Chris@16:     BOOST_MPI_CHECK_RESULT(MPI_Scan,
Chris@16:                            (const_cast<T*>(in_values), out_values, n,
Chris@16:                             boost::mpi::get_mpi_datatype<T>(*in_values),
Chris@16:                             mpi_op.get_mpi_op(), comm));
Chris@16:   }
Chris@16: 
Chris@16:   /**********************************************************************
Chris@16:    * User-defined, tree-based reduction for non-MPI data types          *
Chris@16:    **********************************************************************/
Chris@16: 
Chris@16:   template<typename T, typename Op>
Chris@16:   void
Chris@16:   upper_lower_scan(const communicator& comm, const T* in_values, int n,
Chris@16:                    T* out_values, Op& op, int lower, int upper)
Chris@16:   {
Chris@16:     int tag = environment::collectives_tag();
Chris@16:     int rank = comm.rank();
Chris@16: 
Chris@16:     if (lower + 1 == upper) {
Chris@16:       std::copy(in_values, in_values + n, out_values);
Chris@16:     } else {
Chris@16:       int middle = (lower + upper) / 2;
Chris@16:       
Chris@16:       if (rank < middle) {
Chris@16:         // Lower half
Chris@16:         upper_lower_scan(comm, in_values, n, out_values, op, lower, middle);
Chris@16: 
Chris@16:         // If we're the last process in the lower half, send our values
Chris@16:         // to everyone in the upper half.
Chris@16:         if (rank == middle - 1) {
Chris@16:           packed_oarchive oa(comm);
Chris@16:           for (int i = 0; i < n; ++i)
Chris@16:             oa << out_values[i];
Chris@16: 
Chris@16:           for (int p = middle; p < upper; ++p)
Chris@16:             comm.send(p, tag, oa);
Chris@16:         }
Chris@16:       } else {
Chris@16:         // Upper half
Chris@16:         upper_lower_scan(comm, in_values, n, out_values, op, middle, upper);
Chris@16: 
Chris@16:         // Receive value from the last process in the lower half.
Chris@16:         packed_iarchive ia(comm);
Chris@16:         comm.recv(middle - 1, tag, ia);
Chris@16: 
Chris@16:         // Combine value that came from the left with our value
Chris@16:         T left_value;
Chris@16:         for (int i = 0; i < n; ++i)
Chris@16:           {
Chris@16:             ia >> left_value;
Chris@16:             out_values[i] = op(left_value, out_values[i]);
Chris@16:           }
Chris@16:       }
Chris@16:     }
Chris@16:   }
Chris@16: 
Chris@16:   // We are performing prefix reduction for a type that has no
Chris@16:   // associated MPI datatype and operation, so we'll use a simple
Chris@16:   // upper/lower algorithm.
Chris@16:   template<typename T, typename Op>
Chris@16:   inline void
Chris@16:   scan_impl(const communicator& comm, const T* in_values, int n, T* out_values, 
Chris@16:             Op op, mpl::false_ /*is_mpi_op*/, mpl::false_/*is_mpi_datatype*/)
Chris@16:   {
Chris@16:     upper_lower_scan(comm, in_values, n, out_values, op, 0, comm.size());
Chris@16:   }
Chris@16: } // end namespace detail
Chris@16: 
Chris@16: 
Chris@16: template<typename T, typename Op>
Chris@16: inline void
Chris@16: scan(const communicator& comm, const T& in_value, T& out_value, Op op)
Chris@16: {
Chris@16:   detail::scan_impl(comm, &in_value, 1, &out_value, op, 
Chris@16:                     is_mpi_op<Op, T>(), is_mpi_datatype<T>());
Chris@16: }
Chris@16: 
Chris@16: template<typename T, typename Op>
Chris@16: inline void
Chris@16: scan(const communicator& comm, const T* in_values, int n, T* out_values, Op op)
Chris@16: {
Chris@16:   detail::scan_impl(comm, in_values, n, out_values, op, 
Chris@16:                     is_mpi_op<Op, T>(), is_mpi_datatype<T>());
Chris@16: }
Chris@16: 
Chris@16: template<typename T, typename Op>
Chris@16: inline T
Chris@16: scan(const communicator& comm, const T& in_value, Op op)
Chris@16: {
Chris@16:   T out_value;
Chris@16:   detail::scan_impl(comm, &in_value, 1, &out_value, op, 
Chris@16:                     is_mpi_op<Op, T>(), is_mpi_datatype<T>());
Chris@16:   return out_value;
Chris@16: }
Chris@16: 
Chris@16: } } // end namespace boost::mpi
Chris@16: 
Chris@16: #endif // BOOST_MPI_SCAN_HPP