17 #include <deal.II/meshworker/dof_info.h>
18 #include <deal.II/meshworker/integration_info.h>
19 #include <deal.II/base/quadrature_lib.h>
26 template<
int dim,
int sdim>
32 const unsigned int nqp = fevalv[0]->n_quadrature_points;
34 values.resize(global_data->n_values());
38 for (
unsigned int i=0; i<values.size(); ++i)
40 values[i].resize(n_components);
43 for (
unsigned int j=0; j<values[i].size(); ++j)
45 values[i][j].resize(nqp);
52 gradients.resize(global_data->n_gradients());
55 for (
unsigned int i=0; i<gradients.size(); ++i)
57 gradients[i].resize(n_components);
59 for (
unsigned int j=0; j<gradients[i].size(); ++j)
61 gradients[i][j].resize(nqp);
65 hessians.resize(global_data->n_hessians());
68 for (
unsigned int i=0; i<hessians.size(); ++i)
70 hessians[i].resize(n_components);
72 for (
unsigned int j=0; j<hessians[i].size(); ++j)
74 hessians[i][j].resize(nqp);
80 template<
int dim,
int sdim>
89 template<
int dim,
int sdim>
90 template <
typename number>
96 unsigned int comp = 0;
98 for (
unsigned int b=0; b<info.
block_info->local().size(); ++b)
100 const unsigned int fe_no = info.
block_info->base_element(b);
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);
107 this->global_data->mg_fill(values, gradients, hessians, fe, info.
cell->level(), info.indices,
108 comp, n_comp, block_start, block_size);
110 this->global_data->fill(values, gradients, hessians, fe, info.indices,
111 comp, n_comp, block_start, block_size);
118 const unsigned int n_comp = fe.
get_fe().n_components();
120 this->global_data->mg_fill(values, gradients, hessians, fe, info.
cell->level(), info.indices,
121 0, n_comp, 0, info.indices.size());
123 this->global_data->fill(values, gradients, hessians, fe, info.indices,
124 0, n_comp, 0, info.indices.size());
129 template<
int dim,
int sdim>
133 std::size_t mem =
sizeof(*this)
136 for (
unsigned int i=0; i<fevalv.size(); ++i)
137 mem += fevalv[i]->memory_consumption();
143 template<
int dim,
int sdim>
153 template<
int dim,
int sdim>
159 face_flags |= boundary_flags;
160 neighbor_flags |= neighbor_geometry
164 if (cell_selector.has_values() != 0) cell_flags |=
update_values;
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;
172 if (face_selector.has_values() != 0) face_flags |=
update_values;
176 if (face_selector.has_values() != 0) neighbor_flags |=
update_values;
178 if (face_selector.has_hessians() != 0) neighbor_flags |=
update_hessians;
182 template <
int dim,
int sdim>
191 if (cell) cell_flags |= flags;
192 if (boundary) boundary_flags |= flags;
193 if (face) face_flags |= flags;
194 if (neighbor) neighbor_flags |= flags;
198 template<
int dim,
int sdim>
202 std::size_t mem =
sizeof(*this)
204 -
sizeof (cell_quadrature)
206 -
sizeof (boundary_quadrature)
208 -
sizeof (face_quadrature)
210 -
sizeof (cell_selector)
212 -
sizeof (boundary_selector)
214 -
sizeof (face_selector)
237 DEAL_II_NAMESPACE_CLOSE
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)
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
std::size_t memory_consumption() const
std::size_t memory_consumption() const
std::size_t memory_consumption(const T &t)
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
void initialize_update_flags(bool neighbor_geometry=false)
void initialize_data(const std_cxx1x::shared_ptr< VectorDataBase< dim, spacedim > > &data)
Second derivatives of shape functions.
Shape function gradients.
const FiniteElement< dim, spacedim > & get_fe() const
void fill_local_data(const DoFInfo< dim, spacedim, number > &info, bool split_fevalues)