Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3075-gc235bd4825 2025-04-15 08:10:00+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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
vector_access_internal.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 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
16#ifndef dealii_matrix_free_vector_access_internal_h
17#define dealii_matrix_free_vector_access_internal_h
18
19#include <deal.II/base/config.h>
20
22
26
27#include <boost/algorithm/string/join.hpp>
28
30
31
32namespace internal
33{
34 // below we use type-traits from matrix-free/type_traits.h
35
36
37
38 // access to serial const vectors that have operator[].
39 template <typename VectorType,
40 std::enable_if_t<is_serial_vector_or_array<VectorType>::value,
41 VectorType> * = nullptr>
42 inline typename VectorType::value_type
43 vector_access(const VectorType &vec, const unsigned int entry)
44 {
45 return vec[entry];
46 }
47
48
49
50 // access to serial non-const vectors that have operator[].
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 &
55 vector_access(VectorType &vec, const unsigned int entry)
56 {
57 return vec[entry];
58 }
59
60
61
62 // access to distributed MPI vectors that have a local_element(uint)
63 // method to access data in local index space, which is what we use in
64 // DoFInfo and hence in read_dof_values etc.
65 template <
66 typename VectorType,
67 std::enable_if_t<has_local_element<VectorType>, VectorType> * = nullptr>
68 inline typename VectorType::value_type &
69 vector_access(VectorType &vec, const unsigned int entry)
70 {
71 return vec.local_element(entry);
72 }
73
74
75
76 // same for const access
77 template <
78 typename VectorType,
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)
82 {
83 return vec.local_element(entry);
84 }
85
86
87
88 template <
89 typename VectorType,
90 std::enable_if_t<has_add_local_element<VectorType>, VectorType> * = nullptr>
91 inline void
92 vector_access_add(VectorType &vec,
93 const unsigned int entry,
94 const typename VectorType::value_type &val)
95 {
96 vec.add_local_element(entry, val);
97 }
98
99
100
101 template <typename VectorType,
102 std::enable_if_t<!has_add_local_element<VectorType>, VectorType> * =
103 nullptr>
104 inline void
105 vector_access_add(VectorType &vec,
106 const unsigned int entry,
107 const typename VectorType::value_type &val)
108 {
109 vector_access(vec, entry) += val;
110 }
111
112
113
114 template <
115 typename VectorType,
116 std::enable_if_t<has_add_local_element<VectorType>, VectorType> * = nullptr>
117 inline void
119 const types::global_dof_index entry,
120 const typename VectorType::value_type &val)
121 {
122 vec.add(entry, val);
123 }
124
125
126
127 template <typename VectorType,
128 std::enable_if_t<!has_add_local_element<VectorType>, VectorType> * =
129 nullptr>
130 inline void
131 vector_access_add_global(VectorType &vec,
132 const types::global_dof_index entry,
133 const typename VectorType::value_type &val)
134 {
135 vec[entry] += val;
136 }
137
138
139
140 template <
141 typename VectorType,
142 std::enable_if_t<has_set_local_element<VectorType>, VectorType> * = nullptr>
143 inline void
144 vector_access_set(VectorType &vec,
145 const unsigned int entry,
146 const typename VectorType::value_type &val)
147 {
148 vec.set_local_element(entry, val);
149 }
150
151
152
153 template <typename VectorType,
154 std::enable_if_t<!has_set_local_element<VectorType>, VectorType> * =
155 nullptr>
156 inline void
157 vector_access_set(VectorType &vec,
158 const unsigned int entry,
159 const typename VectorType::value_type &val)
160 {
161 vector_access(vec, entry) = val;
162 }
163
164
165
166 // this is to make sure that the parallel partitioning in VectorType
167 // is really the same as stored in MatrixFree.
168 // version below is when has_partitioners_are_compatible == false
169 // FIXME: this is incorrect for PETSc/Trilinos MPI vectors
170 template <int dim,
171 typename Number,
172 typename VectorizedArrayType,
173 typename VectorType,
174 std::enable_if_t<!has_partitioners_are_compatible<VectorType>,
175 VectorType> * = nullptr>
176 inline void
178 const VectorType &vec,
181 {
182 (void)vec;
183 (void)matrix_free;
184 (void)dof_info;
185
186 AssertDimension(vec.size(), dof_info.vector_partitioner->size());
187 }
188
189
190
191 // same as above for has_partitioners_are_compatible == true
192 template <int dim,
193 typename Number,
194 typename VectorizedArrayType,
195 typename VectorType,
196 std::enable_if_t<has_partitioners_are_compatible<VectorType>,
197 VectorType> * = nullptr>
198 inline void
200 const VectorType &vec,
203 {
204 (void)vec;
205 (void)matrix_free;
206 (void)dof_info;
207
208 if constexpr (running_in_debug_mode())
209 {
210 if (vec.partitioners_are_compatible(*dof_info.vector_partitioner) ==
211 false)
212 {
213 unsigned int dof_index = numbers::invalid_unsigned_int;
214
215 for (unsigned int i = 0; i < matrix_free.n_components(); ++i)
216 if (&matrix_free.get_dof_info(i) == &dof_info)
217 {
218 dof_index = i;
219 break;
220 }
221
224
225 std::vector<std::string> dof_indices_with_compatible_partitioners;
226
227 for (unsigned int i = 0; i < matrix_free.n_components(); ++i)
228 if (vec.partitioners_are_compatible(
229 *matrix_free.get_dof_info(i).vector_partitioner))
230 dof_indices_with_compatible_partitioners.push_back(
231 std::to_string(i));
232
233 if (dof_indices_with_compatible_partitioners.empty())
234 {
235 Assert(false,
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."));
245 }
246 else
247 {
248 Assert(
249 false,
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) +
265 " to get a "
266 "compatible vector?"));
267 }
268 }
269 }
270 }
271
272
273
274 // Below, three classes (VectorReader, VectorSetter,
275 // VectorDistributorLocalToGlobal) implement the same interface and can be
276 // used to to read from vector, set elements of a vector and add to elements
277 // of the vector.
278
279 // 1. A class to read data from vector
280 template <typename Number, typename VectorizedArrayType>
282 {
283 template <typename VectorType>
284 void
285 process_dof(const unsigned int index,
286 const VectorType &vec,
287 Number &res) const
288 {
289 res = vector_access(vec, index);
290 }
291
292
293
294 template <typename VectorNumberType>
295 void
296 process_dof(const VectorNumberType &global, Number &local) const
297 {
298 local = global;
299 }
300
301
302
303 template <typename VectorType>
304 void
305 process_dofs_vectorized(const unsigned int dofs_per_cell,
306 const unsigned int dof_index,
307 VectorType &vec,
308 VectorizedArrayType *dof_values,
309 std::bool_constant<true>) const
310 {
311 if constexpr (running_in_debug_mode())
312 {
313 // in debug mode, run non-vectorized version because this path
314 // has additional checks (e.g., regarding ghosting)
315 process_dofs_vectorized(dofs_per_cell,
316 dof_index,
317 vec,
318 dof_values,
319 std::bool_constant<false>());
320 }
321 else
322 {
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);
327 }
328 }
329
330
331
332 template <typename VectorType>
333 void
334 process_dofs_vectorized(const unsigned int dofs_per_cell,
335 const unsigned int dof_index,
336 const VectorType &vec,
337 VectorizedArrayType *dof_values,
338 std::bool_constant<false>) const
339 {
340 for (unsigned int i = 0; i < dofs_per_cell; ++i)
341 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
342 dof_values[i][v] =
343 vector_access(vec, dof_index + v + i * VectorizedArrayType::size());
344 }
345
346
347
348 template <typename VectorType>
349 void
350 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
351 const unsigned int *dof_indices,
352 VectorType &vec,
353 const unsigned int constant_offset,
354 VectorizedArrayType *dof_values,
355 std::bool_constant<true>) const
356 {
358 vec.begin() + constant_offset,
359 dof_indices,
360 dof_values);
361 }
362
363
364
365 template <typename VectorType>
366 void
367 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
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
373 {
374 for (unsigned int d = 0; d < dofs_per_cell; ++d)
375 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
376 dof_values[d][v] =
377 vector_access(vec, dof_indices[v] + constant_offset + d);
378 }
379
380
381
382 template <typename VectorType>
383 void
384 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
385 const unsigned int *dof_indices,
386 VectorType &vec,
387 VectorizedArrayType *dof_values,
388 std::bool_constant<true> type) const
389 {
391 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
392 }
393
394
395
396 template <typename VectorType>
397 void
398 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
399 const unsigned int *dof_indices,
400 const VectorType &vec,
401 VectorizedArrayType *dof_values,
402 std::bool_constant<false> type) const
403 {
405 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
406 }
407
408
409
410 template <typename Number2>
411 void
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
417 {
419 global_ptr,
420 dof_values);
421 }
422
423
424
425 template <typename Number2>
426 void
428 const unsigned int,
429 const std::array<Number2 *, VectorizedArrayType::size()> &,
430 VectorizedArrayType *,
431 std::bool_constant<false>) const
432 {
434 }
435
436
437
438 // variant where VectorType::value_type is the same as Number -> can call
439 // gather
440 template <typename VectorType>
441 void
442 process_dof_gather(const unsigned int *indices,
443 VectorType &vec,
444 const unsigned int constant_offset,
445 typename VectorType::value_type *vec_ptr,
446 VectorizedArrayType &res,
447 std::bool_constant<true>) const
448 {
449 (void)constant_offset;
450 (void)vec;
451
452 if constexpr (running_in_debug_mode())
453 {
454 // in debug mode, run non-vectorized version because this path
455 // has additional checks (e.g., regarding ghosting)
456 Assert(vec_ptr == vec.begin() + constant_offset, ExcInternalError());
457 process_dof_gather(indices,
458 vec,
459 constant_offset,
460 vec_ptr,
461 res,
462 std::bool_constant<false>());
463 }
464 else
465 {
466 res.gather(vec_ptr, indices);
467 }
468 }
469
470
471
472 // variant where VectorType::value_type is not the same as Number -> must
473 // manually load the data
474 template <typename VectorType>
475 void
476 process_dof_gather(const unsigned int *indices,
477 const VectorType &vec,
478 const unsigned int constant_offset,
479 typename VectorType::value_type *,
480 VectorizedArrayType &res,
481 std::bool_constant<false>) const
482 {
483 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
484 res[v] = vector_access(vec, indices[v] + constant_offset);
485 }
486
487
488
489 template <typename VectorType>
490 void
492 const VectorType &vec,
493 Number &res) const
494 {
495 res = vec[index];
496 }
497
498
499
500 void
501 pre_constraints(const Number &, Number &res) const
502 {
503 res = Number();
504 }
505
506
507
508 template <typename VectorType>
509 void
510 process_constraint(const unsigned int index,
511 const Number weight,
512 const VectorType &vec,
513 Number &res) const
514 {
515 res += weight * vector_access(vec, index);
516 }
517
518
519
520 void
521 post_constraints(const Number &sum, Number &write_pos) const
522 {
523 write_pos = sum;
524 }
525
526
527
528 void
529 process_empty(VectorizedArrayType &res) const
530 {
531 res = VectorizedArrayType();
532 }
533 };
534
535
536
537 // 2. A class to add values to the vector during
538 // FEEvaluation::distribute_local_to_global() call
539 template <typename Number, typename VectorizedArrayType>
541 {
542 template <typename VectorType>
543 void
544 process_dof(const unsigned int index, VectorType &vec, Number &res) const
545 {
546 vector_access_add(vec, index, res);
547 }
548
549
550 template <typename VectorNumberType>
551 void
552 process_dof(VectorNumberType &global, Number &local) const
553 {
554 global += local;
555 }
556
557
558
559 template <typename VectorType>
560 void
561 process_dofs_vectorized(const unsigned int dofs_per_cell,
562 const unsigned int dof_index,
563 VectorType &vec,
564 VectorizedArrayType *dof_values,
565 std::bool_constant<true>) const
566 {
567 Number *vec_ptr = vec.begin() + dof_index;
568 for (unsigned int i = 0; i < dofs_per_cell;
569 ++i, vec_ptr += VectorizedArrayType::size())
570 {
571 VectorizedArrayType tmp;
572 tmp.load(vec_ptr);
573 tmp += dof_values[i];
574 tmp.store(vec_ptr);
575 }
576 }
577
578
579
580 template <typename VectorType>
581 void
582 process_dofs_vectorized(const unsigned int dofs_per_cell,
583 const unsigned int dof_index,
584 VectorType &vec,
585 VectorizedArrayType *dof_values,
586 std::bool_constant<false>) const
587 {
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(),
592 dof_values[i][v]);
593 }
594
595
596
597 template <typename VectorType>
598 void
599 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
600 const unsigned int *dof_indices,
601 VectorType &vec,
602 const unsigned int constant_offset,
603 VectorizedArrayType *dof_values,
604 std::bool_constant<true>) const
605 {
607 dofs_per_cell,
608 dof_values,
609 dof_indices,
610 vec.begin() + constant_offset);
611 }
612
613
614
615 template <typename VectorType>
616 void
617 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
618 const unsigned int *dof_indices,
619 VectorType &vec,
620 const unsigned int constant_offset,
621 VectorizedArrayType *dof_values,
622 std::bool_constant<false>) const
623 {
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,
628 dof_values[d][v]);
629 }
630
631
632
633 template <typename VectorType>
634 void
635 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
636 const unsigned int *dof_indices,
637 VectorType &vec,
638 VectorizedArrayType *dof_values,
639 std::bool_constant<true> type) const
640 {
642 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
643 }
644
645
646
647 template <typename VectorType>
648 void
649 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
650 const unsigned int *dof_indices,
651 VectorType &vec,
652 VectorizedArrayType *dof_values,
653 std::bool_constant<false> type) const
654 {
656 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
657 }
658
659
660
661 template <typename Number2>
662 void
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
668 {
670 dofs_per_cell,
671 dof_values,
672 global_ptr);
673 }
674
675
676
677 template <typename Number2>
678 void
680 const unsigned int,
681 std::array<Number2 *, VectorizedArrayType::size()> &,
682 VectorizedArrayType *,
683 std::bool_constant<false>) const
684 {
686 }
687
688
689
690 // variant where VectorType::value_type is the same as Number -> can call
691 // scatter
692 template <typename VectorType>
693 void
694 process_dof_gather(const unsigned int *indices,
695 VectorType &vec,
696 const unsigned int constant_offset,
697 typename VectorType::value_type *vec_ptr,
698 VectorizedArrayType &res,
699 std::bool_constant<true>) const
700 {
701 (void)constant_offset;
702 (void)vec_ptr;
703 (void)vec;
704
705#if DEAL_II_VECTORIZATION_WIDTH_IN_BITS < 512
706 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
707 vector_access(vec, indices[v] + constant_offset) += res[v];
708#else
709 // only use gather in case there is also scatter.
710 VectorizedArrayType tmp;
711 tmp.gather(vec_ptr, indices);
712 tmp += res;
713 tmp.scatter(indices, vec_ptr);
714#endif
715 }
716
717
718
719 // variant where VectorType::value_type is not the same as Number -> must
720 // manually append all data
721 template <typename VectorType>
722 void
723 process_dof_gather(const unsigned int *indices,
724 VectorType &vec,
725 const unsigned int constant_offset,
726 typename VectorType::value_type *,
727 VectorizedArrayType &res,
728 std::bool_constant<false>) const
729 {
730 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
731 vector_access_add(vec, indices[v] + constant_offset, res[v]);
732 }
733
734
735
736 template <typename VectorType>
737 void
739 VectorType &vec,
740 Number &res) const
741 {
743 }
744
745
746
747 void
748 pre_constraints(const Number &input, Number &res) const
749 {
750 res = input;
751 }
752
753
754
755 template <typename VectorType>
756 void
757 process_constraint(const unsigned int index,
758 const Number weight,
759 VectorType &vec,
760 Number &res) const
761 {
762 vector_access_add(vec, index, weight * res);
763 }
764
765
766
767 void
768 post_constraints(const Number &, Number &) const
769 {}
770
771
772
773 void
774 process_empty(VectorizedArrayType &) const
775 {}
776 };
777
778
779
780 // 3. A class to set elements of the vector
781 template <typename Number, typename VectorizedArrayType>
783 {
784 template <typename VectorType>
785 void
786 process_dof(const unsigned int index, VectorType &vec, Number &res) const
787 {
788 vector_access(vec, index) = res;
789 }
790
791
792
793 template <typename VectorNumberType>
794 void
795 process_dof(VectorNumberType &global, Number &local) const
796 {
797 global = local;
798 }
799
800
801
802 template <typename VectorType>
803 void
804 process_dofs_vectorized(const unsigned int dofs_per_cell,
805 const unsigned int dof_index,
806 VectorType &vec,
807 VectorizedArrayType *dof_values,
808 std::bool_constant<true>) const
809 {
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);
814 }
815
816
817
818 template <typename VectorType>
819 void
820 process_dofs_vectorized(const unsigned int dofs_per_cell,
821 const unsigned int dof_index,
822 VectorType &vec,
823 VectorizedArrayType *dof_values,
824 std::bool_constant<false>) const
825 {
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()) =
829 dof_values[i][v];
830 }
831
832
833
834 template <typename VectorType>
835 void
836 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
837 const unsigned int *dof_indices,
838 VectorType &vec,
839 const unsigned int constant_offset,
840 VectorizedArrayType *dof_values,
841 std::bool_constant<true>) const
842 {
844 dofs_per_cell,
845 dof_values,
846 dof_indices,
847 vec.begin() + constant_offset);
848 }
849
850
851
852 template <typename VectorType, bool booltype>
853 void
854 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
855 const unsigned int *dof_indices,
856 VectorType &vec,
857 const unsigned int constant_offset,
858 VectorizedArrayType *dof_values,
859 std::bool_constant<false>) const
860 {
861 for (unsigned int i = 0; i < dofs_per_cell; ++i)
862 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
863 vector_access(vec, constant_offset + dof_indices[v] + i) =
864 dof_values[i][v];
865 }
866
867
868
869 template <typename VectorType>
870 void
871 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
872 const unsigned int *dof_indices,
873 VectorType &vec,
874 VectorizedArrayType *dof_values,
875 std::bool_constant<true> type) const
876 {
878 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
879 }
880
881
882
883 template <typename VectorType, bool booltype>
884 void
885 process_dofs_vectorized_transpose(const unsigned int dofs_per_cell,
886 const unsigned int *dof_indices,
887 VectorType &vec,
888 VectorizedArrayType *dof_values,
889 std::bool_constant<false> type) const
890 {
892 dofs_per_cell, dof_indices, vec, 0, dof_values, type);
893 }
894
895
896
897 template <typename Number2>
898 void
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
904 {
906 dofs_per_cell,
907 dof_values,
908 global_ptr);
909 }
910
911
912
913 template <typename Number2>
914 void
916 const unsigned int,
917 std::array<Number2 *, VectorizedArrayType::size()> &,
918 VectorizedArrayType *,
919 std::bool_constant<false>) const
920 {
922 }
923
924
925
926 template <typename VectorType>
927 void
928 process_dof_gather(const unsigned int *indices,
929 VectorType &vec,
930 const unsigned int constant_offset,
931 typename VectorType::value_type *vec_ptr,
932 VectorizedArrayType &res,
933 std::bool_constant<true>) const
934 {
935 (void)vec_ptr;
936 (void)constant_offset;
937 (void)vec;
938 Assert(vec_ptr == vec.begin() + constant_offset, ExcInternalError());
939 res.scatter(indices, vec_ptr);
940 }
941
942
943
944 template <typename VectorType>
945 void
946 process_dof_gather(const unsigned int *indices,
947 VectorType &vec,
948 const unsigned int constant_offset,
949 typename VectorType::value_type *,
950 VectorizedArrayType &res,
951 std::bool_constant<false>) const
952 {
953 for (unsigned int v = 0; v < VectorizedArrayType::size(); ++v)
954 vector_access(vec, indices[v] + constant_offset) = res[v];
955 }
956
957
958
959 template <typename VectorType>
960 void
962 VectorType &vec,
963 Number &res) const
964 {
965 vec[index] = res;
966 }
967
968
969
970 void
971 pre_constraints(const Number &, Number &) const
972 {}
973
974
975
976 template <typename VectorType>
977 void
978 process_constraint(const unsigned int,
979 const Number,
980 VectorType &,
981 Number &) const
982 {}
983
984
985
986 void
987 post_constraints(const Number &, Number &) const
988 {}
989
990
991
992 void
993 process_empty(VectorizedArrayType &) const
994 {}
995 };
996} // namespace internal
997
998
1000
1001#endif
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
Definition config.h:35
constexpr bool running_in_debug_mode()
Definition config.h:73
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
#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
Definition types.h:238
std::shared_ptr< const Utilities::MPI::Partitioner > vector_partitioner
Definition dof_info.h:593
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)