Compadre  1.3.3
Compadre_NeighborLists.hpp
Go to the documentation of this file.
1 #ifndef _COMPADRE_NEIGHBORLISTS_HPP_
2 #define _COMPADRE_NEIGHBORLISTS_HPP_
3 
4 #include "Compadre_Typedefs.hpp"
5 #include <Kokkos_Core.hpp>
6 
7 namespace Compadre {
8 
9 //! NeighborLists assists in accessing entries of compressed row neighborhood lists
10 template <typename view_type>
12 public:
13 
14  typedef view_type internal_view_type;
15  typedef Kokkos::View<global_index_type*, typename view_type::array_layout,
16  typename view_type::memory_space, typename view_type::memory_traits> internal_row_offsets_view_type;
17 
18 protected:
19 
20 
24 
26  view_type _cr_neighbor_lists;
28 
29  typename internal_row_offsets_view_type::HostMirror _host_row_offsets;
30  typename view_type::HostMirror _host_cr_neighbor_lists;
31  typename view_type::HostMirror _host_number_of_neighbors_list;
32 
33 
34 public:
35 
36 /** @name Constructors
37  * Ways to initialize a NeighborLists object
38  */
39 ///@{
40 
41  //! \brief Constructor for the purpose of classes who have NeighborLists as a member object
44  _needs_sync_to_host = true;
46  }
47 
48  /*! \brief Constructor for when compressed row `cr_neighbor_lists` is preallocated/populated,
49  * `number_of_neighbors_list` and `neighbor_lists_row_offsets` have already been populated.
50  */
51  NeighborLists(view_type cr_neighbor_lists, view_type number_of_neighbors_list,
52  internal_row_offsets_view_type neighbor_lists_row_offsets, bool compute_max = true) {
53  compadre_assert_release((view_type::rank==1) &&
54  "cr_neighbor_lists and number_neighbors_list and neighbor_lists_row_offsets must be a 1D Kokkos view.");
55 
56  _number_of_targets = number_of_neighbors_list.extent(0);
57  _number_of_neighbors_list = number_of_neighbors_list;
58  _cr_neighbor_lists = cr_neighbor_lists;
59  _row_offsets = neighbor_lists_row_offsets;
60 
61  _host_cr_neighbor_lists = Kokkos::create_mirror_view(_cr_neighbor_lists);
62  _host_number_of_neighbors_list = Kokkos::create_mirror_view(_number_of_neighbors_list);
63  _host_row_offsets = Kokkos::create_mirror_view(_row_offsets);
64 
65  Kokkos::deep_copy(_host_cr_neighbor_lists, _cr_neighbor_lists);
67  Kokkos::deep_copy(_host_row_offsets, _row_offsets);
68  Kokkos::fence();
69 
70  if (compute_max) {
72  } else {
74  }
75 
76  //check neighbor_lists is large enough
77  compadre_assert_release(((size_t)(this->getTotalNeighborsOverAllListsHost())<=cr_neighbor_lists.extent(0))
78  && "neighbor_lists is not large enough to store all neighbors.");
79 
80  _needs_sync_to_host = false;
81  }
82 
83  /*! \brief Constructor for when compressed row `cr_neighbor_lists` is preallocated/populated,
84  * and `number_of_neighbors_list` is already populated, and row offsets still need to be computed.
85  */
86  NeighborLists(view_type cr_neighbor_lists, view_type number_of_neighbors_list) {
87  compadre_assert_release((view_type::rank==1)
88  && "cr_neighbor_lists and number_neighbors_list must be a 1D Kokkos view.");
89 
90  _number_of_targets = number_of_neighbors_list.extent(0);
91 
92  _row_offsets = internal_row_offsets_view_type("row offsets", number_of_neighbors_list.extent(0));
93  _number_of_neighbors_list = number_of_neighbors_list;
94  _cr_neighbor_lists = cr_neighbor_lists;
95 
96  _host_cr_neighbor_lists = Kokkos::create_mirror_view(_cr_neighbor_lists);
97  _host_number_of_neighbors_list = Kokkos::create_mirror_view(_number_of_neighbors_list);
98  _host_row_offsets = Kokkos::create_mirror_view(_row_offsets);
99 
100  Kokkos::deep_copy(_host_cr_neighbor_lists, _cr_neighbor_lists);
102  Kokkos::deep_copy(_host_row_offsets, _row_offsets);
103  Kokkos::fence();
104 
107 
108  //check neighbor_lists is large enough
109  compadre_assert_release(((size_t)(this->getTotalNeighborsOverAllListsHost())<=cr_neighbor_lists.extent(0))
110  && "neighbor_lists is not large enough to store all neighbors.");
111 
112  _needs_sync_to_host = false;
113  }
114 
115  /*! \brief Constructor for when `number_of_neighbors_list` is already populated.
116  * Will allocate space for compressed row neighbor lists data, and will allocate then
117  * populate information for row offsets.
118  */
119  NeighborLists(view_type number_of_neighbors_list) {
120  compadre_assert_release((view_type::rank==1)
121  && "cr_neighbor_lists and number_neighbors_list must be a 1D Kokkos view.");
122 
123  _number_of_targets = number_of_neighbors_list.extent(0);
124 
125  _row_offsets = internal_row_offsets_view_type("row offsets", number_of_neighbors_list.extent(0));
126  _number_of_neighbors_list = number_of_neighbors_list;
127 
128  _host_number_of_neighbors_list = Kokkos::create_mirror_view(_number_of_neighbors_list);
129  _host_row_offsets = Kokkos::create_mirror_view(_row_offsets);
130 
132  Kokkos::deep_copy(_host_row_offsets, _row_offsets);
133  Kokkos::fence();
134 
137 
138  _cr_neighbor_lists = view_type("compressed row neighbor lists data", this->getTotalNeighborsOverAllListsHost());
139  _host_cr_neighbor_lists = Kokkos::create_mirror_view(_cr_neighbor_lists);
140  Kokkos::deep_copy(_host_cr_neighbor_lists, _cr_neighbor_lists);
141  Kokkos::fence();
142 
143  _needs_sync_to_host = false;
144  }
145 ///@}
146 
147 /** @name Public modifiers
148  */
149 ///@{
150 
151  //! Setter function for N(i,j) indexing where N(i,j) is the index of the jth neighbor of i
152  KOKKOS_INLINE_FUNCTION
153  void setNeighborDevice(int target_index, int neighbor_num, int new_value) {
154  _cr_neighbor_lists(_row_offsets(target_index)+neighbor_num) = new_value;
155  // indicate that host view is now out of sync with device
156  // but only in debug mode (notice the next line is both setting the variable and checking it was set)
158  }
159 
160  //! Calculate the maximum number of neighbors of all targets' neighborhoods (host)
162  if (_number_of_neighbors_list.extent(0)==0) {
164  } else {
165  auto number_of_neighbors_list = _number_of_neighbors_list;
166  Kokkos::parallel_reduce("max number of neighbors",
167  Kokkos::RangePolicy<typename view_type::execution_space>(0, _number_of_neighbors_list.extent(0)),
168  KOKKOS_LAMBDA(const int i, int& t_max_num_neighbors) {
169  t_max_num_neighbors = (number_of_neighbors_list(i) > t_max_num_neighbors) ? number_of_neighbors_list(i) : t_max_num_neighbors;
170  }, Kokkos::Max<int>(_max_neighbor_list_row_storage_size));
171  Kokkos::fence();
172  }
173  }
174 
175  //! Calculate the row offsets for each target's neighborhood (host)
177  auto number_of_neighbors_list = _number_of_neighbors_list;
178  auto row_offsets = _row_offsets;
179  Kokkos::parallel_scan("number of neighbors offsets",
180  Kokkos::RangePolicy<typename view_type::execution_space>(0, _number_of_neighbors_list.extent(0)),
181  KOKKOS_LAMBDA(const int i, global_index_type& lsum, bool final) {
182  row_offsets(i) = lsum;
183  lsum += number_of_neighbors_list(i);
184  });
185  Kokkos::deep_copy(_host_row_offsets, _row_offsets);
186  Kokkos::fence();
187  }
188 
189  //! Sync the host from the device (copy device data to host)
191  Kokkos::deep_copy(_host_cr_neighbor_lists, _cr_neighbor_lists);
192  Kokkos::fence();
193  _needs_sync_to_host = false;
194  }
195 
196  //! Device view into neighbor lists data (use with caution)
197  view_type getNeighborLists() {
198  return _cr_neighbor_lists;
199  }
200 
201 ///@}
202 /** @name Public accessors
203  */
204 ///@{
205  //! Get number of total targets having neighborhoods (host/device).
206  KOKKOS_INLINE_FUNCTION
207  int getNumberOfTargets() const {
208  return _number_of_targets;
209  }
210 
211  //! Get number of neighbors for a given target (host)
212  int getNumberOfNeighborsHost(int target_index) const {
213  return _host_number_of_neighbors_list(target_index);
214  }
215 
216  //! Get number of neighbors for a given target (device)
217  KOKKOS_INLINE_FUNCTION
218  int getNumberOfNeighborsDevice(int target_index) const {
219  return _number_of_neighbors_list(target_index);
220  }
221 
222  //! Get offset into compressed row neighbor lists (host)
223  global_index_type getRowOffsetHost(int target_index) const {
224  return _host_row_offsets(target_index);
225  }
226 
227  //! Get offset into compressed row neighbor lists (device)
228  KOKKOS_INLINE_FUNCTION
229  global_index_type getRowOffsetDevice(int target_index) const {
230  return _row_offsets(target_index);
231  }
232 
233  //! Offers N(i,j) indexing where N(i,j) is the index of the jth neighbor of i (host)
234  int getNeighborHost(int target_index, int neighbor_num) const {
236  && "Stale information in host_cr_neighbor_lists. Call CopyDeviceDataToHost() to refresh.");
237  compadre_assert_debug((neighbor_num<_number_of_neighbors_list(target_index))
238  && "neighor_num exceeds number of neighbors for this target_index.");
239  return _host_cr_neighbor_lists(_row_offsets(target_index)+neighbor_num);
240  }
241 
242  //! Offers N(i,j) indexing where N(i,j) is the index of the jth neighbor of i (device)
243  KOKKOS_INLINE_FUNCTION
244  int getNeighborDevice(int target_index, int neighbor_num) const {
245  return _cr_neighbor_lists(_row_offsets(target_index)+neighbor_num);
246  }
247 
248  //! Get the maximum number of neighbors of all targets' neighborhoods (host/device)
249  KOKKOS_INLINE_FUNCTION
250  int getMaxNumNeighbors() const {
251  compadre_kernel_assert_debug((_max_neighbor_list_row_storage_size > -1) && "getMaxNumNeighbors() called but maximum never calculated.");
253  }
254 
255  //! Get the sum of the number of neighbors of all targets' neighborhoods (host)
257  if (this->getNumberOfTargets()==0) {
258  return 0;
259  } else {
261  }
262  }
263 
264  //! Get the sum of the number of neighbors of all targets' neighborhoods (device)
265  KOKKOS_INLINE_FUNCTION
268  }
269 ///@}
270 
271 }; // NeighborLists
272 
273 //! CreateNeighborLists allows for the construction of an object of type NeighborLists with template deduction
274 template <typename view_type>
275 NeighborLists<view_type> CreateNeighborLists(view_type number_of_neighbors_list) {
276  return NeighborLists<view_type>(number_of_neighbors_list);
277 }
278 
279 //! CreateNeighborLists allows for the construction of an object of type NeighborLists with template deduction
280 template <typename view_type>
281 NeighborLists<view_type> CreateNeighborLists(view_type neighbor_lists, view_type number_of_neighbors_list) {
282  return NeighborLists<view_type>(neighbor_lists, number_of_neighbors_list);
283 }
284 
285 //! CreateNeighborLists allows for the construction of an object of type NeighborLists with template deduction
286 template <typename view_type>
287 NeighborLists<view_type> CreateNeighborLists(view_type neighbor_lists, view_type number_of_neighbors_list, view_type neighbor_lists_row_offsets) {
288  return NeighborLists<view_type>(neighbor_lists, number_of_neighbors_list, neighbor_lists_row_offsets);
289 }
290 
291 //! Converts 2D neighbor lists to compressed row neighbor lists
292 template <typename view_type_2d, typename view_type_1d = Kokkos::View<int*, typename view_type_2d::memory_space, typename view_type_2d::memory_traits> >
294 
295  // gets total number of neighbors over all lists
296  // computes calculation where the data resides (device/host)
297  global_index_type total_storage_size = 0;
298  Kokkos::parallel_reduce("total number of neighbors over all lists", Kokkos::RangePolicy<typename view_type_2d::execution_space>(0, neighbor_lists.extent(0)),
299  KOKKOS_LAMBDA(const int i, global_index_type& t_total_num_neighbors) {
300  t_total_num_neighbors += neighbor_lists(i,0);
301  }, Kokkos::Sum<global_index_type>(total_storage_size));
302  Kokkos::fence();
303 
304  // view_type_1d may be on host or device, and view_type_2d may be either as well (could even be opposite)
305  view_type_1d new_cr_neighbor_lists("compressed row neighbor lists", total_storage_size);
306  view_type_1d new_number_of_neighbors_list("number of neighbors list", neighbor_lists.extent(0));
307 
308  // copy number of neighbors list over to view_type_1d
309  // d_neighbor_lists will be accessible from view_type_1d's execution space
310  auto d_neighbor_lists = create_mirror_view(typename view_type_1d::execution_space(), neighbor_lists);
311  Kokkos::deep_copy(d_neighbor_lists, neighbor_lists);
312  Kokkos::fence();
313  Kokkos::parallel_for("copy number of neighbors to compressed row",
314  Kokkos::RangePolicy<typename view_type_1d::execution_space>(0, neighbor_lists.extent(0)),
315  KOKKOS_LAMBDA(const int i) {
316  new_number_of_neighbors_list(i) = d_neighbor_lists(i,0);
317  });
318  Kokkos::fence();
319 
320 
321  // this will calculate row offsets
322  auto nla(CreateNeighborLists(new_cr_neighbor_lists, new_number_of_neighbors_list));
323  auto cr_data = nla.getNeighborLists();
324 
325  // if device_execution_space can access this view, then write directly into the view
326  if (Kokkos::SpaceAccessibility<device_execution_space, typename view_type_1d::memory_space>::accessible==1) {
327  Kokkos::parallel_for("copy neighbor lists to compressed row", Kokkos::RangePolicy<typename view_type_1d::execution_space>(0, neighbor_lists.extent(0)),
328  KOKKOS_LAMBDA(const int i) {
329  for (int j=0; j<d_neighbor_lists(i,0); ++j) {
330  cr_data(nla.getRowOffsetDevice(i)+j) = d_neighbor_lists(i,j+1);
331  }
332  });
333  Kokkos::fence();
334  nla.copyDeviceDataToHost(); // has a fence at the end
335  }
336  // otherwise we are writing to a view that can't be seen from device (must be host space),
337  // and d_neighbor_lists was already made to be a view_type that is accessible from view_type_1d's execution_space
338  // (which we know is host) so we can do a parallel_for over the host_execution_space
339  else {
340  Kokkos::parallel_for("copy neighbor lists to compressed row", Kokkos::RangePolicy<host_execution_space>(0, neighbor_lists.extent(0)),
341  KOKKOS_LAMBDA(const int i) {
342  for (int j=0; j<neighbor_lists(i,0); ++j) {
343  cr_data(nla.getRowOffsetHost(i)+j) = d_neighbor_lists(i,j+1);
344  }
345  });
346  Kokkos::fence();
347  }
348 
349  return nla;
350 }
351 
352 } // Compadre namespace
353 
354 #endif
355 
KOKKOS_INLINE_FUNCTION int getNumberOfNeighborsDevice(int target_index) const
Get number of neighbors for a given target (device)
view_type::HostMirror _host_number_of_neighbors_list
std::size_t global_index_type
KOKKOS_INLINE_FUNCTION int getMaxNumNeighbors() const
Get the maximum number of neighbors of all targets&#39; neighborhoods (host/device)
internal_row_offsets_view_type::HostMirror _host_row_offsets
view_type getNeighborLists()
Device view into neighbor lists data (use with caution)
NeighborLists< view_type > CreateNeighborLists(view_type number_of_neighbors_list)
CreateNeighborLists allows for the construction of an object of type NeighborLists with template dedu...
Kokkos::View< global_index_type *, typename view_type::array_layout, typename view_type::memory_space, typename view_type::memory_traits > internal_row_offsets_view_type
int getNeighborHost(int target_index, int neighbor_num) const
Offers N(i,j) indexing where N(i,j) is the index of the jth neighbor of i (host)
void computeMaxNumNeighbors()
Calculate the maximum number of neighbors of all targets&#39; neighborhoods (host)
#define TO_GLOBAL(variable)
KOKKOS_INLINE_FUNCTION global_index_type getTotalNeighborsOverAllListsDevice() const
Get the sum of the number of neighbors of all targets&#39; neighborhoods (device)
NeighborLists()
Constructor for the purpose of classes who have NeighborLists as a member object. ...
KOKKOS_INLINE_FUNCTION void setNeighborDevice(int target_index, int neighbor_num, int new_value)
Setter function for N(i,j) indexing where N(i,j) is the index of the jth neighbor of i...
NeighborLists(view_type number_of_neighbors_list)
Constructor for when number_of_neighbors_list is already populated. Will allocate space for compresse...
NeighborLists assists in accessing entries of compressed row neighborhood lists.
KOKKOS_INLINE_FUNCTION int getNeighborDevice(int target_index, int neighbor_num) const
Offers N(i,j) indexing where N(i,j) is the index of the jth neighbor of i (device) ...
KOKKOS_INLINE_FUNCTION global_index_type getRowOffsetDevice(int target_index) const
Get offset into compressed row neighbor lists (device)
NeighborLists(view_type cr_neighbor_lists, view_type number_of_neighbors_list)
Constructor for when compressed row cr_neighbor_lists is preallocated/populated, and number_of_neighb...
global_index_type getTotalNeighborsOverAllListsHost() const
Get the sum of the number of neighbors of all targets&#39; neighborhoods (host)
void computeRowOffsets()
Calculate the row offsets for each target&#39;s neighborhood (host)
void copyDeviceDataToHost()
Sync the host from the device (copy device data to host)
#define compadre_assert_debug(condition)
compadre_assert_debug is used for assertions that are checked in loops, as these significantly impact...
int getNumberOfNeighborsHost(int target_index) const
Get number of neighbors for a given target (host)
view_type::HostMirror _host_cr_neighbor_lists
NeighborLists< view_type_1d > Convert2DToCompressedRowNeighborLists(view_type_2d neighbor_lists)
Converts 2D neighbor lists to compressed row neighbor lists.
NeighborLists(view_type cr_neighbor_lists, view_type number_of_neighbors_list, internal_row_offsets_view_type neighbor_lists_row_offsets, bool compute_max=true)
Constructor for when compressed row cr_neighbor_lists is preallocated/populated, number_of_neighbors_...
global_index_type getRowOffsetHost(int target_index) const
Get offset into compressed row neighbor lists (host)
#define compadre_assert_release(condition)
compadre_assert_release is used for assertions that should always be checked, but generally are not e...
internal_row_offsets_view_type _row_offsets
KOKKOS_INLINE_FUNCTION int getNumberOfTargets() const
Get number of total targets having neighborhoods (host/device).
#define compadre_kernel_assert_debug(condition)