Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
integration_info.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: integration_info.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 #include <deal.II/meshworker/dof_info.h>
18 #include <deal.II/meshworker/integration_info.h>
19 #include <deal.II/base/quadrature_lib.h>
20 
22 
23 
24 namespace MeshWorker
25 {
26  template<int dim, int sdim>
27  void
29  const std_cxx1x::shared_ptr<VectorDataBase<dim,sdim> > &data)
30  {
31  global_data = data;
32  const unsigned int nqp = fevalv[0]->n_quadrature_points;
33 
34  values.resize(global_data->n_values());
35 // deallog << "values: " << values.size() << " [";
36  // For all selected finite
37  // element functions
38  for (unsigned int i=0; i<values.size(); ++i)
39  {
40  values[i].resize(n_components);
41 // deallog << ' ' << values[i].size() << " {";
42  // For all components
43  for (unsigned int j=0; j<values[i].size(); ++j)
44  {
45  values[i][j].resize(nqp);
46 // deallog << ' ' << values[i][j].size();
47  }
48 // deallog << " }";
49  }
50 // deallog << " ]" << std::endl;
51 
52  gradients.resize(global_data->n_gradients());
53  // For all selected finite
54  // element functions
55  for (unsigned int i=0; i<gradients.size(); ++i)
56  {
57  gradients[i].resize(n_components);
58  // For all components
59  for (unsigned int j=0; j<gradients[i].size(); ++j)
60  {
61  gradients[i][j].resize(nqp);
62  }
63  }
64 
65  hessians.resize(global_data->n_hessians());
66  // For all selected finite
67  // element functions
68  for (unsigned int i=0; i<hessians.size(); ++i)
69  {
70  hessians[i].resize(n_components);
71  // For all components
72  for (unsigned int j=0; j<hessians[i].size(); ++j)
73  {
74  hessians[i][j].resize(nqp);
75  }
76  }
77  }
78 
79 
80  template<int dim, int sdim>
81  void
83  {
84  fevalv.resize(0);
85  }
86 
87 
88 
89  template<int dim, int sdim>
90  template <typename number>
91  void
93  {
94  if (split_fevalues)
95  {
96  unsigned int comp = 0;
97  // Loop over all blocks
98  for (unsigned int b=0; b<info.block_info->local().size(); ++b)
99  {
100  const unsigned int fe_no = info.block_info->base_element(b);
101  const FEValuesBase<dim,sdim> &fe = this->fe_values(fe_no);
102  const unsigned int n_comp = fe.get_fe().n_components();
103  const unsigned int block_start = info.block_info->local().block_start(b);
104  const unsigned int block_size = info.block_info->local().block_size(b);
105 
106  if (info.level_cell)
107  this->global_data->mg_fill(values, gradients, hessians, fe, info.cell->level(), info.indices,
108  comp, n_comp, block_start, block_size);
109  else
110  this->global_data->fill(values, gradients, hessians, fe, info.indices,
111  comp, n_comp, block_start, block_size);
112  comp += n_comp;
113  }
114  }
115  else
116  {
117  const FEValuesBase<dim,sdim> &fe = this->fe_values(0);
118  const unsigned int n_comp = fe.get_fe().n_components();
119  if (info.level_cell)
120  this->global_data->mg_fill(values, gradients, hessians, fe, info.cell->level(), info.indices,
121  0, n_comp, 0, info.indices.size());
122  else
123  this->global_data->fill(values, gradients, hessians, fe, info.indices,
124  0, n_comp, 0, info.indices.size());
125  }
126  }
127 
128 
129  template<int dim, int sdim>
130  std::size_t
132  {
133  std::size_t mem = sizeof(*this)
135  - sizeof (fevalv);
136  for (unsigned int i=0; i<fevalv.size(); ++i)
137  mem += fevalv[i]->memory_consumption();
138  return mem;
139  }
140 
141 //----------------------------------------------------------------------//
142 
143  template<int dim, int sdim>
145  {
146  cell_flags = update_default;
147  boundary_flags = update_default;
148  face_flags = update_default;
149  neighbor_flags = update_default;
150  }
151 
152 
153  template<int dim, int sdim>
154  void
156  {
157  cell_flags |= update_JxW_values;
159  face_flags |= boundary_flags;
160  neighbor_flags |= neighbor_geometry
161  ? boundary_flags
162  : update_default;
163 
164  if (cell_selector.has_values() != 0) cell_flags |= update_values;
165  if (cell_selector.has_gradients() != 0) cell_flags |= update_gradients;
166  if (cell_selector.has_hessians() != 0) cell_flags |= update_hessians;
167 
168  if (boundary_selector.has_values() != 0) boundary_flags |= update_values;
169  if (boundary_selector.has_gradients() != 0) boundary_flags |= update_gradients;
170  if (boundary_selector.has_hessians() != 0) boundary_flags |= update_hessians;
171 
172  if (face_selector.has_values() != 0) face_flags |= update_values;
173  if (face_selector.has_gradients() != 0) face_flags |= update_gradients;
174  if (face_selector.has_hessians() != 0) face_flags |= update_hessians;
175 
176  if (face_selector.has_values() != 0) neighbor_flags |= update_values;
177  if (face_selector.has_gradients() != 0) neighbor_flags |= update_gradients;
178  if (face_selector.has_hessians() != 0) neighbor_flags |= update_hessians;
179  }
180 
181 
182  template <int dim, int sdim>
183  void
185  const UpdateFlags flags,
186  bool cell,
187  bool boundary,
188  bool face,
189  bool neighbor)
190  {
191  if (cell) cell_flags |= flags;
192  if (boundary) boundary_flags |= flags;
193  if (face) face_flags |= flags;
194  if (neighbor) neighbor_flags |= flags;
195  }
196 
197 
198  template<int dim, int sdim>
199  std::size_t
201  {
202  std::size_t mem = sizeof(*this)
203  + MemoryConsumption::memory_consumption(cell_quadrature)
204  - sizeof (cell_quadrature)
205  + MemoryConsumption::memory_consumption(boundary_quadrature)
206  - sizeof (boundary_quadrature)
207  + MemoryConsumption::memory_consumption(face_quadrature)
208  - sizeof (face_quadrature)
210  -sizeof (cell_selector)
211  + MemoryConsumption::memory_consumption(boundary_selector)
212  -sizeof (boundary_selector)
214  -sizeof (face_selector)
216  - sizeof(cell)
218  - sizeof(boundary)
220  - sizeof(face)
222  - sizeof(subface)
224  - sizeof(neighbor);
225 // if (cell_data != 0)
226 // mem += MemoryConsumption::memory_consumption(*cell_data);
227 // if (boundary_data != 0)
228 // mem += MemoryConsumption::memory_consumption(*boundary_data);
229 // if (face_data != 0)
230 // mem += MemoryConsumption::memory_consumption(*face_data);
231 
232  return mem;
233  }
234 }
235 
236 
237 DEAL_II_NAMESPACE_CLOSE
238 
Transformed quadrature weights.
void add_update_flags(const UpdateFlags flags, const bool cell=true, const bool boundary=true, const bool face=true, const bool neighbor=true)
Shape function values.
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:173
void initialize_data(const std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > &data)
No update.
UpdateFlags
std::size_t memory_consumption(const T &t)
void initialize_update_flags(bool neighbor_geometry=false)
Second derivatives of shape functions.
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:73
const FiniteElement< dim, spacedim > & get_fe() const
Shape function gradients.
Normal vectors.