16#ifndef dealii_cuda_hanging_nodes_internal_h
17#define dealii_cuda_hanging_nodes_internal_h
21#ifdef DEAL_II_COMPILER_CUDA_AWARE
39 template <
unsigned int size>
40 __device__
inline unsigned int
41 index2(
unsigned int i,
unsigned int j)
48 template <
unsigned int size>
49 __device__
inline unsigned int
50 index3(
unsigned int i,
unsigned int j,
unsigned int k)
52 return i + size * j + size * size * k;
57 template <
unsigned int fe_degree,
58 unsigned int direction,
61 __device__
inline void
67 const unsigned int x_idx = threadIdx.x % (fe_degree + 1);
68 const unsigned int y_idx = threadIdx.y;
70 const auto this_type =
72 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_x :
73 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_y;
75 const unsigned int interp_idx = (direction == 0) ? x_idx : y_idx;
80 const bool constrained_face =
83 ::internal::MatrixFreeFunctions::ConstraintKinds::face_y :
87 ::internal::MatrixFreeFunctions::ConstraintKinds::face_x :
90 ::internal::MatrixFreeFunctions::ConstraintKinds::unconstrained;
95 const bool constrained_dof =
97 (((constraint_mask & ::internal::MatrixFreeFunctions::
98 ConstraintKinds::subcell_y) !=
99 ::internal::MatrixFreeFunctions::ConstraintKinds::
102 (y_idx == fe_degree))) ||
104 (((constraint_mask & ::internal::MatrixFreeFunctions::
105 ConstraintKinds::subcell_x) !=
109 (x_idx == fe_degree)));
111 if (constrained_face && constrained_dof)
113 const bool type = (constraint_mask & this_type) !=
114 ::internal::MatrixFreeFunctions::
115 ConstraintKinds::unconstrained;
119 for (
unsigned int i = 0; i <= fe_degree; ++i)
121 const unsigned int real_idx =
122 (direction == 0) ? index2<fe_degree + 1>(i, y_idx) :
123 index2<fe_degree + 1>(x_idx, i);
134 for (
unsigned int i = 0; i <= fe_degree; ++i)
136 const unsigned int real_idx =
137 (direction == 0) ? index2<fe_degree + 1>(i, y_idx) :
138 index2<fe_degree + 1>(x_idx, i);
143 fe_degree - interp_idx] :
155 if (constrained_face && constrained_dof)
156 values[index2<fe_degree + 1>(x_idx, y_idx)] = t;
163 template <
unsigned int fe_degree,
164 unsigned int direction,
167 __device__
inline void
173 const unsigned int x_idx = threadIdx.x % (fe_degree + 1);
174 const unsigned int y_idx = threadIdx.y;
175 const unsigned int z_idx = threadIdx.z;
177 const auto this_type =
179 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_x :
181 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_y :
182 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_z;
183 const auto face1_type =
185 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_y :
187 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_z :
188 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_x;
189 const auto face2_type =
191 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_z :
193 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_x :
194 ::internal::MatrixFreeFunctions::ConstraintKinds::subcell_y;
200 ::internal::MatrixFreeFunctions::ConstraintKinds::face_y :
202 ::internal::MatrixFreeFunctions::ConstraintKinds::face_z :
203 ::internal::MatrixFreeFunctions::ConstraintKinds::face_x;
206 ::internal::MatrixFreeFunctions::ConstraintKinds::face_z :
208 ::internal::MatrixFreeFunctions::ConstraintKinds::face_x :
209 ::internal::MatrixFreeFunctions::ConstraintKinds::face_y;
212 ::internal::MatrixFreeFunctions::ConstraintKinds::edge_x :
214 ::internal::MatrixFreeFunctions::ConstraintKinds::edge_y :
215 ::internal::MatrixFreeFunctions::ConstraintKinds::edge_z;
216 const auto constrained_face = constraint_mask & (face1 | face2 | edge);
218 const unsigned int interp_idx = (direction == 0) ? x_idx :
219 (direction == 1) ? y_idx :
221 const unsigned int face1_idx = (direction == 0) ? y_idx :
222 (direction == 1) ? z_idx :
224 const unsigned int face2_idx = (direction == 0) ? z_idx :
225 (direction == 1) ? x_idx :
229 const bool on_face1 = ((constraint_mask & face1_type) !=
230 ::internal::MatrixFreeFunctions::
231 ConstraintKinds::unconstrained) ?
233 (face1_idx == fe_degree);
234 const bool on_face2 = ((constraint_mask & face2_type) !=
235 ::internal::MatrixFreeFunctions::
236 ConstraintKinds::unconstrained) ?
238 (face2_idx == fe_degree);
239 const bool constrained_dof =
240 ((((constraint_mask & face1) != ::internal::MatrixFreeFunctions::
241 ConstraintKinds::unconstrained) &&
243 (((constraint_mask & face2) != ::internal::MatrixFreeFunctions::
244 ConstraintKinds::unconstrained) &&
246 (((constraint_mask & edge) != ::internal::MatrixFreeFunctions::
247 ConstraintKinds::unconstrained) &&
248 on_face1 && on_face2));
250 if ((constrained_face != ::internal::MatrixFreeFunctions::
251 ConstraintKinds::unconstrained) &&
254 const bool type = (constraint_mask & this_type) !=
255 ::internal::MatrixFreeFunctions::
256 ConstraintKinds::unconstrained;
259 for (
unsigned int i = 0; i <= fe_degree; ++i)
261 const unsigned int real_idx =
262 (direction == 0) ? index3<fe_degree + 1>(i, y_idx, z_idx) :
263 (direction == 1) ? index3<fe_degree + 1>(x_idx, i, z_idx) :
264 index3<fe_degree + 1>(x_idx, y_idx, i);
275 for (
unsigned int i = 0; i <= fe_degree; ++i)
277 const unsigned int real_idx =
278 (direction == 0) ? index3<fe_degree + 1>(i, y_idx, z_idx) :
279 (direction == 1) ? index3<fe_degree + 1>(x_idx, i, z_idx) :
280 index3<fe_degree + 1>(x_idx, y_idx, i);
285 fe_degree - interp_idx] :
298 if ((constrained_face != ::internal::MatrixFreeFunctions::
299 ConstraintKinds::unconstrained) &&
301 values[index3<fe_degree + 1>(x_idx, y_idx, z_idx)] = t;
315 template <
int dim,
int fe_degree,
bool transpose,
typename Number>
324 interpolate_boundary_2d<fe_degree, 0, transpose>(constraint_mask,
327 interpolate_boundary_2d<fe_degree, 1, transpose>(constraint_mask,
333 interpolate_boundary_3d<fe_degree, 0, transpose>(constraint_mask,
336 interpolate_boundary_3d<fe_degree, 1, transpose>(constraint_mask,
339 interpolate_boundary_3d<fe_degree, 2, transpose>(constraint_mask,
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
void interpolate_boundary_2d(const ::internal::MatrixFreeFunctions::ConstraintKinds constraint_mask, Number *values)
void resolve_hanging_nodes(const ::internal::MatrixFreeFunctions::ConstraintKinds constraint_mask, Number *values)
void interpolate_boundary_3d(const ::internal::MatrixFreeFunctions::ConstraintKinds constraint_mask, Number *values)
unsigned int index3(unsigned int i, unsigned int j, unsigned int k)
__constant__ double constraint_weights[(CUDAWrappers::mf_max_elem_degree+1) *(CUDAWrappers::mf_max_elem_degree+1)]
unsigned int index2(unsigned int i, unsigned int j)
constexpr unsigned int mf_max_elem_degree
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)