Reference documentation for deal.II version Git 3d170a2f6b 2020-04-03 17:03:45 -0400
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
cell_id.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2018 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 
21 #include <deal.II/base/exceptions.h>
22 
23 #include <array>
24 #include <cstdint>
25 #include <iostream>
26 #include <vector>
27 
28 #ifdef DEAL_II_WITH_P4EST
29 # include <deal.II/distributed/p4est_wrappers.h>
30 #endif
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 // Forward declarations
35 #ifndef DOXYGEN
36 template <int, int>
37 class Triangulation;
38 #endif
39 
68 class CellId
69 {
70 public:
78  using binary_type = std::array<unsigned int, 4>;
79 
91  const std::vector<std::uint8_t> &child_indices);
92 
105  const unsigned int n_child_indices,
106  const std::uint8_t * child_indices);
107 
112  CellId(const binary_type &binary_representation);
113 
117  CellId();
118 
131  std::string
132  to_string() const;
133 
137  template <int dim>
139  to_binary() const;
140 
144  template <int dim, int spacedim>
146  to_cell(const Triangulation<dim, spacedim> &tria) const;
147 
151  bool
152  operator==(const CellId &other) const;
153 
157  bool
158  operator!=(const CellId &other) const;
159 
165  bool
166  operator<(const CellId &other) const;
167 
171  bool
172  is_parent_of(const CellId &other) const;
173 
177  bool
178  is_ancestor_of(const CellId &other) const;
179 
183  template <class Archive>
184  void
185  serialize(Archive &ar, const unsigned int version);
186 
191  get_coarse_cell_id() const;
192 
193 private:
199 
204  unsigned int n_child_indices;
205 
215 #ifdef DEAL_II_WITH_P4EST
216  std::array<std::uint8_t, internal::p4est::functions<2>::max_level>
218 #else
219  std::array<std::uint8_t, 30> child_indices;
220 #endif
221 
222  friend std::istream &
223  operator>>(std::istream &is, CellId &cid);
224  friend std::ostream &
225  operator<<(std::ostream &os, const CellId &cid);
226 };
227 
228 
229 
233 inline std::ostream &
234 operator<<(std::ostream &os, const CellId &cid)
235 {
236  os << cid.coarse_cell_id << '_' << cid.n_child_indices << ':';
237  for (unsigned int i = 0; i < cid.n_child_indices; ++i)
238  // write the child indices. because they are between 0 and 2^dim-1, they all
239  // just have one digit, so we could write them as one character
240  // objects. it's probably clearer to write them as one-digit characters
241  // starting at '0'
242  os << static_cast<unsigned char>('0' + cid.child_indices[i]);
243  return os;
244 }
245 
246 
247 
251 template <class Archive>
252 void
253 CellId::serialize(Archive &ar, const unsigned int /*version*/)
254 {
255  ar &coarse_cell_id;
256  ar &n_child_indices;
257  ar &child_indices;
258 }
259 
263 inline std::istream &
264 operator>>(std::istream &is, CellId &cid)
265 {
266  unsigned int cellid;
267  is >> cellid;
268  if (is.eof())
269  return is;
270 
271  cid.coarse_cell_id = cellid;
272  char dummy;
273  is >> dummy;
274  Assert(dummy == '_', ExcMessage("invalid CellId"));
275  is >> cid.n_child_indices;
276  is >> dummy;
277  Assert(dummy == ':', ExcMessage("invalid CellId"));
278 
279  unsigned char value;
280  for (unsigned int i = 0; i < cid.n_child_indices; ++i)
281  {
282  // read the one-digit child index (as an integer number) and
283  // convert it back into unsigned integer type
284  is >> value;
285  cid.child_indices[i] = value - '0';
286  }
287  return is;
288 }
289 
290 
291 
292 inline bool
293 CellId::operator==(const CellId &other) const
294 {
295  if (this->coarse_cell_id != other.coarse_cell_id)
296  return false;
297  if (n_child_indices != other.n_child_indices)
298  return false;
299 
300  for (unsigned int i = 0; i < n_child_indices; ++i)
301  if (child_indices[i] != other.child_indices[i])
302  return false;
303 
304  return true;
305 }
306 
307 
308 
309 inline bool
310 CellId::operator!=(const CellId &other) const
311 {
312  return !(*this == other);
313 }
314 
315 
316 
317 inline bool
318 CellId::operator<(const CellId &other) const
319 {
320  if (this->coarse_cell_id != other.coarse_cell_id)
321  return this->coarse_cell_id < other.coarse_cell_id;
322 
323  unsigned int idx = 0;
324  while (idx < n_child_indices)
325  {
326  if (idx >= other.n_child_indices)
327  return false;
328 
329  if (child_indices[idx] != other.child_indices[idx])
330  return child_indices[idx] < other.child_indices[idx];
331 
332  ++idx;
333  }
334 
335  if (n_child_indices == other.n_child_indices)
336  return false;
337  return true; // other.id is longer
338 }
339 
340 
341 
342 inline bool
343 CellId::is_parent_of(const CellId &other) const
344 {
345  if (this->coarse_cell_id != other.coarse_cell_id)
346  return false;
347 
348  if (n_child_indices + 1 != other.n_child_indices)
349  return false;
350 
351  for (unsigned int idx = 0; idx < n_child_indices; ++idx)
352  if (child_indices[idx] != other.child_indices[idx])
353  return false;
354 
355  return true; // other.id is longer
356 }
357 
358 
359 
360 inline bool
361 CellId::is_ancestor_of(const CellId &other) const
362 {
363  if (this->coarse_cell_id != other.coarse_cell_id)
364  return false;
365 
366  if (n_child_indices >= other.n_child_indices)
367  return false;
368 
369  for (unsigned int idx = 0; idx < n_child_indices; ++idx)
370  if (child_indices[idx] != other.child_indices[idx])
371  return false;
372 
373  return true; // other.id is longer
374 }
375 
376 
379 {
380  return coarse_cell_id;
381 }
382 
383 
384 DEAL_II_NAMESPACE_CLOSE
385 
386 #endif
void serialize(Archive &ar, const unsigned int version)
Definition: cell_id.h:253
types::coarse_cell_id get_coarse_cell_id() const
Definition: cell_id.h:378
Triangulation< dim, spacedim >::cell_iterator to_cell(const Triangulation< dim, spacedim > &tria) const
Definition: cell_id.cc:159
binary_type to_binary() const
Definition: cell_id.cc:102
bool is_parent_of(const CellId &other) const
Definition: cell_id.h:343
CellId()
Definition: cell_id.cc:24
std::array< unsigned int, 4 > binary_type
Definition: cell_id.h:78
friend std::ostream & operator<<(std::ostream &os, const CellId &cid)
Definition: cell_id.h:234
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1419
bool operator<(const CellId &other) const
Definition: cell_id.h:318
bool operator!=(const CellId &other) const
Definition: cell_id.h:310
friend std::istream & operator>>(std::istream &is, CellId &cid)
Definition: cell_id.h:264
bool operator==(const CellId &other) const
Definition: cell_id.h:293
std::string to_string() const
Definition: cell_id.cc:148
Definition: cell_id.h:68
bool is_ancestor_of(const CellId &other) const
Definition: cell_id.h:361
unsigned int n_child_indices
Definition: cell_id.h:204
types::coarse_cell_id coarse_cell_id
Definition: cell_id.h:198
global_cell_index coarse_cell_id
Definition: types.h:113
std::array< std::uint8_t, internal::p4est::functions< 2 >::max_level > child_indices
Definition: cell_id.h:217