]> https://gitweb.dealii.org/ - dealii.git/commitdiff
extend new functions to non-tensor product
authorDenis Davydov <davydden@gmail.com>
Thu, 19 May 2016 12:15:22 +0000 (14:15 +0200)
committerDenis Davydov <davydden@gmail.com>
Thu, 19 May 2016 14:55:54 +0000 (16:55 +0200)
include/deal.II/fe/fe_tools.h
source/fe/fe_tools.cc

index 14fb9995e7137806a7fcef03b862810ab2894061..892e8a90e50909e6e0d55b08bc1575636dd4fb33 100644 (file)
@@ -950,11 +950,13 @@ namespace FETools
                     const bool do_tensor_product = true);
 
   /**
-   * build...
-   * @param face_system_to_base_table
-   * @param face_system_to_component_table
-   * @param finite_element
-   * @param do_tensor_product
+   * For a given (composite) @p finite_element build @p face_system_to_base_table,
+   * and @p face_system_to_component_table.
+   *
+   * If @p do_tensor_product is <code>true</code>, the underlying finite element
+   * is assumed to be build using the tensor product rule. That is, the number of
+   * composite components is the sum of components in each finite element times
+   * multiplicity.
    */
   template <int dim, int spacedim>
   void
index 8f29669e06ca37705351b2561bff69809831a136..b0a6e7e2a68314a10a2535bdca7e85c91f497f72 100644 (file)
@@ -597,7 +597,7 @@ namespace FETools
         unsigned int comp_start = 0;
         for (unsigned int base=0; base<fe.n_base_elements(); ++base)
           for (unsigned int m=0; m<fe.element_multiplicity(base);
-               ++m, comp_start+=fe.base_element(base).n_components())
+               ++m, comp_start+=fe.base_element(base).n_components() * do_tensor_product)
             for (unsigned int local_index = 0;
                  local_index < fe.base_element(base).dofs_per_vertex;
                  ++local_index, ++total_index)
@@ -634,7 +634,7 @@ namespace FETools
           unsigned int comp_start = 0;
           for (unsigned int base=0; base<fe.n_base_elements(); ++base)
             for (unsigned int m=0; m<fe.element_multiplicity(base);
-                 ++m, comp_start+=fe.base_element(base).n_components())
+                 ++m, comp_start+=fe.base_element(base).n_components() * do_tensor_product)
               for (unsigned int local_index = 0;
                    local_index < fe.base_element(base).dofs_per_line;
                    ++local_index, ++total_index)
@@ -672,7 +672,7 @@ namespace FETools
           unsigned int comp_start = 0;
           for (unsigned int base=0; base<fe.n_base_elements(); ++base)
             for (unsigned int m=0; m<fe.element_multiplicity(base);
-                 ++m, comp_start += fe.base_element(base).n_components())
+                 ++m, comp_start += fe.base_element(base).n_components() * do_tensor_product)
               for (unsigned int local_index = 0;
                    local_index < fe.base_element(base).dofs_per_quad;
                    ++local_index, ++total_index)
@@ -710,7 +710,7 @@ namespace FETools
           unsigned int comp_start = 0;
           for (unsigned int base=0; base<fe.n_base_elements(); ++base)
             for (unsigned int m=0; m<fe.element_multiplicity(base);
-                 ++m, comp_start+=fe.base_element(base).n_components())
+                 ++m, comp_start+=fe.base_element(base).n_components() * do_tensor_product)
               for (unsigned int local_index = 0;
                    local_index < fe.base_element(base).dofs_per_hex;
                    ++local_index, ++total_index)
@@ -764,7 +764,7 @@ namespace FETools
         unsigned int comp_start = 0;
         for (unsigned int base=0; base<fe.n_base_elements(); ++base)
           for (unsigned int m=0; m<fe.element_multiplicity(base);
-               ++m, comp_start += fe.base_element(base).n_components())
+               ++m, comp_start += fe.base_element(base).n_components() * do_tensor_product)
             for (unsigned int local_index = 0;
                  local_index < fe.base_element(base).dofs_per_vertex;
                  ++local_index, ++total_index)
@@ -813,7 +813,7 @@ namespace FETools
           unsigned int comp_start = 0;
           for (unsigned int base = 0; base < fe.n_base_elements(); ++base)
             for (unsigned int m=0; m<fe.element_multiplicity(base);
-                 ++m, comp_start += fe.base_element(base).n_components())
+                 ++m, comp_start += fe.base_element(base).n_components() * do_tensor_product)
               for (unsigned int local_index = 0;
                    local_index < fe.base_element(base).dofs_per_line;
                    ++local_index, ++total_index)
@@ -857,7 +857,7 @@ namespace FETools
           unsigned int comp_start = 0;
           for (unsigned int base=0; base<fe.n_base_elements(); ++base)
             for (unsigned int m=0; m<fe.element_multiplicity(base);
-                 ++m, comp_start += fe.base_element(base).n_components())
+                 ++m, comp_start += fe.base_element(base).n_components() * do_tensor_product)
               for (unsigned int local_index = 0;
                    local_index < fe.base_element(base).dofs_per_quad;
                    ++local_index, ++total_index)

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.