]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New output_indices and set/get_renumbering functions according to those in the Polyno...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 11 Mar 2004 14:55:33 +0000 (14:55 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 11 Mar 2004 14:55:33 +0000 (14:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@8715 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/tensor_product_polynomials.h
deal.II/base/source/tensor_product_polynomials.cc

index 11405be2d9e18c4a7bd5b8138d0287bca3fabf6c..ffc5c6f7268435bff54ae7cd54fcb31f5f1accca 100644 (file)
  * ..., P<sub>1</sub>(x)P<sub>2</sub>(y),
  * P<sub>2</sub>(x)P<sub>2</sub>(y), P<sub>3</sub>(x)P<sub>2</sub>(y),
  * ...</i> and likewise in 3d.
+ *
+ * The @ref{output_indices} function prints the ordering of the
+ * dim-dimensional polynomials, i.e. for each polynomial in the
+ * polynomial space it gives the indices i,j,k of the one-dimensional
+ * polynomials in x,y and z direction. The ordering of the
+ * dim-dimensional polynomials can be changed by using the
+ * @p{set_renumbering} function. 
  * 
- * @author Ralf Hartmann, Guido Kanschat, 2000, Wolfgang Bangerth 2003
+ * @author Ralf Hartmann, 2000, 2004, Guido Kanschat, 2000, Wolfgang Bangerth 2003
  */
 template <int dim>
 class TensorProductPolynomials
@@ -60,6 +67,33 @@ class TensorProductPolynomials
                                      */
     template <class Pol>
     TensorProductPolynomials (const std::vector<Pol> &pols);
+    
+                                    /**
+                                     * Prints the list of the indices
+                                     * to <tt>out</tt>.
+                                     */
+    void output_indices(std::ostream &out) const;
+
+                                    /**
+                                     * Sets the ordering of the
+                                     * polynomials. Requires
+                                     * <tt>renumber.size()==n()</tt>.
+                                     * Stores a copy of
+                                     * <tt>renumber</tt>.
+                                     */
+    void set_renumbering(const std::vector<unsigned int> &renumber);
+
+                                    /**
+                                     * Gives read access to the
+                                     * renumbering vector.
+                                     */
+    const std::vector<unsigned int> &get_renumbering() const;
+
+                                    /**
+                                     * Gives read access to the
+                                     * inverse renumbering vector.
+                                     */
+    const std::vector<unsigned int> &get_renumbering_inverse() const;
 
                                     /**
                                      * Computes the value and the
@@ -200,6 +234,18 @@ class TensorProductPolynomials
                                      */
     unsigned int n_tensor_pols;
 
+                                    /**
+                                     * Index map for reordering the
+                                     * polynomials.
+                                     */
+    std::vector<unsigned int> index_map;
+
+                                    /**
+                                     * Index map for reordering the
+                                     * polynomials.
+                                     */
+    std::vector<unsigned int> index_map_inverse;
+
                                      /**
                                       * Each tensor product polynomial
                                       * <i>i</i> is a product of
@@ -223,7 +269,27 @@ class TensorProductPolynomials
     unsigned int x_to_the_dim (const unsigned int x);
 };
 
+/// @if NoDoc
+
+template <int dim>
+inline
+const std::vector<unsigned int> &
+TensorProductPolynomials<dim>::get_renumbering() const
+{
+  return index_map;
+}
+
 
+template <int dim>
+inline
+const std::vector<unsigned int> &
+TensorProductPolynomials<dim>::get_renumbering_inverse() const
+{
+  return index_map_inverse;
+}
+
+
+/// @endif
 
 
 /**
@@ -480,8 +546,20 @@ TensorProductPolynomials<dim>::
 TensorProductPolynomials(const std::vector<Pol> &pols)
                :
                polynomials (pols.begin(), pols.end()),
-               n_tensor_pols(x_to_the_dim(pols.size()))
-{}
+               n_tensor_pols(x_to_the_dim(pols.size())),
+               index_map(n_tensor_pols),
+               index_map_inverse(n_tensor_pols)
+{
+                                  // per default set this index map
+                                  // to identity. This map can be
+                                  // changed by the user through the
+                                  // set_renumbering function
+  for (unsigned int i=0; i<n_tensor_pols; ++i)
+    {
+      index_map[i]=i;
+      index_map_inverse[i]=i;
+    }  
+}
 
 
 
index 83e6dbc8c6b11039e702b781add7b175f42dccb5..78f608ae24ab99fe4103f67bfe0665105e021372 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -27,7 +27,7 @@ compute_index (const unsigned int i,
                unsigned int       (&indices)[1]) const
 {
   Assert (i<polynomials.size(), ExcInternalError());
-  indices[0] = i;
+  indices[0] = index_map[i];
 }
 
 
@@ -40,9 +40,10 @@ compute_index (const unsigned int i,
 {
   const unsigned int n_pols = polynomials.size();
   Assert (i<n_pols*n_pols, ExcInternalError());
+  const unsigned int n=index_map[i];
 
-  indices[0] = i % n_pols;
-  indices[1] = i / n_pols;
+  indices[0] = n % n_pols;
+  indices[1] = n / n_pols;
 }
 
 
@@ -55,10 +56,42 @@ compute_index (const unsigned int i,
 {
   const unsigned int n_pols = polynomials.size();
   Assert (i<n_pols*n_pols*n_pols, ExcInternalError());
+  const unsigned int n=index_map[i];
+
+  indices[0] = n % n_pols;
+  indices[1] = (n/n_pols) % n_pols;
+  indices[2] = n / (n_pols*n_pols);
+}
+
+
+template <int dim>
+void
+TensorProductPolynomials<dim>::output_indices(std::ostream &out) const
+{
+  unsigned int ix[dim];
+  for (unsigned int i=0; i<n_tensor_pols; ++i)
+    {
+      compute_index(i,ix);
+      out << i << "\t";
+      for (unsigned int d=0; d<dim; ++d)
+       out << ix[d] << " ";
+      out << std::endl;
+    }
+}
+
+
+
+template <int dim>
+void
+TensorProductPolynomials<dim>::set_renumbering(
+  const std::vector<unsigned int> &renumber)
+{
+  Assert(renumber.size()==index_map.size(),
+        ExcDimensionMismatch(renumber.size(), index_map.size()));
 
-  indices[0] = i % n_pols;
-  indices[1] = (i/n_pols) % n_pols;
-  indices[2] = i / (n_pols*n_pols);
+  index_map=renumber;
+  for (unsigned int i=0; i<index_map.size(); ++i)
+    index_map_inverse[index_map[i]]=i;
 }
 
 

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.