MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
face_map_utils.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2024, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#ifndef MFEM_FACE_MAP_UTILS_HPP
13#define MFEM_FACE_MAP_UTILS_HPP
14
17#include <utility> // std::pair
18#include <vector>
19
20namespace mfem
21{
22
23namespace internal
24{
25
26/// Each face of a hexahedron is given by a level set x_i = l, where x_i is one
27/// of x, y, or z (corresponding to i=0, i=1, i=2), and l is either 0 or 1.
28/// Returns i and level.
29std::pair<int,int> GetFaceNormal3D(const int face_id);
30
31/// @brief Fills in the entries of the lexicographic face_map.
32///
33/// For use in FiniteElement::GetFaceMap.
34///
35/// n_face_dofs_per_component is the number of DOFs for each vector component
36/// on the face (there is only one vector component in all cases except for 3D
37/// Nedelec elements, where the face DOFs have two components to span the
38/// tangent space).
39///
40/// The DOFs for the i-th vector component begin at offsets[i] (i.e. the number
41/// of vector components is given by offsets.size()).
42///
43/// The DOFs for each vector component are arranged in a Cartesian grid defined
44/// by strides and n_dofs_per_dim.
45void FillFaceMap(const int n_face_dofs_per_component,
46 const std::vector<int> &offsets,
47 const std::vector<int> &strides,
48 const std::vector<int> &n_dofs_per_dim,
49 Array<int> &face_map);
50
51/// Return the face map for nodal tensor elements (H1, L2, and Bernstein basis).
52void GetTensorFaceMap(const int dim, const int order, const int face_id,
53 Array<int> &face_map);
54
55/// @brief Given a face DOF index in native (counter-clockwise) ordering, return
56/// the corresponding DOF index in lexicographic ordering (for a quadrilateral
57/// element).
58MFEM_HOST_DEVICE
59inline int ToLexOrdering2D(const int face_id, const int size1d, const int i)
60{
61 if (face_id==2 || face_id==3)
62 {
63 return size1d-1-i;
64 }
65 else
66 {
67 return i;
68 }
69}
70
71/// @brief Given a face DOF index on a shared face, ordered lexicographically
72/// relative to element 1, return the corresponding face DOF index ordered
73/// lexicographically relative to element 2.
74MFEM_HOST_DEVICE
75inline int PermuteFace2D(const int face_id1, const int face_id2,
76 const int orientation, const int size1d,
77 const int index)
78{
79 int new_index;
80 // Convert from element 1 lex ordering to native ordering
81 if (face_id1 == 2 || face_id1 == 3)
82 {
83 new_index = size1d-1-index;
84 }
85 else
86 {
87 new_index = index;
88 }
89 // Permute based on face orientations
90 if (orientation == 1)
91 {
92 new_index = size1d-1-new_index;
93 }
94 // Covert to element 2 lex ordering
95 return ToLexOrdering2D(face_id2, size1d, new_index);
96}
97
98/// @brief Given a face DOF index in native (counter-clockwise) ordering, return
99/// the corresponding DOF index in lexicographic ordering (for a hexahedral
100/// element).
101MFEM_HOST_DEVICE
102inline int ToLexOrdering3D(const int face_id, const int size1d, const int i,
103 const int j)
104{
105 if (face_id==2 || face_id==1 || face_id==5)
106 {
107 return i + j*size1d;
108 }
109 else if (face_id==3 || face_id==4)
110 {
111 return (size1d-1-i) + j*size1d;
112 }
113 else // face_id==0
114 {
115 return i + (size1d-1-j)*size1d;
116 }
117}
118
119/// @brief Given the index of a face DOF in lexicographic ordering relative
120/// element 1, permute the index so that it is lexicographically ordered
121/// relative to element 2.
122///
123/// The given face corresponds to local face index @a face_id1 relative to
124/// element 1, and @a face_id2 (with @a orientation) relative to element 2.
125MFEM_HOST_DEVICE
126inline int PermuteFace3D(const int face_id1, const int face_id2,
127 const int orientation,
128 const int size1d, const int index)
129{
130 int i=0, j=0, new_i=0, new_j=0;
131 i = index%size1d;
132 j = index/size1d;
133 // Convert from lex ordering
134 if (face_id1==3 || face_id1==4)
135 {
136 i = size1d-1-i;
137 }
138 else if (face_id1==0)
139 {
140 j = size1d-1-j;
141 }
142 // Permute based on face orientations
143 switch (orientation)
144 {
145 case 0:
146 new_i = i;
147 new_j = j;
148 break;
149 case 1:
150 new_i = j;
151 new_j = i;
152 break;
153 case 2:
154 new_i = j;
155 new_j = (size1d-1-i);
156 break;
157 case 3:
158 new_i = (size1d-1-i);
159 new_j = j;
160 break;
161 case 4:
162 new_i = (size1d-1-i);
163 new_j = (size1d-1-j);
164 break;
165 case 5:
166 new_i = (size1d-1-j);
167 new_j = (size1d-1-i);
168 break;
169 case 6:
170 new_i = (size1d-1-j);
171 new_j = i;
172 break;
173 case 7:
174 new_i = i;
175 new_j = (size1d-1-j);
176 break;
177 }
178 return ToLexOrdering3D(face_id2, size1d, new_i, new_j);
179}
180
181/// @brief Given a face DOF (or quadrature) index ordered lexicographically
182/// relative to element 1, return the associated (i, j) coordinates.
183///
184/// The returned coordinates will be relative to element 1 or element 2
185/// according to the value of side (side == 0 corresponds element 1).
186MFEM_HOST_DEVICE
187inline void FaceIdxToVolIdx2D(const int qi, const int nq, const int face_id0,
188 const int face_id1, const int side, int &i, int &j)
189{
190 // Note: in 2D, a consistently ordered mesh will always have the element 2
191 // face reversed relative to element 1, so orientation is determined entirely
192 // by side. (In 3D, separate orientation information is needed).
193 const int orientation = side;
194
195 const int face_id = (side == 0) ? face_id0 : face_id1;
196 const int edge_idx = (side == 0) ? qi : PermuteFace2D(face_id0, face_id1,
197 orientation, nq, qi);
198
199 const int level = (face_id == 0 || face_id == 3) ? 0 : (nq-1);
200 const bool x_axis = (face_id == 0 || face_id == 2);
201
202 i = x_axis ? edge_idx : level;
203 j = x_axis ? level : edge_idx;
204}
205
206/// @brief Given a face DOF (or quadrature) index ordered lexicographically
207/// relative to element 1, return the associated (i, j, k) coordinates.
208///
209/// The returned coordinates will be relative to element 1 or element 2
210/// according to the value of side (side == 0 corresponds element 1).
211MFEM_HOST_DEVICE
212inline void FaceIdxToVolIdx3D(const int index, const int size1d,
213 const int face_id0, const int face_id1,
214 const int side, const int orientation,
215 int& i, int& j, int& k)
216{
217 MFEM_VERIFY_KERNEL(face_id1 >= 0 || side == 0,
218 "Accessing second side but face_id1 is not valid.");
219
220 const int face_id = (side == 0) ? face_id0 : face_id1;
221 const int fidx = (side == 0) ? index
222 : PermuteFace3D(face_id0, face_id1, orientation, size1d, index);
223
224 const bool xy_plane = (face_id == 0 || face_id == 5);
225 const bool yz_plane = (face_id == 2 || face_id == 4);
226
227 const int level = (face_id == 0 || face_id == 1 || face_id == 4)
228 ? 0 : (size1d-1);
229
230 const int _i = fidx % size1d;
231 const int _j = fidx / size1d;
232
233 k = xy_plane ? level : _j;
234 j = yz_plane ? _i : xy_plane ? _j : level;
235 i = yz_plane ? level : _i;
236}
237
238} // namespace internal
239
240} // namespace mfem
241
242#endif
int dim
Definition ex24.cpp:53
int index(int i, int j, int nx, int ny)
Definition life.cpp:235