Reference documentation for deal.II version Git 1aa49a030b 2021-12-07 05:12:24 -0500
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
cell_id.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_cell_id_h
17 #define dealii_cell_id_h
18 
19 #include <deal.II/base/config.h>
20 
23 
24 #include <array>
25 #include <cstdint>
26 #include <iostream>
27 #include <vector>
28 
29 #ifdef DEAL_II_WITH_P4EST
31 #endif
32 
34 
35 // Forward declarations
36 #ifndef DOXYGEN
37 template <int, int>
38 class Triangulation;
39 #endif
40 
70 class CellId
71 {
72 public:
80  using binary_type = std::array<unsigned int, 4>;
81 
93  const std::vector<std::uint8_t> &child_indices);
94 
107  const unsigned int n_child_indices,
108  const std::uint8_t * child_indices);
109 
114  CellId(const binary_type &binary_representation);
115 
120  explicit CellId(const std::string &string_representation);
121 
125  CellId();
126 
139  std::string
140  to_string() const;
141 
145  template <int dim>
147  to_binary() const;
148 
154  template <int dim, int spacedim>
157 
161  bool
162  operator==(const CellId &other) const;
163 
167  bool
168  operator!=(const CellId &other) const;
169 
175  bool
176  operator<(const CellId &other) const;
177 
181  bool
182  is_parent_of(const CellId &other) const;
183 
187  bool
188  is_ancestor_of(const CellId &other) const;
189 
195  template <class Archive>
196  void
197  serialize(Archive &ar, const unsigned int version);
198 
203  get_coarse_cell_id() const;
204 
214  get_child_indices() const;
215 
216 private:
222 
227  unsigned int n_child_indices;
228 
238 #ifdef DEAL_II_WITH_P4EST
239  std::array<std::uint8_t, internal::p4est::functions<2>::max_level>
241 #else
242  std::array<std::uint8_t, 30> child_indices;
243 #endif
244 
245  friend std::istream &
246  operator>>(std::istream &is, CellId &cid);
247  friend std::ostream &
248  operator<<(std::ostream &os, const CellId &cid);
249 };
250 
251 
252 
256 inline std::ostream &
257 operator<<(std::ostream &os, const CellId &cid)
258 {
259  os << cid.coarse_cell_id << '_' << cid.n_child_indices << ':';
260  for (unsigned int i = 0; i < cid.n_child_indices; ++i)
261  // write the child indices. because they are between 0 and 2^dim-1, they all
262  // just have one digit, so we could write them as one character
263  // objects. it's probably clearer to write them as one-digit characters
264  // starting at '0'
265  os << static_cast<unsigned char>('0' + cid.child_indices[i]);
266  return os;
267 }
268 
269 
270 
274 template <class Archive>
275 void
276 CellId::serialize(Archive &ar, const unsigned int /*version*/)
277 {
278  ar &coarse_cell_id;
279  ar &n_child_indices;
280  ar &child_indices;
281 }
282 
286 inline std::istream &
287 operator>>(std::istream &is, CellId &cid)
288 {
289  unsigned int cellid;
290  is >> cellid;
291  if (is.eof())
292  return is;
293 
294  cid.coarse_cell_id = cellid;
295  char dummy;
296  is >> dummy;
297  Assert(dummy == '_', ExcMessage("invalid CellId"));
298  is >> cid.n_child_indices;
299  is >> dummy;
300  Assert(dummy == ':', ExcMessage("invalid CellId"));
301 
302  unsigned char value;
303  for (unsigned int i = 0; i < cid.n_child_indices; ++i)
304  {
305  // read the one-digit child index (as an integer number) and
306  // convert it back into unsigned integer type
307  is >> value;
308  cid.child_indices[i] = value - '0';
309  }
310  return is;
311 }
312 
313 
314 
315 inline bool
316 CellId::operator==(const CellId &other) const
317 {
318  if (this->coarse_cell_id != other.coarse_cell_id)
319  return false;
320  if (n_child_indices != other.n_child_indices)
321  return false;
322 
323  for (unsigned int i = 0; i < n_child_indices; ++i)
324  if (child_indices[i] != other.child_indices[i])
325  return false;
326 
327  return true;
328 }
329 
330 
331 
332 inline bool
333 CellId::operator!=(const CellId &other) const
334 {
335  return !(*this == other);
336 }
337 
338 
339 
340 inline bool
341 CellId::operator<(const CellId &other) const
342 {
343  if (this->coarse_cell_id != other.coarse_cell_id)
344  return this->coarse_cell_id < other.coarse_cell_id;
345 
346  unsigned int idx = 0;
347  while (idx < n_child_indices)
348  {
349  if (idx >= other.n_child_indices)
350  return false;
351 
352  if (child_indices[idx] != other.child_indices[idx])
353  return child_indices[idx] < other.child_indices[idx];
354 
355  ++idx;
356  }
357 
358  if (n_child_indices == other.n_child_indices)
359  return false;
360  return true; // other.id is longer
361 }
362 
363 
364 
365 inline bool
366 CellId::is_parent_of(const CellId &other) const
367 {
368  if (this->coarse_cell_id != other.coarse_cell_id)
369  return false;
370 
371  if (n_child_indices + 1 != other.n_child_indices)
372  return false;
373 
374  for (unsigned int idx = 0; idx < n_child_indices; ++idx)
375  if (child_indices[idx] != other.child_indices[idx])
376  return false;
377 
378  return true; // other.id is longer
379 }
380 
381 
382 
383 inline bool
384 CellId::is_ancestor_of(const CellId &other) const
385 {
386  if (this->coarse_cell_id != other.coarse_cell_id)
387  return false;
388 
389  if (n_child_indices >= other.n_child_indices)
390  return false;
391 
392  for (unsigned int idx = 0; idx < n_child_indices; ++idx)
393  if (child_indices[idx] != other.child_indices[idx])
394  return false;
395 
396  return true; // other.id is longer
397 }
398 
399 
400 
403 {
404  return coarse_cell_id;
405 }
406 
407 
408 
411 {
412  return {child_indices.data(), n_child_indices};
413 }
414 
415 
417 
418 #endif
void serialize(Archive &ar, const unsigned int version)
Definition: cell_id.h:276
types::coarse_cell_id get_coarse_cell_id() const
Definition: cell_id.h:402
Triangulation< dim, spacedim >::cell_iterator to_cell(const Triangulation< dim, spacedim > &tria) const
Definition: cell_id.cc:167
binary_type to_binary() const
Definition: cell_id.cc:110
const ::Triangulation< dim, spacedim > & tria
bool is_parent_of(const CellId &other) const
Definition: cell_id.h:366
CellId()
Definition: cell_id.cc:24
std::array< unsigned int, 4 > binary_type
Definition: cell_id.h:80
friend std::ostream & operator<<(std::ostream &os, const CellId &cid)
Definition: cell_id.h:257
static ::ExceptionBase & ExcMessage(std::string arg1)
ArrayView< const std::uint8_t > get_child_indices() const
Definition: cell_id.h:410
#define Assert(cond, exc)
Definition: exceptions.h:1461
bool operator<(const CellId &other) const
Definition: cell_id.h:341
bool operator!=(const CellId &other) const
Definition: cell_id.h:333
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:404
friend std::istream & operator>>(std::istream &is, CellId &cid)
Definition: cell_id.h:287
bool operator==(const CellId &other) const
Definition: cell_id.h:316
std::string to_string() const
Definition: cell_id.cc:156
Definition: cell_id.h:70
bool is_ancestor_of(const CellId &other) const
Definition: cell_id.h:384
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:403
unsigned int n_child_indices
Definition: cell_id.h:227
types::coarse_cell_id coarse_cell_id
Definition: cell_id.h:221
#define DEAL_II_DEPRECATED
Definition: config.h:160
std::array< std::uint8_t, internal::p4est::functions< 2 >::max_level > child_indices
Definition: cell_id.h:240