From: grahambenharper Date: Tue, 9 Jul 2019 10:14:47 +0000 (-0600) Subject: Add two convenience functions to FE_PolyTensor and improve based on suggestions X-Git-Tag: v9.2.0-rc1~1381^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f10f8ee3141dac20f74da0128d8478418d5739ba;p=dealii.git Add two convenience functions to FE_PolyTensor and improve based on suggestions --- diff --git a/include/deal.II/fe/fe_poly_tensor.h b/include/deal.II/fe/fe_poly_tensor.h index 89034a6f45..4bb34d4dc8 100644 --- a/include/deal.II/fe/fe_poly_tensor.h +++ b/include/deal.II/fe/fe_poly_tensor.h @@ -216,6 +216,18 @@ protected: */ std::vector mapping_type; + /** + * Returns a boolean that is true when the finite element uses a single + * mapping and false when the finite element uses multiple mappings. + */ + bool + single_mapping() const; + + /** + * Returns MappingType @p i for the finite element. + */ + MappingType + get_mapping_type(unsigned int i) const; /* NOTE: The following function has its definition inlined into the class declaration because we otherwise run into a compiler error with MS Visual @@ -254,19 +266,17 @@ protected: // allocate memory const bool update_transformed_shape_values = - (std::find_if(this->mapping_type.begin(), - this->mapping_type.end(), - [](const MappingType t) { return t != mapping_none; }) != - this->mapping_type.end()); + std::any_of(this->mapping_type.begin(), + this->mapping_type.end(), + [](const MappingType t) { return t != mapping_none; }); const bool update_transformed_shape_grads = - (std::find_if(this->mapping_type.begin(), - this->mapping_type.end(), - [](const MappingType t) { - return (t == mapping_raviart_thomas || - t == mapping_piola || t == mapping_nedelec || - t == mapping_contravariant); - }) != this->mapping_type.end()); + std::any_of(this->mapping_type.begin(), + this->mapping_type.end(), + [](const MappingType t) { + return (t == mapping_raviart_thomas || t == mapping_piola || + t == mapping_nedelec || t == mapping_contravariant); + }); const bool update_transformed_shape_hessian_tensors = update_transformed_shape_values; diff --git a/source/fe/fe_poly_tensor.cc b/source/fe/fe_poly_tensor.cc index bc1d0844a3..28bc8e5a52 100644 --- a/source/fe/fe_poly_tensor.cc +++ b/source/fe/fe_poly_tensor.cc @@ -158,7 +158,8 @@ namespace internal // TODO: This is probably only going to work for those elements for // which all dofs are face dofs for (unsigned int l = 0; l < GeometryInfo::lines_per_cell; ++l) - if (!(cell->line_orientation(l)) & mapping_type[0] == mapping_nedelec) + if (!(cell->line_orientation(l)) && + mapping_type[0] == mapping_nedelec) face_sign[l] = -1.0; } } // namespace @@ -193,6 +194,28 @@ FE_PolyTensor::FE_PolyTensor( +template +bool +FE_PolyTensor::single_mapping() const +{ + return mapping_type.size() == 1; +} + + + +template +MappingType +FE_PolyTensor::get_mapping_type( + unsigned int i) const +{ + if (single_mapping()) + return mapping_type[0]; + Assert(i < mapping_type.size(), ExcIndexRange(i, 0, mapping_type.size())); + return mapping_type[i]; +} + + + template double FE_PolyTensor::shape_value( @@ -392,9 +415,7 @@ FE_PolyTensor::fill_fe_values( for (unsigned int i = 0; i < this->dofs_per_cell; ++i) { - MappingType mapping_type = - (this->mapping_type.size() > 1 ? this->mapping_type[i] : - this->mapping_type[0]); + const MappingType mapping_type = get_mapping_type(i); const unsigned int first = output_data.shape_function_to_row_table[i * this->n_components() + @@ -972,9 +993,7 @@ FE_PolyTensor::fill_fe_face_values( for (unsigned int i = 0; i < this->dofs_per_cell; ++i) { - const MappingType mapping_type = - (this->mapping_type.size() > 1 ? this->mapping_type[i] : - this->mapping_type[0]); + const MappingType mapping_type = get_mapping_type(i); const unsigned int first = output_data.shape_function_to_row_table[i * this->n_components() + @@ -1604,9 +1623,7 @@ FE_PolyTensor::fill_fe_subface_values( for (unsigned int i = 0; i < this->dofs_per_cell; ++i) { - MappingType mapping_type = - (this->mapping_type.size() > 1 ? this->mapping_type[i] : - this->mapping_type[0]); + const MappingType mapping_type = get_mapping_type(i); const unsigned int first = output_data.shape_function_to_row_table[i * this->n_components() + @@ -2170,9 +2187,7 @@ FE_PolyTensor::requires_update_flags( for (unsigned int i = 0; i < this->dofs_per_cell; ++i) { - MappingType mapping_type = - (this->mapping_type.size() > 1 ? this->mapping_type[i] : - this->mapping_type[0]); + const MappingType mapping_type = get_mapping_type(i); switch (mapping_type) {