if (cell->is_locally_owned())
{
const unsigned int fe_index = cell->active_fe_index();
- if (fe[fe_index].degree != 0)
- {
- // for each cell:
- // get location of finite element
- // support_points
- fe_values.reinit(cell);
- const std::vector<Point<spacedim> > &support_points =
- fe_values.get_present_fe_values().get_quadrature_points();
-
- // pick out the representative
- // support points
- rep_points.resize (dofs_of_rep_points[fe_index].size());
- for (unsigned int j=0; j<dofs_of_rep_points[fe_index].size(); ++j)
- rep_points[j] = support_points[dofs_of_rep_points[fe_index][j]];
-
- // get indices of the dofs on this cell
- dofs_on_cell.resize (fe[fe_index].dofs_per_cell);
- cell->get_dof_indices (dofs_on_cell);
-
-
- if (fe_is_system)
- {
- // get function values at
- // these points. Here: get
- // all components
- function_values_system[fe_index]
- .resize (n_rep_points[fe_index],
- Vector<double> (fe[fe_index].n_components()));
- function.vector_value_list (rep_points,
- function_values_system[fe_index]);
- // distribute the function
- // values to the global
- // vector
- for (unsigned int i=0; i<fe[fe_index].dofs_per_cell; ++i)
- {
- const unsigned int component
- = fe[fe_index].system_to_component_index(i).first;
- const unsigned int rep_dof=dof_to_rep_index_table[fe_index][i];
- vec(dofs_on_cell[i])
- = function_values_system[fe_index][rep_dof](component);
- }
- }
- else
- {
- // get first component only,
- // which is the only component
- // in the function anyway
- function_values_scalar[fe_index].resize (n_rep_points[fe_index]);
- function.value_list (rep_points,
- function_values_scalar[fe_index],
- 0);
- // distribute the function
- // values to the global
- // vector
- for (unsigned int i=0; i<fe[fe_index].dofs_per_cell; ++i)
- vec(dofs_on_cell[i])
- = function_values_scalar[fe_index][dof_to_rep_index_table[fe_index][i]];
- }
- }
- }
+ if (fe[fe_index].degree != 0)
+ {
+ // for each cell:
+ // get location of finite element
+ // support_points
+ fe_values.reinit(cell);
+ const std::vector<Point<spacedim> > &support_points =
+ fe_values.get_present_fe_values().get_quadrature_points();
+
+ // pick out the representative
+ // support points
+ rep_points.resize (dofs_of_rep_points[fe_index].size());
+ for (unsigned int j=0; j<dofs_of_rep_points[fe_index].size(); ++j)
+ rep_points[j] = support_points[dofs_of_rep_points[fe_index][j]];
+
+ // get indices of the dofs on this cell
+ dofs_on_cell.resize (fe[fe_index].dofs_per_cell);
+ cell->get_dof_indices (dofs_on_cell);
+
+
+ if (fe_is_system)
+ {
+ // get function values at
+ // these points. Here: get
+ // all components
+ function_values_system[fe_index]
+ .resize (n_rep_points[fe_index],
+ Vector<double> (fe[fe_index].n_components()));
+ function.vector_value_list (rep_points,
+ function_values_system[fe_index]);
+ // distribute the function
+ // values to the global
+ // vector
+ for (unsigned int i=0; i<fe[fe_index].dofs_per_cell; ++i)
+ {
+ const unsigned int component
+ = fe[fe_index].system_to_component_index(i).first;
+ const unsigned int rep_dof=dof_to_rep_index_table[fe_index][i];
+ vec(dofs_on_cell[i])
+ = function_values_system[fe_index][rep_dof](component);
+ }
+ }
+ else
+ {
+ // get first component only,
+ // which is the only component
+ // in the function anyway
+ function_values_scalar[fe_index].resize (n_rep_points[fe_index]);
+ function.value_list (rep_points,
+ function_values_scalar[fe_index],
+ 0);
+ // distribute the function
+ // values to the global
+ // vector
+ for (unsigned int i=0; i<fe[fe_index].dofs_per_cell; ++i)
+ vec(dofs_on_cell[i])
+ = function_values_scalar[fe_index][dof_to_rep_index_table[fe_index][i]];
+ }
+ }
+ }
vec.compress(VectorOperation::insert);
}
template <int dim, int spacedim,
- template <int,int> class DH,
- class Vector>
+ template <int,int> class DH,
+ class Vector>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
const Vector &u1,
template <int dim, int spacedim,
- template <int,int> class DH,
- class Vector>
+ template <int,int> class DH,
+ class Vector>
void
interpolate_to_different_mesh (const DH<dim, spacedim> &dof1,
const Vector &u1,
template <int dim, int spacedim,
- template <int,int> class DH,
- class Vector>
+ template <int,int> class DH,
+ class Vector>
void
interpolate_to_different_mesh (const InterGridMap<DH<dim, spacedim> > &intergridmap,
const Vector &u1,
* Compute the boundary values to be used in the project() functions.
*/
template <int dim, int spacedim,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
+ template <int,int> class DH,
+ template <int,int> class M_or_MC,
+ template <int> class Q_or_QC>
void project_compute_b_v (const M_or_MC<dim, spacedim> &mapping,
const DH<dim,spacedim> &dof,
const Function<spacedim> &function,
* Generic implementation of the project() function
*/
template <int dim, int spacedim,
- class Vector,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
+ class Vector,
+ template <int,int> class DH,
+ template <int,int> class M_or_MC,
+ template <int> class Q_or_QC>
void do_project (const M_or_MC<dim, spacedim> &mapping,
const DH<dim,spacedim> &dof,
const ConstraintMatrix &constraints,
// easier to simply work on
// individual vertices
template <class DH,
- template <int,int> class M_or_MC>
+ template <int,int> class M_or_MC>
static inline
void do_interpolate_boundary_values (const M_or_MC<DH::dimension, DH::space_dimension> &,
const DH &dof,
// whenever possible (i.e., if dim==1), the function template above
// will be used
template <class DH,
- template <int,int> class M_or_MC,
- int dim_>
+ template <int,int> class M_or_MC,
+ int dim_>
static inline
void
do_interpolate_boundary_values (const M_or_MC<DH::dimension, DH::space_dimension> &mapping,
:
(dim == 2 ?
(i<2*fe.dofs_per_vertex ? i : i+2*fe.dofs_per_vertex)
- :
- (dim == 3 ?
- (i<4*fe.dofs_per_vertex ?
- i
:
- (i<4*fe.dofs_per_vertex+4*fe.dofs_per_line ?
- i+4*fe.dofs_per_vertex
+ (dim == 3 ?
+ (i<4*fe.dofs_per_vertex ?
+ i
+ :
+ (i<4*fe.dofs_per_vertex+4*fe.dofs_per_line ?
+ i+4*fe.dofs_per_vertex
+ :
+ i+4*fe.dofs_per_vertex+8*fe.dofs_per_line))
:
- i+4*fe.dofs_per_vertex+8*fe.dofs_per_line))
- :
- numbers::invalid_unsigned_int)));
+ numbers::invalid_unsigned_int)));
Assert (cell_i < fe.dofs_per_cell, ExcInternalError());
// make sure that if this is not a primitive
}
else
// fe has only one component, so save some computations
- {
+ {
// get only the one component that this function has
dof_values_scalar.resize (fe.dofs_per_face);
function_map.find(boundary_component)->second
namespace
{
template <int dim, int spacedim,
- template <int,int> class DH,
- template <int,int> class M_or_MC,
- template <int> class Q_or_QC>
+ template <int,int> class DH,
+ template <int,int> class M_or_MC,
+ template <int> class Q_or_QC>
void
do_project_boundary_values (const M_or_MC<dim, spacedim> &mapping,
const DH<dim, spacedim> &dof,
switch (max_n_contributions_per_cell)
{
- // first deal with the case that a number of cells all have
- // registered that they have a normal vector defined at the
- // location of a given vector dof, and that each of them have
- // encountered this vector dof exactly once while looping over all
- // their faces. as stated in the documentation, this is the case
- // where we want to simply average over all normal vectors
- //
- // the typical case is in 2d where multiple cells meet at one
- // vertex sitting on the boundary. same in 3d for a vertex that
- // is associated with only one of the boundary indicators passed
- // to this function
+ // first deal with the case that a number of cells all have
+ // registered that they have a normal vector defined at the
+ // location of a given vector dof, and that each of them have
+ // encountered this vector dof exactly once while looping over all
+ // their faces. as stated in the documentation, this is the case
+ // where we want to simply average over all normal vectors
+ //
+ // the typical case is in 2d where multiple cells meet at one
+ // vertex sitting on the boundary. same in 3d for a vertex that
+ // is associated with only one of the boundary indicators passed
+ // to this function
case 1:
{
// compute the average normal vector from all the ones that have
update_flags |= UpdateFlags (update_gradients);
if (spacedim == dim+1)
update_flags |= UpdateFlags (update_normal_vectors);
- // no break!
+ // no break!
default:
update_flags |= UpdateFlags (update_values);