HavoqGT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
vertex_data.hpp
Go to the documentation of this file.
1 
2 /*
3  * Copyright (c) 2013, Lawrence Livermore National Security, LLC.
4  * Produced at the Lawrence Livermore National Laboratory.
5  * Re-written by Steven Feldman <feldman12@llnl.gov>.
6  * LLNL-CODE-644630.
7  * All rights reserved.
8  *
9  * This file is part of HavoqGT, Version 0.1.
10  * For details, see https://computation.llnl.gov/casc/dcca-pub/dcca/Downloads.html
11  *
12  * Please also read this link – Our Notice and GNU Lesser General Public License.
13  * http://www.gnu.org/licenses/old-licenses/lgpl-2.1.html
14  *
15  * This program is free software; you can redistribute it and/or modify it under
16  * the terms of the GNU Lesser General Public License (as published by the Free
17  * Software Foundation) version 2.1 dated February 1999.
18  *
19  * This program is distributed in the hope that it will be useful, but WITHOUT ANY
20  * WARRANTY; without even the IMPLIED WARRANTY OF MERCHANTABILITY or FITNESS FOR A
21  * PARTICULAR PURPOSE. See the terms and conditions of the GNU General Public
22  * License for more details.
23  *
24  * You should have received a copy of the GNU Lesser General Public License along
25  * with this program; if not, write to the Free Software Foundation, Inc.,
26  * 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27  *
28  * OUR NOTICE AND TERMS AND CONDITIONS OF THE GNU GENERAL PUBLIC LICENSE
29  *
30  * Our Preamble Notice
31  *
32  * A. This notice is required to be provided under our contract with the
33  * U.S. Department of Energy (DOE). This work was produced at the Lawrence
34  * Livermore National Laboratory under Contract No. DE-AC52-07NA27344 with the DOE.
35  *
36  * B. Neither the United States Government nor Lawrence Livermore National
37  * Security, LLC nor any of their employees, makes any warranty, express or
38  * implied, or assumes any liability or responsibility for the accuracy,
39  * completeness, or usefulness of any information, apparatus, product, or process
40  * disclosed, or represents that its use would not infringe privately-owned rights.
41  *
42  * C. Also, reference herein to any specific commercial products, process, or
43  * services by trade name, trademark, manufacturer or otherwise does not
44  * necessarily constitute or imply its endorsement, recommendation, or favoring by
45  * the United States Government or Lawrence Livermore National Security, LLC. The
46  * views and opinions of authors expressed herein do not necessarily state or
47  * reflect those of the United States Government or Lawrence Livermore National
48  * Security, LLC, and shall not be used for advertising or product endorsement
49  * purposes.
50  *
51  */
52 
53 #ifndef HAVOQGT_MPI_IMPL_VERTEX_DATA_HPP_
54 #define HAVOQGT_MPI_IMPL_VERTEX_DATA_HPP_
55 
57 
58 namespace havoqgt {
59 namespace mpi {
60 
61 template <typename SegementManager>
62 template <typename T, typename Allocator >
63 class delegate_partitioned_graph<SegementManager>::vertex_data {
64  public:
65  typedef T value_type;
67 
68  T& operator[] (const vertex_locator& locator)
69  {
70  if(locator.is_delegate()) {
71  assert(locator.local_id() < m_delegate_data.size());
72  return m_delegate_data[locator.local_id()];
73  }
74  assert(locator.local_id() < m_owned_vert_data.size());
75  return m_owned_vert_data[locator.local_id()];
76  }
77 
78 
79  const T& operator[] (const vertex_locator& locator) const
80  {
81  if(locator.is_delegate()) {
82  assert(locator.local_id() < m_delegate_data.size());
83  return m_delegate_data[locator.local_id()];
84  }
85  assert(locator.local_id() < m_owned_vert_data.size());
86  return m_owned_vert_data[locator.local_id()];
87  }
88 
89  void reset(const T& r) {
90  for(size_t i=0; i<m_owned_vert_data.size(); ++i) {
91  m_owned_vert_data[i] = r;
92  }
93  for(size_t i=0; i<m_delegate_data.size(); ++i) {
94  m_delegate_data[i] = r;
95  }
96  }
97 
98  void all_reduce() {
99  std::vector<T> tmp_in(m_delegate_data.begin(), m_delegate_data.end());
100  std::vector<T> tmp_out(tmp_in.size(), 0);
101  mpi_all_reduce(tmp_in, tmp_out, std::plus<T>(), MPI_COMM_WORLD);
102  std::copy(tmp_out.begin(), tmp_out.end(), m_delegate_data.begin());
103  }
104 
105  vertex_data(const delegate_partitioned_graph& dpg, Allocator allocate = Allocator() )
106  : m_owned_vert_data(allocate)
107  , m_delegate_data(allocate) {
108  m_owned_vert_data.resize(dpg.m_owned_info.size());
109  m_delegate_data.resize(dpg.m_delegate_info.size());
110  }
111 
112  T local_accumulate() const {
113  T to_return = T();
114  to_return = std::accumulate(m_owned_vert_data.begin(), m_owned_vert_data.end(), to_return);
115  to_return = std::accumulate(m_delegate_data.begin(), m_delegate_data.end(), to_return);
116  return to_return;
117  }
118 
119  T global_accumulate() const {
120  T local = local_accumulate();
121  return mpi_all_reduce(local,std::plus<T>(), MPI_COMM_WORLD);
122  }
123 
124  private:
125  bip::vector<T, Allocator > m_owned_vert_data;
126  bip::vector<T, Allocator > m_delegate_data;
127 };
128 
129 } // mpi
130 } // namespace havoqgt
131 #endif // HAVOQGT_MPI_IMPL_VERTEX_DATA_HPP_
T mpi_all_reduce(T in_d, Op in_op, MPI_Comm mpi_comm)
Definition: mpi.hpp:176
vertex_data(const delegate_partitioned_graph &dpg, Allocator allocate=Allocator())