18 #ifndef __deal2__mg_transfer_block_templates_h
19 #define __deal2__mg_transfer_block_templates_h
21 #include <deal.II/lac/sparse_matrix.h>
22 #include <deal.II/lac/constraint_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_block.h>
32 DEAL_II_NAMESPACE_OPEN
37 typedef std::vector<std::pair<unsigned int, unsigned int> >::const_iterator IT;
40 template <
typename number>
41 template <
int dim,
typename number2,
int spacedim>
48 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
49 for (IT i= copy_indices[selected_block][level].begin();
50 i != copy_indices[selected_block][level].end(); ++i)
51 dst.
block(selected_block)(i->first) = src[level](i->second);
56 template <
typename number>
57 template <
int dim,
typename number2,
int spacedim>
64 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
65 for (IT i= copy_indices[selected_block][level].begin();
66 i != copy_indices[selected_block][level].end(); ++i)
67 dst(i->first) = src[level](i->second);
72 template <
typename number>
73 template <
int dim,
typename number2,
int spacedim>
80 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
81 for (IT i= copy_indices[selected_block][level].begin();
82 i != copy_indices[selected_block][level].end(); ++i)
83 dst.
block(selected_block)(i->first) += src[level](i->second);
88 template <
typename number>
89 template <
int dim,
typename number2,
int spacedim>
96 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
97 for (IT i= copy_indices[selected_block][level].begin();
98 i != copy_indices[selected_block][level].end(); ++i)
99 dst(i->first) += src[level](i->second);
104 template <
typename number>
116 template <
typename number>
117 template <
int dim,
typename number2,
int spacedim>
124 for (
unsigned int block=0; block<selected.size(); ++block)
126 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
127 for (IT i= copy_indices[block][level].begin();
128 i != copy_indices[block][level].end(); ++i)
129 dst.
block(block)(i->first) = src[level].block(mg_block[block])(i->second);
134 template <
typename number>
135 template <
int dim,
typename number2,
int spacedim>
142 for (
unsigned int block=0; block<selected.size(); ++block)
144 for (
unsigned int level=0; level<mg_dof_handler.
get_tria().n_levels(); ++level)
145 for (IT i= copy_indices[block][level].begin();
146 i != copy_indices[block][level].end(); ++i)
147 dst.
block(block)(i->first) += src[level].block(mg_block[block])(i->second);
150 DEAL_II_NAMESPACE_CLOSE
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
void copy_from_mg_add(const DoFHandler< dim, spacedim > &mg_dof, BlockVector< number2 > &dst, const MGLevelObject< BlockVector< number > > &src) const
void copy_from_mg(const DoFHandler< dim, spacedim > &mg_dof, BlockVector< number2 > &dst, const MGLevelObject< BlockVector< number > > &src) const
BlockType & block(const unsigned int i)
std::size_t memory_consumption() const
std::size_t memory_consumption() const
const Triangulation< dim, spacedim > & get_tria() const