VTK  9.3.0
vtkCellGridSidesQuery.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-License-Identifier: BSD-3-Clause
8 #ifndef vtkCellGridSidesQuery_h
9 #define vtkCellGridSidesQuery_h
10 
11 #include "vtkCellGridQuery.h"
12 
13 #include "vtkStringToken.h" // For API.
14 
15 #include <functional>
16 #include <map>
17 #include <set>
18 #include <unordered_map>
19 #include <vector>
20 
21 VTK_ABI_NAMESPACE_BEGIN
22 
32 class VTKCOMMONDATAMODEL_EXPORT vtkCellGridSidesQuery : public vtkCellGridQuery
33 {
34 public:
37  void PrintSelf(ostream& os, vtkIndent indent) override;
38 
39  void Initialize() override;
40  void Finalize() override;
41 
42  std::map<vtkStringToken,
43  std::unordered_map<vtkStringToken, std::unordered_map<vtkIdType, std::set<int>>>>&
45  {
46  return this->Sides;
47  }
48  // std::map<vtkStringToken, std::unordered_map<vtkIdType, std::set<int>>>& GetSides() { return
49  // this->Sides; }
50 
71 
72  template <typename T, std::size_t NN>
73  void AddSide(vtkStringToken cellType, vtkIdType cell, int side, vtkStringToken shape,
74  const std::array<T, NN>& conn)
75  {
76  // Find where we should start hashing and in what direction.
77  std::size_t ss = 0;
78  T smin = conn[0];
79  for (std::size_t jj = 1; jj < NN; ++jj)
80  {
81  if (conn[jj] < smin)
82  {
83  smin = conn[jj];
84  ss = jj;
85  }
86  }
87  bool forward = conn[(ss + 1) % NN] > conn[(ss + NN - 1) % NN];
88 
89  std::size_t hashedValue = std::hash<std::size_t>{}(NN);
90  vtkCellGridSidesQuery::HashCombiner()(hashedValue, shape.GetId());
91  if (forward)
92  {
93  for (std::size_t ii = 0; ii < NN; ++ii)
94  {
95  std::size_t hashedToken = std::hash<T>{}(conn[(ss + ii) % NN]);
96  vtkCellGridSidesQuery::HashCombiner()(hashedValue, hashedToken);
97  }
98  }
99  else // backward
100  {
101  for (std::size_t ii = 0; ii < NN; ++ii)
102  {
103  std::size_t hashedToken = std::hash<T>{}(conn[(ss + NN - ii) % NN]);
104  // hashedValue = hashedValue ^ (hashedToken << (ii + 1));
105  vtkCellGridSidesQuery::HashCombiner()(hashedValue, hashedToken);
106  }
107  }
108  this->Hashes[hashedValue].Sides.insert(Side{ cellType, shape, cell, side });
109  }
110 
111  template <typename T>
112  void AddSide(vtkStringToken cellType, vtkIdType cell, int side, vtkStringToken shape,
113  const std::vector<T>& conn)
114  {
115  std::size_t ss = 0;
116  std::size_t NN = conn.size();
117  if (NN == 0)
118  {
119  return;
120  }
121 
122  T smin = conn[0];
123  for (std::size_t jj = 1; jj < NN; ++jj)
124  {
125  if (conn[jj] < smin)
126  {
127  smin = conn[jj];
128  ss = jj;
129  }
130  }
131  bool forward = conn[(ss + 1) % NN] > conn[(ss + NN - 1) % NN];
132 
133  std::size_t hashedValue = std::hash<std::size_t>{}(NN);
134  vtkCellGridSidesQuery::HashCombiner()(hashedValue, shape.GetId());
135  // std::cout << "Hash(" << (forward ? "F" : "R") << ")";
136  if (forward)
137  {
138  for (std::size_t ii = 0; ii < NN; ++ii)
139  {
140  std::size_t hashedToken = std::hash<T>{}(conn[(ss + ii) % NN]);
141  vtkCellGridSidesQuery::HashCombiner()(hashedValue, hashedToken);
142  // std::cout << " " << conn[(ss + ii) % NN];
143  }
144  }
145  else // backward
146  {
147  for (std::size_t ii = 0; ii < NN; ++ii)
148  {
149  std::size_t hashedToken = std::hash<T>{}(conn[(ss + NN - ii) % NN]);
150  // hashedValue = hashedValue ^ (hashedToken << (ii + 1));
151  vtkCellGridSidesQuery::HashCombiner()(hashedValue, hashedToken);
152  // std::cout << " " << conn[(ss + NN - ii) % NN];
153  }
154  }
155  // std::cout << " = " << std::hex << hashedValue << std::dec << "\n";
156  this->Hashes[hashedValue].Sides.insert(Side{ cellType, shape, cell, side });
157  }
159 
160 protected:
162  ~vtkCellGridSidesQuery() override = default;
163 
164  // Hash combiner adapted from boost::hash_combine
166  {
167  public:
168  template <typename T>
169  void operator()(T& h, typename std::enable_if<sizeof(T) == 8, std::size_t>::type k)
170  {
171  constexpr T m = 0xc6a4a7935bd1e995ull;
172  constexpr int r = 47;
173 
174  k *= m;
175  k ^= k >> r;
176  k *= m;
177 
178  h ^= k;
179  h *= m;
180 
181  // Completely arbitrary number, to prevent 0's
182  // from hashing to 0.
183  h += 0xe6546b64;
184  }
185 
186  template <typename T>
187  void operator()(T& h, typename std::enable_if<sizeof(T) == 4, std::size_t>::type k)
188  {
189  constexpr std::uint32_t c1 = 0xcc9e2d51;
190  constexpr std::uint32_t c2 = 0x1b873593;
191  constexpr std::uint32_t r1 = 15;
192  constexpr std::uint32_t r2 = 13;
193 
194  k *= c1;
195  k = (k << r1) | (k >> (32 - r1));
196  k *= c2;
197 
198  h ^= k;
199  h = (h << r2) | (h >> (32 - r2));
200  h = h * 5 + 0xe6546b64;
201  }
202  };
203 
204  struct Side
205  {
209  int SideId;
210  bool operator<(const Side& other) const
211  {
212  return (this->CellType < other.CellType) ||
213  (this->CellType == other.CellType &&
214  ((this->DOF < other.DOF) || (this->DOF == other.DOF && this->SideId < other.SideId)));
215  }
216  };
217  struct Entry
218  {
219  std::set<Side> Sides;
220  };
221  std::unordered_map<std::size_t, Entry> Hashes;
222  std::map<vtkStringToken,
223  std::unordered_map<vtkStringToken, std::unordered_map<vtkIdType, std::set<int>>>>
225 
226 private:
228  void operator=(const vtkCellGridSidesQuery&) = delete;
229 };
230 
231 VTK_ABI_NAMESPACE_END
232 #endif // vtkCellGridSidesQuery_h
Perform an operation on cells in a vtkCellMetadata instance.
void operator()(T &h, typename std::enable_if< sizeof(T)==4, std::size_t >::type k)
void operator()(T &h, typename std::enable_if< sizeof(T)==8, std::size_t >::type k)
A cell-grid query for enumerating sides of cells.
~vtkCellGridSidesQuery() override=default
std::map< vtkStringToken, std::unordered_map< vtkStringToken, std::unordered_map< vtkIdType, std::set< int > > > > & GetSides()
void Finalize() override
Override this if your query-result state requires finalization.
void AddSide(vtkStringToken cellType, vtkIdType cell, int side, vtkStringToken shape, const std::vector< T > &conn)
vtkCellGridSidesQuery()=default
std::map< vtkStringToken, std::unordered_map< vtkStringToken, std::unordered_map< vtkIdType, std::set< int > > > > Sides
std::unordered_map< std::size_t, Entry > Hashes
void Initialize() override
Override this if your query-result state requires initialization.
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
static vtkCellGridSidesQuery * New()
void AddSide(vtkStringToken cellType, vtkIdType cell, int side, vtkStringToken shape, const std::array< T, NN > &conn)
Add a side with the given shape and connectivity to the request's state.
a simple class to control print indentation
Definition: vtkIndent.h:29
Computes the portion of a dataset which is inside a selection.
Represent a string by its integer hash.
Hash GetId() const
Return the token's ID (usually its hash but possibly not in the case of collisions).
@ type
Definition: vtkX3D.h:516
@ side
Definition: vtkX3D.h:475
bool operator<(const Side &other) const
int vtkIdType
Definition: vtkType.h:315