libMesh
sparsity_pattern.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 #ifndef LIBMESH_SPARSITY_PATTERN_H
20 #define LIBMESH_SPARSITY_PATTERN_H
21 
22 // Local Includes
23 #include "libmesh/elem_range.h"
24 #include "libmesh/threads_allocators.h"
25 #include "libmesh/parallel_object.h"
26 
27 // C++ includes
28 #include <vector>
29 
30 namespace libMesh
31 {
32 
33 // Forward declarations
34 class MeshBase;
35 class DofMap;
36 class CouplingMatrix;
37 
48 namespace SparsityPattern // use a namespace so member classes can be forward-declared.
49 {
50 typedef std::vector<dof_id_type, Threads::scalable_allocator<dof_id_type>> Row;
51 class Graph : public std::vector<Row> {};
52 
53 class NonlocalGraph : public std::map<dof_id_type, Row> {};
54 
64 template<typename BidirectionalIterator>
65 static void sort_row (const BidirectionalIterator begin,
66  BidirectionalIterator middle,
67  const BidirectionalIterator end);
68 
80 class Build : public ParallelObject
81 {
82 private:
83  const MeshBase & mesh;
84  const DofMap & dof_map;
86  const std::set<GhostingFunctor *> & coupling_functors;
89 
90  void handle_vi_vj(const Elem * partner,
91  const std::vector<dof_id_type> & element_dofs_i,
92  unsigned int vj);
93 
94 public:
95 
98 
99  std::vector<dof_id_type> n_nz;
100  std::vector<dof_id_type> n_oz;
101 
102  Build (const MeshBase & mesh_in,
103  const DofMap & dof_map_in,
104  const CouplingMatrix * dof_coupling_in,
105  std::set<GhostingFunctor *> coupling_functors_in,
106  const bool implicit_neighbor_dofs_in,
107  const bool need_full_sparsity_pattern_in);
108 
109  Build (Build & other, Threads::split);
110 
111  void operator()(const ConstElemRange & range);
112 
113  void join (const Build & other);
114 
115  void parallel_sync ();
116 };
117 
118 #if defined(__GNUC__) && (__GNUC__ < 4) && !defined(__INTEL_COMPILER)
119 
124 void _dummy_function(void);
125 #endif
126 
127 }
128 
129 
130 
131 // ------------------------------------------------------------
132 // SparsityPattern inline member functions
133 template<typename BidirectionalIterator>
134 inline
135 void SparsityPattern::sort_row (const BidirectionalIterator begin,
136  BidirectionalIterator middle,
137  const BidirectionalIterator end)
138 {
139  if ((begin == middle) || (middle == end)) return;
140 
141  libmesh_assert_greater (std::distance (begin, middle), 0);
142  libmesh_assert_greater (std::distance (middle, end), 0);
143  libmesh_assert (std::unique (begin, middle) == middle);
144  libmesh_assert (std::unique (middle, end) == end);
145 
146  while (middle != end)
147  {
148  BidirectionalIterator
149  b = middle,
150  a = b-1;
151 
152  // Bubble-sort the middle value downward
153  while (!(*a < *b)) // *a & *b are less-than comparable, so use <
154  {
155  std::swap (*a, *b);
156 
157 #if defined(__GNUC__) && (__GNUC__ < 4) && !defined(__INTEL_COMPILER)
158  /* Prohibit optimization at this point since gcc 3.3.5 seems
159  to have a bug. */
161 #endif
162 
163  if (a == begin) break;
164 
165  b=a;
166  --a;
167  }
168 
169  ++middle;
170  }
171 
172  // Assure the algorithm worked if we are in DEBUG mode
173 #ifdef DEBUG
174  {
175  // SGI STL extension!
176  // libmesh_assert (std::is_sorted(begin,end));
177 
178  BidirectionalIterator
179  prev = begin,
180  first = begin;
181 
182  for (++first; first != end; prev=first, ++first)
183  if (*first < *prev)
184  libmesh_assert(false);
185  }
186 #endif
187 
188  // Make sure the two ranges did not contain any common elements
189  libmesh_assert (std::unique (begin, end) == end);
190 }
191 
192 } // namespace libMesh
193 
194 #endif // LIBMESH_SPARSITY_PATTERN_H
SparsityPattern::Graph sparsity_pattern
This helper class can be called on multiple threads to compute the sparsity pattern (or graph) of the...
const std::set< GhostingFunctor * > & coupling_functors
Dummy "splitting object" used to distinguish splitting constructors from copy constructors.
Definition: threads_none.h:63
static const unsigned int prev[3]
A lookup table for the decrement modulo 3 operation, for iterating through the three nodes per elemen...
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
The StoredRange class defines a contiguous, divisible set of objects.
Definition: stored_range.h:52
The libMesh namespace provides an interface to certain functionality in the library.
Real distance(const Point &p)
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
std::vector< dof_id_type, Threads::scalable_allocator< dof_id_type > > Row
const CouplingMatrix * dof_coupling
void _dummy_function(void)
Dummy function that does nothing but can be used to prohibit compiler optimization in some situations...
std::vector< dof_id_type > n_nz
This class forms the base class for all other classes that are expected to be implemented in parallel...
SparsityPattern::NonlocalGraph nonlocal_pattern
void swap(Iterator &lhs, Iterator &rhs)
swap, used to implement op=
This class defines a coupling matrix.
static void sort_row(const BidirectionalIterator begin, BidirectionalIterator middle, const BidirectionalIterator end)
Splices the two sorted ranges [begin,middle) and [middle,end) into one sorted range [begin...
std::vector< dof_id_type > n_oz