15#ifndef dealii_portable_fe_evaluation_h
16#define dealii_portable_fe_evaluation_h
28#include <deal.II/matrix_free/portable_matrix_free.templates.h>
31#include <Kokkos_Core.hpp>
65 int n_q_points_1d = fe_degree + 1,
66 int n_components_ = 1,
67 typename Number =
double>
74 using value_type = std::conditional_t<(n_components_ == 1),
84 std::conditional_t<n_components_ == dim,
133 const unsigned int dof_index = 0);
244 template <
typename Functor>
266 , precomputed_data(
data->precomputed_data)
267 , shared_data(
data->shared_data)
268 , cell_id(
data->team_member.league_rank())
318 Kokkos::parallel_for(Kokkos::TeamThreadRange(
data->team_member,
319 tensor_dofs_per_component),
321 for (unsigned int c = 0; c < n_components_; ++c)
322 shared_data->values(i, c) =
323 src[precomputed_data->local_to_global(
324 i + tensor_dofs_per_component * c, cell_id)];
326 data->team_member.team_barrier();
328 for (
unsigned int c = 0; c < n_components_; ++c)
330 if (precomputed_data->constraint_mask(cell_id * n_components + c) !=
333 internal::resolve_hanging_nodes<dim, fe_degree, false, Number>(
335 precomputed_data->constraint_weights,
336 precomputed_data->constraint_mask(cell_id * n_components + c),
337 Kokkos::subview(shared_data->values, Kokkos::ALL, c));
352 for (
unsigned int c = 0; c < n_components_; ++c)
354 if (precomputed_data->constraint_mask(cell_id * n_components + c) !=
357 internal::resolve_hanging_nodes<dim, fe_degree, true, Number>(
359 precomputed_data->constraint_weights,
360 precomputed_data->constraint_mask(cell_id * n_components + c),
361 Kokkos::subview(shared_data->values, Kokkos::ALL, c));
364 if (precomputed_data->use_coloring)
366 Kokkos::parallel_for(
367 Kokkos::TeamThreadRange(
data->team_member, tensor_dofs_per_component),
369 for (unsigned int c = 0; c < n_components_; ++c)
370 dst[precomputed_data->local_to_global(
371 i + tensor_dofs_per_component * c, cell_id)] +=
372 shared_data->values(i, c);
377 Kokkos::parallel_for(
378 Kokkos::TeamThreadRange(
data->team_member, tensor_dofs_per_component),
380 for (unsigned int c = 0; c < n_components_; ++c)
381 Kokkos::atomic_add(&dst[precomputed_data->local_to_global(
382 i + (tensor_dofs_per_component)*c, cell_id)],
383 shared_data->values(i, c));
401 if (fe_degree >= 0 && fe_degree + 1 == n_q_points_1d &&
402 precomputed_data->element_type ==
403 ElementType::tensor_symmetric_collocation)
406 n_components, evaluation_flag,
data);
410 else if (fe_degree >= 0 &&
412 precomputed_data->element_type <= ElementType::tensor_symmetric)
418 Number>::evaluate(n_components, evaluation_flag,
data);
420 else if (fe_degree >= 0 && precomputed_data->element_type <=
421 ElementType::tensor_symmetric_no_collocation)
428 Kokkos::abort(
"The element type is not yet supported by the portable "
429 "matrix-free module.");
446 if (fe_degree >= 0 && fe_degree + 1 == n_q_points_1d &&
447 precomputed_data->element_type ==
448 ElementType::tensor_symmetric_collocation)
455 else if (fe_degree >= 0 &&
457 precomputed_data->element_type <= ElementType::tensor_symmetric)
463 Number>::integrate(n_components, integration_flag,
data);
465 else if (fe_degree >= 0 && precomputed_data->element_type <=
466 ElementType::tensor_symmetric_no_collocation)
473 Kokkos::abort(
"The element type is not yet supported by the portable "
474 "matrix-free module.");
494 if constexpr (n_components_ == 1)
496 return shared_data->values(q_point, 0);
501 for (
unsigned int c = 0; c < n_components; ++c)
502 result[c] = shared_data->values(q_point, c);
523 if constexpr (n_components_ == 1)
525 return shared_data->values(dof, 0);
530 for (
unsigned int c = 0; c < n_components; ++c)
531 result[c] = shared_data->values(dof, c);
548 if constexpr (n_components_ == 1)
550 shared_data->values(q_point, 0) =
551 val_in * precomputed_data->JxW(q_point, cell_id);
555 for (
unsigned int c = 0; c < n_components; ++c)
556 shared_data->values(q_point, c) =
557 val_in[c] * precomputed_data->JxW(q_point, cell_id);
573 if constexpr (n_components_ == 1)
575 shared_data->values(dof, 0) = val_in;
579 for (
unsigned int c = 0; c < n_components; ++c)
580 shared_data->values(dof, c) = val_in[c];
595 Number>::gradient_type
602 if constexpr (n_components_ == 1)
604 for (
unsigned int d_1 = 0; d_1 < dim; ++d_1)
607 for (
unsigned int d_2 = 0; d_2 < dim; ++d_2)
609 precomputed_data->inv_jacobian(q_point, cell_id, d_2, d_1) *
610 shared_data->gradients(q_point, d_2, 0);
616 for (
unsigned int c = 0; c < n_components; ++c)
617 for (
unsigned int d_1 = 0; d_1 < dim; ++d_1)
620 for (
unsigned int d_2 = 0; d_2 < dim; ++d_2)
622 precomputed_data->inv_jacobian(q_point, cell_id, d_2, d_1) *
623 shared_data->gradients(q_point, d_2, c);
643 if constexpr (n_components_ == 1)
645 for (
unsigned int d_1 = 0; d_1 < dim; ++d_1)
648 for (
unsigned int d_2 = 0; d_2 < dim; ++d_2)
650 precomputed_data->inv_jacobian(q_point, cell_id, d_1, d_2) *
652 shared_data->gradients(q_point, d_1, 0) =
653 tmp * precomputed_data->JxW(q_point, cell_id);
658 for (
unsigned int c = 0; c < n_components; ++c)
659 for (
unsigned int d_1 = 0; d_1 < dim; ++d_1)
662 for (
unsigned int d_2 = 0; d_2 < dim; ++d_2)
664 precomputed_data->inv_jacobian(q_point, cell_id, d_1, d_2) *
666 shared_data->gradients(q_point, d_1, c) =
667 tmp * precomputed_data->JxW(q_point, cell_id);
679 template <
typename Functor>
684 Kokkos::parallel_for(Kokkos::TeamThreadRange(
data->team_member, n_q_points),
685 [&](
const int &i) { func(this, i); });
686 data->team_member.team_barrier();
697 constexpr unsigned int
friend class FEEvaluation
int get_current_cell_index()
static constexpr unsigned int n_q_points
SharedData< dim, Number > * shared_data
void integrate(const EvaluationFlags::EvaluationFlags integration_flag)
static constexpr unsigned int dimension
const data_type * get_matrix_free_data()
const MatrixFree< dim, Number >::Data * data
static constexpr unsigned int n_components
void evaluate(const EvaluationFlags::EvaluationFlags evaluate_flag)
void read_dof_values(const DeviceVector< Number > &src)
value_type get_value(int q_point) const
void distribute_local_to_global(DeviceVector< Number > &dst) const
static constexpr unsigned int tensor_dofs_per_cell
value_type get_dof_value(int q_point) const
const MatrixFree< dim, Number >::PrecomputedData * precomputed_data
std::conditional_t< n_components_==1, Tensor< 1, dim, Number >, std::conditional_t< n_components_==dim, Tensor< 2, dim, Number >, Tensor< 1, n_components_, Tensor< 1, dim, Number > > > > gradient_type
std::conditional_t<(n_components_==1), Number, Tensor< 1, n_components_, Number > > value_type
gradient_type get_gradient(int q_point) const
static constexpr unsigned int tensor_dofs_per_component
typename MatrixFree< dim, Number >::Data data_type
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
#define DEAL_II_HOST_DEVICE
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
std::vector< index_type > data
EvaluationFlags
The EvaluationFlags enum.
constexpr bool use_collocation_evaluation(const unsigned int fe_degree, const unsigned int n_q_points_1d)
Kokkos::View< Number *, MemorySpace::Default::kokkos_space > DeviceVector
constexpr T pow(const T base, const int iexp)
static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, const typename MatrixFree< dim, Number >::Data *data)
static void evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const typename MatrixFree< dim, Number >::Data *data)
static void evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const typename MatrixFree< dim, Number >::Data *data)
static void integrate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags integration_flag, const typename MatrixFree< dim, Number >::Data *data)