]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
remove n_quadrature_points member variable after being deprecated for nearly 3 years
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 11 Nov 2010 20:46:11 +0000 (20:46 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 11 Nov 2010 20:46:11 +0000 (20:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@22694 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/base/quadrature.h
deal.II/source/base/quadrature.cc
deal.II/source/fe/fe_nedelec.cc
deal.II/source/fe/fe_q.cc
deal.II/source/fe/mapping_q_eulerian.cc
deal.II/source/numerics/data_out.cc
deal.II/source/numerics/data_out_faces.cc
deal.II/source/numerics/data_out_rotation.cc
deal.II/source/numerics/data_out_stack.cc
deal.II/source/numerics/error_estimator.cc

index fa3829aa8efdbdeb5f08176df8aa437323be9f2e..d99b5a7c421bf6046c01a369a80c9dbd0b3671e2 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -90,26 +90,6 @@ class Quadrature : public Subscriptor
                                      */
     typedef Quadrature<dim-1> SubQuadrature;
 
-                                    /**
-                                     * @deprecated Use size()
-                                     * instead.
-                                     *
-                                     * Number of quadrature points.
-                                     *
-                                     * @warning After introduction of
-                                     * the assignment operator, this
-                                     * number is not constant anymore
-                                     * and erroneous assignment can
-                                     * compromise integrity of the
-                                     * Quadrature object. Since
-                                     * direct data access should be
-                                     * considered a design flaw
-                                     * anyway, it is strongly
-                                     * suggested to use size()
-                                     * instead.
-                                     */
-    unsigned int n_quadrature_points;
-
                                     /**
                                      * Constructor.
                                      *
index 2d7de09456f32fa6d26de89e8b2b630786e9685c..764e145b1654907dcac7a794cf8276d32d7380ea 100644 (file)
@@ -40,11 +40,12 @@ namespace
 }
 
 
-
+//TODO: It would be desirable to have a Tensor<rank,0>
 template <>
 Quadrature<0>::Quadrature (const unsigned int)
-               : n_quadrature_points(1),
-                 weights (1, 1.)
+               :
+//             quadrature_points(1),
+               weights (1, 1.)
 {}
 
 
@@ -57,7 +58,6 @@ Quadrature<0>::~Quadrature ()
 
 template <int dim>
 Quadrature<dim>::Quadrature (const unsigned int n_q) :
-               n_quadrature_points(n_q),
                quadrature_points (n_q, Point<dim>()),
                weights (n_q, 0)
 {}
@@ -68,7 +68,6 @@ template <int dim>
 Quadrature<dim>::Quadrature (const std::vector<Point<dim> > &points,
                             const std::vector<double>      &weights)
                :
-               n_quadrature_points(points.size()),
                quadrature_points(points),
                weights(weights)
 {
@@ -81,7 +80,6 @@ Quadrature<dim>::Quadrature (const std::vector<Point<dim> > &points,
 template <int dim>
 Quadrature<dim>::Quadrature (const std::vector<Point<dim> > &points)
                :
-               n_quadrature_points(points.size()),
                quadrature_points(points),
                weights(points.size(), std::atof("Inf"))
 {
@@ -94,19 +92,17 @@ Quadrature<dim>::Quadrature (const std::vector<Point<dim> > &points)
 template <int dim>
 Quadrature<dim>::Quadrature (const Point<dim> &point)
                 :
-               n_quadrature_points(1),
                quadrature_points(std::vector<Point<dim> > (1, point)),
                weights(std::vector<double> (1, 1.))
 {}
 
 
 template <>
-Quadrature<0>::Quadrature (const Quadrature<-1> &,
-                          const Quadrature<1> &)
-               :
-                n_quadrature_points (1),
-               weights(1, 1.)
-{}
+Quadrature<0>::Quadrature (const SubQuadrature&,
+                          const Quadrature<1>&)
+{
+  Assert(false, ExcImpossibleInDim(0));
+}
 
 
 
@@ -114,10 +110,8 @@ template <int dim>
 Quadrature<dim>::Quadrature (const SubQuadrature &q1,
                             const Quadrature<1> &q2)
                :
-               n_quadrature_points (q1.size() *
-                                    q2.size()),
-               quadrature_points (n_quadrature_points),
-               weights (n_quadrature_points, 0)
+               quadrature_points (q1.size() * q2.size()),
+               weights (q1.size() * q2.size())
 {
   unsigned int present_index = 0;
   for (unsigned int i2=0; i2<q2.size(); ++i2)
@@ -157,9 +151,8 @@ template <>
 Quadrature<1>::Quadrature (const SubQuadrature&,
                           const Quadrature<1>& q2)
                :
-               n_quadrature_points (q2.size()),
-               quadrature_points (n_quadrature_points),
-               weights (n_quadrature_points, 0)
+               quadrature_points (q2.size()),
+               weights (q2.size())
 {
   unsigned int present_index = 0;
   for (unsigned int i2=0; i2<q2.size(); ++i2)
@@ -194,22 +187,21 @@ Quadrature<1>::Quadrature (const SubQuadrature&,
 template <>
 Quadrature<0>::Quadrature (const Quadrature<1> &)
                :
-               n_quadrature_points (1),
-               weights (1, 1.)
+               Subscriptor(),
+//             quadrature_points(1),
+               weights(1,1.)
 {}
 
 
 template <>
 Quadrature<1>::Quadrature (const Quadrature<0> &)
                :
-               n_quadrature_points (numbers::invalid_unsigned_int),
-               quadrature_points (),
-               weights ()
+               Subscriptor()
 {
                                    // this function should never be
                                    // called -- this should be the
                                    // copy constructor in 1d...
-  Assert (false, ExcInternalError());
+  Assert (false, ExcImpossibleInDim(1));
 }
 
 
@@ -218,9 +210,8 @@ template <int dim>
 Quadrature<dim>::Quadrature (const Quadrature<dim != 1 ? 1 : 0> &q)
                :
                Subscriptor(),
-               n_quadrature_points (dimpow<dim>(q.size())),
-               quadrature_points (n_quadrature_points),
-               weights (n_quadrature_points, 0.)
+               quadrature_points (dimpow<dim>(q.size())),
+               weights (dimpow<dim>(q.size()))
 {
   Assert (dim <= 3, ExcNotImplemented());
   
@@ -253,7 +244,6 @@ template <int dim>
 Quadrature<dim>::Quadrature (const Quadrature<dim> &q)
                :
                Subscriptor(),
-               n_quadrature_points (q.size()),
                quadrature_points (q.quadrature_points),
                weights (q.weights)
 {}
@@ -265,7 +255,6 @@ Quadrature<dim>::operator= (const Quadrature<dim>& q)
 {
   weights = q.weights;
   quadrature_points = q.quadrature_points;
-  n_quadrature_points = q.size();
   return *this;
 }
 
index c024ab8133841566ac604ddfe67f7a0ee7a440ed..f744500769673bb0cf5a2f659339e7319057e200 100644 (file)
@@ -23,6 +23,7 @@
 
 DEAL_II_NAMESPACE_OPEN
 
+//#define DEBUG_NEDELEC
 
 template <int dim>
 FE_Nedelec<dim>::FE_Nedelec (const unsigned int p) :
@@ -36,6 +37,10 @@ FE_PolyTensor<PolynomialsNedelec<dim>, dim>
   std::vector<bool> (dim, true))),
 deg (p)
 {
+#ifdef DEBUG_NEDELEC
+  deallog << get_name() << std::endl;
+#endif
+  
   Assert (dim >= 2, ExcImpossibleInDim(dim));
 
   const unsigned int n_dofs = this->dofs_per_cell;
@@ -58,11 +63,20 @@ deg (p)
                                   // matrices to the right sizes.
                                   // Restriction only for isotropic
                                   // refinement
+#ifdef DEBUG_NEDELEC
+  deallog << "Embedding" << std::endl;
+#endif
   this->reinit_restriction_and_prolongation_matrices ();
                                   // Fill prolongation matrices with embedding operators
   FETools::compute_embedding_matrices (*this, this->prolongation, true);
+#ifdef DEBUG_NEDELEC
+  deallog << "Restriction" << std::endl;
+#endif
   initialize_restriction ();
-
+  
+#ifdef DEBUG_NEDELEC
+  deallog << "Face Embedding" << std::endl;
+#endif
   FullMatrix<double> face_embeddings[GeometryInfo<dim>::max_children_per_face];
 
   for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
index c1d4f0ddfad3b42941254215aca54946535c4a08..8bdede7b5183dd007027394ddbf66145686f23e6 100644 (file)
@@ -588,14 +588,14 @@ FE_Q<dim,spacedim>::FE_Q (const Quadrature<1> &points)
                :
                FE_Poly<TensorProductPolynomials<dim>, dim, spacedim> (
                  TensorProductPolynomials<dim>(Polynomials::generate_complete_Lagrange_basis(points.get_points())),
-                 FiniteElementData<dim>(get_dpo_vector(points.n_quadrature_points-1),
-                                        1, points.n_quadrature_points-1,
+                 FiniteElementData<dim>(get_dpo_vector(points.size()-1),
+                                        1, points.size()-1,
                                         FiniteElementData<dim>::H1),
                  std::vector<bool> (1, false),
                  std::vector<std::vector<bool> >(1, std::vector<bool>(1,true))),
-               face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (points.n_quadrature_points-1)))
+               face_index_map(FE_Q_Helper::invert_numbering(face_lexicographic_to_hierarchic_numbering (points.size()-1)))
 {
-  const unsigned int degree = points.n_quadrature_points-1;
+  const unsigned int degree = points.size()-1;
 
   Assert (degree > 0,
           ExcMessage ("This element can only be used for polynomial degrees "
index 6a1e9010cdb2a86189451510584798733c732f59..6f46d41c1666081eb923ca07b64760b746321abf 100644 (file)
@@ -72,7 +72,7 @@ SupportQuadrature (const unsigned int map_degree)
                                   // trapezoidal quadrature rule.
   const QTrapez<1> q1d;
   const QIterated<dim> q_iterated(q1d,map_degree);
-  const unsigned int n_q_points = q_iterated.n_quadrature_points;
+  const unsigned int n_q_points = q_iterated.size();
 
                                   // we then need to define a renumbering
                                   // vector that allows us to go from a
@@ -150,7 +150,7 @@ compute_mapping_support_points
                                   // components appropriately, or create a
                                   // separate dof handler for the displacements.
 
-  const unsigned int n_support_pts = support_quadrature.n_quadrature_points;
+  const unsigned int n_support_pts = support_quadrature.size();
   const unsigned int n_components  = euler_dof_handler->get_fe().n_components();
 
   Assert (n_components >= spacedim, ExcDimensionMismatch(n_components, spacedim) );
index 9a346cbf8255fe6f0ff196b4c8a72552aa79af87..d72ff30abb18038e1a5acfd12abf792ddb6ffeae 100644 (file)
@@ -57,7 +57,7 @@ namespace internal
                    ParallelDataBase<dim,spacedim> (n_components,
                                                    n_datasets,
                                                    n_subdivisions,
-                                                   quadrature.n_quadrature_points,
+                                                   quadrature.size(),
                                                    n_postprocessor_outputs,
                                                    finite_elements),
                    q_collection (quadrature),
@@ -1077,7 +1077,7 @@ void DataOut<dim,DH>::build_patches (const Mapping<DH::dimension,DH::space_dimen
 
   DataOutBase::Patch<DH::dimension, DH::space_dimension> sample_patch;
   sample_patch.n_subdivisions = n_subdivisions;
-  sample_patch.data.reinit (n_datasets, patch_points.n_quadrature_points);
+  sample_patch.data.reinit (n_datasets, patch_points.size());
 
 
 
index 8760e59d5686d077dd94ab1a1accaa7cc49f55a5..cad21ea6b09acc0f692b069383ee3f6a104cd3f1 100644 (file)
@@ -49,7 +49,7 @@ namespace internal
                  ParallelDataBase<dim,spacedim> (n_components,
                                                  n_datasets,
                                                  n_subdivisions,
-                                                 quadrature.n_quadrature_points,
+                                                 quadrature.size(),
                                                  n_postprocessor_outputs,
                                                  finite_elements),
                  q_collection (quadrature),
@@ -357,7 +357,7 @@ void DataOutFaces<dim,DH>::build_patches (const Mapping<DH::dimension> &mapping,
   DataOutBase::Patch<DH::dimension-1,DH::space_dimension> sample_patch;
   sample_patch.n_subdivisions = n_subdivisions;
   sample_patch.data.reinit (n_datasets,
-                           patch_points.n_quadrature_points);
+                           patch_points.size());
 
                                   // now build the patches in parallel
   WorkStream::run (&all_faces[0],
index 94a0e8f58715e5ca709c481996263fa22e82cccf..4d14e22165762e00ec732eaf30d102a0548b0f6f 100644 (file)
@@ -62,7 +62,7 @@ namespace internal
                    ParallelDataBase<dim,spacedim> (n_components,
                                                    n_datasets,
                                                    n_subdivisions,
-                                                   quadrature.n_quadrature_points,
+                                                   quadrature.size(),
                                                    n_postprocessor_outputs,
                                                    finite_elements),
                    n_patches_per_circle (n_patches_per_circle),
@@ -507,7 +507,7 @@ void DataOutRotation<dim,DH>::build_patches (const unsigned int n_patches_per_ci
     {
       new_patches[i].n_subdivisions = n_subdivisions;
       new_patches[i].data.reinit (n_datasets,
-                                 patch_points.n_quadrature_points
+                                 patch_points.size()
                                  * (n_subdivisions+1));
     }
 
index 8c24e7c3a9312ab5b263ee0a9d71370f246bec78..59ed4fcc33758b1d715ce09d87650a9547c97034 100644 (file)
@@ -279,7 +279,7 @@ void DataOutStack<dim,spacedim,DH>::build_patches (const unsigned int nnnn_subdi
   hp::FEValues<dim> x_fe_patch_values (fe_collection, q_collection,
                                        update_values);
 
-  const unsigned int n_q_points = patch_points.n_quadrature_points;
+  const unsigned int n_q_points = patch_points.size();
   std::vector<double>          patch_values (n_q_points);
   std::vector<Vector<double> > patch_values_system (n_q_points,
                                                    Vector<double>(n_components));
index e2efdb979559fac1a91a7966a4202b806e2669df..55106b8c3a28bd80be50a68fd147b085608e36f9 100644 (file)
@@ -330,7 +330,7 @@ namespace internal
     void
     ParallelData<DH>::resize (const unsigned int active_fe_index)
     {
-      const unsigned int n_q_points   = face_quadratures[active_fe_index].n_quadrature_points;
+      const unsigned int n_q_points   = face_quadratures[active_fe_index].size();
       const unsigned int n_components = finite_element.n_components();
 
       normal_vectors.resize(n_q_points);
@@ -420,8 +420,7 @@ namespace internal
       const unsigned int dim = DH::dimension;
 
       const typename DH::face_iterator face = cell->face(face_no);
-      const unsigned int n_q_points         = parallel_data.face_quadratures[cell->active_fe_index()]
-                                             .n_quadrature_points,
+      const unsigned int n_q_points         = parallel_data.face_quadratures[cell->active_fe_index()].size(),
                         n_components       = parallel_data.finite_element.n_components(),
                         n_solution_vectors = solutions.size();
 
@@ -632,8 +631,7 @@ namespace internal
       const unsigned int dim = DH::dimension;
 
       const typename DH::cell_iterator neighbor = cell->neighbor(face_no);
-      const unsigned int n_q_points         = parallel_data.face_quadratures[cell->active_fe_index()]
-                                             .n_quadrature_points,
+      const unsigned int n_q_points         = parallel_data.face_quadratures[cell->active_fe_index()].size(),
                         n_components       = parallel_data.finite_element.n_components(),
                         n_solution_vectors = solutions.size();
       const typename DH::face_iterator

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.