Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
mg_constrained_dofs.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_constrained_dofs.h 31075 2013-10-02 17:51:27Z heister @f$
3 //
4 // Copyright (C) 2010 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__mg_constraints_h
18 #define __deal2__mg_constraints_h
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/subscriptor.h>
22 
23 #include <deal.II/multigrid/mg_tools.h>
24 
25 #include <vector>
26 #include <set>
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 template <int dim, int spacedim> class DoFHandler;
31 template <int dim> struct FunctionMap;
32 
33 
41 {
42 public:
43 
44  typedef std::vector<std::set<types::global_dof_index> >::size_type size_dof;
56  template <int dim, int spacedim>
57  void initialize(const DoFHandler<dim,spacedim> &dof);
58 
66  template <int dim, int spacedim>
67  void initialize(const DoFHandler<dim,spacedim> &dof,
68  const typename FunctionMap<dim>::type &function_map,
69  const ComponentMask &component_mask = ComponentMask());
70 
74  void clear();
75 
80  bool is_boundary_index (const unsigned int level,
81  const types::global_dof_index index) const;
82 
87 // TODO: remove
88  bool non_refinement_edge_index (const unsigned int level,
89  const types::global_dof_index index) const;
90 
94  bool at_refinement_edge (const unsigned int level,
95  const types::global_dof_index index) const;
96 
102  bool at_refinement_edge_boundary (const unsigned int level,
103  const types::global_dof_index index) const;
104 
109 // TODO: remove
110  const std::vector<std::set<types::global_dof_index> > &
111  get_boundary_indices () const;
112 
117  // TODO: remove
118  const std::vector<std::set<types::global_dof_index> > &
120 
126  // TODO: remove
127  const std::vector<std::vector<bool> > &
129 
135  // TODO: remove
136  const std::vector<std::vector<bool> > &
138 
143  bool set_boundary_values () const;
144 
150 private:
151 
155  std::vector<std::set<types::global_dof_index> > boundary_indices;
156 
161  std::vector<std::set<types::global_dof_index> > non_refinement_edge_indices;
162 
167  std::vector<std::vector<bool> > refinement_edge_indices;
168 
175  std::vector<std::vector<bool> > refinement_edge_boundary_indices;
176 };
177 
178 
179 template <int dim, int spacedim>
180 inline
181 void
183 {
184  const unsigned int nlevels = dof.get_tria().n_global_levels();
185  refinement_edge_indices.resize(nlevels);
186  refinement_edge_boundary_indices.resize(nlevels);
187  non_refinement_edge_indices.resize(nlevels);
188  for (unsigned int l=0; l<nlevels; ++l)
189  {
190  refinement_edge_indices[l].resize(dof.n_dofs(l));
191  refinement_edge_boundary_indices[l].resize(dof.n_dofs(l));
192  }
195 
196  MGTools::extract_non_interface_dofs (dof, non_refinement_edge_indices);
197 }
198 
199 
200 template <int dim, int spacedim>
201 inline
202 void
204  const DoFHandler<dim,spacedim> &dof,
205  const typename FunctionMap<dim>::type &function_map,
206  const ComponentMask &component_mask)
207 {
208  const unsigned int nlevels = dof.get_tria().n_global_levels();
209  boundary_indices.resize(nlevels);
210  refinement_edge_indices.resize(nlevels);
211  refinement_edge_boundary_indices.resize(nlevels);
212  non_refinement_edge_indices.resize(nlevels);
213 
214  for (unsigned int l=0; l<nlevels; ++l)
215  {
216  boundary_indices[l].clear();
217  refinement_edge_indices[l].resize(dof.n_dofs(l));
218  refinement_edge_boundary_indices[l].resize(dof.n_dofs(l));
219  }
220 
221  MGTools::make_boundary_list (dof, function_map, boundary_indices, component_mask);
224  MGTools::extract_non_interface_dofs (dof, non_refinement_edge_indices);
225 }
226 
227 
228 inline
229 void
231 {
232  for (unsigned int l=0; l<boundary_indices.size(); ++l)
233  boundary_indices[l].clear();
234 
235  for (unsigned int l=0; l<refinement_edge_indices.size(); ++l)
237 
238  for (unsigned int l=0; l<refinement_edge_boundary_indices.size(); ++l)
240 
241  for (unsigned int l=0; l<non_refinement_edge_indices.size(); ++l)
243 }
244 
245 
246 inline
247 bool
248 MGConstrainedDoFs::is_boundary_index (const unsigned int level,
249  const types::global_dof_index index) const
250 {
251  AssertIndexRange(level, boundary_indices.size());
252  if (boundary_indices[level].find(index) != boundary_indices[level].end())
253  return true;
254  else
255  return false;
256 }
257 
258 inline
259 bool
261  const types::global_dof_index index) const
262 {
264 
265  if (non_refinement_edge_indices[level].find(index) != non_refinement_edge_indices[level].end())
266  return true;
267  else
268  return false;
269 }
270 
271 inline
272 bool
273 MGConstrainedDoFs::at_refinement_edge (const unsigned int level,
274  const types::global_dof_index index) const
275 {
277  AssertIndexRange(index, refinement_edge_indices[level].size());
278 
279  return refinement_edge_indices[level][index];
280 }
281 
282 
283 inline
284 bool
286  const types::global_dof_index index) const
287 {
290 
291  return refinement_edge_boundary_indices[level][index];
292 }
293 
294 inline
295 const std::vector<std::set<types::global_dof_index> > &
297 {
298  return boundary_indices;
299 }
300 
301 inline
302 const std::vector<std::set<types::global_dof_index> > &
304 {
306 }
307 
308 inline
309 const std::vector<std::vector<bool> > &
311 {
313 }
314 
315 inline
316 const std::vector<std::vector<bool> > &
318 {
320 }
321 
322 inline
323 bool
325 {
326  const bool boundary_values_need_to_be_set
327  = boundary_indices.size()!=0;
328  return boundary_values_need_to_be_set;
329 }
330 
331 inline
332 bool
334 {
335  bool is_continuous = false;
336  for (unsigned int l=0; l<refinement_edge_indices.size(); ++l)
337  for (unsigned int i=0; i<refinement_edge_indices[l].size(); ++i)
338  if (refinement_edge_indices[l][i])
339  {
340  is_continuous = true;
341  return is_continuous;
342  }
343  return is_continuous;
344 }
345 
346 DEAL_II_NAMESPACE_CLOSE
347 
348 #endif
bool is_boundary_index(const unsigned int level, const types::global_dof_index index) const
bool at_refinement_edge(const unsigned int level, const types::global_dof_index index) const
std::vector< std::vector< bool > > refinement_edge_indices
std::map< types::boundary_id, const Function< dim > * > type
Definition: function_map.h:56
#define AssertIndexRange(index, range)
Definition: exceptions.h:888
const std::vector< std::set< types::global_dof_index > > & get_non_refinement_edge_indices() const
const std::vector< std::vector< bool > > & get_refinement_edge_boundary_indices() const
void make_boundary_list(const DoFHandler< dim, spacedim > &mg_dof, const typename FunctionMap< dim >::type &function_map, std::vector< std::set< types::global_dof_index > > &boundary_indices, const ComponentMask &component_mask=ComponentMask())
void extract_inner_interface_dofs(const DoFHandler< dim, spacedim > &mg_dof_handler, std::vector< std::vector< bool > > &interface_dofs, std::vector< std::vector< bool > > &boundary_interface_dofs)
unsigned int global_dof_index
Definition: types.h:100
types::global_dof_index n_dofs() const
const std::vector< std::set< types::global_dof_index > > & get_boundary_indices() const
std::vector< std::set< types::global_dof_index > > non_refinement_edge_indices
bool non_refinement_edge_index(const unsigned int level, const types::global_dof_index index) const
std::vector< std::set< types::global_dof_index > > boundary_indices
std::vector< std::vector< bool > > refinement_edge_boundary_indices
bool continuity_across_refinement_edges() const
const Triangulation< dim, spacedim > & get_tria() const
bool set_boundary_values() const
const std::vector< std::vector< bool > > & get_refinement_edge_indices() const
bool at_refinement_edge_boundary(const unsigned int level, const types::global_dof_index index) const
void initialize(const DoFHandler< dim, spacedim > &dof)