Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
mg_transfer.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_transfer.templates.h 30040 2013-07-18 17:06:48Z 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_templates_h
19 #define __deal2__mg_transfer_templates_h
20 
21 #include <deal.II/lac/trilinos_vector.h>
22 #include <deal.II/lac/sparse_matrix.h>
23 #include <deal.II/grid/tria_iterator.h>
24 #include <deal.II/fe/fe.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.h>
29 
30 #include <algorithm>
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 
35 namespace
36 {
49  template <int dim, typename number, int spacedim>
50  void
51  reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
52  const std::vector<unsigned int> &,
54  {
55  for (unsigned int level=v.min_level();
56  level<=v.max_level(); ++level)
57  {
58  unsigned int n = mg_dof.n_dofs (level);
59  v[level].reinit(n);
60  }
61 
62  }
63 
74  template <int dim, typename number, int spacedim>
75  void
76  reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
77  std::vector<unsigned int> target_component,
79  {
80  const unsigned int n_blocks = mg_dof.get_fe().n_blocks();
81  if (target_component.size()==0)
82  {
83  target_component.resize(n_blocks);
84  for (unsigned int i=0; i<n_blocks; ++i)
85  target_component[i] = i;
86  }
87  Assert(target_component.size()==n_blocks,
88  ExcDimensionMismatch(target_component.size(),n_blocks));
89  const unsigned int max_block
90  = *std::max_element (target_component.begin(),
91  target_component.end());
92  const unsigned int n_target_blocks = max_block + 1;
93 
94  std::vector<std::vector<types::global_dof_index> >
95  ndofs(mg_dof.get_tria().n_levels(),
96  std::vector<types::global_dof_index>(n_target_blocks));
97  MGTools::count_dofs_per_block (mg_dof, ndofs, target_component);
98 
99  for (unsigned int level=v.min_level();
100  level<=v.max_level(); ++level)
101  {
102  v[level].reinit(n_target_blocks);
103  for (unsigned int b=0; b<n_target_blocks; ++b)
104  v[level].block(b).reinit(ndofs[level][b]);
105  v[level].collect_sizes();
106  }
107  }
108 
109 #ifdef DEAL_II_WITH_TRILINOS
110 
118  template <int dim, int spacedim>
119  void
120  reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
121  const std::vector<unsigned int> &,
123  {
124  const ::parallel::distributed::Triangulation<dim,spacedim> *tria =
126  (&mg_dof.get_tria()));
127  AssertThrow(tria!=NULL, ExcMessage("multigrid with Trilinos vectors only works with distributed Triangulation!"));
128 
129 #ifdef DEAL_II_WITH_P4EST
130  for (unsigned int level=v.min_level();
131  level<=v.max_level(); ++level)
132  {
133  v[level].reinit(mg_dof.locally_owned_mg_dofs(level), tria->get_communicator());
134  }
135 #else
136  (void)v;
137 #endif
138  }
139 #endif
140 }
141 
142 
143 
144 /* --------------------- MGTransferPrebuilt -------------- */
145 
146 
147 
148 
149 template <class VECTOR>
150 template <int dim, class InVector, int spacedim>
151 void
153  const DoFHandler<dim,spacedim> &mg_dof_handler,
155  const InVector &src) const
156 {
157  reinit_vector(mg_dof_handler, component_to_block_map, dst);
158  bool first = true;
159  for (unsigned int level=mg_dof_handler.get_tria().n_global_levels(); level != 0;)
160  {
161  --level;
162  VECTOR &dst_level = dst[level];
163 
164  typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
165  for (IT i= copy_indices[level].begin();
166  i != copy_indices[level].end(); ++i)
167  dst_level(i->second) = src(i->first);
168 
169  // for (IT i= copy_indices_to_me[level].begin();
170  // i != copy_indices_to_me[level].end(); ++i)
171  // dst_level(i->second) = src(i->first);
172 
173  dst_level.compress(VectorOperation::insert);
174 
175  if (!first)
176  restrict_and_add (level+1, dst[level], dst[level+1]);
177 
178  first = false;
179  }
180 }
181 
182 
183 
184 template <class VECTOR>
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> &src) const
191 {
192  // For non-DG: degrees of
193  // freedom in the refinement
194  // face may need special
195  // attention, since they belong
196  // to the coarse level, but
197  // have fine level basis
198  // functions
199  dst = 0;
200  for (unsigned int level=0; level<mg_dof_handler.get_tria().n_global_levels(); ++level)
201  {
202  typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
203 
204  // First copy all indices local to this process
205  if (constraints==0)
206  for (IT i= copy_indices[level].begin();
207  i != copy_indices[level].end(); ++i)
208  dst(i->first) = src[level](i->second);
209  else
210  for (IT i= copy_indices[level].begin();
211  i != copy_indices[level].end(); ++i)
212  constraints->distribute_local_to_global(i->first, src[level](i->second), dst);
213 
214  // Do the same for the indices where the level index is local,
215  // but the global index is not
216  if (constraints==0)
217  for (IT i= copy_indices_from_me[level].begin();
218  i != copy_indices_from_me[level].end(); ++i)
219  dst(i->first) = src[level](i->second);
220  else
221  for (IT i= copy_indices_from_me[level].begin();
222  i != copy_indices_from_me[level].end(); ++i)
223  constraints->distribute_local_to_global(i->first, src[level](i->second), dst);
224  }
225 }
226 
227 
228 
229 template <class VECTOR>
230 template <int dim, class OutVector, int spacedim>
231 void
233  const DoFHandler<dim,spacedim> &mg_dof_handler,
234  OutVector &dst,
235  const MGLevelObject<VECTOR> &src) const
236 {
237  // For non-DG: degrees of
238  // freedom in the refinement
239  // face may need special
240  // attention, since they belong
241  // to the coarse level, but
242  // have fine level basis
243  // functions
244  for (unsigned int level=0; level<mg_dof_handler.get_tria().n_global_levels(); ++level)
245  {
246  typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
247  if (constraints==0)
248  for (IT i= copy_indices[level].begin();
249  i != copy_indices[level].end(); ++i)
250  dst(i->first) += src[level](i->second);
251  else
252  for (IT i= copy_indices[level].begin();
253  i != copy_indices[level].end(); ++i)
254  constraints->distribute_local_to_global(i->first, src[level](i->second), dst);
255 
256  // Do the same for the indices where the level index is local,
257  // but the global index is not
258  if (constraints==0)
259  for (IT i= copy_indices_from_me[level].begin();
260  i != copy_indices_from_me[level].end(); ++i)
261  dst(i->first) += src[level](i->second);
262  else
263  for (IT i= copy_indices_from_me[level].begin();
264  i != copy_indices_from_me[level].end(); ++i)
265  constraints->distribute_local_to_global(i->first, src[level](i->second), dst);
266  }
267 }
268 
269 
270 
271 template <class VECTOR>
272 void
274 set_component_to_block_map (const std::vector<unsigned int> &map)
275 {
276  component_to_block_map = map;
277 }
278 
279 template <class VECTOR>
280 std::size_t
282 {
283  std::size_t result = sizeof(*this);
284  result += sizeof(unsigned int) * sizes.size();
285 
286  for (unsigned int i=0; i<prolongation_matrices.size(); ++i)
287  result += prolongation_matrices[i]->memory_consumption()
288  + prolongation_sparsities[i]->memory_consumption();
289 
290  return result;
291 }
292 
293 
294 DEAL_II_NAMESPACE_CLOSE
295 
296 #endif
::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:362
#define Assert(cond, exc)
Definition: exceptions.h:299
void count_dofs_per_block(const DH &dof_handler, std::vector< std::vector< types::global_dof_index > > &dofs_per_block, std::vector< unsigned int > target_block=std::vector< unsigned int >())
unsigned int max_level() const
void set_component_to_block_map(const std::vector< unsigned int > &map)
void copy_from_mg(const DoFHandler< dim, spacedim > &mg_dof, OutVector &dst, const MGLevelObject< VECTOR > &src) const
void copy_from_mg_add(const DoFHandler< dim, spacedim > &mg_dof, OutVector &dst, const MGLevelObject< VECTOR > &src) const
std::size_t memory_consumption() const
void copy_to_mg(const DoFHandler< dim, spacedim > &mg_dof, MGLevelObject< VECTOR > &dst, const InVector &src) const
unsigned int min_level() const
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
const Triangulation< dim, spacedim > & get_tria() const