1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2011 - 2020 by the deal.II authors 4 // 5 // This file is part of the deal.II library. 6 // 7 // The deal.II library is free software; you can use it, redistribute 8 // it, and/or modify it under the terms of the GNU Lesser General 9 // Public License as published by the Free Software Foundation; either 10 // version 2.1 of the License, or (at your option) any later version. 11 // The full text of the license can be found in the file LICENSE.md at 12 // the top level directory of deal.II. 13 // 14 // --------------------------------------------------------------------- 15 16 17 #ifndef dealii_matrix_free_dof_info_h 18 #define dealii_matrix_free_dof_info_h 19 20 21 #include <deal.II/base/config.h> 22 23 #include <deal.II/base/exceptions.h> 24 #include <deal.II/base/partitioner.h> 25 #include <deal.II/base/vectorization.h> 26 27 #include <deal.II/dofs/dof_handler.h> 28 29 #include <deal.II/lac/affine_constraints.h> 30 #include <deal.II/lac/dynamic_sparsity_pattern.h> 31 32 #include <deal.II/matrix_free/face_info.h> 33 #include <deal.II/matrix_free/mapping_info.h> 34 #include <deal.II/matrix_free/shape_info.h> 35 #include <deal.II/matrix_free/task_info.h> 36 37 #include <array> 38 #include <memory> 39 40 41 DEAL_II_NAMESPACE_OPEN 42 43 namespace internal 44 { 45 namespace MatrixFreeFunctions 46 { 47 /** 48 * A struct that takes entries describing a constraint and puts them into 49 * a sorted list where duplicates are filtered out 50 */ 51 template <typename Number> 52 struct ConstraintValues 53 { 54 ConstraintValues(); 55 56 /** 57 * This function inserts some constrained entries to the collection of 58 * all values. It stores the (reordered) numbering of the dofs 59 * (according to the ordering that matches with the function) in 60 * new_indices, and returns the storage position the double array for 61 * access later on. 62 */ 63 template <typename number2> 64 unsigned short 65 insert_entries( 66 const std::vector<std::pair<types::global_dof_index, number2>> 67 &entries); 68 69 std::vector<std::pair<types::global_dof_index, double>> 70 constraint_entries; 71 std::vector<types::global_dof_index> constraint_indices; 72 73 std::pair<std::vector<Number>, types::global_dof_index> next_constraint; 74 std::map<std::vector<Number>, 75 types::global_dof_index, 76 FPArrayComparator<Number>> 77 constraints; 78 }; 79 80 /** 81 * The class that stores the indices of the degrees of freedom for all the 82 * cells. Essentially, this is a smart number cache in the style of a 83 * DoFHandler that also embeds the description of constraints directly on 84 * the cell level without the need to refer to the external 85 * AffineConstraints object. 86 * 87 * This class only stores index relations. The weights for hanging node 88 * constraints are stored in a different field. This is because a 89 * different field allows for the same compressed weight data on different 90 * DoFHandlers for vector-valued problems. There, the indices might be 91 * constrained differently on different components (e.g. Dirichlet 92 * conditions only on selected components), whereas the weights from 93 * hanging nodes are the same and need to be stored only once. The 94 * combination will be handled in the MatrixFree class. 95 * 96 * @ingroup matrixfree 97 */ 98 struct DoFInfo 99 { 100 /** 101 * This value is used to define subranges in the vectors which we can 102 * zero inside the MatrixFree::loop() call. The goal is to only clear a 103 * part of the vector at a time to keep the values that are zeroed in 104 * caches, saving one global vector access for the case where this is 105 * applied rather than `vector = 0.;`. 106 * 107 * We set the granularity to 64 - that is a number sufficiently large 108 * to minimize loop peel overhead within the work (and compatible with 109 * vectorization lengths of up to 16) and small enough to not waste on 110 * the size of the individual chunks. 111 */ 112 static const unsigned int chunk_size_zero_vector = 64; 113 114 /** 115 * Default empty constructor. 116 */ 117 DoFInfo(); 118 119 /** 120 * Copy constructor. 121 */ 122 DoFInfo(const DoFInfo &) = default; 123 124 /** 125 * Clear all data fields in this class. 126 */ 127 void 128 clear(); 129 130 /** 131 * Return the FE index for a given finite element degree. If not in hp 132 * mode, this function always returns index 0. If an index is not found 133 * in hp mode, it returns numbers::invalid_unsigned_int. 134 */ 135 unsigned int 136 fe_index_from_degree(const unsigned int first_selected_component, 137 const unsigned int fe_degree) const; 138 139 /** 140 * Populate the vector @p locall_indices with locally owned degrees of freedom 141 * stored on the cell block @p cell. 142 * If @p with_constraints is `true`, then the returned vector will contain indices 143 * required to resolve constraints. 144 * 145 * The image below illustrates the output of this function for cell blocks 146 * zero and one with zero Dirichlet boundary conditions at the bottom of 147 * the domain. Note that due to the presence of constraints, the DoFs 148 * returned by this function for the case `with_constraints = true` are 149 * not a simple union 150 * of per cell DoFs on the cell block @p cell. 151 * 152 * @image html dofinfo_get_dof_indices.png 153 * 154 * @note The returned indices may contain duplicates. The unique set can be 155 * obtain using `std::sort()` followed by `std::unique()` and 156 * `std::vector::erase()`. 157 */ 158 void 159 get_dof_indices_on_cell_batch(std::vector<unsigned int> &locall_indices, 160 const unsigned int cell, 161 const bool with_constraints = true) const; 162 163 /** 164 * This internal method takes the local indices on a cell and fills them 165 * into this class. It resolves the constraints and distributes the 166 * results. Ghost indices, i.e., indices that are located on another 167 * processor, get a temporary number by this function, and will later be 168 * assigned the correct index after all the ghost indices have been 169 * collected by the call to @p assign_ghosts. 170 */ 171 template <typename number> 172 void 173 read_dof_indices( 174 const std::vector<types::global_dof_index> &local_indices, 175 const std::vector<unsigned int> & lexicographic_inv, 176 const dealii::AffineConstraints<number> & constraints, 177 const unsigned int cell_number, 178 ConstraintValues<double> & constraint_values, 179 bool & cell_at_boundary); 180 181 /** 182 * This method assigns the correct indices to ghost indices from the 183 * temporary numbering employed by the @p read_dof_indices function. The 184 * numbers are localized with respect to the MPI process, and ghosts 185 * start at the end of the locally owned range. This way, we get direct 186 * access to all vector entries. 187 */ 188 void 189 assign_ghosts(const std::vector<unsigned int> &boundary_cells); 190 191 /** 192 * This method reorders the way cells are gone through based on a given 193 * renumbering of the cells. It also takes @p vectorization_length cells 194 * together and interprets them as one cell only, as is needed for 195 * vectorization. 196 */ 197 void 198 reorder_cells(const TaskInfo & task_info, 199 const std::vector<unsigned int> & renumbering, 200 const std::vector<unsigned int> & constraint_pool_row_index, 201 const std::vector<unsigned char> &irregular_cells); 202 203 /** 204 * Finds possible compression for the cell indices that we can apply for 205 * increased efficiency. Run at the end of reorder_cells. 206 */ 207 void 208 compute_cell_index_compression( 209 const std::vector<unsigned char> &irregular_cells); 210 211 /** 212 * Finds possible compression for the face indices that we can apply for 213 * increased efficiency. Run at the end of reorder_cells. 214 */ 215 template <int length> 216 void 217 compute_face_index_compression( 218 const std::vector<FaceToCellTopology<length>> &faces); 219 220 /** 221 * This function computes the connectivity of the currently stored 222 * indices in terms of connections between the individual cells and 223 * fills the structure into a sparsity pattern. 224 */ 225 void 226 make_connectivity_graph(const TaskInfo & task_info, 227 const std::vector<unsigned int> &renumbering, 228 DynamicSparsityPattern &connectivity) const; 229 230 /** 231 * In case face integrals are enabled, find out whether certain loops 232 * over the unknowns only access a subset of all the ghost dofs we keep 233 * in the main partitioner. 234 */ 235 void 236 compute_tight_partitioners( 237 const Table<2, ShapeInfo<double>> & shape_info, 238 const unsigned int n_owned_cells, 239 const unsigned int n_lanes, 240 const std::vector<FaceToCellTopology<1>> &inner_faces, 241 const std::vector<FaceToCellTopology<1>> &ghosted_faces, 242 const bool fill_cell_centric); 243 244 /** 245 * Compute a renumbering of the degrees of freedom to improve the data 246 * access patterns for this class that can be utilized by the categories 247 * in the IndexStorageVariants enum. For example, the index ordering can 248 * be improved for typical DG elements by interleaving the degrees of 249 * freedom from batches of cells, which avoids the explicit data 250 * transposition in IndexStorageVariants::contiguous. Currently, these 251 * more advanced features are not implemented, so there is only limited 252 * value of this function. 253 */ 254 void 255 compute_dof_renumbering( 256 std::vector<types::global_dof_index> &renumbering); 257 258 /** 259 * Fills the array that defines how to zero selected ranges in the result 260 * vector within the cell loop, filling the two member variables @p 261 * vector_zero_range_list_index and @p vector_zero_range_list. 262 * 263 * The intent of this pattern is to zero the vector entries in close 264 * temporal proximity to the first access and thus keeping the vector 265 * entries in cache. 266 */ 267 template <int length> 268 void 269 compute_vector_zero_access_pattern( 270 const TaskInfo & task_info, 271 const std::vector<FaceToCellTopology<length>> &faces); 272 273 /** 274 * Return the memory consumption in bytes of this class. 275 */ 276 std::size_t 277 memory_consumption() const; 278 279 /** 280 * Prints a detailed summary of memory consumption in the different 281 * structures of this class to the given output stream. 282 */ 283 template <typename StreamType> 284 void 285 print_memory_consumption(StreamType & out, 286 const TaskInfo &size_info) const; 287 288 /** 289 * Prints a representation of the indices in the class to the given 290 * output stream. 291 */ 292 template <typename Number> 293 void 294 print(const std::vector<Number> & constraint_pool_data, 295 const std::vector<unsigned int> &constraint_pool_row_index, 296 std::ostream & out) const; 297 298 /** 299 * Enum for various storage variants of the indices. This storage format 300 * is used to implement more efficient indexing schemes in case the 301 * underlying data structures allow for them, and to inform the access 302 * functions in FEEvaluationBase::read_write_operation() on which array 303 * to get the data from. One example of more efficient storage is the 304 * enum value IndexStorageVariants::contiguous, which means that one can 305 * get the indices to all degrees of freedom of a cell by reading only 306 * the first index for each cell, whereas all subsequent indices are 307 * merely an offset from the first index. 308 */ 309 enum class IndexStorageVariants : unsigned char 310 { 311 /** 312 * This value indicates that no index compression was found and the 313 * only valid storage is to access all indices present on the cell, 314 * possibly including constraints. For a cell/face of this index type, 315 * the data access in FEEvaluationBase is directed to the array @p 316 * dof_indices with the index 317 * `row_starts[cell_index*n_vectorization*n_components].first`. 318 */ 319 full, 320 /** 321 * This value indicates that the indices are interleaved for access 322 * with vectorized gather and scatter operation. This storage variant 323 * is possible in case there are no constraints on the cell and the 324 * indices in the batch of cells are not pointing to the same global 325 * index in different slots of a vectorized array (in order to support 326 * scatter operations). For a cell/face of this index type, the data 327 * access in FEEvaluationBase is directed to the array 328 * `dof_indices_interleaved` with the index 329 * `row_starts[cell_index*n_vectorization*n_components].first`. 330 */ 331 interleaved, 332 /** 333 * This value indicates that the indices within a cell are all 334 * contiguous, and one can get the index to the cell by reading that 335 * single value for each of the cells in the cell batch. For a 336 * cell/face of this index type, the data access in FEEvaluationBase 337 * is directed to the array `dof_indices_contiguous` with the index 338 * `cell_index*n_vectorization*n_components`. 339 */ 340 contiguous, 341 /** 342 * This value indicates that the indices with a cell are contiguous and 343 * interleaved for vectorization, i.e., the first DoF index on a cell 344 * to the four or eight cells in the vectorization batch come first, 345 * than the second DoF index, and so on. Furthermore, the interleaving 346 * between cells implies that only the batches for vectorization can be 347 * accessed efficiently, whereas there is a strided access for getting 348 * only some of the entries. 349 * 350 * The two additional categories `interleaved_contiguous_strided` and 351 * `interleaved_contiguous_mixed_strides` are a consequence of this 352 * storage type. The former is for faces where at least one of the two 353 * adjacent sides will break with the interleaved storage. We then have 354 * to make a strided access as described in the next category. The last 355 * category `interleaved_contiguous_mixed_strides` appears in the ghost 356 * layer, see the more detailed description of that category below. 357 * Again, this is something that cannot be avoided in general once we 358 * interleave the indices between cells. 359 * 360 * For a cell/face of this index type, the data access in 361 * FEEvaluationBase is directed to the array `dof_indices_contiguous` 362 * with the index `cell_index*n_vectorization*n_components`. 363 */ 364 interleaved_contiguous, 365 /** 366 * Similar to interleaved_contiguous storage, but for the case when the 367 * interleaved indices are only contiguous within the degrees of 368 * freedom, but not also over the components of a vectorized array. 369 * This happens typically on faces with DG where the cells have 370 * `interleaved_contiguous` storage but the faces' numbering is not the 371 * same as the cell's numbering. For a 372 * cell/face of this index type, the data access in FEEvaluationBase 373 * is directed to the array `dof_indices_contiguous` with the index 374 * `cell_index*n_vectorization*n_components`. 375 */ 376 interleaved_contiguous_strided, 377 /** 378 * Similar to interleaved_contiguous_separate storage, but for the case 379 * when the interleaved indices are not `n_vectorization apart`. This 380 * happens typically within the ghost layer of DG where the remote 381 * owner has applied an interleaved storage and the current processor 382 * only sees some of the cells. For a 383 * cell/face of this index type, the data access in FEEvaluationBase 384 * is directed to the array `dof_indices_contiguous` with the index 385 * `cell_index*n_vectorization*n_components`, including the array 386 * `dof_indices_interleave_strides` for the information about the 387 * actual stride. 388 */ 389 interleaved_contiguous_mixed_strides 390 }; 391 392 /** 393 * Enum used to distinguish the data arrays for the vectorization type 394 * in cells and faces. 395 */ 396 enum DoFAccessIndex : unsigned char 397 { 398 /** 399 * The data index for the faces designated as interior 400 */ 401 dof_access_face_interior = 0, 402 /** 403 * The data index for the faces designated as exterior 404 */ 405 dof_access_face_exterior = 1, 406 /** 407 * The data index for the cells 408 */ 409 dof_access_cell = 2 410 }; 411 412 /** 413 * Stores the dimension of the underlying DoFHandler. Since the indices 414 * are not templated, this is the variable that makes the dimension 415 * accessible in the (rare) cases it is needed inside this class. 416 */ 417 unsigned int dimension; 418 419 /** 420 * For efficiency reasons, always keep a fixed number of cells with 421 * similar properties together. This variable controls the number of 422 * cells batched together. As opposed to the other classes which are 423 * templated on the number type, this class as a pure index container is 424 * not templated, so we need to keep the information otherwise contained 425 * in VectorizedArrayType::size(). 426 */ 427 unsigned int vectorization_length; 428 429 /** 430 * Stores the index storage variant of all cell and face batches. 431 * 432 * The three arrays given here address the types for the faces decorated 433 * as interior (0), the faces decorated with as exterior (1), and the 434 * cells (2) according to CellOrFaceAccess. 435 */ 436 std::array<std::vector<IndexStorageVariants>, 3> index_storage_variants; 437 438 /** 439 * Stores the rowstart indices of the compressed row storage in the @p 440 * dof_indices and @p constraint_indicator fields. These two fields are 441 * always accessed together, so it is simpler to keep just one variable 442 * for them. This also obviates keeping two rowstart vectors in sync. 443 */ 444 std::vector<std::pair<unsigned int, unsigned int>> row_starts; 445 446 /** 447 * Stores the indices of the degrees of freedom for each cell. These 448 * indices are computed in MPI-local index space, i.e., each processor 449 * stores the locally owned indices as numbers between <tt>0</tt> and 450 * <tt>n_locally_owned_dofs-1</tt> and ghost indices in the range 451 * <tt>n_locally_owned_dofs</tt> to 452 * <tt>n_locally_owned_dofs+n_ghost_dofs</tt>. The translation between 453 * this MPI-local index space and the global numbering of degrees of 454 * freedom is stored in the @p vector_partitioner data structure. This 455 * array also includes the indirect contributions from constraints, 456 * which are described by the @p constraint_indicator field. Because of 457 * variable lengths of rows, this would be a vector of a vector. 458 * However, we use one contiguous memory region and store the rowstart 459 * in the variable @p row_starts. 460 */ 461 std::vector<unsigned int> dof_indices; 462 463 /** 464 * This variable describes the position of constraints in terms of the 465 * local numbering of degrees of freedom on a cell. The first number 466 * stores the distance from one constrained degree of freedom to the 467 * next. This allows to identify the position of constrained DoFs as we 468 * loop through the local degrees of freedom of the cell when reading 469 * from or writing to a vector. The second number stores the index of 470 * the constraint weights, stored in the variable constraint_pool_data. 471 */ 472 std::vector<std::pair<unsigned short, unsigned short>> 473 constraint_indicator; 474 475 /** 476 * Reordered index storage for `IndexStorageVariants::interleaved`. 477 */ 478 std::vector<unsigned int> dof_indices_interleaved; 479 480 /** 481 * Compressed index storage for faster access than through @p 482 * dof_indices used according to the description in IndexStorageVariants. 483 * 484 * The three arrays given here address the types for the faces decorated 485 * as interior (0), the faces decorated with as exterior (1), and the 486 * cells (2) according to CellOrFaceAccess. 487 */ 488 std::array<std::vector<unsigned int>, 3> dof_indices_contiguous; 489 490 /** 491 * Compressed index storage for faster access than through @p 492 * dof_indices used according to the description in IndexStorageVariants. 493 * 494 * The three arrays given here address the types for the faces decorated 495 * as minus (0), the faces decorated with as plus (1), and the cells 496 * (2). 497 */ 498 std::array<std::vector<unsigned int>, 3> dof_indices_interleave_strides; 499 500 /** 501 * Caches the number of indices filled when vectorizing. This 502 * information can implicitly deduced from the row_starts data fields, 503 * but this field allows for faster access. 504 * 505 * The three arrays given here address the types for the faces decorated 506 * as interior (0), the faces decorated with as exterior (1), and the 507 * cells (2) according to CellOrFaceAccess. 508 */ 509 std::array<std::vector<unsigned char>, 3> n_vectorization_lanes_filled; 510 511 /** 512 * This stores the parallel partitioning that can be used to set up 513 * vectors. The partitioner includes the description of the local range 514 * in the vector, and also includes how the ghosts look like. This 515 * enables initialization of vectors based on the DoFInfo field. 516 */ 517 std::shared_ptr<const Utilities::MPI::Partitioner> vector_partitioner; 518 519 /** 520 * This partitioning selects a subset of ghost indices to the full 521 * vector partitioner stored in @p vector_partitioner. These 522 * partitioners are used in specialized loops that only import parts of 523 * the ghosted region for reducing the amount of communication. There 524 * are five variants of the partitioner initialized: 525 * - one that queries only the cell values, 526 * - one that additionally describes the indices for 527 * evaluating the function values on relevant faces, 528 * - one that describes the indices for evaluation both the function 529 * values and the gradients on relevant faces adjacent to the locally 530 * owned cells, 531 * - one that additionally describes the indices for 532 * evaluating the function values on all faces, and 533 * - one that describes the indices for evaluation both the function 534 * values and the gradients on all faces adjacent to the locally owned 535 * cells. 536 */ 537 std::array<std::shared_ptr<const Utilities::MPI::Partitioner>, 5> 538 vector_partitioner_face_variants; 539 540 /** 541 * This stores a (sorted) list of all locally owned degrees of freedom 542 * that are constrained. 543 */ 544 std::vector<unsigned int> constrained_dofs; 545 546 /** 547 * Stores the rowstart indices of the compressed row storage in the @p 548 * plain_dof_indices fields. 549 */ 550 std::vector<unsigned int> row_starts_plain_indices; 551 552 /** 553 * Stores the indices of the degrees of freedom for each cell. This 554 * array does not include the indirect contributions from constraints, 555 * which are included in @p dof_indices. Because of variable lengths of 556 * rows, this would be a vector of a vector. However, we use one 557 * contiguous memory region and store the rowstart in the variable @p 558 * row_starts_plain_indices. 559 */ 560 std::vector<unsigned int> plain_dof_indices; 561 562 /** 563 * Stores the offset in terms of the number of base elements over all 564 * DoFInfo objects. 565 */ 566 unsigned int global_base_element_offset; 567 568 /** 569 * Stores the number of base elements in the DoFHandler where the 570 * indices have been read from. 571 */ 572 unsigned int n_base_elements; 573 574 /** 575 * Stores the number of components of each base element in the finite 576 * element where the indices have been read from. 577 */ 578 std::vector<unsigned int> n_components; 579 580 /** 581 * The ith entry of this vector stores the component number of the given 582 * base element. 583 */ 584 std::vector<unsigned int> start_components; 585 586 /** 587 * For a given component in an FESystem, this variable tells which base 588 * element the index belongs to. 589 */ 590 std::vector<unsigned int> component_to_base_index; 591 592 /** 593 * For a vector-valued element, this gives the constant offset in the 594 * number of degrees of freedom starting at the given component, as the 595 * degrees are numbered by degrees of freedom. This data structure does 596 * not take possible constraints and thus, shorter or longer lists, into 597 * account. This information is encoded in the row_starts variables 598 * directly. 599 * 600 * The outer vector goes through the various fe indices in the hp case, 601 * similarly to the @p dofs_per_cell variable. 602 */ 603 std::vector<std::vector<unsigned int>> component_dof_indices_offset; 604 605 /** 606 * Stores the number of degrees of freedom per cell. 607 */ 608 std::vector<unsigned int> dofs_per_cell; 609 610 /** 611 * Stores the number of degrees of freedom per face. 612 */ 613 std::vector<unsigned int> dofs_per_face; 614 615 /** 616 * Informs on whether plain indices are cached. 617 */ 618 bool store_plain_indices; 619 620 /** 621 * Stores the index of the active finite element in the hp case. 622 */ 623 std::vector<unsigned int> cell_active_fe_index; 624 625 /** 626 * Stores the maximum degree of different finite elements for the hp 627 * case. 628 */ 629 unsigned int max_fe_index; 630 631 /** 632 * To each of the slots in an hp adaptive case, the inner vector stores 633 * the corresponding element degree. This is used by the constructor of 634 * FEEvaluationBase to identify the correct data slot in the hp case. 635 */ 636 std::vector<std::vector<unsigned int>> fe_index_conversion; 637 638 /** 639 * Temporarily stores the numbers of ghosts during setup. Cleared when 640 * calling @p assign_ghosts. Then, all information is collected by the 641 * partitioner. 642 */ 643 std::vector<types::global_dof_index> ghost_dofs; 644 645 /** 646 * Stores an integer to each partition in TaskInfo that indicates 647 * whether to clear certain parts in the result vector if the user 648 * requested it with the respective argument in the MatrixFree::loop. 649 */ 650 std::vector<unsigned int> vector_zero_range_list_index; 651 652 /** 653 * Stores the actual ranges in the vector to be cleared. 654 */ 655 std::vector<std::pair<unsigned int, unsigned int>> vector_zero_range_list; 656 657 /** 658 * Stores an integer to each partition in TaskInfo that indicates when 659 * to schedule operations that will be done before any access to vector 660 * entries. 661 */ 662 std::vector<unsigned int> cell_loop_pre_list_index; 663 664 /** 665 * Stores the actual ranges of the operation before any access to vector 666 * entries. 667 */ 668 std::vector<std::pair<unsigned int, unsigned int>> cell_loop_pre_list; 669 670 /** 671 * Stores an integer to each partition in TaskInfo that indicates when 672 * to schedule operations that will be done after all access to vector 673 * entries. 674 */ 675 std::vector<unsigned int> cell_loop_post_list_index; 676 677 /** 678 * Stores the actual ranges of the operation after all access to vector 679 * entries. 680 */ 681 std::vector<std::pair<unsigned int, unsigned int>> cell_loop_post_list; 682 }; 683 684 685 /*-------------------------- Inline functions ---------------------------*/ 686 687 #ifndef DOXYGEN 688 689 690 inline unsigned int fe_index_from_degree(const unsigned int first_selected_component,const unsigned int fe_degree)691 DoFInfo::fe_index_from_degree(const unsigned int first_selected_component, 692 const unsigned int fe_degree) const 693 { 694 const unsigned int n_indices = fe_index_conversion.size(); 695 if (n_indices <= 1) 696 return 0; 697 for (unsigned int i = 0; i < n_indices; ++i) 698 if (fe_index_conversion[i][first_selected_component] == fe_degree) 699 return i; 700 return numbers::invalid_unsigned_int; 701 } 702 703 #endif // ifndef DOXYGEN 704 705 } // end of namespace MatrixFreeFunctions 706 } // end of namespace internal 707 708 DEAL_II_NAMESPACE_CLOSE 709 710 #endif 711