From: guido Date: Mon, 15 Feb 1999 16:24:03 +0000 (+0000) Subject: system_to_component_table is correct X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4d6d92d2fc4908dfcb584aafa694a3b74a01f4b0;p=dealii-svn.git system_to_component_table is correct git-svn-id: https://svn.dealii.org/trunk@804 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index db56eae77a..7b91913ddd 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -27,39 +27,112 @@ void FESystem::initialize () // Multi-component base elements have to be thought of. // 1. Vertices unsigned total_index = 0; - for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + for (unsigned vertex_number= 0 ; vertex_number < GeometryInfo<2>::vertices_per_cell ; + ++vertex_number) { - for (unsigned local_index = 0 ; - local_index < base_element(comp).dofs_per_vertex ; - ++local_index) + unsigned comp_start = 0; + for(unsigned comp = 0; comp < n_component_elements() ; + ++comp) { - system_to_component_table[total_index++] = pair - (comp, 0); + for (unsigned m = 0; m < element_multiplicity(comp); ++m) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_vertex ; + ++local_index) + { + cerr << "V(sys, comp, m, loc) = (" + << total_index << "," << comp << "," << m << "," << local_index << ")" << endl; + system_to_component_table[total_index++] + = pair + (comp_start+m, + vertex_number*base_element(comp).dofs_per_vertex+local_index); + } + } + comp_start += element_multiplicity(comp); } } + // 2. Lines - for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + for (unsigned line_number= 0 ; line_number < GeometryInfo<2>::lines_per_cell ; + ++line_number) { - for (unsigned local_index = 0 ; - local_index < base_element(comp).dofs_per_vertex ; - ++local_index) + unsigned comp_start = 0; + for(unsigned comp = 0; comp < n_component_elements() ; + ++comp) { - system_to_component_table[total_index++] = pair - (comp, 0); + for (unsigned m = 0; m < element_multiplicity(comp); ++m) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_line ; + ++local_index) + { + cerr << "L(sys, comp, m, loc) = (" + << total_index << "," << comp << "," << m << "," << local_index << ")" << endl; + system_to_component_table[total_index++] + = pair + (comp_start+m, + line_number*base_element(comp).dofs_per_line + +local_index+base_element(comp).first_line_index); + } + } + comp_start += element_multiplicity(comp); } } + // 3. Quads - for(unsigned comp = 0; comp < n_component_elements() ; ++comp) + for (unsigned quad_number= 0 ; quad_number < GeometryInfo<2>::quads_per_cell ; + ++quad_number) { - for (unsigned local_index = 0 ; - local_index < base_element(comp).dofs_per_vertex ; - ++local_index) + unsigned comp_start = 0; + for(unsigned comp = 0; comp < n_component_elements() ; + ++comp) { - system_to_component_table[total_index++] = pair - (comp, 0); + for (unsigned m = 0; m < element_multiplicity(comp); ++m) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_quad ; + ++local_index) + { + cerr << "Q(sys, comp, m, loc) = (" + << total_index << "," << comp << "," << m << "," << local_index << ")" << endl; + system_to_component_table[total_index++] + = pair + (comp_start+m, + quad_number*base_element(comp).dofs_per_quad + +local_index+base_element(comp).first_quad_index); + } + } + comp_start += element_multiplicity(comp); } } + // 4. Hex + for (unsigned hex_number= 0 ; hex_number < GeometryInfo<2>::hexes_per_cell ; + ++hex_number) + { + unsigned comp_start = 0; + for(unsigned comp = 0; comp < n_component_elements() ; + ++comp, comp_start += element_multiplicity(comp)) + { + for (unsigned m = 0; m < element_multiplicity(comp); ++m) + { + for (unsigned local_index = 0 ; + local_index < base_element(comp).dofs_per_hex ; + ++local_index) + { + cerr << "Q(sys, comp, m, loc) = (" + << total_index << "," << comp << "," << m << "," << local_index << ")" << endl; + system_to_component_table[total_index++] + = pair + (comp_start+m, + hex_number*base_element(comp).dofs_per_hex + +local_index+base_element(comp).first_hex_index); + } + } + } + } + + // distribute the matrices of the base // finite element to the matrices of @@ -93,13 +166,27 @@ void FESystem::initialize () template <> FiniteElementData<1> FESystem<1>::multiply_dof_numbers (const FiniteElementData<1> &fe_data, - const unsigned int N) { + const unsigned int N) +{ return FiniteElementData<1> (fe_data.dofs_per_vertex * N, fe_data.dofs_per_line * N, fe_data.n_transform_functions, fe_data.n_components * N); }; +template <> +FiniteElementData<1> +FESystem<1>::multiply_dof_numbers (const FiniteElementData<1> &fe1, + const unsigned int N1, + const FiniteElementData<1> &fe2, + const unsigned int N2) +{ + return FiniteElementData<1> (fe1.dofs_per_vertex * N1 + fe2.dofs_per_vertex * N2 , + fe1.dofs_per_line * N1 + fe2.dofs_per_line * N2 , + fe1.n_transform_functions + fe2.n_transform_functions , + fe1.n_components * N1 + fe2.n_components * N2 ); +}; + #endif @@ -117,17 +204,33 @@ FESystem<2>::multiply_dof_numbers (const FiniteElementData<2> &fe_data, fe_data.n_components * N); }; -#endif +template <> +FiniteElementData<2> +FESystem<2>::multiply_dof_numbers (const FiniteElementData<2> &fe1, + const unsigned int N1, + const FiniteElementData<2> &fe2, + const unsigned int N2) +{ + return FiniteElementData<2> (fe1.dofs_per_vertex * N1 + fe2.dofs_per_vertex * N2 , + fe1.dofs_per_line * N1 + fe2.dofs_per_line * N2 , + fe1.dofs_per_quad * N1 + fe2.dofs_per_quad * N2 , + fe1.n_transform_functions + fe2.n_transform_functions , + fe1.n_components * N1 + fe2.n_components * N2 ); +}; +#endif -/* template double FESystem::shape_value (const unsigned int i, const Point &p) const { + Assert(false, ExcNotImplemented()); + return 0.; + Assert((ishape_value (i / n_sub_elements, p); }; @@ -136,7 +239,11 @@ double FESystem::shape_value (const unsigned int i, template Tensor<1,dim> FESystem::shape_grad (const unsigned int i, - const Point &p) const { + const Point &p) const +{ + Assert(false, ExcNotImplemented()); + return Tensor<1,dim>(); + Assert((ishape_grad (i / n_sub_elements, p); @@ -147,7 +254,11 @@ FESystem::shape_grad (const unsigned int i, template Tensor<2,dim> FESystem::shape_grad_grad (const unsigned int i, - const Point &p) const { + const Point &p) const + { + Assert(false, ExcNotImplemented()); + return Tensor<2,dim>(); + Assert((ishape_grad_grad (i / n_sub_elements, p); @@ -156,8 +267,9 @@ FESystem::shape_grad_grad (const unsigned int i, template -void FESystem::get_unit_support_points (vector > &support_points) const { - Assert (support_points.size() == total_dofs, +void FESystem::get_unit_support_points (vector > &support_points) const +{ +/* Assert (support_points.size() == total_dofs, ExcWrongFieldDimension (support_points.size(), total_dofs)); vector > base_support_points (base_element->total_dofs); @@ -166,6 +278,7 @@ void FESystem::get_unit_support_points (vector > &support_points for (unsigned int i=0; itotal_dofs; ++i) for (unsigned int n=0; n::get_unit_support_points (vector > &support_points template void FESystem::get_support_points (const DoFHandler::cell_iterator &cell, const Boundary &boundary, - vector > &support_points) const { + vector > &support_points) const +{ + Assert(false, ExcNotImplemented()); +/* Assert (support_points.size() == total_dofs, ExcWrongFieldDimension (support_points.size(), total_dofs)); @@ -183,6 +299,7 @@ void FESystem::get_support_points (const DoFHandler::cell_iterator &ce for (unsigned int i=0; itotal_dofs; ++i) for (unsigned int n=0; n::get_support_points (const DoFHandler::cell_iterator &ce template void FESystem::get_face_support_points (const DoFHandler::face_iterator &face, const Boundary &boundary, - vector > &support_points) const { + vector > &support_points) const +{ + Assert(false, ExcNotImplemented()); +/* Assert (support_points.size() == dofs_per_face, ExcWrongFieldDimension (support_points.size(), dofs_per_face)); @@ -200,6 +320,7 @@ void FESystem::get_face_support_points (const DoFHandler::face_iterato for (unsigned int i=0; idofs_per_face; ++i) for (unsigned int n=0; n::get_face_support_points (const DoFHandler::face_iterato template void FESystem::get_local_mass_matrix (const DoFHandler::cell_iterator &cell, const Boundary &boundary, - dFMatrix &local_mass_matrix) const { + dFMatrix &local_mass_matrix) const +{ + Assert(false, ExcNotImplemented()); +/* Assert (local_mass_matrix.n() == total_dofs, ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs)); Assert (local_mass_matrix.m() == total_dofs, @@ -228,22 +352,30 @@ void FESystem::get_local_mass_matrix (const DoFHandler::cell_iterator // only fill diagonals of the blocks local_mass_matrix (i*n_sub_elements + n, j*n_sub_elements + n) = base_mass_matrix (i,j); +*/ }; template double FESystem::shape_value_transform (const unsigned int i, - const Point &p) const { - return base_element->shape_value_transform (i, p); + const Point &p) const +{ + Assert(false, ExcNotImplemented()); + return 0.; + +// return base_element->shape_value_transform (i, p); }; - template Tensor<1,dim> FESystem::shape_grad_transform (const unsigned int i, - const Point &p) const { - return base_element->shape_grad_transform (i, p); + const Point &p) const +{ + Assert(false, ExcNotImplemented()); + return Tensor<1,dim>(); + +// return base_element->shape_grad_transform (i, p); }; @@ -252,8 +384,11 @@ template void FESystem::get_face_jacobians (const DoFHandler::face_iterator &face, const Boundary &boundary, const vector > &unit_points, - vector &face_jacobi_determinants) const { - base_element->get_face_jacobians (face, boundary, unit_points, face_jacobi_determinants); + vector &face_jacobi_determinants) const +{ + Assert(false, ExcNotImplemented()); + +// base_element->get_face_jacobians (face, boundary, unit_points, face_jacobi_determinants); }; @@ -262,8 +397,11 @@ template void FESystem::get_subface_jacobians (const DoFHandler::face_iterator &face, const unsigned int subface_no, const vector > &unit_points, - vector &face_jacobi_determinants) const { - base_element->get_subface_jacobians (face, subface_no, unit_points, face_jacobi_determinants); + vector &face_jacobi_determinants) const +{ + Assert(false, ExcNotImplemented()); + +// base_element->get_subface_jacobians (face, subface_no, unit_points, face_jacobi_determinants); }; @@ -274,8 +412,11 @@ void FESystem::get_normal_vectors (const DoFHandler::cell_iterator &ce const unsigned int face_no, const Boundary &boundary, const vector > &unit_points, - vector > &normal_vectors) const { - base_element->get_normal_vectors (cell, face_no, boundary, unit_points, normal_vectors); + vector > &normal_vectors) const +{ + Assert(false, ExcNotImplemented()); + +// base_element->get_normal_vectors (cell, face_no, boundary, unit_points, normal_vectors); }; @@ -285,11 +426,14 @@ void FESystem::get_normal_vectors (const DoFHandler::cell_iterator &ce const unsigned int face_no, const unsigned int subface_no, const vector > &unit_points, - vector > &normal_vectors) const { - base_element->get_normal_vectors (cell, face_no, subface_no, unit_points, normal_vectors); + vector > &normal_vectors) const +{ + Assert(false, ExcNotImplemented()); + +// base_element->get_normal_vectors (cell, face_no, subface_no, unit_points, normal_vectors); }; -*/ +