Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
dof_info.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_info.h 30040 2013-07-18 17:06:48Z maier @f$
3 //
4 // Copyright (C) 2006 - 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 
18 #ifndef __deal2__mesh_worker_dof_info_h
19 #define __deal2__mesh_worker_dof_info_h
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/quadrature_lib.h>
23 #include <deal.II/base/std_cxx1x/shared_ptr.h>
24 #include <deal.II/dofs/block_info.h>
25 #include <deal.II/fe/fe_values.h>
26 #include <deal.II/meshworker/local_results.h>
27 #include <deal.II/meshworker/vector_selector.h>
28 
30 
31 namespace MeshWorker
32 {
33  template <int dim, class DOFINFO> class DoFInfoBox;
68  template<int dim, int spacedim = dim, typename number = double>
69  class DoFInfo : public LocalResults<number>
70  {
71  public:
74 
77 
88  unsigned int face_number;
99  unsigned int sub_number;
100 
101  /*
102  * The DoF indices of the
103  * current cell
104  */
105  std::vector<types::global_dof_index> indices;
106 
111  std::vector<std::vector<types::global_dof_index> > indices_by_block;
112 
117  DoFInfo(const BlockInfo &block_info);
118 
125  DoFInfo (const DoFHandler<dim, spacedim> &dof_handler);
126 
131  template <class DHCellIterator>
132  void reinit(const DHCellIterator &c);
133 
137  template <class DHCellIterator, class DHFaceIterator>
138  void reinit(const DHCellIterator &c,
139  const DHFaceIterator &f,
140  const unsigned int face_no);
141 
146  template <class DHCellIterator, class DHFaceIterator>
147  void reinit(const DHCellIterator &c,
148  const DHFaceIterator &f,
149  const unsigned int face_no,
150  const unsigned int subface_no);
151 
156  template <class DHFaceIterator>
157  void set_face (const DHFaceIterator &f,
158  const unsigned int face_no);
159 
164  template <class DHFaceIterator>
165  void set_subface (const DHFaceIterator &f,
166  const unsigned int face_no,
167  const unsigned int subface_no);
168 
169  const BlockIndices &local_indices() const;
170 
171 
174 
175  bool level_cell;
176  private:
182  DoFInfo ();
183 
185  void set_block_indices ();
187  template <class DHCellIterator>
188  void get_indices(const DHCellIterator &c);
189 
191  std::vector<types::global_dof_index> indices_org;
192 
202 
203  friend class DoFInfoBox<dim, DoFInfo<dim, spacedim, number> >;
204  };
205 
206 
218  template <int dim, class DOFINFO>
219  class DoFInfoBox
220  {
221  public:
226  DoFInfoBox(const DOFINFO &seed);
227 
234 
238  void reset();
239 
249  template <class ASSEMBLER>
250  void assemble(ASSEMBLER &ass) const;
251 
255  std::size_t memory_consumption () const;
256 
257 
261  DOFINFO cell;
270 
283  };
284 
285 //----------------------------------------------------------------------//
286 
287  template <int dim, int spacedim, typename number>
289  {
290  std::vector<types::global_dof_index> aux(1);
291  aux[0] = dof_handler.get_fe().dofs_per_cell;
292  aux_local_indices.reinit(aux);
293  }
294 
295 
296  template <int dim, int spacedim, typename number>
297  template <class DHCellIterator>
298  inline void
300  {
301  indices.resize(c->get_fe().dofs_per_cell);
302  if (block_info == 0 || block_info->local().size() == 0)
303  c->get_active_or_mg_dof_indices(indices);
304  else
305  {
306  indices_org.resize(c->get_fe().dofs_per_cell);
307  c->get_active_or_mg_dof_indices(indices_org);
308  set_block_indices();
309  }
310  }
311 
312 
313  template <int dim, int spacedim, typename number>
314  template <class DHCellIterator>
315  inline void
316  DoFInfo<dim,spacedim,number>::reinit(const DHCellIterator &c)
317  {
318  get_indices(c);
319  level_cell = c->is_level_cell();
320 
321  cell = typename Triangulation<dim,spacedim>::cell_iterator(*c);
322  face_number = deal_II_numbers::invalid_unsigned_int;
323  sub_number = deal_II_numbers::invalid_unsigned_int;
324  if (block_info)
325  LocalResults<number>::reinit(block_info->local());
326  else
327  LocalResults<number>::reinit(aux_local_indices);
328  }
329 
330 
331  template<int dim, int spacedim, typename number>
332  template <class DHFaceIterator>
333  inline void
335  const DHFaceIterator &f,
336  const unsigned int face_no)
337  {
338  face = static_cast<typename Triangulation<dim>::face_iterator> (f);
339  face_number = face_no;
340  sub_number = deal_II_numbers::invalid_unsigned_int;
341  }
342 
343 
344  template<int dim, int spacedim, typename number>
345  template <class DHCellIterator, class DHFaceIterator>
346  inline void
348  const DHCellIterator &c,
349  const DHFaceIterator &f,
350  const unsigned int face_no)
351  {
352  if ((cell.state() != IteratorState::valid)
353  || cell != typename Triangulation<dim>::cell_iterator(*c))
354  get_indices(c);
355  level_cell = c->is_level_cell();
356 
357  cell = typename Triangulation<dim>::cell_iterator(*c);
358  set_face(f,face_no);
359 
360  if (block_info)
361  LocalResults<number>::reinit(block_info->local());
362  else
363  LocalResults<number>::reinit(aux_local_indices);
364  }
365 
366 
367  template<int dim, int spacedim, typename number>
368  template <class DHFaceIterator>
369  inline void
371  const DHFaceIterator &f,
372  const unsigned int face_no,
373  const unsigned int subface_no)
374  {
375  face = static_cast<typename Triangulation<dim>::face_iterator> (f);
376  face_number = face_no;
377  sub_number = subface_no;
378  }
379 
380 
381  template<int dim, int spacedim, typename number>
382  template <class DHCellIterator, class DHFaceIterator>
383  inline void
385  const DHCellIterator &c,
386  const DHFaceIterator &f,
387  const unsigned int face_no,
388  const unsigned int subface_no)
389  {
390  if (cell.state() != IteratorState::valid
391  || cell != static_cast<typename Triangulation<dim>::cell_iterator> (c))
392  get_indices(c);
393  level_cell = c->is_level_cell();
394 
395  cell = static_cast<typename Triangulation<dim>::cell_iterator> (c);
396  set_subface(f, face_no, subface_no);
397 
398  if (block_info)
399  LocalResults<number>::reinit(block_info->local());
400  else
401  LocalResults<number>::reinit(aux_local_indices);
402  }
403 
404 
405  template<int dim, int spacedim, typename number>
406  inline const BlockIndices &
408  {
409  if (block_info)
410  return block_info->local();
411  return aux_local_indices;
412  }
413 
414 //----------------------------------------------------------------------//
415 
416  template <int dim, class DOFINFO>
417  inline
419  :
420  cell(seed)
421  {
422  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
423  {
424  exterior[i] = seed;
425  interior[i] = seed;
426  interior_face_available[i] = false;
427  exterior_face_available[i] = false;
428  }
429  }
430 
431 
432  template <int dim, class DOFINFO>
433  inline
435  :
436  cell(other.cell)
437  {
438  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
439  {
440  exterior[i] = other.exterior[i];
441  interior[i] = other.interior[i];
442  interior_face_available[i] = false;
443  exterior_face_available[i] = false;
444  }
445  }
446 
447 
448  template <int dim, class DOFINFO>
449  inline void
451  {
452  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
453  {
454  interior_face_available[i] = false;
455  exterior_face_available[i] = false;
456  }
457  }
458 
459 
460  template <int dim, class DOFINFO>
461  template <class ASSEMBLER>
462  inline void
463  DoFInfoBox<dim, DOFINFO>::assemble (ASSEMBLER &assembler) const
464  {
465  assembler.assemble(cell);
466  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
467  {
468  // Only do something if data available
469  if (interior_face_available[i])
470  {
471  // If both data
472  // available, it is an
473  // interior face
474  if (exterior_face_available[i])
475  assembler.assemble(interior[i], exterior[i]);
476  else
477  assembler.assemble(interior[i]);
478  }
479  }
480  }
481 }
482 
483 DEAL_II_NAMESPACE_CLOSE
484 
485 #endif
Iterator points to a valid object.
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:173
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:269
void set_subface(const DHFaceIterator &f, const unsigned int face_no, const unsigned int subface_no)
Definition: dof_info.h:370
Auxiliary class aiding in the handling of block structures like in BlockVector or FESystem...
Definition: block_indices.h:55
void assemble(ASSEMBLER &ass) const
Definition: dof_info.h:463
void get_indices(const DHCellIterator &c)
Fill index vector with active indices.
Definition: dof_info.h:299
bool exterior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:282
void set_face(const DHFaceIterator &f, const unsigned int face_no)
Definition: dof_info.h:334
unsigned int face_number
Definition: dof_info.h:88
std::vector< types::global_dof_index > indices_org
Auxiliary vector.
Definition: dof_info.h:191
BlockIndices aux_local_indices
Definition: dof_info.h:201
DoFInfoBox(const DOFINFO &seed)
Definition: dof_info.h:418
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:265
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:73
void set_block_indices()
Set up local block indices.
std::size_t memory_consumption() const
unsigned int sub_number
Definition: dof_info.h:99
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:276
void reinit(const BlockIndices &local_sizes)
void reinit(const DHCellIterator &c)
Definition: dof_info.h:316
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:93
const FiniteElement< dim, spacedim > & get_fe() const
std::vector< std::vector< types::global_dof_index > > indices_by_block
Definition: dof_info.h:111
Triangulation< dim, spacedim >::face_iterator face
The current face.
Definition: dof_info.h:76