From: leicht Date: Thu, 18 Jan 2007 10:43:46 +0000 (+0000) Subject: Reimplement the changes from patch 14330: Do not permute the shape functions on faces... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b812dc366f8a2658a6fabd5be2233e7348eb9616;p=dealii-svn.git Reimplement the changes from patch 14330: Do not permute the shape functions on faces in 3D, but reorder the corresponding dofs in cell->get_dof_indices(). git-svn-id: https://svn.dealii.org/trunk@14339 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_accessor.templates.h b/deal.II/deal.II/include/dofs/dof_accessor.templates.h index 4a732ec914..124e2f9065 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.templates.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.templates.h @@ -589,8 +589,13 @@ DoFObjectAccessor<3,DH>::get_dof_indices (std::vector &dof_indices for (unsigned int d=0; dline(line)->dof_index(d,fe_index); for (unsigned int quad=0; quad<6; ++quad) - for (unsigned int d=0; dquad(quad)->dof_index(d,fe_index); + if (this->face_orientation(quad)) + for (unsigned int d=0; dquad(quad)->dof_index(d,fe_index); + else + for (unsigned int d=0; dquad(quad)->dof_index(this->dof_handler->get_fe()[fe_index]. + adjust_quad_dof_index_for_face_orientation(d),fe_index); for (unsigned int d=0; ddof_index(d,fe_index); } diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 76a768dbb3..25a3752841 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -1195,6 +1195,20 @@ class FiniteElement : public Subscriptor, std::pair face_system_to_component_index (const unsigned int index) const; + /** + * For faces with non-standard + * face_orientation in 3D, the dofs on + * faces (quads) have to be permuted in + * order to be combined with the correct + * shape functions. Given a local dof @p + * index on a quad, return the local index, + * if the face has non-standard + * face_orientation. In 2D and 1D there is + * no need for permutation so the identity + * is returned. + */ + unsigned int adjust_quad_dof_index_for_face_orientation (const unsigned int index) const; + /** * Return in which of the vector * components of this finite @@ -1459,31 +1473,6 @@ class FiniteElement : public Subscriptor, unsigned int component_to_block_index (const unsigned int component) const; - /** - * For faces with non-standard - * face_orientation in 3D, the shape - * functions on faces have to be permutat - * in order to be combined with the correct - * dofs. This function returns a vector of - * integer values, that have to be added to - * the index of a shape function in order - * to result in the permuted index. Prior - * content of the vector @p shifts is - * erased. In 3D a vector of length @p - * dofs_per_quad is returned, in 2D and 1D - * there is no need for permutation and a - * vector of length 0 is returned, the same - * is true for elements which have no dofs - * on quads. The general implementation - * returns a vector of zeros, resulting in - * no permutation at all. This has to be - * overloaded by derived finite element - * classes. - */ - virtual - void - get_face_shape_function_shifts (std::vector &shifts) const; - //@} /** @@ -2078,6 +2067,29 @@ class FiniteElement : public Subscriptor, */ std::vector > generalized_face_support_points; + /** + * For faces with non-standard + * face_orientation in 3D, the dofs on + * faces (quads) have to be permuted in + * order to be combined with the correct + * shape functions. Given a local dof @p + * index on a quad, return the shift in the + * local index, if the face has + * non-standard face_orientation, + * i.e. old_index + shift = + * new_index. In 2D and 1D there is + * no need for permutation so the vector is + * empty. In 3D it has the size of @p + * dofs_per_quad. + * + * The standard implementation fills this + * with zeros, i.e. no permuatation at + * all. Derived finite element classes have + * to fill this vector with the correct + * values. + */ + std::vector adjust_quad_dof_index_for_face_orientation_table; + private: /** * Store what diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index cd57fcf2b6..075439ce24 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -438,26 +438,6 @@ class FE_Q : public FE_Poly,dim> virtual FiniteElementDomination::Domination compare_for_face_domination (const FiniteElement &fe_other) const; - - /** - * For faces with non-standard - * face_orientation in 3D, the shape - * functions on faces have to be permuted - * in order to be combined with the correct - * dofs. This function returns a vector of - * integer values, that have to be added to - * the index of a shape function in order - * to result in the permuted index. Prior - * content of the vector @p shifts is - * erased. In 3D a vector of length @p - * dofs_per_quad is returned, in 2D and 1D - * there is no need for permutation and a - * vector of length 0 is returned. - */ - virtual - void - get_face_shape_function_shifts (std::vector &shifts) const; - //@} /** @@ -546,6 +526,15 @@ class FE_Q : public FE_Poly,dim> * constructor. */ void initialize_unit_face_support_points (); + + /** + * Initialize the + * @p adjust_quad_dof_index_for_face_orientation_table field + * of the FiniteElement + * class. Called from the + * constructor. + */ + void initialize_quad_dof_index_permutation (); /** * Mapping from hierarchic to diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index ff0eb1d84b..901cda2194 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -539,25 +539,6 @@ class FESystem : public FiniteElement virtual FiniteElementDomination::Domination compare_for_face_domination (const FiniteElement &fe_other) const; - - /** - * For faces with non-standard - * face_orientation in 3D, the shape - * functions on faces have to be permuted - * in order to be combined with the correct - * dofs. This function returns a vector of - * integer values, that have to be added to - * the index of a shape function in order - * to result in the permuted index. Prior - * content of the vector @p shifts is - * erased. In 3D a vector of length @p - * dofs_per_quad is returned, in 2D and 1D - * there is no need for permutation and a - * vector of length 0 is returned. - */ - virtual - void - get_face_shape_function_shifts (std::vector &shifts) const; //@} /** @@ -743,6 +724,15 @@ class FESystem : public FiniteElement */ void initialize_unit_face_support_points (); + /** + * Initialize the + * @p adjust_quad_dof_index_for_face_orientation_table field + * of the FiniteElement + * class. Called from the + * constructor. + */ + void initialize_quad_dof_index_permutation (); + /** * Helper function used in the constructor: * take a @p FiniteElementData object diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index a0fcde8046..addb01b923 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -263,38 +263,7 @@ class FEValuesData * non-zero components. */ std::vector shape_function_to_row_table; - - /** - * Vector containing the permutation of - * shape functions necessary if the faces - * of a cell have the wrong - * face_orientation. This is computed only - * once. Actually, this does not contain - * the permutation itself but rather the - * shift of indices needed to calculate the - * permutation. - */ - std::vector shift_in_face_shape_functions; - - /** - * Vector containing the permutation of - * shape functions due to faces with - * non-standard face_orientation on a given - * cell (recomputed on each cell.) If all - * faces are oriented according to the - * standard, this is the identity mapping. - */ - std::vector permuted_shape_functions; - - /** - * Bool flag indicating the need to update - * the @p permuted_shape_functions vector - * on each cell. This is only necessary in - * 3d and if the finite element has - * shape_functions on the face. - */ - bool update_shape_function_permutation; - + /** * Original update flags handed * to the constructor of @@ -657,19 +626,7 @@ class FEValuesBase : protected FEValuesData, const unsigned int point_no, const unsigned int component) const; - /** - * If shape functions belong to a face in - * 3D, they have to be permuted, if the - * face has non-standard face - * orientation. This functuion takes an - * index of a shape function (on a standard - * cell) and returns the corresponding - * shape function on the real cell. - */ - unsigned int - shift_shape_function_index (const unsigned int i) const; - - + //@} /// @name FunctionAccess Access to values of global finite element functions //@{ @@ -1666,12 +1623,6 @@ class FEValuesBase : protected FEValuesData, */ UpdateFlags compute_update_flags (const UpdateFlags update_flags) const; - /** - * Reinit the permutation of (face) shape - * functions to match the present cell. - */ - void reinit(); - private: /** * Copy constructor. Since @@ -2371,17 +2322,15 @@ double FEValuesBase::shape_value (const unsigned int i, const unsigned int j) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_values, ExcAccessToUninitializedField()); - Assert (fe->is_primitive (I), - ExcShapeFunctionNotPrimitive(I)); + Assert (fe->is_primitive (i), + ExcShapeFunctionNotPrimitive(i)); // if the entire FE is primitive, // then we can take a short-cut: if (fe->is_primitive()) - return this->shape_values(I,j); + return this->shape_values(i,j); else // otherwise, use the mapping // between shape function numbers @@ -2392,7 +2341,7 @@ FEValuesBase::shape_value (const unsigned int i, // question to which vector // component the call of this // function refers - return this->shape_values(this->shape_function_to_row_table[I], j); + return this->shape_values(this->shape_function_to_row_table[i], j); } @@ -2404,8 +2353,6 @@ FEValuesBase::shape_value_component (const unsigned int i, const unsigned int j, const unsigned int component) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_values, ExcAccessToUninitializedField()); Assert (component < fe->n_components(), @@ -2420,10 +2367,10 @@ FEValuesBase::shape_value_component (const unsigned int i, // system_to_component_table only // works if the shape function is // primitive): - if (fe->is_primitive(I)) + if (fe->is_primitive(i)) { - if (component == fe->system_to_component_index(I).first) - return this->shape_values(this->shape_function_to_row_table[I],j); + if (component == fe->system_to_component_index(i).first) + return this->shape_values(this->shape_function_to_row_table[i],j); else return 0; } @@ -2438,7 +2385,7 @@ FEValuesBase::shape_value_component (const unsigned int i, // whether the shape function // is non-zero at all within // this component: - if (fe->get_nonzero_components(I)[component] == false) + if (fe->get_nonzero_components(i)[component] == false) return 0.; // count how many non-zero @@ -2450,10 +2397,10 @@ FEValuesBase::shape_value_component (const unsigned int i, // shape function in the arrays // we index presently: const unsigned int - row = (this->shape_function_to_row_table[I] + row = (this->shape_function_to_row_table[i] + - std::count (fe->get_nonzero_components(I).begin(), - fe->get_nonzero_components(I).begin()+component, + std::count (fe->get_nonzero_components(i).begin(), + fe->get_nonzero_components(i).begin()+component, true)); return this->shape_values(row, j); }; @@ -2467,21 +2414,19 @@ const Tensor<1,dim> & FEValuesBase::shape_grad (const unsigned int i, const unsigned int j) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_gradients, ExcAccessToUninitializedField()); - Assert (fe->is_primitive (I), - ExcShapeFunctionNotPrimitive(I)); + Assert (fe->is_primitive (i), + ExcShapeFunctionNotPrimitive(i)); Assert (ishape_gradients.size(), - ExcIndexRange (I, 0, this->shape_gradients.size())); + ExcIndexRange (i, 0, this->shape_gradients.size())); Assert (jshape_gradients[0].size(), ExcIndexRange (j, 0, this->shape_gradients[0].size())); // if the entire FE is primitive, // then we can take a short-cut: if (fe->is_primitive()) - return this->shape_gradients[I][j]; + return this->shape_gradients[i][j]; else // otherwise, use the mapping // between shape function numbers @@ -2492,7 +2437,7 @@ FEValuesBase::shape_grad (const unsigned int i, // question to which vector // component the call of this // function refers - return this->shape_gradients[this->shape_function_to_row_table[I]][j]; + return this->shape_gradients[this->shape_function_to_row_table[i]][j]; } @@ -2504,8 +2449,6 @@ FEValuesBase::shape_grad_component (const unsigned int i, const unsigned int j, const unsigned int component) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_gradients, ExcAccessToUninitializedField()); Assert (component < fe->n_components(), @@ -2520,10 +2463,10 @@ FEValuesBase::shape_grad_component (const unsigned int i, // system_to_component_table only // works if the shape function is // primitive): - if (fe->is_primitive(I)) + if (fe->is_primitive(i)) { - if (component == fe->system_to_component_index(I).first) - return this->shape_gradients[this->shape_function_to_row_table[I]][j]; + if (component == fe->system_to_component_index(i).first) + return this->shape_gradients[this->shape_function_to_row_table[i]][j]; else return Tensor<1,dim>(); } @@ -2538,7 +2481,7 @@ FEValuesBase::shape_grad_component (const unsigned int i, // whether the shape function // is non-zero at all within // this component: - if (fe->get_nonzero_components(I)[component] == false) + if (fe->get_nonzero_components(i)[component] == false) return Tensor<1,dim>(); // count how many non-zero @@ -2550,10 +2493,10 @@ FEValuesBase::shape_grad_component (const unsigned int i, // shape function in the arrays // we index presently: const unsigned int - row = (this->shape_function_to_row_table[I] + row = (this->shape_function_to_row_table[i] + - std::count (fe->get_nonzero_components(I).begin(), - fe->get_nonzero_components(I).begin()+component, + std::count (fe->get_nonzero_components(i).begin(), + fe->get_nonzero_components(i).begin()+component, true)); return this->shape_gradients[row][j]; }; @@ -2567,21 +2510,19 @@ const Tensor<2,dim> & FEValuesBase::shape_2nd_derivative (const unsigned int i, const unsigned int j) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField()); - Assert (fe->is_primitive (I), + Assert (fe->is_primitive (i), ExcShapeFunctionNotPrimitive(i)); - Assert (Ishape_2nd_derivatives.size(), - ExcIndexRange (I, 0, this->shape_2nd_derivatives.size())); + Assert (ishape_2nd_derivatives.size(), + ExcIndexRange (i, 0, this->shape_2nd_derivatives.size())); Assert (jshape_2nd_derivatives[0].size(), ExcIndexRange (j, 0, this->shape_2nd_derivatives[0].size())); // if the entire FE is primitive, // then we can take a short-cut: if (fe->is_primitive()) - return this->shape_2nd_derivatives[I][j]; + return this->shape_2nd_derivatives[i][j]; else // otherwise, use the mapping // between shape function numbers @@ -2592,7 +2533,7 @@ FEValuesBase::shape_2nd_derivative (const unsigned int i, // question to which vector // component the call of this // function refers - return this->shape_2nd_derivatives[this->shape_function_to_row_table[I]][j]; + return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j]; } @@ -2604,8 +2545,6 @@ FEValuesBase::shape_2nd_derivative_component (const unsigned int i, const unsigned int j, const unsigned int component) const { - const unsigned int I=shift_shape_function_index(i); - Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField()); Assert (component < fe->n_components(), @@ -2620,10 +2559,10 @@ FEValuesBase::shape_2nd_derivative_component (const unsigned int i, // system_to_component_table only // works if the shape function is // primitive): - if (fe->is_primitive(I)) + if (fe->is_primitive(i)) { - if (component == fe->system_to_component_index(I).first) - return this->shape_2nd_derivatives[this->shape_function_to_row_table[I]][j]; + if (component == fe->system_to_component_index(i).first) + return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j]; else return Tensor<2,dim>(); } @@ -2650,10 +2589,10 @@ FEValuesBase::shape_2nd_derivative_component (const unsigned int i, // shape function in the arrays // we index presently: const unsigned int - row = (this->shape_function_to_row_table[I] + row = (this->shape_function_to_row_table[i] + - std::count (fe->get_nonzero_components(I).begin(), - fe->get_nonzero_components(I).begin()+component, + std::count (fe->get_nonzero_components(i).begin(), + fe->get_nonzero_components(i).begin()+component, true)); return this->shape_2nd_derivatives[row][j]; }; @@ -2737,28 +2676,6 @@ FEValuesBase::JxW (const unsigned int i) const } - -template -inline -unsigned int -FEValuesBase::shift_shape_function_index (const unsigned int i) const -{ - // standard implementation for 1D and 2D - Assert(idofs_per_cell, ExcInternalError()); - return i; -} - -template <> -inline -unsigned int -FEValuesBase<3>::shift_shape_function_index (const unsigned int i) const -{ - Assert(idofs_per_cell, ExcInternalError()); - return this->permuted_shape_functions[i]; -} - - - /*------------------------ Inline functions: FEValues ----------------------------*/ @@ -2844,8 +2761,6 @@ FEFaceValuesBase::boundary_form (const unsigned int i) const return this->boundary_forms[i]; } - - #endif // DOXYGEN DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index b0f24874c5..13043e4642 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -253,8 +253,13 @@ DoFCellAccessor >::update_cell_dof_indices_cache () const for (unsigned int d=0; dline(line)->dof_index(d); for (unsigned int quad=0; quad<6; ++quad) - for (unsigned int d=0; dquad(quad)->dof_index(d); + if (this->face_orientation(quad)) + for (unsigned int d=0; dquad(quad)->dof_index(d); + else + for (unsigned int d=0; dquad(quad)->dof_index(this->dof_handler->get_fe(). + adjust_quad_dof_index_for_face_orientation(d)); for (unsigned int d=0; d::FiniteElement ( : FiniteElementData (fe_data), cached_primitivity(false), + adjust_quad_dof_index_for_face_orientation_table (this->dofs_per_quad, 0), system_to_base_table(this->dofs_per_cell), face_system_to_base_table(this->dofs_per_face), component_to_base_table (this->components, @@ -333,25 +334,25 @@ FiniteElement::component_to_block_index (const unsigned int index) const } -#if deal_II_dimension < 3 - template -void -FiniteElement::get_face_shape_function_shifts (std::vector &shifts) const +unsigned int +FiniteElement::adjust_quad_dof_index_for_face_orientation (const unsigned int) const { - // general template for 1D and 2D, return an - // empty vector - shifts.clear(); + // general template for 1D and 2D: not implemented + Assert (false, ExcNotImplemented()); + return deal_II_numbers::invalid_unsigned_int; } -#else +#if deal_II_dimension == 3 template <> -void -FiniteElement<3>::get_face_shape_function_shifts (std::vector &shifts) const +unsigned int +FiniteElement<3>::adjust_quad_dof_index_for_face_orientation (const unsigned int index) const { - shifts.clear(); - shifts.resize(this->dofs_per_quad,0); + Assert (indexdofs_per_quad, ExcIndexRange(index,0,this->dofs_per_quad)); + Assert (adjust_quad_dof_index_for_face_orientation_table.size()==this->dofs_per_quad, + ExcInternalError()); + return index+adjust_quad_dof_index_for_face_orientation_table[index]; } #endif diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index b584801c38..199ecb5e5c 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -210,6 +210,8 @@ FE_Q::FE_Q (const unsigned int degree) initialize_constraints (); initialize_embedding (); initialize_restriction (); + + initialize_quad_dof_index_permutation(); } @@ -709,38 +711,6 @@ compare_for_face_domination (const FiniteElement &fe_other) const } - -template -void -FE_Q::get_face_shape_function_shifts (std::vector &shifts) const -{ - // general template for 1D and 2D, return an - // empty vector - shifts.clear(); -} - - - -#if deal_II_dimension == 3 - -template <> -void -FE_Q<3>::get_face_shape_function_shifts (std::vector &shifts) const -{ - shifts.resize(this->dofs_per_quad); - - unsigned int points=this->degree-1; - Assert(points*points==this->dofs_per_quad, ExcInternalError()); - - for (unsigned int local=0; localdofs_per_quad; ++local) - // face support points are in lexicographic - // ordering with x running fastest. invert - // that (y running fastest) - shifts[local] = (local%points)*points + local/points - local; -} - -#endif - //--------------------------------------------------------------------------- // Auxiliary functions //--------------------------------------------------------------------------- @@ -825,6 +795,40 @@ void FE_Q::initialize_unit_face_support_points () +template +void +FE_Q::initialize_quad_dof_index_permutation () +{ + // general template for 1D and 2D, do nothing +} + + + +#if deal_II_dimension == 3 + +template <> +void +FE_Q<3>::initialize_quad_dof_index_permutation () +{ + + Assert (adjust_quad_dof_index_for_face_orientation_table.size()==this->dofs_per_quad, + ExcInternalError()); + + unsigned int points=this->degree-1; + Assert(points*points==this->dofs_per_quad, ExcInternalError()); + + for (unsigned int local=0; localdofs_per_quad; ++local) + // face support points are in lexicographic + // ordering with x running fastest. invert + // that (y running fastest) + this->adjust_quad_dof_index_for_face_orientation_table[local] + =(local%points)*points + local/points - local; +} + +#endif + + + template std::vector FE_Q::get_dpo_vector(const unsigned int deg) diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 48f2e4b7a8..b6db8a9000 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -1729,6 +1729,8 @@ void FESystem::initialize () // on cell and face initialize_unit_support_points (); initialize_unit_face_support_points (); + + initialize_quad_dof_index_permutation (); } @@ -1847,6 +1849,40 @@ initialize_unit_face_support_points () +template +void +FESystem::initialize_quad_dof_index_permutation () +{ + // general template for 1D and 2D, do nothing +} + + + +#if deal_II_dimension == 3 + +template <> +void +FESystem<3>::initialize_quad_dof_index_permutation () +{ + // to obtain the shifts for this composed + // element, concatenate the shift vectors of + // the base elements + for (unsigned int b=0; b &temp=this->base_element(b).adjust_quad_dof_index_for_face_orientation_table; + for (unsigned int c=0; cdofs_per_quad, + ExcInternalError()); +} + +#endif + + + template bool FESystem:: @@ -2926,41 +2962,6 @@ FESystem::unit_face_support_point (const unsigned index) const -template -void -FESystem::get_face_shape_function_shifts (std::vector &shifts) const -{ - // general template for 1D and 2D, return an - // empty vector - shifts.clear(); -} - - - -#if deal_II_dimension == 3 - -template <> -void -FESystem<3>::get_face_shape_function_shifts (std::vector &shifts) const -{ - shifts.clear(); - std::vector temp; - // to obtain the shifts for this composed - // element, concatenate the shift vectors of - // the base elements - for (unsigned int b=0; bbase_element(b).get_face_shape_function_shifts(temp); - for (unsigned int c=0; cdofs_per_quad, ExcInternalError()); -} - -#endif - - - template unsigned int FESystem::memory_consumption () const diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 7d5e1eadb6..84aeac4bf5 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -316,34 +316,6 @@ FEValuesData::initialize (const unsigned int n_quadrature_points, if (flags & update_cell_JxW_values) this->cell_JxW_values.resize(n_quadrature_points); - - // initialize the permutation fields, if they - // are needed - if (dim==3) - { - permuted_shape_functions.resize(fe.dofs_per_cell); - // initialize cell permutation mapping - // with identity - for (unsigned int i=0; ipermuted_shape_functions[i]=i; - - if (fe.dofs_per_quad>0 && - (flags & update_values || - flags & update_gradients || - flags & update_second_derivatives)) - { - // ask the fe to fill the vector of - // shifts - fe.get_face_shape_function_shifts(shift_in_face_shape_functions); - Assert (shift_in_face_shape_functions.size()==fe.dofs_per_quad, - ExcInternalError()); - update_shape_function_permutation=true; - } - else - update_shape_function_permutation=false; - } - else - update_shape_function_permutation=false; } @@ -393,43 +365,6 @@ FEValuesBase::~FEValuesBase () } -#if deal_II_dimension <3 -template -void FEValuesBase::reinit () -{ - // do nothing in 1D and 2D -} -#else - -template <> -void FEValuesBase<3>::reinit () -{ - // in 3D reinit the permutation of face shape - // functions, if necessary - if (update_shape_function_permutation) - { - Assert (this->shift_in_face_shape_functions.size()==this->fe->dofs_per_quad, - ExcInternalError()); - Triangulation<3>::cell_iterator thiscell = *(this->present_cell); - unsigned int offset=fe->first_quad_index; - for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) - { - if (thiscell->face_orientation(face_no)) - for (unsigned int i=0; idofs_per_quad; ++i) - this->permuted_shape_functions[offset+i]=offset+i; - else - for (unsigned int i=0; idofs_per_quad; ++i) - this->permuted_shape_functions[offset+i]=offset+i+this->shift_in_face_shape_functions[i]; - - offset+=this->fe->dofs_per_quad; - } - Assert (offset-fe->first_quad_index==GeometryInfo<3>::faces_per_cell*this->fe->dofs_per_quad, - ExcInternalError()); - } -} -#endif - - template template @@ -1076,15 +1011,12 @@ FEValuesBase::memory_consumption () const return (MemoryConsumption::memory_consumption (this->shape_values) + MemoryConsumption::memory_consumption (this->shape_gradients) + MemoryConsumption::memory_consumption (this->shape_2nd_derivatives) + - MemoryConsumption::memory_consumption (this->permuted_shape_functions) + - MemoryConsumption::memory_consumption (this->shift_in_face_shape_functions) + MemoryConsumption::memory_consumption (this->JxW_values) + MemoryConsumption::memory_consumption (this->quadrature_points) + MemoryConsumption::memory_consumption (this->normal_vectors) + MemoryConsumption::memory_consumption (this->boundary_forms) + MemoryConsumption::memory_consumption (this->cell_JxW_values) + sizeof(this->update_flags) + - sizeof(this->update_shape_function_permutation) + MemoryConsumption::memory_consumption (n_quadrature_points) + MemoryConsumption::memory_consumption (dofs_per_cell) + MemoryConsumption::memory_consumption (mapping) + @@ -1328,8 +1260,6 @@ void FEValues::do_reinit () this->fe_data->clear_first_cell (); this->mapping_data->clear_first_cell (); - - FEValuesBase::reinit(); } @@ -1607,7 +1537,6 @@ void FEFaceValues::do_reinit (const unsigned int face_no) this->fe_data->clear_first_cell (); this->mapping_data->clear_first_cell (); - FEValuesBase::reinit(); } @@ -1860,7 +1789,6 @@ void FESubfaceValues::do_reinit (const unsigned int face_no, this->fe_data->clear_first_cell (); this->mapping_data->clear_first_cell (); - FEValuesBase::reinit(); } diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 77293dbeaa..cce5704b4f 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -959,9 +959,9 @@ inconvenience this causes.

deal.II

    -
  1. Fixed: On faces with wrong face_orientation the shape - functions have to reordered in order to be combined with the correct - dofs. This is only relevant for continuous elements in 3D. At least for +

  2. Fixed: On faces with wrong face_orientation the dofs + have to reordered in order to be combined with the correct shape + functions. This is only relevant for continuous elements in 3D. At least for FE_Q and systems of FE_Q this works now, for other finite elements the reordering vector still has to be implemented.