16#ifndef dealii_matrix_free_vector_access_internal_h
17#define dealii_matrix_free_vector_access_internal_h
27#include <boost/algorithm/string/join.hpp>
40 std::enable_if_t<is_serial_vector_or_array<VectorType>::value,
42 inline typename VectorType::value_type
51 template <
typename VectorType,
52 std::enable_if_t<is_serial_vector_or_array<VectorType>::value,
53 VectorType> * =
nullptr>
54 inline typename VectorType::value_type &
67 std::enable_if_t<has_local_element<VectorType>, VectorType> * =
nullptr>
68 inline typename VectorType::value_type &
71 return vec.local_element(entry);
79 std::enable_if_t<has_local_element<VectorType>,
VectorType> * =
nullptr>
80 inline typename VectorType::value_type
81 vector_access(
const VectorType &vec,
const unsigned int entry)
83 return vec.local_element(entry);
90 std::enable_if_t<has_add_local_element<VectorType>,
VectorType> * =
nullptr>
93 const unsigned int entry,
94 const typename VectorType::value_type &val)
96 vec.add_local_element(entry, val);
101 template <
typename VectorType,
102 std::enable_if_t<!has_add_local_element<VectorType>, VectorType> * =
106 const unsigned int entry,
107 const typename VectorType::value_type &val)
116 std::enable_if_t<has_add_local_element<VectorType>,
VectorType> * =
nullptr>
120 const typename VectorType::value_type &val)
127 template <
typename VectorType,
128 std::enable_if_t<!has_add_local_element<VectorType>, VectorType> * =
133 const typename VectorType::value_type &val)
142 std::enable_if_t<has_set_local_element<VectorType>,
VectorType> * =
nullptr>
145 const unsigned int entry,
146 const typename VectorType::value_type &val)
148 vec.set_local_element(entry, val);
153 template <
typename VectorType,
154 std::enable_if_t<!has_set_local_element<VectorType>, VectorType> * =
158 const unsigned int entry,
159 const typename VectorType::value_type &val)
172 typename VectorizedArrayType,
174 std::enable_if_t<!has_partitioners_are_compatible<VectorType>,
178 const VectorType &vec,
194 typename VectorizedArrayType,
196 std::enable_if_t<has_partitioners_are_compatible<VectorType>,
197 VectorType> * =
nullptr>
200 const VectorType &vec,
215 for (
unsigned int i = 0; i < matrix_free.
n_components(); ++i)
225 std::vector<std::string> dof_indices_with_compatible_partitioners;
227 for (
unsigned int i = 0; i < matrix_free.
n_components(); ++i)
228 if (vec.partitioners_are_compatible(
230 dof_indices_with_compatible_partitioners.push_back(
233 if (dof_indices_with_compatible_partitioners.empty())
237 "The parallel layout of the given vector is "
238 "compatible neither with the Partitioner of the "
239 "current FEEvaluation with dof_handler_index=" +
240 std::to_string(dof_index) +
241 " nor with any Partitioner in MatrixFree. A "
242 "potential reason is that you did not use "
243 "MatrixFree::initialize_dof_vector() to get a "
244 "compatible vector."));
251 "The parallel layout of the given vector is "
252 "not compatible with the Partitioner of the "
253 "current FEEvaluation with dof_handler_index=" +
254 std::to_string(dof_index) +
255 ". However, the underlying "
256 "MatrixFree contains Partitioner objects that are compatible. "
257 "They have the following dof_handler_index values: " +
258 boost::algorithm::join(
259 dof_indices_with_compatible_partitioners,
", ") +
260 ". Did you want to pass any of these values to the "
261 "constructor of the current FEEvaluation object or "
262 "did you not use MatrixFree::initialize_dof_vector() "
263 "with dof_handler_index=" +
264 std::to_string(dof_index) +
266 "compatible vector?"));
280 template <
typename Number,
typename VectorizedArrayType>
283 template <
typename VectorType>
286 const VectorType &vec,
294 template <
typename VectorNumberType>
303 template <
typename VectorType>
306 const unsigned int dof_index,
308 VectorizedArrayType *dof_values,
309 std::bool_constant<true>)
const
319 std::bool_constant<false>());
323 const Number *vec_ptr = vec.begin() + dof_index;
324 for (
unsigned int i = 0; i < dofs_per_cell;
325 ++i, vec_ptr += VectorizedArrayType::size())
326 dof_values[i].load(vec_ptr);
332 template <
typename VectorType>
335 const unsigned int dof_index,
336 const VectorType &vec,
337 VectorizedArrayType *dof_values,
338 std::bool_constant<false>)
const
340 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
341 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
343 vector_access(vec, dof_index + v + i * VectorizedArrayType::size());
348 template <
typename VectorType>
351 const unsigned int *dof_indices,
353 const unsigned int constant_offset,
354 VectorizedArrayType *dof_values,
355 std::bool_constant<true>)
const
358 vec.begin() + constant_offset,
365 template <
typename VectorType>
368 const unsigned int *dof_indices,
369 const VectorType &vec,
370 const unsigned int constant_offset,
371 VectorizedArrayType *dof_values,
372 std::bool_constant<false>)
const
374 for (
unsigned int d = 0; d < dofs_per_cell; ++d)
375 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
382 template <
typename VectorType>
385 const unsigned int *dof_indices,
387 VectorizedArrayType *dof_values,
388 std::bool_constant<true> type)
const
391 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
396 template <
typename VectorType>
399 const unsigned int *dof_indices,
400 const VectorType &vec,
401 VectorizedArrayType *dof_values,
402 std::bool_constant<false> type)
const
405 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
410 template <
typename Number2>
413 const unsigned int dofs_per_cell,
414 const std::array<Number2 *, VectorizedArrayType::size()> &global_ptr,
415 VectorizedArrayType *dof_values,
416 std::bool_constant<true>)
const
425 template <
typename Number2>
429 const std::array<Number2 *, VectorizedArrayType::size()> &,
430 VectorizedArrayType *,
431 std::bool_constant<false>)
const
440 template <
typename VectorType>
444 const unsigned int constant_offset,
445 typename VectorType::value_type *vec_ptr,
446 VectorizedArrayType &res,
447 std::bool_constant<true>)
const
449 (void)constant_offset;
462 std::bool_constant<false>());
466 res.gather(vec_ptr, indices);
474 template <
typename VectorType>
477 const VectorType &vec,
478 const unsigned int constant_offset,
479 typename VectorType::value_type *,
480 VectorizedArrayType &res,
481 std::bool_constant<false>)
const
483 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
489 template <
typename VectorType>
492 const VectorType &vec,
508 template <
typename VectorType>
512 const VectorType &vec,
531 res = VectorizedArrayType();
539 template <
typename Number,
typename VectorizedArrayType>
542 template <
typename VectorType>
550 template <
typename VectorNumberType>
559 template <
typename VectorType>
562 const unsigned int dof_index,
564 VectorizedArrayType *dof_values,
565 std::bool_constant<true>)
const
567 Number *vec_ptr = vec.begin() + dof_index;
568 for (
unsigned int i = 0; i < dofs_per_cell;
569 ++i, vec_ptr += VectorizedArrayType::size())
571 VectorizedArrayType tmp;
573 tmp += dof_values[i];
580 template <
typename VectorType>
583 const unsigned int dof_index,
585 VectorizedArrayType *dof_values,
586 std::bool_constant<false>)
const
588 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
589 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
591 dof_index + v + i * VectorizedArrayType::size(),
597 template <
typename VectorType>
600 const unsigned int *dof_indices,
602 const unsigned int constant_offset,
603 VectorizedArrayType *dof_values,
604 std::bool_constant<true>)
const
610 vec.begin() + constant_offset);
615 template <
typename VectorType>
618 const unsigned int *dof_indices,
620 const unsigned int constant_offset,
621 VectorizedArrayType *dof_values,
622 std::bool_constant<false>)
const
624 for (
unsigned int d = 0; d < dofs_per_cell; ++d)
625 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
627 dof_indices[v] + constant_offset + d,
633 template <
typename VectorType>
636 const unsigned int *dof_indices,
638 VectorizedArrayType *dof_values,
639 std::bool_constant<true> type)
const
642 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
647 template <
typename VectorType>
650 const unsigned int *dof_indices,
652 VectorizedArrayType *dof_values,
653 std::bool_constant<false> type)
const
656 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
661 template <
typename Number2>
664 const unsigned int dofs_per_cell,
665 std::array<Number2 *, VectorizedArrayType::size()> &global_ptr,
666 VectorizedArrayType *dof_values,
667 std::bool_constant<true>)
const
677 template <
typename Number2>
681 std::array<Number2 *, VectorizedArrayType::size()> &,
682 VectorizedArrayType *,
683 std::bool_constant<false>)
const
692 template <
typename VectorType>
696 const unsigned int constant_offset,
697 typename VectorType::value_type *vec_ptr,
698 VectorizedArrayType &res,
699 std::bool_constant<true>)
const
701 (void)constant_offset;
705#if DEAL_II_VECTORIZATION_WIDTH_IN_BITS < 512
706 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
710 VectorizedArrayType tmp;
711 tmp.gather(vec_ptr, indices);
713 tmp.scatter(indices, vec_ptr);
721 template <
typename VectorType>
725 const unsigned int constant_offset,
726 typename VectorType::value_type *,
727 VectorizedArrayType &res,
728 std::bool_constant<false>)
const
730 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
736 template <
typename VectorType>
755 template <
typename VectorType>
781 template <
typename Number,
typename VectorizedArrayType>
784 template <
typename VectorType>
793 template <
typename VectorNumberType>
802 template <
typename VectorType>
805 const unsigned int dof_index,
807 VectorizedArrayType *dof_values,
808 std::bool_constant<true>)
const
810 Number *vec_ptr = vec.begin() + dof_index;
811 for (
unsigned int i = 0; i < dofs_per_cell;
812 ++i, vec_ptr += VectorizedArrayType::size())
813 dof_values[i].store(vec_ptr);
818 template <
typename VectorType>
821 const unsigned int dof_index,
823 VectorizedArrayType *dof_values,
824 std::bool_constant<false>)
const
826 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
827 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
828 vector_access(vec, dof_index + v + i * VectorizedArrayType::size()) =
834 template <
typename VectorType>
837 const unsigned int *dof_indices,
839 const unsigned int constant_offset,
840 VectorizedArrayType *dof_values,
841 std::bool_constant<true>)
const
847 vec.begin() + constant_offset);
852 template <
typename VectorType,
bool booltype>
855 const unsigned int *dof_indices,
857 const unsigned int constant_offset,
858 VectorizedArrayType *dof_values,
859 std::bool_constant<false>)
const
861 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
862 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
869 template <
typename VectorType>
872 const unsigned int *dof_indices,
874 VectorizedArrayType *dof_values,
875 std::bool_constant<true> type)
const
878 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
883 template <
typename VectorType,
bool booltype>
886 const unsigned int *dof_indices,
888 VectorizedArrayType *dof_values,
889 std::bool_constant<false> type)
const
892 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
897 template <
typename Number2>
900 const unsigned int dofs_per_cell,
901 std::array<Number2 *, VectorizedArrayType::size()> &global_ptr,
902 VectorizedArrayType *dof_values,
903 std::bool_constant<true>)
const
913 template <
typename Number2>
917 std::array<Number2 *, VectorizedArrayType::size()> &,
918 VectorizedArrayType *,
919 std::bool_constant<false>)
const
926 template <
typename VectorType>
930 const unsigned int constant_offset,
931 typename VectorType::value_type *vec_ptr,
932 VectorizedArrayType &res,
933 std::bool_constant<true>)
const
936 (void)constant_offset;
939 res.scatter(indices, vec_ptr);
944 template <
typename VectorType>
948 const unsigned int constant_offset,
949 typename VectorType::value_type *,
950 VectorizedArrayType &res,
951 std::bool_constant<false>)
const
953 for (
unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
959 template <
typename VectorType>
976 template <
typename VectorType>
const internal::MatrixFreeFunctions::DoFInfo & get_dof_info(const unsigned int dof_handler_index_component=0) const
unsigned int n_components() const
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
void vector_access_add(VectorType &vec, const unsigned int entry, const typename VectorType::value_type &val)
void check_vector_compatibility(const VectorType &vec, const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free, const internal::MatrixFreeFunctions::DoFInfo &dof_info)
void vector_access_add_global(VectorType &vec, const types::global_dof_index entry, const typename VectorType::value_type &val)
VectorType::value_type vector_access(const VectorType &vec, const unsigned int entry)
void vector_access_set(VectorType &vec, const unsigned int entry, const typename VectorType::value_type &val)
constexpr unsigned int invalid_unsigned_int
std::shared_ptr< const Utilities::MPI::Partitioner > vector_partitioner
void process_dof_gather(const unsigned int *indices, VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *vec_ptr, VectorizedArrayType &res, std::bool_constant< true >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false > type) const
void pre_constraints(const Number &input, Number &res) const
void process_empty(VectorizedArrayType &) const
void process_dof(const unsigned int index, VectorType &vec, Number &res) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void post_constraints(const Number &, Number &) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true > type) const
void process_dofs_vectorized_transpose(const unsigned int, std::array< Number2 *, VectorizedArrayType::size()> &, VectorizedArrayType *, std::bool_constant< false >) const
void process_dof_gather(const unsigned int *indices, VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *, VectorizedArrayType &res, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, std::array< Number2 *, VectorizedArrayType::size()> &global_ptr, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dof(VectorNumberType &global, Number &local) const
void process_constraint(const unsigned int index, const Number weight, VectorType &vec, Number &res) const
void process_dof_global(const types::global_dof_index index, VectorType &vec, Number &res) const
void pre_constraints(const Number &, Number &res) const
void post_constraints(const Number &sum, Number &write_pos) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dof_global(const types::global_dof_index index, const VectorType &vec, Number &res) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, const VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, const VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void process_constraint(const unsigned int index, const Number weight, const VectorType &vec, Number &res) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true > type) const
void process_dof_gather(const unsigned int *indices, const VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *, VectorizedArrayType &res, std::bool_constant< false >) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dof(const VectorNumberType &global, Number &local) const
void process_empty(VectorizedArrayType &res) const
void process_dof(const unsigned int index, const VectorType &vec, Number &res) const
void process_dofs_vectorized_transpose(const unsigned int, const std::array< Number2 *, VectorizedArrayType::size()> &, VectorizedArrayType *, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, const VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false > type) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const std::array< Number2 *, VectorizedArrayType::size()> &global_ptr, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dof_gather(const unsigned int *indices, VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *vec_ptr, VectorizedArrayType &res, std::bool_constant< true >) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, std::array< Number2 *, VectorizedArrayType::size()> &global_ptr, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_empty(VectorizedArrayType &) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true > type) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< false >) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< false > type) const
void process_dofs_vectorized_transpose(const unsigned int dofs_per_cell, const unsigned int *dof_indices, VectorType &vec, const unsigned int constant_offset, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dof_global(const types::global_dof_index index, VectorType &vec, Number &res) const
void pre_constraints(const Number &, Number &) const
void process_dof_gather(const unsigned int *indices, VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *vec_ptr, VectorizedArrayType &res, std::bool_constant< true >) const
void process_dofs_vectorized(const unsigned int dofs_per_cell, const unsigned int dof_index, VectorType &vec, VectorizedArrayType *dof_values, std::bool_constant< true >) const
void process_dofs_vectorized_transpose(const unsigned int, std::array< Number2 *, VectorizedArrayType::size()> &, VectorizedArrayType *, std::bool_constant< false >) const
void process_constraint(const unsigned int, const Number, VectorType &, Number &) const
void process_dof(const unsigned int index, VectorType &vec, Number &res) const
void post_constraints(const Number &, Number &) const
void process_dof(VectorNumberType &global, Number &local) const
void process_dof_gather(const unsigned int *indices, VectorType &vec, const unsigned int constant_offset, typename VectorType::value_type *, VectorizedArrayType &res, std::bool_constant< false >) const
void vectorized_load_and_transpose(const unsigned int n_entries, const Number *in, const unsigned int *offsets, VectorizedArray< Number, width > *out)
void vectorized_transpose_and_store(const bool add_into, const unsigned int n_entries, const VectorizedArray< Number, width > *in, const unsigned int *offsets, Number *out)