Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
mg_transfer_block.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_transfer_block.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_block_templates_h
19 #define __deal2__mg_transfer_block_templates_h
20 
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>
29 
30 #include <algorithm>
31 
33 
34 /* --------------------- MGTransferBlockSelect -------------- */
35 
36 // Simplify some things below
37 typedef std::vector<std::pair<unsigned int, unsigned int> >::const_iterator IT;
38 
39 
40 template <typename number>
41 template <int dim, typename number2, int spacedim>
42 void
44  const DoFHandler<dim,spacedim> &mg_dof_handler,
46  const MGLevelObject<Vector<number> > &src) const
47 {
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);
52 }
53 
54 
55 
56 template <typename number>
57 template <int dim, typename number2, int spacedim>
58 void
60  const DoFHandler<dim,spacedim> &mg_dof_handler,
61  Vector<number2> &dst,
62  const MGLevelObject<Vector<number> > &src) const
63 {
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);
68 }
69 
70 
71 
72 template <typename number>
73 template <int dim, typename number2, int spacedim>
74 void
76  const DoFHandler<dim,spacedim> &mg_dof_handler,
78  const MGLevelObject<Vector<number> > &src) const
79 {
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);
84 }
85 
86 
87 
88 template <typename number>
89 template <int dim, typename number2, int spacedim>
90 void
92  const DoFHandler<dim,spacedim> &mg_dof_handler,
93  Vector<number2> &dst,
94  const MGLevelObject<Vector<number> > &src) const
95 {
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);
100 }
101 
102 
103 
104 template <typename number>
105 std::size_t
107 {
108  return sizeof(int) + MGTransferBlockBase::memory_consumption();
109 }
110 
111 
112 /* --------------------- MGTransferBlock -------------- */
113 
114 
115 
116 template <typename number>
117 template <int dim, typename number2, int spacedim>
118 void
120  const DoFHandler<dim,spacedim> &mg_dof_handler,
122  const MGLevelObject<BlockVector<number> > &src) const
123 {
124  for (unsigned int block=0; block<selected.size(); ++block)
125  if (selected[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);
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,
140  const MGLevelObject<BlockVector<number> > &src) const
141 {
142  for (unsigned int block=0; block<selected.size(); ++block)
143  if (selected[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);
148 }
149 
150 DEAL_II_NAMESPACE_CLOSE
151 
152 #endif
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
std::size_t memory_consumption() const
std::size_t memory_consumption() const
BlockType & block(const unsigned int i)
const Triangulation< dim, spacedim > & get_tria() const