18 #ifndef __deal2__mg_transfer_templates_h
19 #define __deal2__mg_transfer_templates_h
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>
49 template <
int dim,
typename number,
int spacedim>
51 reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
52 const std::vector<unsigned int> &,
55 for (
unsigned int level=v.min_level();
56 level<=v.max_level(); ++level)
58 unsigned int n = mg_dof.n_dofs (level);
74 template <
int dim,
typename number,
int spacedim>
76 reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
77 std::vector<unsigned int> target_component,
80 const unsigned int n_blocks = mg_dof.get_fe().n_blocks();
81 if (target_component.size()==0)
83 target_component.resize(n_blocks);
84 for (
unsigned int i=0; i<n_blocks; ++i)
85 target_component[i] = i;
87 Assert(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;
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));
99 for (
unsigned int level=v.min_level();
100 level<=v.max_level(); ++level)
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();
109 #ifdef DEAL_II_WITH_TRILINOS
118 template <
int dim,
int spacedim>
120 reinit_vector (const ::DoFHandler<dim,spacedim> &mg_dof,
121 const std::vector<unsigned int> &,
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!"));
129 #ifdef DEAL_II_WITH_P4EST
133 v[level].reinit(mg_dof.locally_owned_mg_dofs(level), tria->get_communicator());
149 template <
class VECTOR>
150 template <
int dim,
class InVector,
int spacedim>
155 const InVector &src)
const
157 reinit_vector(mg_dof_handler, component_to_block_map, dst);
159 for (
unsigned int level=mg_dof_handler.
get_tria().n_global_levels(); level != 0;)
162 VECTOR &dst_level = dst[level];
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);
173 dst_level.compress(VectorOperation::insert);
176 restrict_and_add (level+1, dst[level], dst[level+1]);
184 template <
class VECTOR>
185 template <
int dim,
class OutVector,
int spacedim>
200 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_global_levels(); ++level)
202 typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
206 for (IT i= copy_indices[level].begin();
207 i != copy_indices[level].end(); ++i)
208 dst(i->first) = src[level](i->second);
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);
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);
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);
229 template <
class VECTOR>
230 template <
int dim,
class OutVector,
int spacedim>
244 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_global_levels(); ++level)
246 typedef std::vector<std::pair<types::global_dof_index, unsigned int> >::const_iterator IT;
248 for (IT i= copy_indices[level].begin();
249 i != copy_indices[level].end(); ++i)
250 dst(i->first) += src[level](i->second);
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);
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);
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);
271 template <
class VECTOR>
276 component_to_block_map = map;
279 template <
class VECTOR>
283 std::size_t result =
sizeof(*this);
284 result +=
sizeof(
unsigned int) * sizes.size();
286 for (
unsigned int i=0; i<prolongation_matrices.size(); ++i)
287 result += prolongation_matrices[i]->memory_consumption()
288 + prolongation_sparsities[i]->memory_consumption();
294 DEAL_II_NAMESPACE_CLOSE
::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
#define Assert(cond, exc)
void set_component_to_block_map(const std::vector< unsigned int > &map)
const Triangulation< dim, spacedim > & get_tria() const
unsigned int min_level() const
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
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
unsigned int max_level() const