Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
grid_reordering.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2000 - 2023 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
17#include <deal.II/base/timer.h>
19
22
23#include <algorithm>
24#include <fstream>
25#include <functional>
26#include <iostream>
27
29
30
31// anonymous namespace for internal helper functions
32namespace
33{
43 constexpr std::array<unsigned int, 8> local_vertex_numbering{
44 {0, 1, 5, 4, 2, 3, 7, 6}};
45
52 void
53 reorder_new_to_old_style(std::vector<CellData<1>> &)
54 {}
55
56
57 void
58 reorder_new_to_old_style(std::vector<CellData<2>> &cells)
59 {
60 for (auto &cell : cells)
61 std::swap(cell.vertices[2], cell.vertices[3]);
62 }
63
64
65 void
66 reorder_new_to_old_style(std::vector<CellData<3>> &cells)
67 {
68 unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
69 for (auto &cell : cells)
70 {
71 for (const unsigned int i : GeometryInfo<3>::vertex_indices())
72 tmp[i] = cell.vertices[i];
73 for (const unsigned int i : GeometryInfo<3>::vertex_indices())
74 cell.vertices[i] = tmp[local_vertex_numbering[i]];
75 }
76 }
77
78
82 void
83 reorder_old_to_new_style(std::vector<CellData<1>> &)
84 {}
85
86
87 void
88 reorder_old_to_new_style(std::vector<CellData<2>> &cells)
89 {
90 // just invert the permutation:
91 reorder_new_to_old_style(cells);
92 }
93
94
95 void
96 reorder_old_to_new_style(std::vector<CellData<3>> &cells)
97 {
98 // undo the ordering above
99 unsigned int tmp[GeometryInfo<3>::vertices_per_cell];
100 for (auto &cell : cells)
101 {
102 for (const unsigned int i : GeometryInfo<3>::vertex_indices())
103 tmp[i] = cell.vertices[i];
104 for (const unsigned int i : GeometryInfo<3>::vertex_indices())
105 cell.vertices[local_vertex_numbering[i]] = tmp[i];
106 }
107 }
108} // namespace
109
110
111
112template <int dim, int spacedim>
113void
115 const bool use_new_style_ordering)
116{
117 Assert(cells.size() != 0,
118 ExcMessage("List of elements to orient must have at least one cell"));
119
120 // there is nothing for us to do in 1d
121 if (dim == 1)
122 return;
123
124 // if necessary, convert to new-style format
125 if (use_new_style_ordering == false)
126 reorder_old_to_new_style(cells);
127
129
130 // and convert back if necessary
131 if (use_new_style_ordering == false)
132 reorder_new_to_old_style(cells);
133}
134
135
136
137template <>
138void
140 const std::vector<Point<1>> &,
141 std::vector<CellData<1>> &,
142 const bool)
143{
144 // nothing to be done in 1d
145}
146
147
148
149template <>
150void
152 const std::vector<Point<2>> &,
153 std::vector<CellData<1>> &,
154 const bool)
155{
156 // nothing to be done in 1d
157}
158
159
160
161template <>
162void
164 const std::vector<Point<3>> &,
165 std::vector<CellData<1>> &,
166 const bool)
167{
168 // nothing to be done in 1d
169}
170
171
172template <>
173void
175 const std::vector<Point<2>> &all_vertices,
176 std::vector<CellData<2>> & cells,
177 const bool use_new_style_ordering)
178{
179 if (!use_new_style_ordering)
180 reorder_old_to_new_style(cells);
181
183
184 if (!use_new_style_ordering)
185 reorder_new_to_old_style(cells);
186}
187
188
189
190template <>
191void
193 const std::vector<Point<3>> &,
194 std::vector<CellData<2>> &,
195 const bool)
196{
197 Assert(false, ExcNotImplemented());
198}
199
200
201
202template <>
203void
205 const std::vector<Point<3>> &all_vertices,
206 std::vector<CellData<3>> & cells,
207 const bool use_new_style_ordering)
208{
209 if (!use_new_style_ordering)
210 reorder_old_to_new_style(cells);
211
213
214 if (!use_new_style_ordering)
215 reorder_new_to_old_style(cells);
216}
217
218
219
220/* ------------------------ explicit instantiations ------------------- */
221template class GridReordering<1, 1>;
222template class GridReordering<1, 2>;
223template class GridReordering<1, 3>;
224template class GridReordering<2, 2>;
225template class GridReordering<2, 3>;
226template class GridReordering<3, 3>;
227
static void reorder_cells(std::vector< CellData< dim > > &original_cells, const bool use_new_style_ordering=false)
static void invert_all_cells_of_negative_grid(const std::vector< Point< spacedim > > &all_vertices, std::vector< CellData< dim > > &original_cells, const bool use_new_style_ordering=false)
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
unsigned int vertex_indices[2]
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
void consistently_order_cells(std::vector< CellData< dim > > &cells)
void invert_all_negative_measure_cells(const std::vector< Point< spacedim > > &all_vertices, std::vector< CellData< dim > > &cells)
STL namespace.
void swap(SmartPointer< T, P > &t1, SmartPointer< T, Q > &t2)