Reference documentation for deal.II version 8.1.0
vector_selector.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: vector_selector.templates.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2009 - 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 #include <deal.II/meshworker/vector_selector.h>
19 #include <deal.II/base/vector_slice.h>
20 #include <deal.II/fe/fe_values.h>
21 
23 
24 namespace MeshWorker
25 {
26  template <int dim, int spacedim>
28  {}
29 
30 
31  template <int dim, int spacedim>
33  :
35  {}
36 
37 
38  template <int dim, int spacedim>
40  {}
41 
42 
43  template <int dim, int spacedim>
44  void
46  std::vector<std::vector<std::vector<double> > > &,
47  std::vector<std::vector<std::vector<Tensor<1,dim> > > > &,
48  std::vector<std::vector<std::vector<Tensor<2,dim> > > > &,
50  const std::vector<types::global_dof_index> &,
51  const unsigned int,
52  const unsigned int,
53  const unsigned int,
54  const unsigned int) const
55  {
56  Assert(false, ExcNotImplemented());
57  }
58 
59 
60  template <int dim, int spacedim>
61  void
63  std::vector<std::vector<std::vector<double> > > &,
64  std::vector<std::vector<std::vector<Tensor<1,dim> > > > &,
65  std::vector<std::vector<std::vector<Tensor<2,dim> > > > &,
67  const unsigned int,
68  const std::vector<types::global_dof_index> &,
69  const unsigned int,
70  const unsigned int,
71  const unsigned int,
72  const unsigned int) const
73  {
74  Assert(false, ExcNotImplemented());
75  }
76 
77 
78 //----------------------------------------------------------------------//
79 
80  template <class VECTOR, int dim, int spacedim>
82  {}
83 
84 
85  template <class VECTOR, int dim, int spacedim>
87  :
88  VectorDataBase<dim, spacedim>(s)
89  {}
90 
91 
92  template <class VECTOR, int dim, int spacedim>
93  void
95  {
96  data = d;
98  }
99 
100 
101  template <class VECTOR, int dim, int spacedim>
102  void
103  VectorData<VECTOR, dim, spacedim>::initialize(const VECTOR *v, const std::string &name)
104  {
106  data.add(p, name);
108  }
109 
110 
111  template <class VECTOR, int dim, int spacedim>
112  void
114  std::vector<std::vector<std::vector<double> > > &values,
115  std::vector<std::vector<std::vector<Tensor<1,dim> > > > &gradients,
116  std::vector<std::vector<std::vector<Tensor<2,dim> > > > &hessians,
117  const FEValuesBase<dim,spacedim> &fe,
118  const std::vector<types::global_dof_index> &index,
119  const unsigned int component,
120  const unsigned int n_comp,
121  const unsigned int start,
122  const unsigned int size) const
123  {
124  AssertDimension(values.size(), this->n_values());
125  AssertDimension(gradients.size(), this->n_gradients());
126  AssertDimension(hessians.size(), this->n_hessians());
127 
128  for (unsigned int i=0; i<this->n_values(); ++i)
129  {
130  const VECTOR &src = *data(this->value_index(i));
131  VectorSlice<std::vector<std::vector<double> > > dst(values[i], component, n_comp);
132  fe.get_function_values(src, make_slice(index, start, size), dst, true);
133  }
134 
135  for (unsigned int i=0; i<this->n_gradients(); ++i)
136  {
137  const VECTOR &src = *data(this->gradient_index(i));
138  VectorSlice<std::vector<std::vector<Tensor<1,dim> > > > dst(gradients[i], component, n_comp);
139  fe.get_function_gradients(src, make_slice(index, start, size), dst, true);
140  }
141 
142  for (unsigned int i=0; i<this->n_hessians(); ++i)
143  {
144  const VECTOR &src = *data(this->hessian_index(i));
145  VectorSlice<std::vector<std::vector<Tensor<2,dim> > > > dst(hessians[i], component, n_comp);
146  fe.get_function_hessians(src, make_slice(index, start, size), dst, true);
147  }
148  }
149 
150 
151  template <class VECTOR, int dim, int spacedim>
152  std::size_t
154  {
155  std::size_t mem = VectorSelector::memory_consumption();
156  mem += sizeof (data);
157  return mem;
158  }
159 
160 //----------------------------------------------------------------------//
161 
162  template <class VECTOR, int dim, int spacedim>
164  {}
165 
166 
167  template <class VECTOR, int dim, int spacedim>
169  :
170  VectorDataBase<dim, spacedim>(s)
171  {}
172 
173 
174  template <class VECTOR, int dim, int spacedim>
175  void
177  {
178  data = d;
180  }
181 
182 
183  template <class VECTOR, int dim, int spacedim>
184  void
186  {
188  p = v;
189  data.add(p, name);
191  }
192 
193 
194  template <class VECTOR, int dim, int spacedim>
195  void
197  std::vector<std::vector<std::vector<double> > > &values,
198  std::vector<std::vector<std::vector<Tensor<1,dim> > > > &gradients,
199  std::vector<std::vector<std::vector<Tensor<2,dim> > > > &hessians,
200  const FEValuesBase<dim,spacedim> &fe,
201  const unsigned int level,
202  const std::vector<types::global_dof_index> &index,
203  const unsigned int component,
204  const unsigned int n_comp,
205  const unsigned int start,
206  const unsigned int size) const
207  {
208  AssertDimension(values.size(), this->n_values());
209  AssertDimension(gradients.size(), this->n_gradients());
210  AssertDimension(hessians.size(), this->n_hessians());
211 
212  for (unsigned int i=0; i<this->n_values(); ++i)
213  {
214  const VECTOR &src = (*data(this->value_index(i)))[level];
215  VectorSlice<std::vector<std::vector<double> > > dst(values[i], component, n_comp);
216  fe.get_function_values(src, make_slice(index, start, size), dst, true);
217  }
218 
219  for (unsigned int i=0; i<this->n_gradients(); ++i)
220  {
221  const VECTOR &src = (*data(this->value_index(i)))[level];
222  VectorSlice<std::vector<std::vector<Tensor<1,dim> > > > dst(gradients[i], component, n_comp);
223  fe.get_function_gradients(src, make_slice(index, start, size), dst, true);
224  }
225 
226  for (unsigned int i=0; i<this->n_hessians(); ++i)
227  {
228  const VECTOR &src = (*data(this->value_index(i)))[level];
229  VectorSlice<std::vector<std::vector<Tensor<2,dim> > > > dst(hessians[i], component, n_comp);
230  fe.get_function_hessians(src, make_slice(index, start, size), dst, true);
231  }
232  }
233 
234 
235  template <class VECTOR, int dim, int spacedim>
236  std::size_t
238  {
239  std::size_t mem = VectorSelector::memory_consumption();
240  mem += sizeof (data);
241  return mem;
242  }
243 }
244 
245 DEAL_II_NAMESPACE_CLOSE
void initialize(const NamedData< VECTOR * > &)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
std::size_t memory_consumption() const
void initialize(const NamedData< DATA > &)
virtual void fill(std::vector< std::vector< std::vector< double > > > &values, std::vector< std::vector< std::vector< Tensor< 1, dim > > > > &gradients, std::vector< std::vector< std::vector< Tensor< 2, dim > > > > &hessians, const FEValuesBase< dim, spacedim > &fe, const std::vector< types::global_dof_index > &index, const unsigned int component, const unsigned int n_comp, const unsigned int start, const unsigned int size) const
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim > > &hessians) const
void initialize(const NamedData< MGLevelObject< VECTOR > * > &)
void get_function_values(const InputVector &fe_function, std::vector< number > &values) const
#define Assert(cond, exc)
Definition: exceptions.h:299
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim > > &gradients) const
std::size_t memory_consumption() const
std::size_t memory_consumption() const
virtual void mg_fill(std::vector< std::vector< std::vector< double > > > &values, std::vector< std::vector< std::vector< Tensor< 1, dim > > > > &gradients, std::vector< std::vector< std::vector< Tensor< 2, dim > > > > &hessians, const FEValuesBase< dim, spacedim > &fe, const unsigned int level, const std::vector< types::global_dof_index > &index, const unsigned int component, const unsigned int n_comp, const unsigned int start, const unsigned int size) const
virtual void fill(std::vector< std::vector< std::vector< double > > > &values, std::vector< std::vector< std::vector< Tensor< 1, dim > > > > &gradients, std::vector< std::vector< std::vector< Tensor< 2, dim > > > > &hessians, const FEValuesBase< dim, spacedim > &fe, const std::vector< types::global_dof_index > &index, const unsigned int component, const unsigned int n_comp, const unsigned int start, const unsigned int size) const
::ExceptionBase & ExcNotImplemented()
virtual void mg_fill(std::vector< std::vector< std::vector< double > > > &values, std::vector< std::vector< std::vector< Tensor< 1, dim > > > > &gradients, std::vector< std::vector< std::vector< Tensor< 2, dim > > > > &hessians, const FEValuesBase< dim, spacedim > &fe, const unsigned int level, const std::vector< types::global_dof_index > &index, const unsigned int component, const unsigned int n_comp, const unsigned int start, const unsigned int size) const