]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
system_to_component_table is correct
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Feb 1999 16:24:03 +0000 (16:24 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Feb 1999 16:24:03 +0000 (16:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@804 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/source/fe/fe_system.cc

index db56eae77afc990a0940f506c840e8c8a58e2540..7b91913ddd485d7c3e4f39df9f7fe518ca3b5a90 100644 (file)
@@ -27,39 +27,112 @@ void FESystem<dim>::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<unsigned,unsigned>
-                                                  (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<unsigned,unsigned>
+                   (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<unsigned,unsigned>
-                                                  (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<unsigned,unsigned>
+                   (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<unsigned,unsigned>
-                                                  (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<unsigned,unsigned>
+                   (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<unsigned,unsigned>
+                   (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<dim>::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 <int dim>
 double FESystem<dim>::shape_value (const unsigned int i,
                                   const Point<dim>  &p) const
 {
+  Assert(false, ExcNotImplemented());
+  return 0.;
+  
   Assert((i<total_dofs), ExcInvalidIndex(i));
 
+
 //  return base_element->shape_value (i / n_sub_elements, p);
 };
 
@@ -136,7 +239,11 @@ double FESystem<dim>::shape_value (const unsigned int i,
 template <int dim>
 Tensor<1,dim>
 FESystem<dim>::shape_grad (const unsigned int  i,
-                          const Point<dim>   &p) const {
+                          const Point<dim>   &p) const
+{
+  Assert(false, ExcNotImplemented());
+  return Tensor<1,dim>();
+
   Assert((i<total_dofs), ExcInvalidIndex(i));
 
 //  return base_element->shape_grad (i / n_sub_elements, p);
@@ -147,7 +254,11 @@ FESystem<dim>::shape_grad (const unsigned int  i,
 template <int dim>
 Tensor<2,dim>
 FESystem<dim>::shape_grad_grad (const unsigned int  i,
-                               const Point<dim>   &p) const {
+                               const Point<dim>   &p) const
+ {
+  Assert(false, ExcNotImplemented());
+  return Tensor<2,dim>();
+
   Assert((i<total_dofs), ExcInvalidIndex(i));
 
 //  return base_element->shape_grad_grad (i / n_sub_elements, p);
@@ -156,8 +267,9 @@ FESystem<dim>::shape_grad_grad (const unsigned int  i,
 
 
 template <int dim>
-void FESystem<dim>::get_unit_support_points (vector<Point<dim> > &support_points) const {
-  Assert (support_points.size() == total_dofs,
+void FESystem<dim>::get_unit_support_points (vector<Point<dim> > &support_points) const
+{
+/*  Assert (support_points.size() == total_dofs,
          ExcWrongFieldDimension (support_points.size(), total_dofs));
 
   vector<Point<dim> > base_support_points (base_element->total_dofs);
@@ -166,6 +278,7 @@ void FESystem<dim>::get_unit_support_points (vector<Point<dim> > &support_points
   for (unsigned int i=0; i<base_element->total_dofs; ++i) 
     for (unsigned int n=0; n<n_sub_elements; ++n)
       support_points[i*n_sub_elements+n] = base_support_points[i];
+*/
 };
 
 
@@ -173,7 +286,10 @@ void FESystem<dim>::get_unit_support_points (vector<Point<dim> > &support_points
 template <int dim>
 void FESystem<dim>::get_support_points (const DoFHandler<dim>::cell_iterator &cell,
                                        const Boundary<dim> &boundary,
-                                       vector<Point<dim> > &support_points) const {
+                                       vector<Point<dim> > &support_points) const
+{
+  Assert(false, ExcNotImplemented());
+/*
   Assert (support_points.size() == total_dofs,
          ExcWrongFieldDimension (support_points.size(), total_dofs));
 
@@ -183,6 +299,7 @@ void FESystem<dim>::get_support_points (const DoFHandler<dim>::cell_iterator &ce
   for (unsigned int i=0; i<base_element->total_dofs; ++i) 
     for (unsigned int n=0; n<n_sub_elements; ++n)
       support_points[i*n_sub_elements+n] = base_support_points[i];
+*/
 };
 
 
@@ -190,7 +307,10 @@ void FESystem<dim>::get_support_points (const DoFHandler<dim>::cell_iterator &ce
 template <int dim>
 void FESystem<dim>::get_face_support_points (const DoFHandler<dim>::face_iterator &face,
                                             const Boundary<dim> &boundary,
-                                            vector<Point<dim> > &support_points) const {
+                                            vector<Point<dim> > &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<dim>::get_face_support_points (const DoFHandler<dim>::face_iterato
   for (unsigned int i=0; i<base_element->dofs_per_face; ++i) 
     for (unsigned int n=0; n<n_sub_elements; ++n)
       support_points[i*n_sub_elements+n] = base_support_points[i];
+*/
 };
 
 
@@ -207,7 +328,10 @@ void FESystem<dim>::get_face_support_points (const DoFHandler<dim>::face_iterato
 template <int dim>
 void FESystem<dim>::get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
                                           const Boundary<dim> &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<dim>::get_local_mass_matrix (const DoFHandler<dim>::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 <int dim>
 double FESystem<dim>::shape_value_transform (const unsigned int i,
-                                            const Point<dim>  &p) const {
-  return base_element->shape_value_transform (i, p);
+                                            const Point<dim>  &p) const
+{
+  Assert(false, ExcNotImplemented());
+  return 0.;
+  
+//  return base_element->shape_value_transform (i, p);
 };
 
 
-
 template <int dim>
 Tensor<1,dim> FESystem<dim>::shape_grad_transform (const unsigned int i,
-                                                  const Point<dim>  &p) const {
-  return base_element->shape_grad_transform (i, p);
+                                                  const Point<dim>  &p) const
+{
+  Assert(false, ExcNotImplemented());
+  return Tensor<1,dim>();
+
+//  return base_element->shape_grad_transform (i, p);
 };
 
 
@@ -252,8 +384,11 @@ template <int dim>
 void FESystem<dim>::get_face_jacobians (const DoFHandler<dim>::face_iterator &face,
                                        const Boundary<dim>         &boundary,
                                        const vector<Point<dim-1> > &unit_points,
-                                       vector<double>      &face_jacobi_determinants) const {
-  base_element->get_face_jacobians (face, boundary, unit_points, face_jacobi_determinants);
+                                       vector<double>      &face_jacobi_determinants) const
+{
+  Assert(false, ExcNotImplemented());
+
+//  base_element->get_face_jacobians (face, boundary, unit_points, face_jacobi_determinants);
 };
 
 
@@ -262,8 +397,11 @@ template <int dim>
 void FESystem<dim>::get_subface_jacobians (const DoFHandler<dim>::face_iterator &face,
                                           const unsigned int           subface_no,
                                           const vector<Point<dim-1> > &unit_points,
-                                          vector<double>      &face_jacobi_determinants) const {
-  base_element->get_subface_jacobians (face, subface_no, unit_points, face_jacobi_determinants);
+                                          vector<double>      &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<dim>::get_normal_vectors (const DoFHandler<dim>::cell_iterator &ce
                                        const unsigned int          face_no,
                                        const Boundary<dim>         &boundary,
                                        const vector<Point<dim-1> > &unit_points,
-                                       vector<Point<dim> >         &normal_vectors) const {
-  base_element->get_normal_vectors (cell, face_no, boundary, unit_points, normal_vectors);
+                                       vector<Point<dim> >         &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<dim>::get_normal_vectors (const DoFHandler<dim>::cell_iterator &ce
                                        const unsigned int          face_no,
                                        const unsigned int          subface_no,
                                        const vector<Point<dim-1> > &unit_points,
-                                       vector<Point<dim> >         &normal_vectors) const {
-  base_element->get_normal_vectors (cell, face_no, subface_no, unit_points, normal_vectors);
+                                       vector<Point<dim> >         &normal_vectors) const
+{
+  Assert(false, ExcNotImplemented());
+
+//  base_element->get_normal_vectors (cell, face_no, subface_no, unit_points, normal_vectors);
 };
 
-*/
+
 
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.