Reference documentation for deal.II version GIT relicensing-1321-g19b0133728 2024-07-26 07:50:01+00:00
\(\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
portable_fe_evaluation.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_portable_fe_evaluation_h
16#define dealii_portable_fe_evaluation_h
17
18#include <deal.II/base/config.h>
19
21#include <deal.II/base/tensor.h>
23
27#include <deal.II/matrix_free/portable_matrix_free.templates.h>
29
30#include <Kokkos_Core.hpp>
31
33
37namespace Portable
38{
64 template <int dim,
65 int fe_degree,
66 int n_q_points_1d = fe_degree + 1,
67 int n_components_ = 1,
68 typename Number = double>
70 {
71 public:
75 using value_type = Number;
76
81
86
90 static constexpr unsigned int dimension = dim;
91
95 static constexpr unsigned int n_components = n_components_;
96
100 static constexpr unsigned int n_q_points =
101 Utilities::pow(n_q_points_1d, dim);
102
106 static constexpr unsigned int tensor_dofs_per_cell =
107 Utilities::pow(fe_degree + 1, dim);
108
114
124 read_dof_values(const Number *src);
125
133 distribute_local_to_global(Number *dst) const;
134
143 evaluate(const EvaluationFlags::EvaluationFlags evaluate_flag);
144
152 DEAL_II_DEPRECATED_WITH_COMMENT("Use the version taking EvaluationFlags.")
154 void
155 evaluate(const bool evaluate_val, const bool evaluate_grad);
156
165 integrate(const EvaluationFlags::EvaluationFlags integration_flag);
166
174 DEAL_II_DEPRECATED_WITH_COMMENT("Use the version taking EvaluationFlags.")
176 void
177 integrate(const bool integrate_val, const bool integrate_grad);
178
184 get_value(int q_point) const;
185
191 get_dof_value(int q_point) const;
192
198 submit_value(const value_type &val_in, int q_point);
199
205 submit_dof_value(const value_type &val_in, int q_point);
206
212 get_gradient(int q_point) const;
213
219 submit_gradient(const gradient_type &grad_in, int q_point);
220
221 // clang-format off
232 // clang-format on
233 template <typename Functor>
235 apply_for_each_quad_point(const Functor &func);
236
237 private:
239 SharedData<dim, Number> *shared_data;
241 };
242
243
244
245 template <int dim,
246 int fe_degree,
247 int n_q_points_1d,
248 int n_components_,
249 typename Number>
251 FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
252 FEEvaluation(const data_type *data, SharedData<dim, Number> *shdata)
253 : data(data)
254 , shared_data(shdata)
255 , cell_id(shared_data->team_member.league_rank())
256 {}
257
258
259
260 template <int dim,
261 int fe_degree,
262 int n_q_points_1d,
263 int n_components_,
264 typename Number>
267 read_dof_values(const Number *src)
268 {
269 static_assert(n_components_ == 1, "This function only supports FE with one \
270 components");
271 // Populate the scratch memory
272 Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
273 n_q_points),
274 [&](const int &i) {
275 shared_data->values(i) =
276 src[data->local_to_global(cell_id, i)];
277 });
278 shared_data->team_member.team_barrier();
279
280 internal::resolve_hanging_nodes<dim, fe_degree, false>(
281 shared_data->team_member,
282 data->constraint_weights,
283 data->constraint_mask(cell_id),
284 shared_data->values);
285 }
286
287
288
289 template <int dim,
290 int fe_degree,
291 int n_q_points_1d,
292 int n_components_,
293 typename Number>
296 distribute_local_to_global(Number *dst) const
297 {
298 static_assert(n_components_ == 1, "This function only supports FE with one \
299 components");
300
301 internal::resolve_hanging_nodes<dim, fe_degree, true>(
302 shared_data->team_member,
303 data->constraint_weights,
304 data->constraint_mask(cell_id),
305 shared_data->values);
306
307 if (data->use_coloring)
308 {
309 Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
310 n_q_points),
311 [&](const int &i) {
312 dst[data->local_to_global(cell_id, i)] +=
313 shared_data->values(i);
314 });
315 }
316 else
317 {
318 Kokkos::parallel_for(
319 Kokkos::TeamThreadRange(shared_data->team_member, n_q_points),
320 [&](const int &i) {
321 Kokkos::atomic_add(&dst[data->local_to_global(cell_id, i)],
322 shared_data->values(i));
323 });
324 }
325 }
326
327
328
329 template <int dim,
330 int fe_degree,
331 int n_q_points_1d,
332 int n_components_,
333 typename Number>
336 const EvaluationFlags::EvaluationFlags evaluate_flag)
337 {
338 // First evaluate the gradients because it requires values that will be
339 // changed if evaluate_val is true
342 dim,
343 fe_degree,
344 n_q_points_1d,
345 Number>
346 evaluator_tensor_product(shared_data->team_member,
347 data->shape_values,
348 data->shape_gradients,
349 data->co_shape_gradients);
350
351 if ((evaluate_flag & EvaluationFlags::values) &&
352 (evaluate_flag & EvaluationFlags::gradients))
353 {
354 evaluator_tensor_product.evaluate_values_and_gradients(
355 shared_data->values, shared_data->gradients);
356 shared_data->team_member.team_barrier();
357 }
358 else if (evaluate_flag & EvaluationFlags::gradients)
359 {
360 evaluator_tensor_product.evaluate_gradients(shared_data->values,
361 shared_data->gradients);
362 shared_data->team_member.team_barrier();
363 }
364 else if (evaluate_flag & EvaluationFlags::values)
365 {
366 evaluator_tensor_product.evaluate_values(shared_data->values);
367 shared_data->team_member.team_barrier();
368 }
369 }
370
371
372
373 template <int dim,
374 int fe_degree,
375 int n_q_points_1d,
376 int n_components_,
377 typename Number>
380 const bool evaluate_val,
381 const bool evaluate_grad)
382 {
383 evaluate(
386 }
387
388
389
390 template <int dim,
391 int fe_degree,
392 int n_q_points_1d,
393 int n_components_,
394 typename Number>
397 const EvaluationFlags::EvaluationFlags integration_flag)
398 {
401 dim,
402 fe_degree,
403 n_q_points_1d,
404 Number>
405 evaluator_tensor_product(shared_data->team_member,
406 data->shape_values,
407 data->shape_gradients,
408 data->co_shape_gradients);
409
410 if ((integration_flag & EvaluationFlags::values) &&
411 (integration_flag & EvaluationFlags::gradients))
412 {
413 evaluator_tensor_product.integrate_values_and_gradients(
414 shared_data->values, shared_data->gradients);
415 }
416 else if (integration_flag & EvaluationFlags::values)
417 {
418 evaluator_tensor_product.integrate_values(shared_data->values);
419 shared_data->team_member.team_barrier();
420 }
421 else if (integration_flag & EvaluationFlags::gradients)
422 {
423 evaluator_tensor_product.template integrate_gradients<false>(
424 shared_data->values, shared_data->gradients);
425 shared_data->team_member.team_barrier();
426 }
427 }
428
429
430
431 template <int dim,
432 int fe_degree,
433 int n_q_points_1d,
434 int n_components_,
435 typename Number>
438 const bool integrate_val,
439 const bool integrate_grad)
440 {
441 integrate(
444 }
445
446
447
448 template <int dim,
449 int fe_degree,
450 int n_q_points_1d,
451 int n_components_,
452 typename Number>
454 fe_degree,
455 n_q_points_1d,
456 n_components_,
457 Number>::value_type
459 int q_point) const
460 {
461 return shared_data->values(q_point);
462 }
463
464
465
466 template <int dim,
467 int fe_degree,
468 int n_q_points_1d,
469 int n_components_,
470 typename Number>
472 fe_degree,
473 n_q_points_1d,
474 n_components_,
475 Number>::value_type
477 get_dof_value(int q_point) const
478 {
479 return shared_data->values(q_point);
480 }
481
482
483
484 template <int dim,
485 int fe_degree,
486 int n_q_points_1d,
487 int n_components_,
488 typename Number>
491 submit_value(const value_type &val_in, int q_point)
492 {
493 shared_data->values(q_point) = val_in * data->JxW(cell_id, q_point);
494 }
495
496
497
498 template <int dim,
499 int fe_degree,
500 int n_q_points_1d,
501 int n_components_,
502 typename Number>
505 submit_dof_value(const value_type &val_in, int q_point)
506 {
507 shared_data->values(q_point) = val_in;
508 }
509
510
511
512 template <int dim,
513 int fe_degree,
514 int n_q_points_1d,
515 int n_components_,
516 typename Number>
518 fe_degree,
519 n_q_points_1d,
520 n_components_,
521 Number>::gradient_type
523 get_gradient(int q_point) const
524 {
525 static_assert(n_components_ == 1, "This function only supports FE with one \
526 components");
527
528 gradient_type grad;
529 for (unsigned int d_1 = 0; d_1 < dim; ++d_1)
530 {
531 Number tmp = 0.;
532 for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
533 tmp += data->inv_jacobian(cell_id, q_point, d_2, d_1) *
534 shared_data->gradients(q_point, d_2);
535 grad[d_1] = tmp;
536 }
537
538 return grad;
539 }
540
541
542
543 template <int dim,
544 int fe_degree,
545 int n_q_points_1d,
546 int n_components_,
547 typename Number>
550 submit_gradient(const gradient_type &grad_in, int q_point)
551 {
552 for (unsigned int d_1 = 0; d_1 < dim; ++d_1)
553 {
554 Number tmp = 0.;
555 for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
556 tmp += data->inv_jacobian(cell_id, q_point, d_1, d_2) * grad_in[d_2];
557 shared_data->gradients(q_point, d_1) =
558 tmp * data->JxW(cell_id, q_point);
559 }
560 }
561
562
563
564 template <int dim,
565 int fe_degree,
566 int n_q_points_1d,
567 int n_components_,
568 typename Number>
569 template <typename Functor>
572 apply_for_each_quad_point(const Functor &func)
573 {
574 Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
575 n_q_points),
576 [&](const int &i) { func(this, i); });
577 shared_data->team_member.team_barrier();
578 }
579
580
581
582#ifndef DOXYGEN
583 template <int dim,
584 int fe_degree,
585 int n_q_points_1d,
586 int n_components_,
587 typename Number>
588 constexpr unsigned int
591#endif
592} // namespace Portable
593
595
596#endif
static constexpr unsigned int n_q_points
SharedData< dim, Number > * shared_data
void integrate(const EvaluationFlags::EvaluationFlags integration_flag)
static constexpr unsigned int dimension
static constexpr unsigned int n_components
value_type get_value(int q_point) const
static constexpr unsigned int tensor_dofs_per_cell
value_type get_dof_value(int q_point) const
Tensor< 1, dim, Number > gradient_type
void evaluate(const EvaluationFlags::EvaluationFlags evaluate_flag)
void distribute_local_to_global(Number *dst) const
gradient_type get_gradient(int q_point) const
typename MatrixFree< dim, Number >::Data data_type
void read_dof_values(const Number *src)
void submit_dof_value(const value_type &val_in, int q_point)
void submit_value(const value_type &val_in, int q_point)
void submit_gradient(const gradient_type &grad_in, int q_point)
void apply_for_each_quad_point(const Functor &func)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_DEPRECATED_WITH_COMMENT(comment)
Definition config.h:208
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
The namespace for the EvaluationFlags enum.
EvaluationFlags
The EvaluationFlags enum.
constexpr T pow(const T base, const int iexp)
Definition utilities.h:966
#define DEAL_II_HOST_DEVICE
Definition numbers.h:34