40 template<
unsigned int spacedim>
50 template<
unsigned int spacedim>
53 unsigned int first_component_idx,
57 for (
unsigned int ip=0; ip<
n_points; ip++)
58 for (
unsigned int id=0;
id<dof_indices.size();
id++)
67 template<
unsigned int spacedim>
68 template<
unsigned int DIM>
88 ASSERT(fe_sys !=
nullptr).error(
"Mixed system must be represented by FESystem.");
93 unsigned int comp_offset = 0;
94 for (
auto fe : fe_sys->
fe())
110 ASSERT(
false).error(
"Not implemented.");
122 template<
unsigned int spacedim>
130 template<
unsigned int spacedim>
135 template<
unsigned int spacedim>
136 template<
unsigned int DIM>
152 ASSERT(fe !=
nullptr).error(
"Mixed system must be represented by FESystem.");
155 for (
unsigned int f=0; f<fe->
fe().size(); f++)
160 if ( q.
dim() == DIM )
164 else if ( q.
dim() + 1 == DIM )
167 for (
unsigned int sid = 0; sid < RefElement<DIM>::n_sides; sid++)
173 ASSERT(
false)(q.
dim())(DIM).error(
"Dimension mismatch in FEValues::initialize().");
178 template<
unsigned int spacedim>
179 template<
unsigned int DIM>
188 ASSERT(_fe.
n_components() == spacedim).error(
"FEVector must have spacedim components.");
190 ASSERT(_fe.
n_components() == spacedim*spacedim).error(
"FETensor must have spacedim*spacedim components.");
203 if (fe_sys !=
nullptr)
205 for (
unsigned int f=0; f<fe_sys->
fe().size(); f++)
227 template<
unsigned int spacedim>
228 template<
unsigned int DIM>
233 std::shared_ptr<FEInternalData> data = std::make_shared<FEInternalData>(q.
size(),
n_dofs_);
236 for (
unsigned int i=0; i<q.
size(); i++)
238 for (
unsigned int j=0; j<
n_dofs_; j++)
243 data->ref_shape_values[i][j] = trans(
shape_values.row(j));
248 for (
unsigned int i=0; i<q.
size(); i++)
250 for (
unsigned int j=0; j<
n_dofs_; j++)
256 data->ref_shape_grads[i][j] = grad;
276 template<
unsigned int spacedim>
278 const unsigned int point_no,
279 const unsigned int comp)
const
493 template<
unsigned int spacedim>
522 template<
unsigned int spacedim>
523 template<
class MapType>
536 template<
unsigned int spacedim>
551 template<
unsigned int spacedim>
576 fv[0].initialize(quadrature[0], *fe[0_d], flags);
577 fv[1].initialize(quadrature[1], *fe[1_d], flags);
578 fv[2].initialize(quadrature[2], *fe[2_d], flags);
579 fv[3].initialize(quadrature[3], *fe[3_d], flags);
#define ASSERT_PERMANENT(expr)
Allow use shorter versions of macro names if these names is not used with external library.
#define ASSERT_LT(a, b)
Definition of comparative assert macro (Less Than) only for debug mode.
#define ASSERT_EQ(a, b)
Definition of comparative assert macro (EQual) only for debug mode.
Class for computation of data on cell and side.
Compound finite element on dim dimensional simplex.
const std::vector< std::shared_ptr< FiniteElement< dim > > > & fe() const
std::vector< unsigned int > fe_dofs(unsigned int fe_index)
Return dof indices belonging to given sub-FE.
Structure for storing the precomputed finite element data.
std::vector< std::vector< arma::vec > > ref_shape_values
Precomputed values of basis functions at the quadrature points.
FEInternalData(unsigned int np, unsigned int nd)
std::vector< std::vector< arma::mat > > ref_shape_grads
Precomputed gradients of basis functions at the quadrature points.
unsigned int n_points
Number of quadrature points.
Calculates finite element data on the actual cell.
~FEValues()
Correct deallocation of objects created by 'initialize' methods.
std::shared_ptr< ElementValues< spacedim > > elm_values
Auxiliary object for calculation of element-dependent data.
std::vector< unsigned int > fe_sys_n_space_components_
Numbers of components of FESystem sub-elements in real space.
std::shared_ptr< FEInternalData > init_fe_data(const FiniteElement< DIM > &fe, const Quadrature &q)
Precompute finite element data on reference element.
std::vector< std::vector< double > > shape_values
Shape functions evaluated at the quadrature points.
FEValues()
Default constructor with postponed initialization.
arma::vec::fixed< spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int comp) const
Return the gradient of the function_no-th shape function at the point_no-th quadrature point.
FEType fe_type_
Type of finite element (scalar, vector, tensor).
std::vector< unsigned int > fe_sys_n_components_
Numbers of components of FESystem sub-elements in reference space.
std::vector< shared_ptr< FEInternalData > > side_fe_data
Precomputed FE data (shape functions on reference element) for all side quadrature points.
ViewsCache views_cache_
Auxiliary storage of FEValuesViews accessors.
unsigned int n_dofs_
Number of finite element dofs.
std::vector< std::vector< unsigned int > > fe_sys_dofs_
Dof indices of FESystem sub-elements.
unsigned int dim_
Dimension of reference space.
void initialize(Quadrature &_quadrature, FiniteElement< DIM > &_fe, UpdateFlags _flags)
Initialize structures and calculates cell-independent data.
std::vector< std::vector< arma::vec::fixed< spacedim > > > shape_gradients
void reinit(const ElementAccessor< spacedim > &cell)
Update cell-dependent data (gradients, Jacobians etc.)
unsigned int n_points() const
Returns the number of quadrature points.
UpdateFlags update_flags
Flags that indicate which finite element quantities are to be computed.
void fill_data(const ElementValues< spacedim > &elm_values, const FEInternalData &fe_data)
Computes the shape function values and gradients on the actual cell and fills the FEValues structure.
unsigned int n_components_
Number of components of the FE.
unsigned int n_points_
Number of integration points.
void fill_data_specialized(const ElementValues< spacedim > &elm_values, const FEInternalData &fe_data)
Computes the shape function values and gradients on the actual cell and fills the FEValues structure....
std::shared_ptr< FEInternalData > fe_data
Precomputed finite element data.
void allocate(unsigned int n_points, FiniteElement< DIM > &_fe, UpdateFlags flags)
Allocates space for computed data.
std::vector< FEValues< spacedim > > fe_values_vec
Vector of FEValues for sub-elements of FESystem.
Abstract class for the description of a general finite element on a reference simplex in dim dimensio...
unsigned int n_space_components(unsigned int spacedim)
Number of components of FE in a mapped space with dimension spacedim.
virtual UpdateFlags update_each(UpdateFlags flags)
Decides which additional quantities have to be computed for each cell.
unsigned int n_dofs() const
Returns the number of degrees of freedom needed by the finite element.
FEType type_
Type of FiniteElement.
double shape_value(const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
Calculates the value of the comp-th component of the i-th shape function at the point p on the refere...
arma::vec::fixed< dim > shape_grad(const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
Calculates the comp-th component of the gradient of the i-th shape function at the point p on the ref...
unsigned int n_components() const
Returns numer of components of the basis function.
Helper class allows update values and gradients of FEValues of FEVectorContravariant type.
Helper class allows update values and gradients of FEValues of FEVectorPiola type.
Helper class allows update values and gradients of FEValues of FEScalar type.
Helper class allows update values and gradients of FEValues of FEMixedSystem type.
Helper class allows update values and gradients of FEValues of FETensor type.
Helper class allows update values and gradients of FEValues of FEVector type.
static UpdateFlags update_each(UpdateFlags flags)
Determines which additional quantities have to be computed.
std::array< QGauss, 4 > array
Base class for quadrature rules on simplices in arbitrary dimensions.
Quadrature make_from_side(unsigned int sid) const
unsigned int size() const
Returns number of quadrature points.
Armor::ArmaVec< double, point_dim > point(unsigned int i) const
Returns the ith quadrature point.
unsigned int side_idx() const
Returns local index of the side on the element.
unsigned int dim() const
Returns dimension of the side, that is dimension of the element minus one.
Class ElementValues calculates data related to transformation of reference cell to actual cell (Jacob...
Class FESystem for compound finite elements.
std::vector< FEValues< 3 > > mixed_fe_values(QGauss::array &quadrature, MixedPtr< FiniteElement > fe, UpdateFlags flags)
Class FEValues calculates finite element data on the actual cells such as shape function values,...
Abstract class for description of finite elements.
int LongIdx
Define type that represents indices of large arrays (elements, nodes, dofs etc.)
Class MappingP1 implements the affine transformation of the unit cell onto the actual cell.
ArmaMat< double, N, M > mat
Basic definitions of numerical quadrature rules.
void initialize(const FEValues &fv, const FiniteElement< DIM > &fe)
UpdateFlags
Enum type UpdateFlags indicates which quantities are to be recomputed on each finite element cell.
@ update_values
Shape function values.
@ update_gradients
Shape function gradients.