XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
vertex_list.hpp
Go to the documentation of this file.
1 #ifndef VERTEX_LIST_HPP
2 #define VERTEX_LIST_HPP
3 
4 #include "space_settings.hpp"
5 
6 template<typename Device, typename F>
7 int get_starts_and_ends(int n, F condition, View<bool*, CLayout, HostType>& starts_h, View<bool*, CLayout, HostType>& ends_h){
8  View<bool*, CLayout, Device> starts(NoInit("starts"), n);
9  View<bool*, CLayout, Device> ends(NoInit("ends"), n);
10 
11  // First, determine which vertices go on the list
12  View<bool*, CLayout, Device> is_in(NoInit("is_in"), n);
13  Kokkos::parallel_for("check_condition", Kokkos::RangePolicy<typename Device::execution_space>(0, n), KOKKOS_LAMBDA(const int i){
14  is_in(i) = condition(i);
15  });
16 
17  // Next, subtract to get list starts and stops
18  int nranges = 0;
19  Kokkos::parallel_reduce("check_condition", Kokkos::RangePolicy<typename Device::execution_space>(0, n), KOKKOS_LAMBDA(const int i, int& nranges_l){
20  starts(i) = false;
21  ends(i) = false;
22  if(is_in(i)){
23  if(i==n-1){
24  ends(i) = true;
25  nranges_l++;
26  }else{
27  if(!is_in(i+1)){ // List just ended
28  ends(i) = true;
29  nranges_l++;
30  }
31  }
32 
33  if(i==0){
34  starts(i) = true;
35  }else{
36  if(!is_in(i-1)){ // List just started
37  starts(i) = true;
38  }
39  }
40  }
41  }, nranges);
42 
43  ends_h = View<bool*, CLayout, HostType>(NoInit("ends"), n);
44  starts_h = View<bool*, CLayout, HostType>(NoInit("starts"), n);
45  Kokkos::deep_copy(ends_h, ends);
46  Kokkos::deep_copy(starts_h, starts);
47 
48  return nranges;
49 }
50 
51 // A compressed way to represent a list of vertices.
52 // Note: using C/python convention so range {0, 2} means vertices 0 and 1
53 class VertexList{
54  struct IntegerRange{
55  int start;
56  int end; // Note: using C/python convention so range {0, 2} means vertices 0 and 1
57 
58  KOKKOS_INLINE_FUNCTION bool is_in_range(int i) const{
59  return (i>=start && i<end);
60  }
61 
62  KOKKOS_INLINE_FUNCTION int size() const{
63  return (end - start);
64  }
65  };
66  public:
67  View<IntegerRange*,CLayout,HostType> ranges;
68 
69  public:
70 
72  : ranges("ranges", 0)
73  {}
74 
75  // Initialize with single range
76  VertexList(int start, int end);
77 
78  // Initialize with a condition, e.g.:
79  // VertexList(grid.nnode, KOKKOS_LAMBDA( const int i){return grid_h.is_in_psi_range(i)});
80  template<typename F>
81  VertexList(int n, F condition){
82  View<bool*, CLayout, HostType> starts_h;
83  View<bool*, CLayout, HostType> ends_h;
84  int nranges = get_starts_and_ends<DeviceType>(n, condition, starts_h, ends_h);
85  ranges = View<IntegerRange*,CLayout,HostType>(NoInit("ranges"), nranges);
86 
87  // Now populate views
88  int i_range = 0;
89  for(int i=0; i<n; i++){
90  if(starts_h(i)){
91  ranges(i_range).start = i;
92  }
93  if(ends_h(i)){
94  ranges(i_range).end = i+1; // i is included, therefore the list ends at i+1
95  i_range += 1;
96  }
97  }
98  }
99 
100  // Initialize the list by passing in the full uncompressed list of vertices.
101  VertexList(const View<int*,CLayout,HostType>& unordered_full_list, bool one_indexed);
102 
110  KOKKOS_INLINE_FUNCTION bool is_in_list(int i) const {
111  // First check ranges
112  for (int is=0; is<ranges.size(); is++){
113  if(ranges(is).is_in_range(i)) return true;
114  }
115 
116  return false;
117  }
118 
119  // Number of vertices
120  KOKKOS_INLINE_FUNCTION int size() const{
121  int n = 0;
122 
123  // Add ranges
124  for (int is=0; is<ranges.size(); is++){
125  n += ranges(is).size();
126  }
127 
128  return n;
129  }
130 
136  template<typename DeviceExSpace, typename F>
137  inline void for_all(const std::string label, F lambda_func) const {
138  // Loop through ranges
139  for (int is=0; is<ranges.size(); is++){
140  Kokkos::parallel_for(label, Kokkos::RangePolicy<DeviceExSpace>(ranges(is).start, ranges(is).end), lambda_func);
141  }
142  }
143 
144  template<class Device>
145  void pack_contiguous(const View<double*, CLayout, Device>& input, const View<double*, CLayout, Device>& contiguous) const;
146 
147  template<class Device>
148  void unpack_contiguous(const View<double*, CLayout, Device>& contiguous, const View<double*, CLayout, Device>& output) const;
149 
150  VertexList operator&(const VertexList& list2) const;
151 
152  VertexList& operator|=(const VertexList& rhs); // Custom +=
153 
154  void batch(int first_val, int nvals, const View<int*, CLayout, HostType>& view) const;
155 
156  int min() const;
157 
158  bool is_contiguous() const;
159 
160  void print(const std::string& msg) const;
161 
162  bool is_subset(const VertexList& list2) const;
163 
164  bool intersects(const VertexList& list2) const;
165 
166  bool intersects(int start, int end) const;
167 
168  template<typename T>
169  void set_to_value(const T& view, double value) const{
170  for_all<typename T::execution_space>("vertex_list_set_to_value", KOKKOS_LAMBDA(const int i){
171  view(i) = value;
172  });
173  Kokkos::fence();
174  }
175 
176  template<typename T, typename F>
177  void set_to_value(const T& view, F func) const{
178  for_all<typename T::execution_space>("vertex_list_set_to_value", KOKKOS_LAMBDA(const int i){
179  view(i) = func(i);
180  });
181  Kokkos::fence();
182  }
183 
184  template<class Device>
185  View<int*, CLayout, Device> get_view_int_not_in_list(int nnode_in) const{
186  View<int*, CLayout, Device> view(NoInit("is_inside"), nnode_in);
187  // Set all values to 1
188  Kokkos::deep_copy(view, 1);
189 
190  // Set values inside list to 0
191  for_all<typename Device::execution_space>("set_is_inside", KOKKOS_LAMBDA(const int i){
192  view(i) = 0;
193  });
194  Kokkos::fence();
195 
196  return view;
197  }
198 };
199 
200 #endif
bool is_contiguous() const
Definition: vertex_list.cpp:273
void pack_contiguous(const View< double *, CLayout, Device > &input, const View< double *, CLayout, Device > &contiguous) const
Definition: vertex_list.cpp:56
KOKKOS_INLINE_FUNCTION bool is_in_range(int i) const
Definition: vertex_list.hpp:58
int start
Definition: vertex_list.hpp:55
int get_starts_and_ends(int n, F condition, View< bool *, CLayout, HostType > &starts_h, View< bool *, CLayout, HostType > &ends_h)
Definition: vertex_list.hpp:7
KOKKOS_INLINE_FUNCTION int size() const
Definition: vertex_list.hpp:120
void unpack_contiguous(const View< double *, CLayout, Device > &contiguous, const View< double *, CLayout, Device > &output) const
Definition: vertex_list.cpp:72
void for_all(const std::string label, F lambda_func) const
Definition: vertex_list.hpp:137
void print(const std::string &msg) const
Definition: vertex_list.cpp:277
VertexList & operator|=(const VertexList &rhs)
Definition: vertex_list.cpp:171
void set_to_value(const T &view, double value) const
Definition: vertex_list.hpp:169
VertexList()
Definition: vertex_list.hpp:71
VertexList operator&(const VertexList &list2) const
Definition: vertex_list.cpp:123
int min() const
Definition: vertex_list.cpp:268
void batch(int first_val, int nvals, const View< int *, CLayout, HostType > &view) const
Definition: vertex_list.cpp:243
Definition: vertex_list.hpp:54
bool is_subset(const VertexList &list2) const
Definition: vertex_list.cpp:288
View< int *, CLayout, Device > get_view_int_not_in_list(int nnode_in) const
Definition: vertex_list.hpp:185
VertexList(int n, F condition)
Definition: vertex_list.hpp:81
KOKKOS_INLINE_FUNCTION bool is_in_list(int i) const
Definition: vertex_list.hpp:110
View< IntegerRange *, CLayout, HostType > ranges
Definition: vertex_list.hpp:67
int end
Definition: vertex_list.hpp:56
Definition: vertex_list.hpp:53
void parallel_for(const std::string name, int n_ptl, Function func, Option option, HostAoSoA aosoa_h, DeviceAoSoA aosoa_d)
Definition: streamed_parallel_for.hpp:252
void set_to_value(const T &view, F func) const
Definition: vertex_list.hpp:177
KOKKOS_INLINE_FUNCTION int size() const
Definition: vertex_list.hpp:62
Kokkos::ViewAllocateWithoutInitializing NoInit
Definition: space_settings.hpp:68
bool intersects(const VertexList &list2) const
Definition: vertex_list.cpp:116