198 const std::vector<std::pair<types::global_dof_index, number2>> &entries)
200 next_constraint.first.resize(entries.size());
201 if (entries.size() > 0)
203 constraint_indices.resize(entries.size());
206 constraint_entries.assign(entries.begin(), entries.end());
207 std::sort(constraint_entries.begin(),
208 constraint_entries.end(),
209 [](
const std::pair<types::global_dof_index, double> &p1,
210 const std::pair<types::global_dof_index, double> &p2) {
211 return p1.second < p2.second;
217 constraint_indices[j] = constraint_entries[j].first;
220 next_constraint.first[j] = constraint_entries[j].second;
232 const auto position = constraints.find(next_constraint.first);
233 if (position != constraints.end())
234 insert_position = position->second;
237 next_constraint.second = constraints.size();
238 constraints.insert(next_constraint);
239 insert_position = next_constraint.second;
244 Assert(insert_position < (1 << (8 *
sizeof(
unsigned short))),
246 return static_cast<unsigned short>(insert_position);
255 const unsigned int n_cells,
256 const bool use_fast_hanging_node_algorithm)
258 this->dof_indices_per_cell.resize(n_cells);
259 this->plain_dof_indices_per_cell.resize(n_cells);
260 this->constraint_indicator_per_cell.resize(n_cells);
263 const bool has_hanging_nodes =
266 if (use_fast_hanging_node_algorithm && has_hanging_nodes)
268 hanging_nodes = std::make_unique<HangingNodes<dim>>(
271 hanging_node_constraint_masks.resize(n_cells);
275 lexicographic_numbering.resize(fes.size());
276 shape_infos.resize(fes.size());
278 for (
unsigned int i = 0; i < fes.size(); ++i)
280 if (fes[i].reference_cell().is_hyper_cube())
284 shape_infos[i].reinit(dummy_quadrature, fes[i], 0);
288 const auto dummy_quadrature =
289 fes[i].reference_cell().template get_gauss_type_quadrature<dim>(
291 shape_infos[i].reinit(dummy_quadrature, fes[i], 0);
294 lexicographic_numbering[i] = shape_infos[i].lexicographic_numbering;
296 active_fe_indices.resize(n_cells);
315 const unsigned int cell_no,
316 const unsigned int mg_level,
318 const ::AffineConstraints<typename Number::value_type> &constraints,
319 const std::shared_ptr<const Utilities::MPI::Partitioner> & partitioner)
321 std::vector<types::global_dof_index> local_dof_indices(
322 cell->get_fe().n_dofs_per_cell());
323 std::vector<types::global_dof_index> local_dof_indices_lex(
324 cell->get_fe().n_dofs_per_cell());
327 cell->get_dof_indices(local_dof_indices);
329 cell->get_mg_dof_indices(local_dof_indices);
334 const auto &lexicographic_numbering =
335 shape_infos[cell->active_fe_index()].lexicographic_numbering;
338 local_dof_indices.size());
340 for (
unsigned int i = 0; i < cell->get_fe().n_dofs_per_cell(); ++i)
341 local_dof_indices_lex[i] =
342 local_dof_indices[lexicographic_numbering[i]];
345 std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
352 auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
353 auto &dof_indices = this->dof_indices_per_cell[cell_no];
354 auto &plain_dof_indices = this->plain_dof_indices_per_cell[cell_no];
360 const auto global_to_local =
363 return partitioner->global_to_local(global_index);
369 plain_dof_indices.resize(local_dof_indices_lex.size());
370 for (
unsigned int i = 0; i < local_dof_indices_lex.size(); ++i)
371 plain_dof_indices[i] = global_to_local(local_dof_indices_lex[i]);
378 std::vector<ConstraintKinds>
mask(cell->get_fe().n_components());
379 hanging_nodes->setup_constraints(
380 cell, {}, lexicographic_numbering, local_dof_indices_lex,
mask);
382 hanging_node_constraint_masks[cell_no] =
compress(
mask[0], dim);
383 active_fe_indices[cell_no] = cell->active_fe_index();
386 for (
auto current_dof : local_dof_indices_lex)
388 const auto *entries_ptr =
389 constraints.get_constraint_entries(current_dof);
392 if (entries_ptr !=
nullptr)
394 const auto & entries = *entries_ptr;
396 if (n_entries == 1 &&
398 100 * std::numeric_limits<double>::epsilon())
400 current_dof = entries[0].first;
404 constraint_indicator.push_back(constraint_iterator);
405 constraint_indicator.back().second =
406 constraint_values.insert_entries(entries);
409 constraint_iterator.first = 0;
413 const std::vector<types::global_dof_index>
414 &constraint_indices = constraint_values.constraint_indices;
415 for (
unsigned int j = 0; j < n_entries; ++j)
417 dof_indices.push_back(
418 global_to_local(constraint_indices[j]));
425 dof_indices.push_back(global_to_local(current_dof));
429 Assert(constraint_iterator.first <
430 (1 << (8 *
sizeof(
unsigned short))) - 1,
432 constraint_iterator.first++;
442 const unsigned int cell_no,
443 const std::vector<types::global_dof_index> &local_dof_indices_lex,
444 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
446 const auto global_to_local =
449 return partitioner->global_to_local(global_index);
454 std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
456 auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
457 auto &dof_indices = this->dof_indices_per_cell[cell_no];
459 for (
const auto current_dof : local_dof_indices_lex)
465 std::pair<types::global_dof_index, typename Number::value_type>>
468 constraint_indicator.push_back(constraint_iterator);
469 constraint_indicator.back().second =
470 constraint_values.insert_entries(entries);
472 constraint_iterator.first = 0;
476 dof_indices.push_back(global_to_local(current_dof));
480 Assert(constraint_iterator.first <
481 (1 << (8 *
sizeof(
unsigned short))) - 1,
483 constraint_iterator.first++;
494 this->dof_indices = {};
495 this->plain_dof_indices = {};
496 this->constraint_indicator = {};
498 this->row_starts = {};
499 this->row_starts.emplace_back(0, 0);
501 if (plain_dof_indices_per_cell.empty() ==
false)
503 this->row_starts_plain_indices = {};
504 this->row_starts_plain_indices.emplace_back(0);
507 for (
unsigned int i = 0; i < dof_indices_per_cell.size(); ++i)
509 this->dof_indices.insert(this->dof_indices.end(),
510 dof_indices_per_cell[i].begin(),
511 dof_indices_per_cell[i].end());
512 this->constraint_indicator.insert(
513 this->constraint_indicator.end(),
514 constraint_indicator_per_cell[i].begin(),
515 constraint_indicator_per_cell[i].end());
517 this->row_starts.emplace_back(this->dof_indices.size(),
518 this->constraint_indicator.size());
520 if (plain_dof_indices_per_cell.empty() ==
false)
522 this->plain_dof_indices.insert(
523 this->plain_dof_indices.end(),
524 plain_dof_indices_per_cell[i].begin(),
525 plain_dof_indices_per_cell[i].end());
527 this->row_starts_plain_indices.emplace_back(
528 this->plain_dof_indices.size());
532 std::vector<const std::vector<double> *> constraints(
533 constraint_values.constraints.size());
534 unsigned int length = 0;
535 for (
const auto &it : constraint_values.constraints)
538 constraints[it.second] = &it.first;
539 length += it.first.size();
542 constraint_pool_data.clear();
543 constraint_pool_data.reserve(length);
544 constraint_pool_row_index.reserve(constraint_values.constraints.size() +
546 constraint_pool_row_index.resize(1, 0);
548 for (
const auto &constraint : constraints)
551 constraint_pool_data.insert(constraint_pool_data.end(),
554 constraint_pool_row_index.push_back(constraint_pool_data.size());
559 dof_indices_per_cell.clear();
560 constraint_indicator_per_cell.clear();
563 std::all_of(hanging_node_constraint_masks.begin(),
564 hanging_node_constraint_masks.end(),
566 return i == unconstrained_compressed_constraint_kind;
568 hanging_node_constraint_masks.clear();
578 VectorType & global_vector,
579 Number * local_vector,
580 const unsigned int first_cell,
581 const unsigned int n_cells,
582 const unsigned int n_dofs_per_cell,
583 const bool apply_constraints)
const
585 if ((row_starts_plain_indices.empty() ==
false) &&
586 (apply_constraints ==
false))
588 for (
unsigned int v = 0; v < n_cells; ++v)
590 const unsigned int cell_index = first_cell + v;
591 const unsigned int *dof_indices =
592 this->plain_dof_indices.data() +
595 for (
unsigned int i = 0; i < n_dofs_per_cell; ++dof_indices, ++i)
596 operation.process_dof(*dof_indices,
604 for (
unsigned int v = 0; v < n_cells; ++v)
606 const unsigned int cell_index = first_cell + v;
607 const unsigned int *dof_indices =
608 this->dof_indices.data() + this->row_starts[
cell_index].first;
609 unsigned int index_indicators = this->row_starts[
cell_index].second;
610 unsigned int next_index_indicators =
613 unsigned int ind_local = 0;
614 for (; index_indicators != next_index_indicators; ++index_indicators)
616 const std::pair<unsigned short, unsigned short> indicator =
617 this->constraint_indicator[index_indicators];
620 for (
unsigned int j = 0; j < indicator.first; ++j)
621 operation.process_dof(dof_indices[j],
623 local_vector[ind_local + j][v]);
625 ind_local += indicator.first;
626 dof_indices += indicator.first;
630 typename Number::value_type
value;
631 operation.pre_constraints(local_vector[ind_local][v],
value);
633 const typename Number::value_type *data_val =
634 this->constraint_pool_begin(indicator.second);
635 const typename Number::value_type *end_pool =
636 this->constraint_pool_end(indicator.second);
637 for (; data_val != end_pool; ++data_val, ++dof_indices)
638 operation.process_constraint(*dof_indices,
643 operation.post_constraints(
value, local_vector[ind_local][v]);
649 for (; ind_local < n_dofs_per_cell; ++dof_indices, ++ind_local)
650 operation.process_dof(*dof_indices,
652 local_vector[ind_local][v]);