dof_map.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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/auto_ptr.h"
26 #include "libmesh/enum_order.h"
28 #include "libmesh/libmesh.h" // libMesh::invalid_uint
29 #include "libmesh/variable.h"
30 #include "libmesh/threads.h"
32 #include "libmesh/elem_range.h"
35 
36 // C++ Includes -----------------------------------
37 #include <algorithm>
38 #include <cstddef>
39 #include <iterator>
40 #include <map>
41 #include <string>
42 #include <vector>
43 
44 namespace libMesh
45 {
46 
47 // Forward Declarations
48 class CouplingMatrix;
49 class DirichletBoundary;
50 class DirichletBoundaries;
51 class DofMap;
52 class DofObject;
53 class Elem;
54 class FEType;
55 class MeshBase;
56 class Mesh;
57 class PeriodicBoundaryBase;
58 class PeriodicBoundaries;
59 namespace SparsityPattern { class Build; }
60 class System;
61 template <typename T> class DenseVectorBase;
62 template <typename T> class DenseVector;
63 template <typename T> class DenseMatrix;
64 template <typename T> class SparseMatrix;
65 template <typename T> class NumericVector;
66 
67 
68 
69 // ------------------------------------------------------------
70 // Do we need constraints for anything?
71 
72 #if defined(LIBMESH_ENABLE_AMR) || \
73  defined(LIBMESH_ENABLE_PERIODIC) || \
74  defined(LIBMESH_ENABLE_DIRICHLET)
75 # define LIBMESH_ENABLE_CONSTRAINTS 1
76 #endif
77 
78 // ------------------------------------------------------------
79 // AMR constraint matrix types
80 
81 #ifdef LIBMESH_ENABLE_CONSTRAINTS
82 
85 typedef std::map<dof_id_type, Real,
86  std::less<dof_id_type>,
88 
95 class DofConstraints : public std::map<dof_id_type,
96  DofConstraintRow,
97  std::less<dof_id_type>,
98  Threads::scalable_allocator<std::pair<const dof_id_type, Number> > >
99 {
100 };
101 
108  public std::map<dof_id_type, Number,
109  std::less<dof_id_type>,
110  Threads::scalable_allocator<std::pair<const dof_id_type, Number> > >
111 {
112 };
113 
119  public std::map<unsigned int, DofConstraintValueMap,
120  std::less<unsigned int>,
121  Threads::scalable_allocator
122  <std::pair<const unsigned int, DofConstraintValueMap> > >
123 {
124 };
125 
126 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
127 
133 typedef std::map<const Node *, Real,
134  std::less<const Node *>,
136 
143 class NodeConstraints : public std::map<const Node *,
144  std::pair<NodeConstraintRow,Point>,
145  std::less<const Node *>,
146  Threads::scalable_allocator<std::pair<const Node * const, std::pair<NodeConstraintRow,Point> > > >
147 {
148 };
149 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
150 
151 #endif // LIBMESH_ENABLE_CONSTRAINTS
152 
153 
154 
155 // ------------------------------------------------------------
156 // DofMap class definition
157 
167  class DofMap : public ReferenceCountedObject<DofMap>,
168  public ParallelObject
169 {
170 public:
171 
177  explicit
178  DofMap(const unsigned int sys_number,
179  const ParallelObject &parent_decomp);
180 
184  ~DofMap();
185 
191  {
192  public:
194 
198  virtual void augment_sparsity_pattern (SparsityPattern::Graph & sparsity,
199  std::vector<dof_id_type> & n_nz,
200  std::vector<dof_id_type> & n_oz) = 0;
201  };
202 
208  {
209  public:
210  virtual ~AugmentSendList () {}
211 
215  virtual void augment_send_list (std::vector<dof_id_type> & send_list) = 0;
216  };
217 
223  void attach_matrix (SparseMatrix<Number>& matrix);
224 
229  bool is_attached (SparseMatrix<Number>& matrix);
230 
236  void distribute_dofs (MeshBase&);
237 
243  void compute_sparsity (const MeshBase&);
244 
248  void clear_sparsity();
249 
261  {
263  }
264 
276  std::vector<dof_id_type> & n_nz,
277  std::vector<dof_id_type> & n_oz,
278  void *),
279  void * context = NULL)
281 
290  {
291  _augment_send_list = &asl;
292  }
293 
298  void attach_extra_send_list_function(void (*func)(std::vector<dof_id_type> &, void *), void * context = NULL)
300 
307  void prepare_send_list ();
308 
316  const std::vector<dof_id_type>& get_send_list() const { return _send_list; }
317 
324  const std::vector<dof_id_type>& get_n_nz() const {
326  return *_n_nz;
327  }
328 
335  const std::vector<dof_id_type>& get_n_oz() const {
337  return *_n_oz;
338  }
339 
340  // /**
341  // * Add an unknown of order \p order and finite element type
342  // * \p type to the system of equations.
343  // */
344  // void add_variable (const Variable &var);
345 
350  void add_variable_group (const VariableGroup &var_group);
351 
355  const VariableGroup& variable_group (const unsigned int c) const;
356 
360  const Variable& variable (const unsigned int c) const;
361 
365  Order variable_order (const unsigned int c) const;
366 
370  Order variable_group_order (const unsigned int vg) const;
371 
375  const FEType& variable_type (const unsigned int c) const;
376 
380  const FEType& variable_group_type (const unsigned int vg) const;
381 
387  unsigned int n_variable_groups() const
388  { return libmesh_cast_int<unsigned int>(_variable_groups.size()); }
389 
395  unsigned int n_variables() const
396  { return libmesh_cast_int<unsigned int>(_variables.size()); }
397 
404  {
405 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
406  return ((this->n_variable_groups() == 1) && (this->n_variables() > 1));
407 #else
408  return false;
409 #endif
410  }
411 
416  unsigned int block_size() const
417  {
418 #ifdef LIBMESH_ENABLE_BLOCKED_STORAGE
419  return (this->has_blocked_representation() ? this->n_variables() : 1);
420 #else
421  return 1;
422 #endif
423  }
424 
428  dof_id_type n_dofs() const { return _n_dfs; }
429 
434 
439  { return this->n_dofs_on_processor (this->processor_id()); }
440 
445  libmesh_assert_less (proc, _first_df.size());
446  return libmesh_cast_int<dof_id_type>(_end_df[proc] - _first_df[proc]);
447  }
448 
453  { libmesh_assert_less (proc, _first_df.size()); return _first_df[proc]; }
454 
456  { return this->first_dof(this->processor_id()); }
457 
458 #ifdef LIBMESH_ENABLE_AMR
459 
463  { libmesh_assert_less (proc, _first_old_df.size()); return _first_old_df[proc]; }
464 
466  { return this->first_old_dof(this->processor_id()); }
467 
468 #endif //LIBMESH_ENABLE_AMR
469 
476  libmesh_deprecated();
477  libmesh_assert_less (proc, _end_df.size());
478  return libmesh_cast_int<dof_id_type>(_end_df[proc] - 1);
479  }
480 
482  { return this->last_dof(this->processor_id()); }
483 
489  { libmesh_assert_less (proc, _end_df.size()); return _end_df[proc]; }
490 
492  { return this->end_dof(this->processor_id()); }
493 
494 #ifdef LIBMESH_ENABLE_AMR
495 
500  { libmesh_assert_less (proc, _end_old_df.size()); return _end_old_df[proc]; }
501 
503  { return this->end_old_dof(this->processor_id()); }
504 
505 #endif //LIBMESH_ENABLE_AMR
506 
512  void dof_indices (const Elem* const elem,
513  std::vector<dof_id_type>& di,
514  const unsigned int vn = libMesh::invalid_uint) const;
515 
522  void SCALAR_dof_indices (std::vector<dof_id_type>& di,
523  const unsigned int vn,
524  const bool old_dofs=false) const;
525 
532  bool all_semilocal_indices (const std::vector<dof_id_type>& dof_indices) const;
533 
540  bool use_coupled_neighbor_dofs(const MeshBase& mesh) const;
541 
553  const std::vector<dof_id_type>& dof_indices,
554  DenseVectorBase<Number>& Ue) const;
555 
556 
557 #ifdef LIBMESH_ENABLE_CONSTRAINTS
558 
559  //--------------------------------------------------------------------
560  // Constraint-specific methods
566 
572 
573 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
574 
579  { return libmesh_cast_int<dof_id_type>(_node_constraints.size()); }
580 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
581 
587  void create_dof_constraints (const MeshBase&, Real time=0);
588 
593 
598 
607 
612  void add_constraint_row (const dof_id_type dof_number,
613  const DofConstraintRow& constraint_row,
614  const Number constraint_rhs,
615  const bool forbid_constraint_overwrite);
616 
627  void add_adjoint_constraint_row (const unsigned int qoi_index,
628  const dof_id_type dof_number,
629  const DofConstraintRow& constraint_row,
630  const Number constraint_rhs,
631  const bool forbid_constraint_overwrite);
632 
638  void add_constraint_row (const dof_id_type dof_number,
639  const DofConstraintRow& constraint_row,
640  const bool forbid_constraint_overwrite = true)
641  { add_constraint_row(dof_number, constraint_row, 0., forbid_constraint_overwrite); }
642 
646  DofConstraints::const_iterator constraint_rows_begin() const
647  { return _dof_constraints.begin(); }
648 
652  DofConstraints::const_iterator constraint_rows_end() const
653  { return _dof_constraints.end(); }
654 
655 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
656 
659  NodeConstraints::const_iterator node_constraint_rows_begin() const
660  { return _node_constraints.begin(); }
661 
665  NodeConstraints::const_iterator node_constraint_rows_end() const
666  { return _node_constraints.end(); }
667 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
668 
673  bool is_constrained_dof (const dof_id_type dof) const;
674 
679  bool has_heterogenous_adjoint_constraint (const unsigned int qoi_num,
680  const dof_id_type dof) const;
681 
686  bool is_constrained_node (const Node* node) const;
687 
694  void print_dof_constraints(std::ostream& os=libMesh::out,
695  bool print_nonlocal=false) const;
696 
702  std::string get_local_constraints(bool print_nonlocal=false) const;
703 
704 
713  std::pair<Real, Real> max_constraint_error(const System &system,
714  NumericVector<Number> *v = NULL) const;
715 
716 #endif // LIBMESH_ENABLE_CONSTRAINTS
717 
718  //--------------------------------------------------------------------
719  // Constraint-specific methods
720  // Some of these methods are enabled (but inlined away to nothing)
721  // when constraints are disabled at configure-time. This is to
722  // increase API compatibility of user code with different library
723  // builds.
724 
739  std::vector<dof_id_type>& elem_dofs,
740  bool asymmetric_constraint_rows = true) const;
741 
749  std::vector<dof_id_type>& row_dofs,
750  std::vector<dof_id_type>& col_dofs,
751  bool asymmetric_constraint_rows = true) const;
752 
757  std::vector<dof_id_type>& dofs,
758  bool asymmetric_constraint_rows = true) const;
759 
769  DenseVector<Number>& rhs,
770  std::vector<dof_id_type>& elem_dofs,
771  bool asymmetric_constraint_rows = true) const;
772 
797  DenseVector<Number>& rhs,
798  std::vector<dof_id_type>& elem_dofs,
799  bool asymmetric_constraint_rows = true,
800  int qoi_index = -1) const;
801 
824  DenseVector<Number>& rhs,
825  std::vector<dof_id_type>& elem_dofs,
826  bool asymmetric_constraint_rows = true,
827  int qoi_index = -1) const;
828 
829 
830 
841  std::vector<dof_id_type>& row_dofs,
842  bool asymmetric_constraint_rows = true) const;
843 
850  void constrain_nothing (std::vector<dof_id_type>& dofs) const;
851 
865  void enforce_constraints_exactly (const System &system,
866  NumericVector<Number> *v = NULL,
867  bool homogeneous = false) const;
868 
876  unsigned int q) const;
877 
878 
879 
880 #ifdef LIBMESH_ENABLE_PERIODIC
881 
882  //--------------------------------------------------------------------
883  // PeriodicBoundary-specific methods
884 
888  void add_periodic_boundary (const PeriodicBoundaryBase& periodic_boundary);
889 
896  void add_periodic_boundary (const PeriodicBoundaryBase& boundary, const PeriodicBoundaryBase& inverse_boundary);
897 
902  bool is_periodic_boundary (const boundary_id_type boundaryid) const;
903 
905  {
906  return _periodic_boundaries;
907  }
908 
909 #endif // LIBMESH_ENABLE_PERIODIC
910 
911 
912 #ifdef LIBMESH_ENABLE_DIRICHLET
913 
914  //--------------------------------------------------------------------
915  // DirichletBoundary-specific methods
916 
920  void add_dirichlet_boundary (const DirichletBoundary& dirichlet_boundary);
921 
927  void add_adjoint_dirichlet_boundary (const DirichletBoundary& dirichlet_boundary,
928  unsigned int q);
929 
933  void remove_dirichlet_boundary (const DirichletBoundary& dirichlet_boundary);
934 
939  void remove_adjoint_dirichlet_boundary (const DirichletBoundary& dirichlet_boundary,
940  unsigned int q);
941 
943  {
944  return _dirichlet_boundaries;
945  }
946 
948  {
949  return _dirichlet_boundaries;
950  }
951 
952  bool has_adjoint_dirichlet_boundaries(unsigned int q) const;
953 
954  const DirichletBoundaries *
955  get_adjoint_dirichlet_boundaries(unsigned int q) const;
956 
958  get_adjoint_dirichlet_boundaries(unsigned int q);
959 
960 #endif // LIBMESH_ENABLE_DIRICHLET
961 
962 
963 #ifdef LIBMESH_ENABLE_AMR
964 
965  //--------------------------------------------------------------------
966  // AMR-specific methods
967 
976  // void augment_send_list_for_projection(const MeshBase &);
977 
984  void old_dof_indices (const Elem* const elem,
985  std::vector<dof_id_type>& di,
986  const unsigned int vn = libMesh::invalid_uint) const;
991  dof_id_type n_old_dofs() const { return _n_old_dfs; }
992 
998  void constrain_p_dofs (unsigned int var,
999  const Elem *elem,
1000  unsigned int s,
1001  unsigned int p);
1002 
1003 #endif // LIBMESH_ENABLE_AMR
1004 
1008  void reinit (MeshBase& mesh);
1009 
1013  void clear ();
1014 
1018  void print_info(std::ostream& os=libMesh::out) const;
1019 
1023  std::string get_info() const;
1024 
1039 
1043  unsigned int sys_number() const;
1044 
1045 private:
1046 
1051 
1055  void invalidate_dofs(MeshBase& mesh) const;
1056 
1061 
1066 
1070  typedef DofObject* (DofMap::*dofobject_accessor)
1072 
1076  template<typename iterator_type>
1077  void set_nonlocal_dof_objects(iterator_type objects_begin,
1078  iterator_type objects_end,
1079  MeshBase &mesh,
1080  dofobject_accessor objects);
1081 
1090  void distribute_local_dofs_var_major (dof_id_type& next_free_dof,
1091  MeshBase& mesh);
1092 
1105  void distribute_local_dofs_node_major (dof_id_type& next_free_dof,
1106  MeshBase& mesh);
1107 
1113 
1114 #ifdef LIBMESH_ENABLE_CONSTRAINTS
1115 
1127  std::vector<dof_id_type>& elem_dofs,
1128  const bool called_recursively=false) const;
1129 
1147  std::vector<dof_id_type>& elem_dofs,
1148  int qoi_index = -1,
1149  const bool called_recursively=false) const;
1150 
1155  void find_connected_dofs (std::vector<dof_id_type> &elem_dofs) const;
1156 
1161  void find_connected_dof_objects (std::vector<const DofObject *> &objs) const;
1162 
1169 
1170 #endif // LIBMESH_ENABLE_CONSTRAINTS
1171 
1175  std::vector<Variable> _variables;
1176 
1180  std::vector<VariableGroup> _variable_groups;
1181 
1185  const unsigned int _sys_number;
1186 
1192  std::vector<SparseMatrix<Number>* > _matrices;
1193 
1197  std::vector<dof_id_type> _first_df;
1198 
1202  std::vector<dof_id_type> _end_df;
1203 
1208  std::vector<dof_id_type> _send_list;
1209 
1214 
1219  std::vector<dof_id_type> & n_nz,
1220  std::vector<dof_id_type> & n_oz,
1221  void *);
1226 
1231 
1235  void (*_extra_send_list_function)(std::vector<dof_id_type> &, void *);
1236 
1241 
1247 
1253 
1260  std::vector<dof_id_type>* _n_nz;
1261 
1266  std::vector<dof_id_type>* _n_oz;
1267 
1272 
1278 
1279 #ifdef LIBMESH_ENABLE_AMR
1280 
1285 
1289  std::vector<dof_id_type> _first_old_df;
1290 
1294  std::vector<dof_id_type> _end_old_df;
1295 
1296 #endif
1297 
1298 #ifdef LIBMESH_ENABLE_CONSTRAINTS
1299 
1304 
1306 
1308 #endif
1309 
1310 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
1311 
1315 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
1316 
1317 
1318 #ifdef LIBMESH_ENABLE_PERIODIC
1319 
1324 #endif
1325 
1326 #ifdef LIBMESH_ENABLE_DIRICHLET
1327 
1332 
1337  std::vector<DirichletBoundaries *> _adjoint_dirichlet_boundaries;
1338 #endif
1339 
1341 };
1342 
1343 
1344 // ------------------------------------------------------------
1345 // Dof Map inline member functions
1346 inline
1347 unsigned int DofMap::sys_number() const
1348 {
1349  return _sys_number;
1350 }
1351 
1352 
1353 
1354 inline
1355 const VariableGroup & DofMap::variable_group (const unsigned int g) const
1356 {
1357  libmesh_assert_less (g, _variable_groups.size());
1358 
1359  return _variable_groups[g];
1360 }
1361 
1362 
1363 
1364 inline
1365 const Variable & DofMap::variable (const unsigned int c) const
1366 {
1367  libmesh_assert_less (c, _variables.size());
1368 
1369  return _variables[c];
1370 }
1371 
1372 
1373 
1374 inline
1375 Order DofMap::variable_order (const unsigned int c) const
1376 {
1377  libmesh_assert_less (c, _variables.size());
1378 
1379  return _variables[c].type().order;
1380 }
1381 
1382 
1383 
1384 inline
1385 Order DofMap::variable_group_order (const unsigned int vg) const
1386 {
1387  libmesh_assert_less (vg, _variable_groups.size());
1388 
1389  return _variable_groups[vg].type().order;
1390 }
1391 
1392 
1393 
1394 inline
1395 const FEType& DofMap::variable_type (const unsigned int c) const
1396 {
1397  libmesh_assert_less (c, _variables.size());
1398 
1399  return _variables[c].type();
1400 }
1401 
1402 
1403 
1404 inline
1405 const FEType& DofMap::variable_group_type (const unsigned int vg) const
1406 {
1407  libmesh_assert_less (vg, _variable_groups.size());
1408 
1409  return _variable_groups[vg].type();
1410 }
1411 
1412 
1413 
1414 inline
1416 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
1417 node
1418 #endif
1419 ) const
1420 {
1421 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
1422  if (_node_constraints.count(node))
1423  return true;
1424 #endif
1425 
1426  return false;
1427 }
1428 
1429 #ifdef LIBMESH_ENABLE_CONSTRAINTS
1430 
1431 inline
1433 {
1434  if (_dof_constraints.count(dof))
1435  return true;
1436 
1437  return false;
1438 }
1439 
1440 inline
1441 bool DofMap::has_heterogenous_adjoint_constraint (const unsigned int qoi_num,
1442  const dof_id_type dof) const
1443 {
1444  AdjointDofConstraintValues::const_iterator it =
1445  _adjoint_constraint_values.find(qoi_num);
1446  if (it != _adjoint_constraint_values.end() &&
1447  it->second.count(dof))
1448  return true;
1449 
1450  return false;
1451 }
1452 
1453 
1454 #else
1455 
1456  //--------------------------------------------------------------------
1457  // Constraint-specific methods get inlined into nothing if
1458  // constraints are disabled, so there's no reason for users not to
1459  // use them.
1460 
1462  std::vector<dof_id_type>&,
1463  bool) const {}
1464 
1465 inline void DofMap::constrain_element_matrix (DenseMatrix<Number>&,
1466  std::vector<dof_id_type>&,
1467  std::vector<dof_id_type>&,
1468  bool) const {}
1469 
1470 inline void DofMap::constrain_element_vector (DenseVector<Number>&,
1471  std::vector<dof_id_type>&,
1472  bool) const {}
1473 
1474 inline void DofMap::constrain_element_matrix_and_vector (DenseMatrix<Number>&,
1475  DenseVector<Number>&,
1476  std::vector<dof_id_type>&,
1477  bool) const {}
1478 
1479 inline void DofMap::constrain_element_dyad_matrix (DenseVector<Number>&,
1480  DenseVector<Number>&,
1481  std::vector<dof_id_type>&,
1482  bool) const {}
1483 
1484 inline void DofMap::enforce_constraints_exactly (const System &,
1485  NumericVector<Number> *,
1486  bool = false) const {}
1487 
1488 inline void DofMap::enforce_adjoint_constraints_exactly (NumericVector<Number> *,
1489  unsigned int) const {}
1490 
1491 #endif // LIBMESH_ENABLE_CONSTRAINTS
1492 
1493 } // namespace libMesh
1494 
1495 #endif // LIBMESH_DOF_MAP_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo