]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Minor cleanup: Make a couple of variables const and move initialization appropriately.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 17 Dec 2008 01:24:36 +0000 (01:24 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 17 Dec 2008 01:24:36 +0000 (01:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@17961 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/mapping_q.h
deal.II/deal.II/source/fe/mapping_q.cc

index c9caead15715da7eec64ff792f4fb41d5e739f91..0ae4531ef9118d0e74ac8d7e7a488e83d3687e7d 100644 (file)
@@ -522,7 +522,7 @@ class MappingQ : public MappingQ1<dim,spacedim>
                                      * mapping of cells at the
                                      * boundary.
                                      */
-    TensorProductPolynomials<dim> *tensor_pols;
+    const TensorProductPolynomials<dim> *tensor_pols;
     
                                     /**
                                      * Number of the Qp tensor
@@ -536,7 +536,7 @@ class MappingQ : public MappingQ1<dim,spacedim>
                                      * numbering. Its size is
                                      * @p dofs_per_cell.
                                      */
-    std::vector<unsigned int> renumber;
+    const std::vector<unsigned int> renumber;
 
                                     /**
                                      * If this flag is set @p true
index cb54d08cd3dac090c9915bcbaec565e34a350f9f..880139da6ac2c9f46cca99c528eb46b8415d3763 100644 (file)
@@ -97,6 +97,22 @@ MappingQ<1>::~MappingQ ()
 
 
 
+namespace
+{
+  template <int dim>
+  std::vector<unsigned int>
+  get_dpo_vector (const unsigned int degree)
+  {
+    std::vector<unsigned int> dpo(dim+1, 1U);
+    for (unsigned int i=1; i<dpo.size(); ++i)
+      dpo[i]=dpo[i-1]*(degree-1);
+    return dpo;
+  }
+}
+
+
+
+
 template<int dim, int spacedim>
 MappingQ<dim,spacedim>::MappingQ (const unsigned int p,
                         const bool use_mapping_q_on_all_cells)
@@ -107,7 +123,10 @@ MappingQ<dim,spacedim>::MappingQ (const unsigned int p,
                        :8+12*(degree-1)+6*(degree-1)*(degree-1)),
                tensor_pols(0),
                n_shape_functions(Utilities::fixed_power<dim>(degree+1)),
-               renumber(0),
+               renumber(FETools::
+                        lexicographic_to_hierarchic_numbering (
+                          FiniteElementData<dim> (get_dpo_vector<dim>(degree), 1,
+                                                  degree))),
                use_mapping_q_on_all_cells (use_mapping_q_on_all_cells),
                feq(degree)
 {
@@ -124,25 +143,6 @@ MappingQ<dim,spacedim>::MappingQ (const unsigned int p,
          ExcInternalError());
   Assert(n_inner+n_outer==n_shape_functions, ExcInternalError());
   
-                                  // build the renumbering of the
-                                  // shape functions of the Qp
-                                  // mapping.
-  renumber.resize(n_shape_functions,0);
-                                  // instead of creating a full FE_Q
-                                  // object to be passed to the
-                                  // lexicographic_to_hierarchic_numbering
-                                  // function we only create the
-                                  // underlying FiniteElementData
-                                  // object. We can't access
-                                  // FE_Q::get_dpo_vector as MappingQ
-                                  // is not friend of FE_Q. Create
-                                  // the dpo vector ourselve.
-  std::vector<unsigned int> dpo(dim+1, 1U);
-  for (unsigned int i=1; i<dpo.size(); ++i)
-    dpo[i]=dpo[i-1]*(degree-1);
-  FETools::lexicographic_to_hierarchic_numbering (
-    FiniteElementData<dim> (dpo, 1, degree), renumber);
-
                                   // build laplace_on_quad_vector
   if (degree>1)
     {

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.