Reference documentation for deal.II version 8.1.0
mg_transfer_component.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_transfer_component.templates.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2003 - 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__mg_transfer_component_templates_h
19 #define __deal2__mg_transfer_component_templates_h
20 
21 #include <deal.II/lac/sparse_matrix.h>
22 #include <deal.II/grid/tria_iterator.h>
23 #include <deal.II/fe/fe.h>
24 #include <deal.II/lac/constraint_matrix.h>
25 #include <deal.II/multigrid/mg_base.h>
26 #include <deal.II/dofs/dof_accessor.h>
27 #include <deal.II/multigrid/mg_tools.h>
28 #include <deal.II/multigrid/mg_transfer_component.h>
29 #include <deal.II/numerics/data_out.h>
30 
31 #include <algorithm>
32 #include <sstream>
33 #include <fstream>
34 
36 
37 /* --------------------- MGTransferSelect -------------- */
38 
39 
40 
41 template <typename number>
42 template <int dim, typename number2, int spacedim>
43 void
45  const DoFHandler<dim,spacedim> &mg_dof_handler,
47  const BlockVector<number2> &src) const
48 {
49  do_copy_to_mg (mg_dof_handler, dst, src.block(target_component[selected_component]));
50 }
51 
52 
53 
54 template <typename number>
55 template <int dim, typename number2, int spacedim>
56 void
58  const DoFHandler<dim,spacedim> &mg_dof_handler,
60  const Vector<number2> &src) const
61 {
62  do_copy_to_mg (mg_dof_handler, dst, src);
63 }
64 
65 
66 
67 template <typename number>
68 template <int dim, typename number2, int spacedim>
69 void
71  const DoFHandler<dim,spacedim> &mg_dof_handler,
73  const MGLevelObject<Vector<number> > &src) const
74 {
75  dst = 0;
76  do_copy_from_mg (mg_dof_handler,
77  dst.block(target_component[selected_component]), src);
78  if (constraints != 0)
79  constraints->condense(dst);
80 }
81 
82 
83 
84 template <typename number>
85 template <int dim, typename number2, int spacedim>
86 void
88  const DoFHandler<dim,spacedim> &mg_dof_handler,
89  Vector<number2> &dst,
90  const MGLevelObject<Vector<number> > &src) const
91 {
92  dst = 0;
93  do_copy_from_mg (mg_dof_handler, dst, src);
94  if (constraints != 0)
95  {
96  //If we were given constraints
97  //apply them to the dst that goes
98  //back now to the linear solver.
99  //Since constraints are globally
100  //defined create a global vector here
101  //and copy dst to the right component,
102  //apply the constraints then and copy
103  //the block back to dst.
104  const unsigned int n_blocks =
105  *std::max_element(target_component.begin(), target_component.end()) + 1;
106  std::vector<types::global_dof_index> dofs_per_block (n_blocks);
107  DoFTools::count_dofs_per_block (mg_dof_handler, dofs_per_block, target_component);
109  tmp.reinit(n_blocks);
110  for (unsigned int b=0; b<n_blocks; ++b)
111  tmp.block(b).reinit(dofs_per_block[b]);
112  tmp.collect_sizes ();
113  tmp.block(target_component[selected_component]) = dst;
114  constraints->condense(tmp);
115  dst = tmp.block(target_component[selected_component]);
116  }
117 }
118 
119 
120 
121 template <typename number>
122 template <int dim, typename number2, int spacedim>
123 void
125  const DoFHandler<dim,spacedim> &mg_dof_handler,
127  const MGLevelObject<Vector<number> > &src) const
128 {
129  do_copy_from_mg_add (mg_dof_handler, dst, src);
130 }
131 
132 
133 
134 template <typename number>
135 template <int dim, typename number2, int spacedim>
136 void
138  const DoFHandler<dim,spacedim> &mg_dof_handler,
139  Vector<number2> &dst,
140  const MGLevelObject<Vector<number> > &src) const
141 {
142  do_copy_from_mg_add (mg_dof_handler, dst, src);
143 }
144 
145 
146 
147 template <typename number>
148 template <int dim, class OutVector, int spacedim>
149 void
151  const DoFHandler<dim,spacedim> &mg_dof_handler,
152  OutVector &dst,
153  const MGLevelObject<Vector<number> > &src) const
154 {
155 // const FiniteElement<dim>& fe = mg_dof_handler.get_fe();
156 
157 // const unsigned int dofs_per_cell = fe.dofs_per_cell;
158 // std::vector<unsigned int> global_dof_indices (dofs_per_cell);
159 // std::vector<unsigned int> level_dof_indices (dofs_per_cell);
160 
161  typename DoFHandler<dim,spacedim>::active_cell_iterator
162  level_cell = mg_dof_handler.begin_active();
163  const typename DoFHandler<dim,spacedim>::active_cell_iterator
164  endc = mg_dof_handler.end();
165 
166  // traverse all cells and copy the
167  // data appropriately to the output
168  // vector
169 
170  // Note that the level is
171  // monotonuosly increasing
172  dst = 0;
173  for (; level_cell != endc; ++level_cell)
174  {
175  const unsigned int level = level_cell->level();
176  typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
177  for (IT i=copy_to_and_from_indices[level].begin();
178  i != copy_to_and_from_indices[level].end(); ++i)
179  dst(i->first) = src[level](i->second);
180  }
181 }
182 
183 
184 template <typename number>
185 template <int dim, class OutVector, int spacedim>
186 void
188  const DoFHandler<dim,spacedim> &mg_dof_handler,
189  OutVector &dst,
190  const MGLevelObject<Vector<number> > &src) const
191 {
192  const FiniteElement<dim> &fe = mg_dof_handler.get_fe();
193  const unsigned int dofs_per_cell = fe.dofs_per_cell;
194 
195  std::vector<types::global_dof_index> global_dof_indices (dofs_per_cell);
196  std::vector<types::global_dof_index> level_dof_indices (dofs_per_cell);
197 
198  typename DoFHandler<dim,spacedim>::active_cell_iterator
199  level_cell = mg_dof_handler.begin_active();
200  const typename DoFHandler<dim,spacedim>::active_cell_iterator
201  endc = mg_dof_handler.end();
202 
203  // traverse all cells and copy the
204  // data appropriately to the output
205  // vector
206 
207  // Note that the level is
208  // monotonuosly increasing
209  dst = 0;
210  for (; level_cell != endc; ++level_cell)
211  {
212  const unsigned int level = level_cell->level();
213  typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
214  for (IT i=copy_to_and_from_indices[level].begin();
215  i != copy_to_and_from_indices[level].end(); ++i)
216  dst(i->first) += src[level](i->second);
217  }
218 }
219 
220 
221 template <typename number>
222 std::size_t
224 {
226 }
227 
228 
229 
230 DEAL_II_NAMESPACE_CLOSE
231 
232 #endif
active_cell_iterator begin_active(const unsigned int level=0) const
std::size_t memory_consumption() const
cell_iterator end() const
void collect_sizes()
void do_copy_from_mg(const DoFHandler< dim, spacedim > &mg_dof, OutVector &dst, const MGLevelObject< Vector< number > > &src) const
void reinit(const unsigned int num_blocks, const size_type block_size=0, const bool fast=false)
void copy_to_mg(const DoFHandler< dim, spacedim > &mg_dof, MGLevelObject< Vector< number > > &dst, const Vector< number2 > &src) const
std::size_t memory_consumption() const
const unsigned int dofs_per_cell
Definition: fe_base.h:271
void count_dofs_per_block(const DH &dof, std::vector< types::global_dof_index > &dofs_per_block, const std::vector< unsigned int > &target_block=std::vector< unsigned int >())
BlockType & block(const unsigned int i)
void do_copy_from_mg_add(const DoFHandler< dim, spacedim > &mg_dof, OutVector &dst, const MGLevelObject< Vector< number > > &src) const
void copy_from_mg_add(const DoFHandler< dim, spacedim > &mg_dof, Vector< number2 > &dst, const MGLevelObject< Vector< number > > &src) const
void copy_from_mg(const DoFHandler< dim, spacedim > &mg_dof, Vector< number2 > &dst, const MGLevelObject< Vector< number > > &src) const
const FiniteElement< dim, spacedim > & get_fe() const