libMesh
dof_map.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_DOF_MAP_H
21 #define LIBMESH_DOF_MAP_H
22 
23 // Local Includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/reference_counted_object.h"
26 #include "libmesh/libmesh.h" // libMesh::invalid_uint
27 #include "libmesh/variable.h"
28 #include "libmesh/threads.h"
29 #include "libmesh/threads_allocators.h"
30 #include "libmesh/elem_range.h"
31 #include "libmesh/ghosting_functor.h"
32 #include "libmesh/sparsity_pattern.h"
33 #include "libmesh/parallel_object.h"
34 #include "libmesh/point.h"
35 #include "libmesh/utility.h"
36 #include "libmesh/elem.h"
37 #include "libmesh/fe_interface.h"
38 #include "libmesh/libmesh_logging.h"
39 #include "libmesh/enum_elem_type.h"
40 #include "libmesh/mesh_subdivision_support.h"
41 #include "libmesh/dof_map_base.h"
42 
43 // TIMPI includes
45 #include "timpi/parallel_sync.h"
46 
47 // C++ Includes
48 #include <algorithm>
49 #include <cstddef>
50 #include <iterator>
51 #include <map>
52 #include <string>
53 #include <vector>
54 #include <memory>
55 
56 namespace libMesh
57 {
58 
59 // Forward Declarations
60 class CouplingMatrix;
61 class DefaultCoupling;
62 class DirichletBoundary;
63 class DirichletBoundaries;
64 class DofMap;
65 class DofObject;
66 class FEType;
67 class MeshBase;
68 class PeriodicBoundaryBase;
69 class PeriodicBoundaries;
70 class System;
71 class NonlinearImplicitSystem;
72 class StaticCondensationDofMap;
73 template <typename T> class DenseVectorBase;
74 template <typename T> class DenseVector;
75 template <typename T> class DenseMatrix;
76 template <typename T> class SparseMatrix;
77 template <typename T> class NumericVector;
78 enum Order : int;
79 
80 
81 
82 // ------------------------------------------------------------
83 // Do we need constraints for anything?
84 
85 #if defined(LIBMESH_ENABLE_AMR) || \
86  defined(LIBMESH_ENABLE_PERIODIC) || \
87  defined(LIBMESH_ENABLE_DIRICHLET)
88 # define LIBMESH_ENABLE_CONSTRAINTS 1
89 #endif
90 
91 // ------------------------------------------------------------
92 // AMR constraint matrix types
93 
94 #ifdef LIBMESH_ENABLE_CONSTRAINTS
95 
98 typedef std::map<dof_id_type, Real,
99  std::less<dof_id_type>,
101 
108 class DofConstraints : public std::map<dof_id_type,
109  DofConstraintRow,
110  std::less<dof_id_type>,
111  Threads::scalable_allocator<std::pair<const dof_id_type, DofConstraintRow>>>
112 {
113 };
114 
121  public std::map<dof_id_type, Number,
122  std::less<dof_id_type>,
123  Threads::scalable_allocator<std::pair<const dof_id_type, Number>>>
124 {
125 };
126 
132  public std::map<unsigned int, DofConstraintValueMap,
133  std::less<unsigned int>,
134  Threads::scalable_allocator
135  <std::pair<const unsigned int, DofConstraintValueMap>>>
136 {
137 };
138 
139 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
140 
146 typedef std::map<const Node *, Real,
147  std::less<const Node *>,
149 
156 class NodeConstraints : public std::map<const Node *,
157  std::pair<NodeConstraintRow,Point>,
158  std::less<const Node *>,
159  Threads::scalable_allocator<std::pair<const Node * const, std::pair<NodeConstraintRow,Point>>>>
160 {
161 };
162 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
163 
164 #endif // LIBMESH_ENABLE_CONSTRAINTS
165 
166 
167 
179 class DofMap : public DofMapBase,
180  public ReferenceCountedObject<DofMap>
181 {
182 public:
183 
189  explicit
190  DofMap(const unsigned int sys_number,
191  MeshBase & mesh);
192 
196  ~DofMap();
197 
202  {};
203 
209  {
210  public:
211  virtual ~AugmentSendList () = default;
212 
216  virtual void augment_send_list (std::vector<dof_id_type> & send_list) = 0;
217  };
218 
224  void attach_matrix (SparseMatrix<Number> & matrix);
225 
232  void update_sparsity_pattern(SparseMatrix<Number> & matrix) const;
233 
238  bool is_attached (SparseMatrix<Number> & matrix);
239 
246  std::size_t distribute_dofs (MeshBase &);
247 
253  void compute_sparsity (const MeshBase &);
254 
258  bool computed_sparsity_already () const;
259 
270  void set_constrained_sparsity_construction(bool use_constraints);
271 
276 
283 
287  void clear_sparsity();
288 
303 
308  void add_default_ghosting();
309 
334  void add_coupling_functor(GhostingFunctor & coupling_functor,
335  bool to_mesh = true);
336 
344  void add_coupling_functor(std::shared_ptr<GhostingFunctor> coupling_functor,
345  bool to_mesh = true)
346  { _shared_functors[coupling_functor.get()] = coupling_functor;
347  this->add_coupling_functor(*coupling_functor, to_mesh); }
348 
354  void remove_coupling_functor(GhostingFunctor & coupling_functor);
355 
359  std::set<GhostingFunctor *>::const_iterator coupling_functors_begin() const
360  { return _coupling_functors.begin(); }
361 
365  std::set<GhostingFunctor *>::const_iterator coupling_functors_end() const
366  { return _coupling_functors.end(); }
367 
372 
396  void add_algebraic_ghosting_functor(GhostingFunctor & evaluable_functor,
397  bool to_mesh = true);
398 
406  void add_algebraic_ghosting_functor(std::shared_ptr<GhostingFunctor> evaluable_functor,
407  bool to_mesh = true)
408  { _shared_functors[evaluable_functor.get()] = evaluable_functor;
409  this->add_algebraic_ghosting_functor(*evaluable_functor, to_mesh); }
410 
416  void remove_algebraic_ghosting_functor(GhostingFunctor & evaluable_functor);
417 
421  std::set<GhostingFunctor *>::const_iterator algebraic_ghosting_functors_begin() const
422  { return _algebraic_ghosting_functors.begin(); }
423 
427  std::set<GhostingFunctor *>::const_iterator algebraic_ghosting_functors_end() const
428  { return _algebraic_ghosting_functors.end(); }
429 
434 
446  {
448  }
449 
461  std::vector<dof_id_type> & n_nz,
462  std::vector<dof_id_type> & n_oz,
463  void *),
464  void * context = nullptr)
466 
475  {
476  _augment_send_list = &asl;
477  }
478 
483  void attach_extra_send_list_function(void (*func)(std::vector<dof_id_type> &, void *),
484  void * context = nullptr)
486 
493  void prepare_send_list ();
494 
501  {
502  _send_list.clear();
503  }
504 
514  void reinit_send_list (MeshBase & mesh);
515 
516 
526  const std::vector<dof_id_type> & get_send_list() const { return _send_list; }
527 
535  const std::vector<dof_id_type> & get_n_nz() const
536  {
538  return _sp->get_n_nz();
539  }
540 
548  const std::vector<dof_id_type> & get_n_oz() const
549  {
551  return _sp->get_n_oz();
552  }
553 
554 
564  {
565  return _sp.get();
566  }
567 
568  // /**
569  // * Add an unknown of order \p order and finite element type
570  // * \p type to the system of equations.
571  // */
572  // void add_variable (const Variable & var);
573 
578  void add_variable_group (VariableGroup var_group);
579 
590  void set_error_on_cyclic_constraint(bool error_on_cyclic_constraint);
591  void set_error_on_constraint_loop(bool error_on_constraint_loop);
592 
596  const VariableGroup & variable_group (const unsigned int c) const;
597 
598  const Variable & variable (const unsigned int c) const override;
599 
603  Order variable_order (const unsigned int c) const;
604 
608  Order variable_group_order (const unsigned int vg) const;
609 
613  const FEType & variable_type (const unsigned int c) const;
614 
618  const FEType & variable_group_type (const unsigned int vg) const;
619 
625  unsigned int n_variable_groups() const
626  { return cast_int<unsigned int>(_variable_groups.size()); }
627 
628  unsigned int n_variables() const override
629  { return cast_int<unsigned int>(_variables.size()); }
630 
634  unsigned int var_group_from_var_number(unsigned int var_num) const;
635 
642  {
643  return ((this->n_variable_groups() == 1) && (this->n_variables() > 1));
644  }
645 
658  unsigned int block_size() const
659  {
660  return (this->has_blocked_representation() ? this->n_variables() : 1);
661  }
662 
663  using DofMapBase::n_dofs;
668  dof_id_type n_dofs(const unsigned int vn) const
669  {
670  dof_id_type n = this->n_local_dofs(vn);
671  this->comm().sum(n);
672  return n;
673  }
674 
679 
686  dof_id_type n_local_dofs(const unsigned int vn) const
687  {
688  dof_id_type n;
689  this->local_variable_indices(n, _mesh, vn);
690  return n;
691  }
692 
697  std::vector<dof_id_type> n_dofs_per_processor(const unsigned int vn) const
698  {
699  std::vector<dof_id_type> n_local_dofs(this->n_processors(), 0);
700  this->comm().allgather(this->n_local_dofs(vn), n_local_dofs);
701  return n_local_dofs;
702  }
703 
708  { std::vector<dof_id_type>::const_iterator ub =
709  std::upper_bound(_end_df.begin(), _end_df.end(), dof);
710  libmesh_assert (ub != _end_df.end());
711  return cast_int<processor_id_type>(ub - _end_df.begin());
712  }
713 
714  void dof_indices (const Elem * const elem,
715  std::vector<dof_id_type> & di) const;
716 
722  void dof_indices (const Elem * const elem,
723  std::vector<dof_id_type> & di,
724  const unsigned int vn,
725  int p_level = -12345) const override;
726 
756  template <typename ScalarDofsFunctor, typename FieldDofsFunctor>
757  void dof_indices(const Elem * const elem,
758  std::vector<dof_id_type> & di,
759  const unsigned int vn,
760  ScalarDofsFunctor scalar_dofs_functor,
761  FieldDofsFunctor field_dofs_functor,
762  int p_level = -12345) const;
763 
768  void dof_indices (const Node * const node,
769  std::vector<dof_id_type> & di) const;
770 
775  void dof_indices (const Node * const node,
776  std::vector<dof_id_type> & di,
777  const unsigned int vn) const override;
778 
785  void dof_indices (const Elem & elem,
786  unsigned int n,
787  std::vector<dof_id_type> & di,
788  const unsigned int vn) const;
789 
790 #ifdef LIBMESH_ENABLE_AMR
791 
798  void old_dof_indices (const Elem & elem,
799  unsigned int n,
800  std::vector<dof_id_type> & di,
801  const unsigned int vn) const;
802 
803 #endif // LIBMESH_ENABLE_AMR
804 
813  void SCALAR_dof_indices (std::vector<dof_id_type> & di,
814  const unsigned int vn,
815  const bool old_dofs=false) const;
816 
824  bool semilocal_index (dof_id_type dof_index) const;
825 
833  bool all_semilocal_indices (const std::vector<dof_id_type> & dof_indices) const;
834 
839  bool local_index (dof_id_type dof_index) const
840  { return (dof_index >= this->first_dof()) && (dof_index < this->end_dof()); }
841 
847  template <typename DofObjectSubclass>
848  bool is_evaluable(const DofObjectSubclass & obj,
849  unsigned int var_num = libMesh::invalid_uint) const;
850 
858  void set_implicit_neighbor_dofs(bool implicit_neighbor_dofs);
859 
864 
871  bool use_coupled_neighbor_dofs(const MeshBase & mesh) const;
872 
884  const std::vector<dof_id_type> & dof_indices,
885  DenseVectorBase<Number> & Ue) const;
886 
892  template <typename T, std::enable_if_t<std::is_same_v<T, dof_id_type> ||
893  std::is_same_v<T, std::vector<dof_id_type>>, int> = 0>
894  void local_variable_indices(T & idx,
895  const MeshBase & mesh,
896  unsigned int var_num) const;
897 
903  template <typename T,
904  std::enable_if_t<std::is_same_v<T, dof_id_type> ||
905  std::is_same_v<T, std::vector<dof_id_type>>,
906  int> = 0>
907  void local_variable_indices(T & idx, unsigned int var_num) const
908  { this->local_variable_indices(idx, this->_mesh, var_num); }
909 
910 #ifdef LIBMESH_ENABLE_CONSTRAINTS
911 
912  //--------------------------------------------------------------------
913  // Constraint-specific methods
919 
925 
926 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
927 
932  { return cast_int<dof_id_type>(_node_constraints.size()); }
933 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
934 
943  void create_dof_constraints (const MeshBase &, Real time=0);
944 
949 
953  void scatter_constraints (MeshBase &);
954 
969  std::set<dof_id_type> & unexpanded_dofs,
970  bool look_for_constrainees);
971 
979  void process_constraints (MeshBase &);
980 
997 
1002  void add_constraint_row (const dof_id_type dof_number,
1003  const DofConstraintRow & constraint_row,
1004  const Number constraint_rhs,
1005  const bool forbid_constraint_overwrite);
1006 
1017  void add_adjoint_constraint_row (const unsigned int qoi_index,
1018  const dof_id_type dof_number,
1019  const DofConstraintRow & constraint_row,
1020  const Number constraint_rhs,
1021  const bool forbid_constraint_overwrite);
1022 
1028  void add_constraint_row (const dof_id_type dof_number,
1029  const DofConstraintRow & constraint_row,
1030  const bool forbid_constraint_overwrite = true)
1031  { add_constraint_row(dof_number, constraint_row, 0., forbid_constraint_overwrite); }
1032 
1036  DofConstraints::const_iterator constraint_rows_begin() const
1037  { return _dof_constraints.begin(); }
1038 
1042  DofConstraints::const_iterator constraint_rows_end() const
1043  { return _dof_constraints.end(); }
1044 
1050 
1052  {
1055  }
1056 
1058  {
1061  }
1062 
1076  {
1078  }
1079 
1080 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
1081 
1084  NodeConstraints::const_iterator node_constraint_rows_begin() const
1085  { return _node_constraints.begin(); }
1086 
1090  NodeConstraints::const_iterator node_constraint_rows_end() const
1091  { return _node_constraints.end(); }
1092 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
1093 
1098  bool is_constrained_dof (const dof_id_type dof) const;
1099 
1104  bool has_heterogeneous_adjoint_constraints (const unsigned int qoi_num) const;
1105 
1109  bool has_heterogenous_adjoint_constraints (const unsigned int qoi_num) const
1110  {
1111  return this->has_heterogeneous_adjoint_constraints (qoi_num);
1112  }
1113 
1119  Number has_heterogeneous_adjoint_constraint (const unsigned int qoi_num,
1120  const dof_id_type dof) const;
1121 
1125  Number has_heterogenous_adjoint_constraint (const unsigned int qoi_num,
1126  const dof_id_type dof) const
1127  {
1128  return this->has_heterogeneous_adjoint_constraint (qoi_num, dof);
1129  }
1130 
1136 
1141  bool is_constrained_node (const Node * node) const;
1142 
1149  void print_dof_constraints(std::ostream & os=libMesh::out,
1150  bool print_nonlocal=false) const;
1151 
1157  std::string get_local_constraints(bool print_nonlocal=false) const;
1158 
1159 
1168  std::pair<Real, Real> max_constraint_error(const System & system,
1169  NumericVector<Number> * v = nullptr) const;
1170 
1171 #endif // LIBMESH_ENABLE_CONSTRAINTS
1172 
1173  //--------------------------------------------------------------------
1174  // Constraint-specific methods
1175  // Some of these methods are enabled (but inlined away to nothing)
1176  // when constraints are disabled at configure-time. This is to
1177  // increase API compatibility of user code with different library
1178  // builds.
1179 
1194  std::vector<dof_id_type> & elem_dofs,
1195  bool asymmetric_constraint_rows = true) const;
1196 
1204  std::vector<dof_id_type> & row_dofs,
1205  std::vector<dof_id_type> & col_dofs,
1206  bool asymmetric_constraint_rows = true) const;
1207 
1212  std::vector<dof_id_type> & dofs,
1213  bool asymmetric_constraint_rows = true) const;
1214 
1224  DenseVector<Number> & rhs,
1225  std::vector<dof_id_type> & elem_dofs,
1226  bool asymmetric_constraint_rows = true) const;
1227 
1252  DenseVector<Number> & rhs,
1253  std::vector<dof_id_type> & elem_dofs,
1254  bool asymmetric_constraint_rows = true,
1255  int qoi_index = -1) const;
1256 
1257  /*
1258  * Backwards compatibility with misspelling.
1259  */
1261  DenseVector<Number> & rhs,
1262  std::vector<dof_id_type> & elem_dofs,
1263  bool asymmetric_constraint_rows = true,
1264  int qoi_index = -1) const
1265  {
1267  (matrix, rhs, elem_dofs, asymmetric_constraint_rows, qoi_index);
1268  }
1269 
1296  DenseVector<Number> & rhs,
1297  std::vector<dof_id_type> & elem_dofs,
1298  bool asymmetric_constraint_rows = true,
1299  int qoi_index = -1) const;
1300 
1301  /*
1302  * Backwards compatibility with misspelling.
1303  */
1305  DenseVector<Number> & rhs,
1306  std::vector<dof_id_type> & elem_dofs,
1307  bool asymmetric_constraint_rows = true,
1308  int qoi_index = -1) const
1309  {
1311  (matrix, rhs, elem_dofs, asymmetric_constraint_rows, qoi_index);
1312  }
1313 
1334  DenseVector<Number> & rhs,
1335  std::vector<dof_id_type> & elem_dofs,
1336  NumericVector<Number> & solution_local) const;
1337 
1353  std::vector<dof_id_type> & elem_dofs,
1354  NumericVector<Number> & solution_local) const;
1355 
1356 
1373  std::vector<dof_id_type> & elem_dofs,
1374  NumericVector<Number> & solution_local) const;
1375 
1385  DenseVector<Number> & w,
1386  std::vector<dof_id_type> & row_dofs,
1387  bool asymmetric_constraint_rows = true) const;
1388 
1395  void constrain_nothing (std::vector<dof_id_type> & dofs) const;
1396 
1410  void enforce_constraints_exactly (const System & system,
1411  NumericVector<Number> * v = nullptr,
1412  bool homogeneous = false) const;
1413 
1421  unsigned int q) const;
1422 
1424  NumericVector<Number> * rhs,
1425  NumericVector<Number> const * solution,
1426  bool homogeneous = true) const;
1427 
1429  SparseMatrix<Number> * jac) const;
1430 
1431 #ifdef LIBMESH_ENABLE_PERIODIC
1432 
1433  //--------------------------------------------------------------------
1434  // PeriodicBoundary-specific methods
1435 
1439  void add_periodic_boundary (const PeriodicBoundaryBase & periodic_boundary);
1440 
1447  void add_periodic_boundary (const PeriodicBoundaryBase & boundary, const PeriodicBoundaryBase & inverse_boundary);
1448 
1453  bool is_periodic_boundary (const boundary_id_type boundaryid) const;
1454 
1456  {
1457  return _periodic_boundaries.get();
1458  }
1459 
1461  {
1462  return _periodic_boundaries.get();
1463  }
1464 
1465 #endif // LIBMESH_ENABLE_PERIODIC
1466 
1467 
1468 #ifdef LIBMESH_ENABLE_DIRICHLET
1469 
1470  //--------------------------------------------------------------------
1471  // DirichletBoundary-specific methods
1472 
1486  void add_dirichlet_boundary (const DirichletBoundary & dirichlet_boundary);
1487 
1493  void add_adjoint_dirichlet_boundary (const DirichletBoundary & dirichlet_boundary,
1494  unsigned int q);
1495 
1499  void remove_dirichlet_boundary (const DirichletBoundary & dirichlet_boundary);
1500 
1505  void remove_adjoint_dirichlet_boundary (const DirichletBoundary & dirichlet_boundary,
1506  unsigned int q);
1507 
1509  {
1510  return _dirichlet_boundaries.get();
1511  }
1512 
1514  {
1515  return _dirichlet_boundaries.get();
1516  }
1517 
1518  bool has_adjoint_dirichlet_boundaries(unsigned int q) const;
1519 
1520  const DirichletBoundaries *
1521  get_adjoint_dirichlet_boundaries(unsigned int q) const;
1522 
1524  get_adjoint_dirichlet_boundaries(unsigned int q);
1525 
1531  const DirichletBoundary & boundary) const;
1532 #endif // LIBMESH_ENABLE_DIRICHLET
1533 
1534 
1535 #ifdef LIBMESH_ENABLE_AMR
1536 
1537  //--------------------------------------------------------------------
1538  // AMR-specific methods
1539 
1548  // void augment_send_list_for_projection(const MeshBase &);
1549 
1550 #ifdef LIBMESH_ENABLE_AMR
1551 
1558  void old_dof_indices (const Elem * const elem,
1559  std::vector<dof_id_type> & di,
1560  const unsigned int vn = libMesh::invalid_uint) const;
1561 
1562 #endif // LIBMESH_ENABLE_AMR
1563 
1569  void constrain_p_dofs (unsigned int var,
1570  const Elem * elem,
1571  unsigned int s,
1572  unsigned int p);
1573 
1574 #endif // LIBMESH_ENABLE_AMR
1575 
1579  void reinit
1580  (MeshBase & mesh,
1581  const std::map<const Node *, std::set<subdomain_id_type>> &
1582  constraining_subdomains);
1583 
1588  virtual void clear () override;
1589 
1593  void print_info(std::ostream & os=libMesh::out) const;
1594 
1598  std::string get_info() const;
1599 
1614 
1618  unsigned int sys_number() const;
1619 
1634  std::unique_ptr<SparsityPattern::Build> build_sparsity(const MeshBase & mesh,
1635  bool calculate_constrained = false,
1636  bool use_condensed_system = false) const;
1637 
1642  void should_p_refine(unsigned int g, bool p_refine);
1643 
1647  bool should_p_refine(unsigned int g) const;
1648 
1652  bool should_p_refine_var(unsigned int var) const;
1653 
1654  // Prevent bad user implicit conversions
1655  void should_p_refine(FEFamily, bool) = delete;
1656  void should_p_refine(Order, bool) = delete;
1657  bool should_p_refine(FEFamily) const = delete;
1658  bool should_p_refine(Order) const = delete;
1659 
1663  void create_static_condensation(MeshBase & mesh, System & system);
1664 
1668  bool has_static_condensation() const { return _sc.get(); }
1669 
1675 
1680 
1681 private:
1682 
1693  void _dof_indices (const Elem & elem,
1694  int p_level,
1695  std::vector<dof_id_type> & di,
1696  const unsigned int vg,
1697  const unsigned int vig,
1698  const Node * const * nodes,
1699  unsigned int n_nodes,
1700  const unsigned int v
1701 #ifdef DEBUG
1702  ,
1703  std::size_t & tot_size
1704 #endif
1705  ) const;
1706 
1722  template <typename FieldDofsFunctor>
1723  void _dof_indices (const Elem & elem,
1724  int p_level,
1725  std::vector<dof_id_type> & di,
1726  const unsigned int vg,
1727  const unsigned int vig,
1728  const Node * const * nodes,
1729  unsigned int n_nodes,
1730  const unsigned int v,
1731 #ifdef DEBUG
1732  std::size_t & tot_size,
1733 #endif
1734  FieldDofsFunctor field_dofs_functor) const;
1735 
1740  void _node_dof_indices (const Elem & elem,
1741  unsigned int n,
1742  const DofObject & obj,
1743  std::vector<dof_id_type> & di,
1744  const unsigned int vn) const;
1745 
1749  void invalidate_dofs(MeshBase & mesh) const;
1750 
1754  DofObject * node_ptr(MeshBase & mesh, dof_id_type i) const;
1755 
1759  DofObject * elem_ptr(MeshBase & mesh, dof_id_type i) const;
1760 
1764  typedef DofObject * (DofMap::*dofobject_accessor)
1765  (MeshBase & mesh, dof_id_type i) const;
1766 
1770  template<typename iterator_type>
1771  void set_nonlocal_dof_objects(iterator_type objects_begin,
1772  iterator_type objects_end,
1773  MeshBase & mesh,
1774  dofobject_accessor objects);
1775 
1790  std::map<const Node *, std::set<subdomain_id_type>>
1792 
1806  (dof_id_type & next_free_dof,
1807  MeshBase & mesh,
1808  const std::map<const Node *, std::set<subdomain_id_type>> &
1809  constraining_subdomains);
1810 
1827  (dof_id_type & next_free_dof,
1828  MeshBase & mesh,
1829  const std::map<const Node *, std::set<subdomain_id_type>> &
1830  constraining_subdomains);
1831 
1832  /*
1833  * Helper method for the above two to count + distriubte SCALAR dofs
1834  */
1835  void distribute_scalar_dofs (dof_id_type & next_free_dof);
1836 
1837 #ifdef DEBUG
1838  /*
1839  * Internal assertions for distribute_local_dofs_*
1840  */
1842 #endif
1843 
1844  /*
1845  * A utility method for obtaining a set of elements to ghost along
1846  * with merged coupling matrices.
1847  */
1848  typedef std::set<std::unique_ptr<CouplingMatrix>, Utility::CompareUnderlying> CouplingMatricesSet;
1849  static void
1851  CouplingMatricesSet & temporary_coupling_matrices,
1852  const std::set<GhostingFunctor *>::iterator & gf_begin,
1853  const std::set<GhostingFunctor *>::iterator & gf_end,
1854  const MeshBase::const_element_iterator & elems_begin,
1855  const MeshBase::const_element_iterator & elems_end,
1856  processor_id_type p);
1857 
1863 
1864 #ifdef LIBMESH_ENABLE_CONSTRAINTS
1865 
1877  std::vector<dof_id_type> & elem_dofs,
1878  const bool called_recursively=false) const;
1879 
1896  DenseVector<Number> & H,
1897  std::vector<dof_id_type> & elem_dofs,
1898  int qoi_index = -1,
1899  const bool called_recursively=false) const;
1900 
1905  void find_connected_dofs (std::vector<dof_id_type> & elem_dofs) const;
1906 
1911  void find_connected_dof_objects (std::vector<const DofObject *> & objs) const;
1912 
1919 
1926 
1927 #endif // LIBMESH_ENABLE_CONSTRAINTS
1928 
1935 
1941 
1945  std::vector<Variable> _variables;
1946 
1950  std::vector<VariableGroup> _variable_groups;
1951 
1955  std::vector<unsigned int> _variable_group_numbers;
1956 
1960  std::unordered_map<unsigned int, unsigned int> _var_to_vg;
1961 
1965  const unsigned int _sys_number;
1966 
1971 
1977  std::vector<SparseMatrix<Number> * > _matrices;
1978 
1983  std::vector<dof_id_type> _first_scalar_df;
1984 
1989  std::vector<dof_id_type> _send_list;
1990 
1995 
2000  std::vector<dof_id_type> & n_nz,
2001  std::vector<dof_id_type> & n_oz,
2002  void *);
2007 
2012 
2016  void (*_extra_send_list_function)(std::vector<dof_id_type> &, void *);
2017 
2022 
2029  std::unique_ptr<DefaultCoupling> _default_coupling;
2030 
2037  std::unique_ptr<DefaultCoupling> _default_evaluating;
2038 
2047  std::set<GhostingFunctor *> _algebraic_ghosting_functors;
2048 
2060  std::set<GhostingFunctor *> _coupling_functors;
2061 
2066  std::map<GhostingFunctor *, std::shared_ptr<GhostingFunctor> > _shared_functors;
2067 
2073 
2079  std::unique_ptr<SparsityPattern::Build> _sp;
2080 
2086 
2087 #ifdef LIBMESH_ENABLE_AMR
2088 
2093  std::vector<dof_id_type> _first_old_scalar_df;
2094 
2098  std::unordered_set<unsigned int> _dont_p_refine;
2099 #endif
2100 
2101 #ifdef LIBMESH_ENABLE_CONSTRAINTS
2102 
2107 
2109 
2111 #endif
2112 
2113 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
2114 
2118 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
2119 
2120 
2121 #ifdef LIBMESH_ENABLE_PERIODIC
2122 
2126  std::unique_ptr<PeriodicBoundaries> _periodic_boundaries;
2127 #endif
2128 
2129 #ifdef LIBMESH_ENABLE_DIRICHLET
2130 
2134  std::unique_ptr<DirichletBoundaries> _dirichlet_boundaries;
2135 
2140  std::vector<std::unique_ptr<DirichletBoundaries>> _adjoint_dirichlet_boundaries;
2141 #endif
2142 
2144 
2151 
2163 
2165  std::unique_ptr<StaticCondensationDofMap> _sc;
2166 };
2167 
2168 
2169 // ------------------------------------------------------------
2170 // Dof Map inline member functions
2171 inline
2172 unsigned int DofMap::sys_number() const
2173 {
2174  return _sys_number;
2175 }
2176 
2177 
2178 
2179 inline
2180 const VariableGroup & DofMap::variable_group (const unsigned int g) const
2181 {
2182  libmesh_assert_less (g, _variable_groups.size());
2183 
2184  return _variable_groups[g];
2185 }
2186 
2187 
2188 
2189 inline
2190 const Variable & DofMap::variable (const unsigned int c) const
2191 {
2192  libmesh_assert_less (c, _variables.size());
2193 
2194  return _variables[c];
2195 }
2196 
2197 
2198 
2199 inline
2200 Order DofMap::variable_order (const unsigned int c) const
2201 {
2202  libmesh_assert_less (c, _variables.size());
2203 
2204  return _variables[c].type().order;
2205 }
2206 
2207 
2208 
2209 inline
2210 Order DofMap::variable_group_order (const unsigned int vg) const
2211 {
2212  libmesh_assert_less (vg, _variable_groups.size());
2213 
2214  return _variable_groups[vg].type().order;
2215 }
2216 
2217 
2218 
2219 inline
2220 const FEType & DofMap::variable_type (const unsigned int c) const
2221 {
2222  libmesh_assert_less (c, _variables.size());
2223 
2224  return _variables[c].type();
2225 }
2226 
2227 
2228 
2229 inline
2230 const FEType & DofMap::variable_group_type (const unsigned int vg) const
2231 {
2232  libmesh_assert_less (vg, _variable_groups.size());
2233 
2234  return _variable_groups[vg].type();
2235 }
2236 
2237 
2238 #ifdef LIBMESH_ENABLE_CONSTRAINTS
2239 
2240 
2241 inline
2243 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
2244  node
2245 #endif
2246  ) const
2247 {
2248 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
2249  if (_node_constraints.count(node))
2250  return true;
2251 #endif
2252 
2253  return false;
2254 }
2255 
2256 
2257 inline
2259 {
2260  if (_dof_constraints.count(dof))
2261  return true;
2262 
2263  return false;
2264 }
2265 
2266 
2267 inline
2268 bool DofMap::has_heterogeneous_adjoint_constraints (const unsigned int qoi_num) const
2269 {
2270  AdjointDofConstraintValues::const_iterator it =
2271  _adjoint_constraint_values.find(qoi_num);
2272  if (it == _adjoint_constraint_values.end())
2273  return false;
2274  if (it->second.empty())
2275  return false;
2276 
2277  return true;
2278 }
2279 
2280 
2281 inline
2283  const dof_id_type dof) const
2284 {
2285  AdjointDofConstraintValues::const_iterator it =
2286  _adjoint_constraint_values.find(qoi_num);
2287  if (it != _adjoint_constraint_values.end())
2288  {
2289  DofConstraintValueMap::const_iterator rhsit =
2290  it->second.find(dof);
2291  if (rhsit == it->second.end())
2292  return 0;
2293  else
2294  return rhsit->second;
2295  }
2296 
2297  return 0;
2298 }
2299 
2300 
2301 
2302 inline
2304 {
2306 }
2307 
2308 
2309 
2310 #else
2311 
2312 //--------------------------------------------------------------------
2313 // Constraint-specific methods get inlined into nothing if
2314 // constraints are disabled, so there's no reason for users not to
2315 // use them.
2316 
2318  std::vector<dof_id_type> &,
2319  bool) const {}
2320 
2322  std::vector<dof_id_type> &,
2323  std::vector<dof_id_type> &,
2324  bool) const {}
2325 
2327  std::vector<dof_id_type> &,
2328  bool) const {}
2329 
2332  std::vector<dof_id_type> &,
2333  bool) const {}
2334 
2337  std::vector<dof_id_type> &, bool, int) const {}
2338 
2341  std::vector<dof_id_type> &, bool, int) const {}
2342 
2345  std::vector<dof_id_type> &,
2346  bool) const {}
2347 
2348 inline void DofMap::constrain_nothing (std::vector<dof_id_type> &) const {}
2349 
2352  bool) const {}
2353 
2355  unsigned int) const {}
2356 
2357 
2361  NumericVector<Number> const *,
2362  bool) const {}
2363 
2366  SparseMatrix<Number> *) const {}
2367 
2368 #endif // LIBMESH_ENABLE_CONSTRAINTS
2369 
2370 
2371 
2372 inline
2374 {
2375  // This got only partly finished...
2376  if (use_constraints)
2377  libmesh_not_implemented();
2378 
2379 #ifdef LIBMESH_ENABLE_CONSTRAINTS
2380  _constrained_sparsity_construction = use_constraints;
2381 #endif
2382  libmesh_ignore(use_constraints);
2383 }
2384 
2385 inline
2387 {
2389 }
2390 
2391 inline
2393 {
2394 #ifdef LIBMESH_ENABLE_CONSTRAINTS
2396 #else
2397  return true;
2398 #endif
2399 }
2400 
2401 inline
2402 void DofMap::should_p_refine(const unsigned int g, const bool p_refine)
2403 {
2404 #ifdef LIBMESH_ENABLE_AMR
2405  if (p_refine)
2406  {
2407  auto it = _dont_p_refine.find(g);
2408  if (it != _dont_p_refine.end())
2409  _dont_p_refine.erase(it);
2410  }
2411  else
2412  _dont_p_refine.insert(g);
2413 #else
2414  libmesh_ignore(g, p_refine);
2415 #endif
2416 }
2417 
2418 inline
2419 bool DofMap::should_p_refine(const unsigned int g) const
2420 {
2421 #ifdef LIBMESH_ENABLE_AMR
2422  return !_dont_p_refine.count(g);
2423 #else
2424  libmesh_ignore(g);
2425  return false;
2426 #endif
2427 }
2428 
2429 inline
2430 unsigned int DofMap::var_group_from_var_number(const unsigned int var_num) const
2431 {
2432  libmesh_assert(var_num < n_variables());
2433  return libmesh_map_find(_var_to_vg, var_num);
2434 }
2435 
2436 inline
2437 bool DofMap::should_p_refine_var(const unsigned int var) const
2438 {
2439 #ifdef LIBMESH_ENABLE_AMR
2440  const auto vg = this->var_group_from_var_number(var);
2441  return !_dont_p_refine.count(vg);
2442 #else
2443  libmesh_ignore(var);
2444  return false;
2445 #endif
2446 }
2447 
2448 template <typename FieldDofsFunctor>
2449 void DofMap::_dof_indices (const Elem & elem,
2450  int p_level,
2451  std::vector<dof_id_type> & di,
2452  const unsigned int vg,
2453  const unsigned int vig,
2454  const Node * const * nodes,
2455  unsigned int n_nodes,
2456  const unsigned int v,
2457 #ifdef DEBUG
2458  std::size_t & tot_size,
2459 #endif
2460  FieldDofsFunctor field_dofs_functor) const
2461 {
2462  const VariableGroup & var = this->variable_group(vg);
2463 
2464  if (var.active_on_subdomain(elem.subdomain_id()))
2465  {
2466  const ElemType type = elem.type();
2467  const unsigned int sys_num = this->sys_number();
2468 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
2469  const bool is_inf = elem.infinite();
2470 #endif
2471 
2472  const bool extra_hanging_dofs =
2474 
2475  FEType fe_type = var.type();
2476 
2477  const bool add_p_level =
2478 #ifdef LIBMESH_ENABLE_AMR
2479  !_dont_p_refine.count(vg);
2480 #else
2481  false;
2482 #endif
2483 
2484 #ifdef DEBUG
2485  // The number of dofs per element is non-static for subdivision FE
2486  if (var.type().family == SUBDIVISION)
2487  tot_size += n_nodes;
2488  else
2489  // FIXME: Is the passed-in p_level just elem.p_level()? If so,
2490  // this seems redundant.
2491  tot_size += FEInterface::n_dofs(fe_type, add_p_level*p_level, &elem);
2492 #endif
2493 
2494  // The total Order is not required when getting the function
2495  // pointer, it is only needed when the function is called (see
2496  // below).
2497  const FEInterface::n_dofs_at_node_ptr ndan =
2498  FEInterface::n_dofs_at_node_function(fe_type, &elem);
2499 
2500  // Get the node-based DOF numbers
2501  for (unsigned int n=0; n != n_nodes; n++)
2502  {
2503  const Node & node = *nodes[n];
2504 
2505  // Cache the intermediate lookups that are common to every
2506  // component
2507 #ifdef DEBUG
2508  const std::pair<unsigned int, unsigned int>
2509  vg_and_offset = node.var_to_vg_and_offset(sys_num,v);
2510  libmesh_assert_equal_to (vg, vg_and_offset.first);
2511  libmesh_assert_equal_to (vig, vg_and_offset.second);
2512 #endif
2513  const unsigned int n_comp = node.n_comp_group(sys_num,vg);
2514 
2515  // There is a potential problem with h refinement. Imagine a
2516  // quad9 that has a linear FE on it. Then, on the hanging side,
2517  // it can falsely identify a DOF at the mid-edge node. This is why
2518  // we go through FEInterface instead of node.n_comp() directly.
2519  const unsigned int nc =
2520 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
2521  is_inf ?
2522  FEInterface::n_dofs_at_node(fe_type, add_p_level*p_level, &elem, n) :
2523 #endif
2524  ndan (type, fe_type.order + add_p_level*p_level, n);
2525 
2526  // If this is a non-vertex on a hanging node with extra
2527  // degrees of freedom, we use the non-vertex dofs (which
2528  // come in reverse order starting from the end, to
2529  // simplify p refinement)
2530  if (extra_hanging_dofs && !elem.is_vertex(n))
2531  {
2532  const int dof_offset = n_comp - nc;
2533 
2534  // We should never have fewer dofs than necessary on a
2535  // node unless we're getting indices on a parent element,
2536  // and we should never need the indices on such a node
2537  if (dof_offset < 0)
2538  {
2539  libmesh_assert(!elem.active());
2540  di.resize(di.size() + nc, DofObject::invalid_id);
2541  }
2542  else
2543  for (int i=int(n_comp)-1; i>=dof_offset; i--)
2544  {
2545  const dof_id_type d =
2546  node.dof_number(sys_num, vg, vig, i, n_comp);
2547  libmesh_assert_not_equal_to (d, DofObject::invalid_id);
2548  field_dofs_functor(elem, n, v, di, d);
2549  }
2550  }
2551  // If this is a vertex or an element without extra hanging
2552  // dofs, our dofs come in forward order coming from the
2553  // beginning
2554  else
2555  {
2556  // We have a good component index only if it's being
2557  // used on this FE type (nc) *and* it's available on
2558  // this DofObject (n_comp).
2559  const unsigned int good_nc = std::min(n_comp, nc);
2560  for (unsigned int i=0; i!=good_nc; ++i)
2561  {
2562  const dof_id_type d =
2563  node.dof_number(sys_num, vg, vig, i, n_comp);
2564  libmesh_assert_not_equal_to (d, DofObject::invalid_id);
2565  libmesh_assert_less (d, this->n_dofs());
2566  field_dofs_functor(elem, n, v, di, d);
2567  }
2568 
2569  // With fewer good component indices than we need, e.g.
2570  // due to subdomain expansion, the remaining expected
2571  // indices are marked invalid.
2572  if (n_comp < nc)
2573  for (unsigned int i=n_comp; i!=nc; ++i)
2574  di.push_back(DofObject::invalid_id);
2575  }
2576  }
2577 
2578  // If there are any element-based DOF numbers, get them
2579  const unsigned int nc = FEInterface::n_dofs_per_elem(fe_type, add_p_level*p_level, &elem);
2580 
2581  // We should never have fewer dofs than necessary on an
2582  // element unless we're getting indices on a parent element
2583  // (and we should never need those indices) or off-domain for a
2584  // subdomain-restricted variable (where invalid_id is the
2585  // correct thing to return)
2586  if (nc != 0)
2587  {
2588  const unsigned int n_comp = elem.n_comp_group(sys_num,vg);
2589  if (elem.n_systems() > sys_num && nc <= n_comp)
2590  {
2591  for (unsigned int i=0; i<nc; i++)
2592  {
2593  const dof_id_type d =
2594  elem.dof_number(sys_num, vg, vig, i, n_comp);
2595  libmesh_assert_not_equal_to (d, DofObject::invalid_id);
2596 
2597  field_dofs_functor(elem, invalid_uint, v, di, d);
2598  }
2599  }
2600  else
2601  {
2602  libmesh_assert(!elem.active() || fe_type.family == LAGRANGE || fe_type.family == SUBDIVISION);
2603  di.resize(di.size() + nc, DofObject::invalid_id);
2604  }
2605  }
2606  }
2607 }
2608 
2609 
2610 
2611 template <typename ScalarDofsFunctor, typename FieldDofsFunctor>
2612 void DofMap::dof_indices (const Elem * const elem,
2613  std::vector<dof_id_type> & di,
2614  const unsigned int vn,
2615  ScalarDofsFunctor scalar_dofs_functor,
2616  FieldDofsFunctor field_dofs_functor,
2617  int p_level) const
2618 {
2619  // We now allow elem==nullptr to request just SCALAR dofs
2620  // libmesh_assert(elem);
2621 
2622  LOG_SCOPE("dof_indices()", "DofMap");
2623 
2624  // Clear the DOF indices vector
2625  di.clear();
2626 
2627  // Use the default p refinement level?
2628  if (p_level == -12345)
2629  p_level = elem ? elem->p_level() : 0;
2630 
2631  const unsigned int vg = this->_variable_group_numbers[vn];
2632  const VariableGroup & var = this->variable_group(vg);
2633  const unsigned int vig = vn - var.number();
2634 
2635 #ifdef DEBUG
2636  // Check that sizes match in DEBUG mode
2637  std::size_t tot_size = 0;
2638 #endif
2639 
2640  if (elem && elem->type() == TRI3SUBDIVISION)
2641  {
2642  // Subdivision surface FE require the 1-ring around elem
2643  const Tri3Subdivision * sd_elem = static_cast<const Tri3Subdivision *>(elem);
2644 
2645  // Ghost subdivision elements have no real dofs
2646  if (!sd_elem->is_ghost())
2647  {
2648  // Determine the nodes contributing to element elem
2649  std::vector<const Node *> elem_nodes;
2650  MeshTools::Subdivision::find_one_ring(sd_elem, elem_nodes);
2651 
2652  _dof_indices(*elem, p_level, di, vg, vig, elem_nodes.data(),
2653  cast_int<unsigned int>(elem_nodes.size()), vn,
2654 #ifdef DEBUG
2655  tot_size,
2656 #endif
2657  field_dofs_functor);
2658  }
2659 
2660  return;
2661  }
2662 
2663  // Get the dof numbers
2664  if (var.type().family == SCALAR &&
2665  (!elem ||
2666  var.active_on_subdomain(elem->subdomain_id())))
2667  {
2668 #ifdef DEBUG
2669  tot_size += var.type().order;
2670 #endif
2671  std::vector<dof_id_type> di_new;
2672  this->SCALAR_dof_indices(di_new,vn);
2673  scalar_dofs_functor(*elem, di, di_new);
2674  }
2675  else if (elem)
2676  _dof_indices(*elem, p_level, di, vg, vig, elem->get_nodes(),
2677  elem->n_nodes(), vn,
2678 #ifdef DEBUG
2679  tot_size,
2680 #endif
2681  field_dofs_functor);
2682 
2683 #ifdef DEBUG
2684  libmesh_assert_equal_to (tot_size, di.size());
2685 #endif
2686 }
2687 
2688 inline
2690 {
2692  return *_sc;
2693 }
2694 
2695 } // namespace libMesh
2696 
2697 #endif // LIBMESH_DOF_MAP_H
std::vector< VariableGroup > _variable_groups
The finite element type for each variable group.
Definition: dof_map.h:1950
void remove_adjoint_dirichlet_boundary(const DirichletBoundary &dirichlet_boundary, unsigned int q)
Removes from the system the specified Dirichlet boundary for the adjoint equation defined by Quantity...
void find_connected_dofs(std::vector< dof_id_type > &elem_dofs) const
Finds all the DOFS associated with the element DOFs elem_dofs.
Definition: dof_map.C:2855
std::unique_ptr< SparsityPattern::Build > _sp
The sparsity pattern of the global matrix.
Definition: dof_map.h:2079
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:530
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:196
std::string get_local_constraints(bool print_nonlocal=false) const
Gets a string reporting all DoF and Node constraints local to this processor.
FEFamily family
The type of finite element.
Definition: fe_type.h:221
void constrain_element_dyad_matrix(DenseVector< Number > &v, DenseVector< Number > &w, std::vector< dof_id_type > &row_dofs, bool asymmetric_constraint_rows=true) const
Constrains a dyadic element matrix B = v w&#39;.
Definition: dof_map.h:2343
A class holding degree of freedom information pertinent to static condensation.
bool _implicit_neighbor_dofs_initialized
Bools to indicate if we override the –implicit_neighbor_dofs commandline options.
Definition: dof_map.h:2149
ElemType
Defines an enum for geometric element types.
bool is_constrained_node(const Node *node) const
Definition: dof_map.h:2242
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
unsigned int n_variable_groups() const
Definition: dof_map.h:625
void add_adjoint_dirichlet_boundary(const DirichletBoundary &dirichlet_boundary, unsigned int q)
Adds a copy of the specified Dirichlet boundary to the system, corresponding to the adjoint problem d...
const unsigned int _sys_number
The number of the system we manage DOFs for.
Definition: dof_map.h:1965
bool _implicit_neighbor_dofs
Definition: dof_map.h:2150
DofConstraintValueMap & get_primal_constraint_values()
Definition: dof_map.h:2303
void check_dirichlet_bcid_consistency(const MeshBase &mesh, const DirichletBoundary &boundary) const
Check that all the ids in dirichlet_bcids are actually present in the mesh.
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
Definition: dof_object.h:1032
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
DefaultCoupling & default_coupling()
Default coupling functor.
Definition: dof_map.h:371
const std::vector< dof_id_type > & get_n_oz() const
Definition: dof_map.h:548
A Node is like a Point, but with more information.
Definition: node.h:52
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:355
This abstract base class defines the interface by which library code and user code can report associa...
dof_id_type n_SCALAR_dofs() const
Definition: dof_map.h:678
std::set< GhostingFunctor * >::const_iterator coupling_functors_begin() const
Beginning of range of coupling functors.
Definition: dof_map.h:359
void build_constraint_matrix_and_vector(DenseMatrix< Number > &C, DenseVector< Number > &H, std::vector< dof_id_type > &elem_dofs, int qoi_index=-1, const bool called_recursively=false) const
Build the constraint matrix C and the forcing vector H associated with the element degree of freedom ...
~DofMap()
Destructor.
Definition: dof_map.C:202
This helper class can be called on multiple threads to compute the sparsity pattern (or graph) of the...
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:310
void local_variable_indices(T &idx, unsigned int var_num) const
If T == dof_id_type, counts, if T == std::vector<dof_id_type>, fills an array of, those dof indices w...
Definition: dof_map.h:907
processor_id_type dof_owner(const dof_id_type dof) const
Definition: dof_map.h:707
void scatter_constraints(MeshBase &)
Sends constraint equations to constraining processors.
void reinit(MeshBase &mesh, const std::map< const Node *, std::set< subdomain_id_type >> &constraining_subdomains)
Reinitialize the underlying data structures conformal to the current mesh.
Definition: dof_map.C:503
void add_periodic_boundary(const PeriodicBoundaryBase &periodic_boundary)
Adds a copy of the specified periodic boundary to the system.
bool _error_on_constraint_loop
This flag indicates whether or not we do an opt-mode check for the presence of constraint loops...
Definition: dof_map.h:1934
void add_adjoint_constraint_row(const unsigned int qoi_index, const dof_id_type dof_number, const DofConstraintRow &constraint_row, const Number constraint_rhs, const bool forbid_constraint_overwrite)
Adds a copy of the user-defined row to the constraint matrix, using an inhomogeneous right-hand-side ...
void process_mesh_constraint_rows(const MeshBase &mesh)
Adds any spline constraints from the Mesh to our DoF constraints.
void * _extra_sparsity_context
A pointer associated with the extra sparsity that can optionally be passed in.
Definition: dof_map.h:2006
void extract_local_vector(const NumericVector< Number > &Ug, const std::vector< dof_id_type > &dof_indices, DenseVectorBase< Number > &Ue) const
Builds the local element vector Ue from the global vector Ug, accounting for any constrained degrees ...
Definition: dof_map.C:2082
void constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix and vector.
Definition: dof_map.h:2330
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:2164
dof_id_type n_dofs() const
Definition: dof_map_base.h:105
The definition of the const_element_iterator struct.
Definition: mesh_base.h:2216
std::size_t distribute_dofs(MeshBase &)
Distribute dofs on the current mesh.
Definition: dof_map.C:978
void set_implicit_neighbor_dofs(bool implicit_neighbor_dofs)
Allow the implicit_neighbor_dofs flag to be set programmatically.
Definition: dof_map.C:1920
void add_default_ghosting()
Add the default functor(s) for coupling and algebraic ghosting.
Definition: dof_map.C:2024
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
void local_variable_indices(T &idx, const MeshBase &mesh, unsigned int var_num) const
If T == dof_id_type, counts, if T == std::vector<dof_id_type>, fills an array of, those dof indices w...
Definition: dof_map.C:1151
bool is_periodic_boundary(const boundary_id_type boundaryid) const
Definition: dof_map.C:215
Number has_heterogenous_adjoint_constraint(const unsigned int qoi_num, const dof_id_type dof) const
Backwards compatibility with misspelling.
Definition: dof_map.h:1125
dof_id_type n_local_constrained_dofs() const
void remove_dirichlet_boundary(const DirichletBoundary &dirichlet_boundary)
Removes the specified Dirichlet boundary from the system.
const FEType & variable_group_type(const unsigned int vg) const
Definition: dof_map.h:2230
std::unique_ptr< StaticCondensationDofMap > _sc
Static condensation class.
Definition: dof_map.h:2165
NodeConstraints::const_iterator node_constraint_rows_begin() const
Definition: dof_map.h:1084
std::vector< dof_id_type > _send_list
A list containing all the global DOF indices that affect the solution on my processor.
Definition: dof_map.h:1989
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2220
std::unique_ptr< SparsityPattern::Build > build_sparsity(const MeshBase &mesh, bool calculate_constrained=false, bool use_condensed_system=false) const
Builds a sparsity pattern for matrices using the current degree-of-freedom numbering and coupling...
Definition: dof_map.C:61
void set_verify_dirichlet_bc_consistency(bool val)
Set the _verify_dirichlet_bc_consistency flag.
Definition: dof_map.C:1926
unsigned int block_size() const
Definition: dof_map.h:658
void sum(T &r) const
bool is_attached(SparseMatrix< Number > &matrix)
Matrices should not be attached more than once.
Definition: dof_map.C:329
void attach_matrix(SparseMatrix< Number > &matrix)
Additional matrices may be attached to this DofMap.
Definition: dof_map.C:274
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
void clear_send_list()
Clears the _send_list vector.
Definition: dof_map.h:500
MeshBase & mesh
void gather_constraints(MeshBase &mesh, std::set< dof_id_type > &unexpanded_dofs, bool look_for_constrainees)
Helper function for querying about constraint equations on other processors.
DefaultCoupling & default_algebraic_ghosting()
Default algebraic ghosting functor.
Definition: dof_map.h:433
std::unique_ptr< DirichletBoundaries > _dirichlet_boundaries
Data structure containing Dirichlet functions.
Definition: dof_map.h:2134
This class allows one to associate Dirichlet boundary values with a given set of mesh boundary ids an...
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
Definition: vector_fe_ex5.C:44
const Parallel::Communicator & comm() const
The Node constraint storage format.
Definition: dof_map.h:156
std::set< GhostingFunctor * > _algebraic_ghosting_functors
The list of all GhostingFunctor objects to be used when distributing ghosted vectors.
Definition: dof_map.h:2047
unsigned int p_level() const
Definition: elem.h:3108
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
bool use_coupled_neighbor_dofs(const MeshBase &mesh) const
Tells other library functions whether or not this problem includes coupling between dofs in neighbori...
Definition: dof_map.C:1932
std::map< const Elem *, const CouplingMatrix *, CompareDofObjectsByPIDAndThenID > map_type
What elements do we care about and what variables do we care about on each element?
void reinit_static_condensation()
Calls reinit on the static condensation map if it exists.
Definition: dof_map.C:3080
void should_p_refine(unsigned int g, bool p_refine)
Describe whether the given variable group should be p-refined.
Definition: dof_map.h:2402
dof_id_type n_dofs(const unsigned int vn) const
Definition: dof_map.h:668
virtual void augment_send_list(std::vector< dof_id_type > &send_list)=0
User-defined function to augment the send list.
The libMesh namespace provides an interface to certain functionality in the library.
void add_constraint_row(const dof_id_type dof_number, const DofConstraintRow &constraint_row, const bool forbid_constraint_overwrite=true)
Adds a copy of the user-defined row to the constraint matrix, using a homogeneous right-hand-side for...
Definition: dof_map.h:1028
dof_id_type n_local_dofs(const unsigned int vn) const
Definition: dof_map.h:686
void set_error_on_constraint_loop(bool error_on_constraint_loop)
Definition: dof_map.C:241
void add_coupling_functor(GhostingFunctor &coupling_functor, bool to_mesh=true)
Adds a functor which can specify coupling requirements for creation of sparse matrices.
Definition: dof_map.C:2033
const SparsityPattern::Build * get_sparsity_pattern() const
Definition: dof_map.h:563
std::string get_info() const
Gets summary info about the sparsity bandwidth and constraints.
Definition: dof_map.C:2925
unsigned int sys_number() const
Definition: dof_map.h:2172
void SCALAR_dof_indices(std::vector< dof_id_type > &di, const unsigned int vn, const bool old_dofs=false) const
Fills the vector di with the global degree of freedom indices corresponding to the SCALAR variable vn...
Definition: dof_map.C:2547
uint8_t processor_id_type
Definition: id_types.h:104
This is the MeshBase class.
Definition: mesh_base.h:75
std::vector< dof_id_type > _first_scalar_df
First DOF index for SCALAR variable v, or garbage for non-SCALAR variable v.
Definition: dof_map.h:1983
const Variable & variable(const unsigned int c) const override
Definition: dof_map.h:2190
This class implements the default algebraic coupling in libMesh: elements couple to themselves...
void enforce_adjoint_constraints_exactly(NumericVector< Number > &v, unsigned int q) const
Heterogeneously constrains the numeric vector v, which represents an adjoint solution defined on the ...
Definition: dof_map.h:2354
AdjointDofConstraintValues _adjoint_constraint_values
Definition: dof_map.h:2110
bool _constrained_sparsity_construction
This flag indicates whether or not we explicitly take constraint equations into account when computin...
Definition: dof_map.h:1940
bool should_p_refine_var(unsigned int var) const
Whether the given variable should be p-refined.
Definition: dof_map.h:2437
unsigned int var_group_from_var_number(unsigned int var_num) const
Definition: dof_map.h:2430
dof_id_type end_dof() const
Definition: dof_map_base.h:83
AugmentSendList * _augment_send_list
Function object to call to add extra entries to the send list.
Definition: dof_map.h:2011
void heterogenously_constrain_element_vector(const DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true, int qoi_index=-1) const
Definition: dof_map.h:1304
bool need_full_sparsity_pattern
Default false; set to true if any attached matrix requires a full sparsity pattern.
Definition: dof_map.h:2072
bool constrained_sparsity_construction()
Returns true iff the current policy when constructing sparsity patterns is to explicitly account for ...
Definition: dof_map.h:2392
void distribute_local_dofs_var_major(dof_id_type &next_free_dof, MeshBase &mesh, const std::map< const Node *, std::set< subdomain_id_type >> &constraining_subdomains)
Distributes the global degrees of freedom, for dofs on this processor.
Definition: dof_map.C:1447
Generic sparse matrix.
Definition: vector_fe_ex5.C:46
void enforce_constraints_on_residual(const NonlinearImplicitSystem &system, NumericVector< Number > *rhs, NumericVector< Number > const *solution, bool homogeneous=true) const
Definition: dof_map.h:2359
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:179
virtual void clear() override
Free all new memory associated with the object, but restore its original state, with the mesh pointer...
Definition: dof_map.C:901
void _node_dof_indices(const Elem &elem, unsigned int n, const DofObject &obj, std::vector< dof_id_type > &di, const unsigned int vn) const
Helper function that implements the element-nodal versions of dof_indices and old_dof_indices.
Definition: dof_map.C:2431
void add_algebraic_ghosting_functor(std::shared_ptr< GhostingFunctor > evaluable_functor, bool to_mesh=true)
Adds a functor which can specify algebraic ghosting requirements for use with distributed vectors...
Definition: dof_map.h:406
void find_one_ring(const Tri3Subdivision *elem, std::vector< const Node *> &nodes)
Determines the 1-ring of element elem, and writes it to the nodes vector.
void reinit_send_list(MeshBase &mesh)
Clears the _send_list vector and then rebuilds it.
Definition: dof_map.C:1906
processor_id_type n_processors() const
void libmesh_ignore(const Args &...)
const dof_id_type n_nodes
Definition: tecplot_io.C:67
dof_id_type first_dof() const
Definition: dof_map_base.h:73
std::unordered_map< unsigned int, unsigned int > _var_to_vg
A map from variable number to variable group number.
Definition: dof_map.h:1960
This class defines the notion of a variable in the system.
Definition: variable.h:50
bool has_static_condensation() const
Checks whether we have static condensation.
Definition: dof_map.h:1668
static void merge_ghost_functor_outputs(GhostingFunctor::map_type &elements_to_ghost, CouplingMatricesSet &temporary_coupling_matrices, const std::set< GhostingFunctor *>::iterator &gf_begin, const std::set< GhostingFunctor *>::iterator &gf_end, const MeshBase::const_element_iterator &elems_begin, const MeshBase::const_element_iterator &elems_end, processor_id_type p)
Definition: dof_map.C:1614
void add_neighbors_to_send_list(MeshBase &mesh)
Adds entries to the _send_list vector corresponding to DoFs on elements neighboring the current proce...
Definition: dof_map.C:1716
std::vector< dof_id_type > _first_old_scalar_df
First old DOF index for SCALAR variable v, or garbage for non-SCALAR variable v.
Definition: dof_map.h:2093
int8_t boundary_id_type
Definition: id_types.h:51
bool has_adjoint_dirichlet_boundaries(unsigned int q) const
void(* _extra_sparsity_function)(SparsityPattern::Graph &, std::vector< dof_id_type > &n_nz, std::vector< dof_id_type > &n_oz, void *)
A function pointer to a function to call to add extra entries to the sparsity pattern.
Definition: dof_map.h:1999
dof_id_type _n_SCALAR_dofs
The total number of SCALAR dofs associated to all SCALAR variables.
Definition: dof_map.h:2085
void assert_no_nodes_missed(MeshBase &mesh)
Definition: dof_map.C:1592
PeriodicBoundaries * get_periodic_boundaries()
Definition: dof_map.h:1455
StaticCondensationDofMap & get_static_condensation()
Definition: dof_map.h:2689
void set_nonlocal_dof_objects(iterator_type objects_begin, iterator_type objects_end, MeshBase &mesh, dofobject_accessor objects)
Helper function for distributing dofs in parallel.
Definition: dof_map.C:352
void add_constraint_row(const dof_id_type dof_number, const DofConstraintRow &constraint_row, const Number constraint_rhs, const bool forbid_constraint_overwrite)
Adds a copy of the user-defined row to the constraint matrix, using an inhomogeneous right-hand-side ...
static bool extra_hanging_dofs(const FEType &fe_t)
virtual unsigned int n_nodes() const =0
std::set< std::unique_ptr< CouplingMatrix >, Utility::CompareUnderlying > CouplingMatricesSet
Definition: dof_map.h:1848
DofConstraints::const_iterator constraint_rows_begin() const
Definition: dof_map.h:1036
void heterogeneously_constrain_element_vector(const DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true, int qoi_index=-1) const
Constrains the element vector.
Definition: dof_map.h:2340
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:96
unsigned int n_variables() const override
Definition: dof_map.h:628
unsigned int n_systems() const
Definition: dof_object.h:937
void add_variable_group(VariableGroup var_group)
Add an unknown of order order and finite element type type to the system of equations.
Definition: dof_map.C:248
void swap_dof_constraints()
Similar to the stash/unstash_dof_constraints() API, but swaps _dof_constraints and _stashed_dof_const...
Definition: dof_map.h:1075
NodeConstraints::const_iterator node_constraint_rows_end() const
Definition: dof_map.h:1090
This base class provides a minimal set of interfaces for satisfying user requests for...
Definition: dof_map_base.h:51
dof_id_type n_constrained_nodes() const
Definition: dof_map.h:931
const Node *const * get_nodes() const
Definition: elem.h:2499
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
void print_dof_constraints(std::ostream &os=libMesh::out, bool print_nonlocal=false) const
Prints (from processor 0) all DoF and Node constraints.
libmesh_assert(ctx)
void attach_extra_sparsity_object(SparsityPattern::AugmentSparsityPattern &asp)
Attach an object to use to populate the sparsity pattern with extra entries.
Definition: dof_map.h:445
DirichletBoundaries * get_dirichlet_boundaries()
Definition: dof_map.h:1513
bool is_constrained_dof(const dof_id_type dof) const
Definition: dof_map.h:2258
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and non-linear solv...
The Tri3Subdivision element is a three-noded subdivision surface shell element used in mechanics calc...
static const dof_id_type invalid_id
An invalid id to distinguish an uninitialized DofObject.
Definition: dof_object.h:482
void allgather_recursive_constraints(MeshBase &)
Gathers constraint equation dependencies from other processors.
const VariableGroup & variable_group(const unsigned int c) const
Definition: dof_map.h:2180
DofConstraints _dof_constraints
Data structure containing DOF constraints.
Definition: dof_map.h:2106
void add_coupling_functor(std::shared_ptr< GhostingFunctor > coupling_functor, bool to_mesh=true)
Adds a functor which can specify coupling requirements for creation of sparse matrices.
Definition: dof_map.h:344
DofMap(const unsigned int sys_number, MeshBase &mesh)
Constructor.
Definition: dof_map.C:136
CouplingMatrix * _dof_coupling
Degree of freedom coupling.
Definition: dof_map.h:1613
bool active_on_subdomain(subdomain_id_type sid) const
Definition: variable.h:167
void create_static_condensation(MeshBase &mesh, System &system)
Add a static condensation class.
Definition: dof_map.C:3075
void print_info(std::ostream &os=libMesh::out) const
Prints summary info about the sparsity bandwidth and constraints.
Definition: dof_map.C:2918
void * _extra_send_list_context
A pointer associated with the extra send list that can optionally be passed in.
Definition: dof_map.h:2021
Number has_heterogeneous_adjoint_constraint(const unsigned int qoi_num, const dof_id_type dof) const
Definition: dof_map.h:2282
This class implements reference counting.
const DirichletBoundaries * get_dirichlet_boundaries() const
Definition: dof_map.h:1508
bool computed_sparsity_already() const
Returns true iff a sparsity pattern has already been computed.
Definition: dof_map.C:292
This class defines a logically grouped set of variables in the system.
Definition: variable.h:203
void update_sparsity_pattern(SparseMatrix< Number > &matrix) const
Additional matrices may be be temporarily initialized by this DofMap.
Definition: dof_map.C:303
std::pair< Real, Real > max_constraint_error(const System &system, NumericVector< Number > *v=nullptr) const
Tests the constrained degrees of freedom on the numeric vector v, which represents a solution defined...
static n_dofs_at_node_ptr n_dofs_at_node_function(const unsigned int dim, const FEType &fe_t)
Definition: fe_interface.C:458
Storage for DofConstraint right hand sides for a particular problem.
Definition: dof_map.h:120
std::vector< SparseMatrix< Number > *> _matrices
Additional matrices handled by this object.
Definition: dof_map.h:1977
DofObject * elem_ptr(MeshBase &mesh, dof_id_type i) const
Definition: dof_map.C:344
Storage for DofConstraint right hand sides for all adjoint problems.
Definition: dof_map.h:131
std::vector< std::unique_ptr< DirichletBoundaries > > _adjoint_dirichlet_boundaries
Data structure containing Dirichlet functions.
Definition: dof_map.h:2140
std::set< GhostingFunctor * >::const_iterator algebraic_ghosting_functors_end() const
End of range of algebraic ghosting functors.
Definition: dof_map.h:427
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
Definition: fe_interface.C:436
void distribute_scalar_dofs(dof_id_type &next_free_dof)
Definition: dof_map.C:1568
DofConstraints::const_iterator constraint_rows_end() const
Definition: dof_map.h:1042
std::unique_ptr< PeriodicBoundaries > _periodic_boundaries
Data structure containing periodic boundaries.
Definition: dof_map.h:2126
void heterogenously_constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true, int qoi_index=-1) const
Definition: dof_map.h:1260
bool has_heterogenous_adjoint_constraints(const unsigned int qoi_num) const
Backwards compatibility with misspelling.
Definition: dof_map.h:1109
void unstash_dof_constraints()
Definition: dof_map.h:1057
void remove_coupling_functor(GhostingFunctor &coupling_functor)
Removes a functor which was previously added to the set of coupling functors, from both this DofMap a...
Definition: dof_map.C:2045
unsigned int(* n_dofs_at_node_ptr)(const ElemType, const Order, const unsigned int)
Definition: fe_interface.h:151
DofObject *(DofMap::* dofobject_accessor)(MeshBase &mesh, dof_id_type i) const
A member function type like node_ptr() or elem_ptr().
Definition: dof_map.h:1765
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Struct which defines a custom comparison object that can be used with std::sets of std::unique_ptrs...
Definition: utility.h:503
void process_constraints(MeshBase &)
Postprocesses any constrained degrees of freedom to be constrained only in terms of unconstrained dof...
subdomain_id_type subdomain_id() const
Definition: elem.h:2582
void stash_dof_constraints()
Definition: dof_map.h:1051
dof_id_type n_constrained_dofs() const
void constrain_element_vector(DenseVector< Number > &rhs, std::vector< dof_id_type > &dofs, bool asymmetric_constraint_rows=true) const
Constrains the element vector.
Definition: dof_map.h:2326
const DofConstraints & get_dof_constraints() const
Provide a const accessor to the DofConstraints map.
Definition: dof_map.h:1049
const std::vector< dof_id_type > & get_n_nz() const
Definition: dof_map.h:535
void _dof_indices(const Elem &elem, int p_level, std::vector< dof_id_type > &di, const unsigned int vg, const unsigned int vig, const Node *const *nodes, unsigned int n_nodes, const unsigned int v #ifdef DEBUG, std::size_t &tot_size #endif) const
Helper function that gets the dof indices on the current element for a non-SCALAR type variable...
Definition: dof_map.C:2515
void invalidate_dofs(MeshBase &mesh) const
Invalidates all active DofObject dofs for this system.
Definition: dof_map.C:886
const DirichletBoundaries * get_adjoint_dirichlet_boundaries(unsigned int q) const
void attach_extra_send_list_function(void(*func)(std::vector< dof_id_type > &, void *), void *context=nullptr)
Attach a function pointer to use as a callback to populate the send_list with extra entries...
Definition: dof_map.h:483
virtual bool is_vertex(const unsigned int i) const =0
bool _verify_dirichlet_bc_consistency
Flag which determines whether we should do some additional checking of the consistency of the Dirichl...
Definition: dof_map.h:2162
std::unique_ptr< DefaultCoupling > _default_coupling
The default coupling GhostingFunctor, used to implement standard libMesh sparsity pattern constructio...
Definition: dof_map.h:2029
void add_constraints_to_send_list()
Adds entries to the _send_list vector corresponding to DoFs which are dependencies for constraint equ...
OStreamProxy out
DofConstraintValueMap _primal_constraint_values
Definition: dof_map.h:2108
void heterogeneously_constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true, int qoi_index=-1) const
Constrains the element matrix and vector.
Definition: dof_map.h:2336
std::vector< Variable > _variables
The finite element type for each variable.
Definition: dof_map.h:1945
DofConstraints _stashed_dof_constraints
Definition: dof_map.h:2106
void remove_default_ghosting()
Remove any default ghosting functor(s).
Definition: dof_map.C:2016
void attach_extra_sparsity_function(void(*func)(SparsityPattern::Graph &sparsity, std::vector< dof_id_type > &n_nz, std::vector< dof_id_type > &n_oz, void *), void *context=nullptr)
Attach a function pointer to use as a callback to populate the sparsity pattern with extra entries...
Definition: dof_map.h:460
SparsityPattern::AugmentSparsityPattern * _augment_sparsity_pattern
Function object to call to add extra entries to the sparsity pattern.
Definition: dof_map.h:1994
void find_connected_dof_objects(std::vector< const DofObject *> &objs) const
Finds all the DofObjects associated with the set in objs.
std::map< dof_id_type, Real, std::less< dof_id_type >, Threads::scalable_allocator< std::pair< const dof_id_type, Real > > > DofConstraintRow
A row of the Dof constraint matrix.
Definition: dof_map.h:100
The DofObject defines an abstract base class for objects that have degrees of freedom associated with...
Definition: dof_object.h:54
void clear_sparsity()
Clears the sparsity pattern.
Definition: dof_map.C:2009
void set_constrained_sparsity_construction(bool use_constraints)
Sets the current policy for constructing sparsity patterns: if use_constraints is true (for robustnes...
Definition: dof_map.h:2373
Order variable_order(const unsigned int c) const
Definition: dof_map.h:2200
dof_id_type n_local_dofs() const
Definition: dof_map_base.h:115
std::set< GhostingFunctor * > _coupling_functors
The list of all GhostingFunctor objects to be used when coupling degrees of freedom in matrix sparsit...
Definition: dof_map.h:2060
std::unordered_set< unsigned int > _dont_p_refine
A container of variable groups that we should not p-refine.
Definition: dof_map.h:2098
void add_dirichlet_boundary(const DirichletBoundary &dirichlet_boundary)
Adds a copy of the specified Dirichlet boundary to the system.
unsigned int number(unsigned int v) const
Definition: variable.h:304
void build_constraint_matrix(DenseMatrix< Number > &C, std::vector< dof_id_type > &elem_dofs, const bool called_recursively=false) const
Build the constraint matrix C associated with the element degree of freedom indices elem_dofs...
std::pair< unsigned int, unsigned int > var_to_vg_and_offset(const unsigned int s, const unsigned int var) const
Definition: dof_object.h:1200
void(* _extra_send_list_function)(std::vector< dof_id_type > &, void *)
A function pointer to a function to call to add extra entries to the send list.
Definition: dof_map.h:2016
Scalable allocator to be used in multithreaded code chunks which allocate a lot of dynamic memory...
DofObject * node_ptr(MeshBase &mesh, dof_id_type i) const
Definition: dof_map.C:337
void heterogeneously_constrain_element_residual(DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, NumericVector< Number > &solution_local) const
Constrains the element residual.
The base class for defining periodic boundaries.
void create_dof_constraints(const MeshBase &, Real time=0)
Rebuilds the raw degree of freedom and DofObject constraints, based on attached DirichletBoundary obj...
std::set< GhostingFunctor * >::const_iterator coupling_functors_end() const
End of range of coupling functors.
Definition: dof_map.h:365
Defines an abstract dense vector base class for use in Finite Element-type computations.
Definition: dof_map.h:73
std::map< const Node *, Real, std::less< const Node * >, Threads::scalable_allocator< std::pair< const Node *const, Real > > > NodeConstraintRow
A row of the Node constraint mapping.
Definition: dof_map.h:148
MeshBase & _mesh
The mesh that system uses.
Definition: dof_map.h:1970
void prepare_send_list()
Takes the _send_list vector (which may have duplicate entries) and sorts it.
Definition: dof_map.C:1864
bool all_semilocal_indices(const std::vector< dof_id_type > &dof_indices) const
Definition: dof_map.C:2599
Defines a dense vector for use in Finite Element-type computations.
Definition: dof_map.h:74
virtual bool infinite() const =0
std::unique_ptr< DefaultCoupling > _default_evaluating
The default algebraic GhostingFunctor, used to implement standard libMesh send_list construction...
Definition: dof_map.h:2037
void full_sparsity_pattern_needed()
Sets need_full_sparsity_pattern to true regardless of the requirements by matrices.
Definition: dof_map.h:2386
std::vector< unsigned int > _variable_group_numbers
The variable group number for each variable.
Definition: dof_map.h:1955
std::map< GhostingFunctor *, std::shared_ptr< GhostingFunctor > > _shared_functors
Hang on to references to any GhostingFunctor objects we were passed in shared_ptr form...
Definition: dof_map.h:2066
void heterogeneously_constrain_element_jacobian_and_residual(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, NumericVector< Number > &solution_local) const
Constrains the element Jacobian and residual.
FEFamily
defines an enum for finite element families.
void remove_algebraic_ghosting_functor(GhostingFunctor &evaluable_functor)
Removes a functor which was previously added to the set of algebraic ghosting functors, from both this DofMap and from the underlying mesh.
Definition: dof_map.C:2070
virtual ~AugmentSendList()=default
Backwards compatibility for prior AugmentSparsityPattern users.
Definition: dof_map.h:201
void attach_extra_send_list_object(DofMap::AugmentSendList &asl)
Attach an object to populate the send_list with extra entries.
Definition: dof_map.h:474
Order variable_group_order(const unsigned int vg) const
Definition: dof_map.h:2210
std::set< GhostingFunctor * >::const_iterator algebraic_ghosting_functors_begin() const
Beginning of range of algebraic ghosting functors.
Definition: dof_map.h:421
void compute_sparsity(const MeshBase &)
Computes the sparsity pattern for the matrices corresponding to proc_id and sends that data to Linear...
Definition: dof_map.C:1988
void add_algebraic_ghosting_functor(GhostingFunctor &evaluable_functor, bool to_mesh=true)
Adds a functor which can specify algebraic ghosting requirements for use with distributed vectors...
Definition: dof_map.C:2058
unsigned int n_comp_group(const unsigned int s, const unsigned int vg) const
Definition: dof_object.h:1015
bool has_heterogeneous_adjoint_constraints(const unsigned int qoi_num) const
Definition: dof_map.h:2268
void constrain_element_residual(DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, NumericVector< Number > &solution_local) const
Constrains the element residual.
const PeriodicBoundaries * get_periodic_boundaries() const
Definition: dof_map.h:1460
std::map< const Node *, std::set< subdomain_id_type > > calculate_constraining_subdomains()
We may have mesh constraint rows with dependent nodes in one subdomain but dependency nodes in anothe...
Definition: dof_map.C:1285
Defines a dense matrix for use in Finite Element-type computations.
Definition: dof_map.h:75
bool active() const
Definition: elem.h:2941
void set_error_on_cyclic_constraint(bool error_on_cyclic_constraint)
Specify whether or not we perform an extra (opt-mode enabled) check for constraint loops...
Definition: dof_map.C:234
virtual ElemType type() const =0
Abstract base class to be used to add user-defined implicit degree of freedom couplings.
bool has_blocked_representation() const
Definition: dof_map.h:641
The constraint matrix storage format.
Definition: dof_map.h:108
void ErrorVector unsigned int
Definition: adjoints_ex3.C:360
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:526
std::vector< dof_id_type > n_dofs_per_processor(const unsigned int vn) const
Definition: dof_map.h:697
Abstract base class to be used to add user-defined parallel degree of freedom couplings.
Definition: dof_map.h:208
void constrain_nothing(std::vector< dof_id_type > &dofs) const
Does not actually constrain anything, but modifies dofs in the same way as any of the constrain funct...
Definition: dof_map.h:2348
bool is_evaluable(const DofObjectSubclass &obj, unsigned int var_num=libMesh::invalid_uint) const
Definition: dof_map.C:2612
bool semilocal_index(dof_id_type dof_index) const
Definition: dof_map.C:2583
void distribute_local_dofs_node_major(dof_id_type &next_free_dof, MeshBase &mesh, const std::map< const Node *, std::set< subdomain_id_type >> &constraining_subdomains)
Distributes the global degrees of freedom for dofs on this processor.
Definition: dof_map.C:1319
void constrain_p_dofs(unsigned int var, const Elem *elem, unsigned int s, unsigned int p)
Constrains degrees of freedom on side s of element elem which correspond to variable number var and t...
void check_for_cyclic_constraints()
Throw an error if we detect any constraint loops, i.e.
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)
A useful inline function which replaces the macros used previously.
void old_dof_indices(const Elem &elem, unsigned int n, std::vector< dof_id_type > &di, const unsigned int vn) const
Appends to the vector di the old global degree of freedom indices for elem.node_ref(n), for one variable vn.
Definition: dof_map.C:2418
uint8_t dof_id_type
Definition: id_types.h:67
std::vector< dof_id_type > _end_df
Last DOF index (plus 1) on processor p.
Definition: dof_map_base.h:159
bool local_index(dof_id_type dof_index) const
Definition: dof_map.h:839
const FEType & type() const
Definition: variable.h:144
NodeConstraints _node_constraints
Data structure containing DofObject constraints.
Definition: dof_map.h:2117
void enforce_constraints_exactly(const System &system, NumericVector< Number > *v=nullptr, bool homogeneous=false) const
Constrains the numeric vector v, which represents a solution defined on the mesh. ...
Definition: dof_map.h:2350
This class defines a coupling matrix.
void constrain_element_matrix(DenseMatrix< Number > &matrix, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix.
Definition: dof_map.h:2317
void enforce_constraints_on_jacobian(const NonlinearImplicitSystem &system, SparseMatrix< Number > *jac) const
Definition: dof_map.h:2365