]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
use SmartPointer instead of a reference in member variable
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 24 Oct 2010 17:32:21 +0000 (17:32 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 24 Oct 2010 17:32:21 +0000 (17:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@22456 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/mapping_q1_eulerian.h
deal.II/deal.II/include/fe/mapping_q_eulerian.h
deal.II/deal.II/source/fe/mapping_q1_eulerian.cc
deal.II/deal.II/source/fe/mapping_q_eulerian.cc

index f0ebbb09aa197f3eeb22ad1ae6feda9194548c27..bf64fa0b13b68d596b5da4f00b9a6f3661dda407 100644 (file)
@@ -130,20 +130,6 @@ class MappingQ1Eulerian : public MappingQ1<dim,spacedim>
                                      */
     bool preserves_vertex_locations () const;
 
-
-                                    /**
-                                     * Exception
-                                     */
-    DeclException0 (ExcWrongNoOfComponents);
-
-                                    /**
-                                     * Exception.
-                                     */
-    DeclException2 (ExcWrongVectorSize,
-                   int, int,
-                   << "Vector has wrong size " << arg1
-                   << ", expected size " << arg2);
-
                                     /**
                                      * Exception.
                                      */
@@ -174,14 +160,14 @@ class MappingQ1Eulerian : public MappingQ1<dim,spacedim>
                                      * Reference to the vector of
                                      * shifts.
                                      */
-    const VECTOR &euler_transform_vectors;
+    SmartPointer<const VECTOR, MappingQ1Eulerian<dim,VECTOR,spacedim> > euler_transform_vectors;
 
                                      /**
                                      * Pointer to the DoFHandler to
                                      * which the mapping vector is
                                      * associated.
                                      */
-    const SmartPointer<const DoFHandler<dim,spacedim>,MappingQ1Eulerian<dim,VECTOR,spacedim> > shiftmap_dof_handler;
+    SmartPointer<const DoFHandler<dim,spacedim>,MappingQ1Eulerian<dim,VECTOR,spacedim> > shiftmap_dof_handler;
 
 
   private:
index c0b23637b04ac1b9e1c46bd3eca3e29a4e21af67..9fd2999ff3d8e71fec6f50fba5b4a22385011f7e 100644 (file)
@@ -133,27 +133,12 @@ class MappingQEulerian : public MappingQ<dim, spacedim>
                                      * locations).
                                      */
     bool preserves_vertex_locations () const;
-
-                                     /**
-                                      * Exception
-                                      */
-
-    DeclException0 (ExcWrongNoOfComponents);
-
+    
                                     /**
                                       * Exception
                                       */
-
     DeclException0 (ExcInactiveCell);
-
-                                    /**
-                                      * Exception
-                                      */
-
-    DeclException2 (ExcWrongVectorSize, int, int,
-                   << "Vector has wrong size " << arg1
-                   << "-- expected size " << arg2);
-
+    
   protected:
                                     /**
                                      * Implementation of the interface in
@@ -178,7 +163,7 @@ class MappingQEulerian : public MappingQ<dim, spacedim>
                                       * shifts.
                                       */
 
-    const VECTOR &euler_vector;
+    SmartPointer<const VECTOR, MappingQEulerian<dim,VECTOR,spacedim> > euler_vector;
 
                                      /**
                                       * Pointer to the DoFHandler to
@@ -186,7 +171,7 @@ class MappingQEulerian : public MappingQ<dim, spacedim>
                                       * associated.
                                       */
 
-    const SmartPointer<const DoFHandler<dim,spacedim>,MappingQEulerian<dim,VECTOR,spacedim> > euler_dof_handler;
+    SmartPointer<const DoFHandler<dim,spacedim>,MappingQEulerian<dim,VECTOR,spacedim> > euler_dof_handler;
 
 
   private:
index 1f85b2c50c7487012132e105055c0f2189055ab7..3dd1292162955ba033cb0f205343ebe16dea9206 100644 (file)
@@ -28,8 +28,8 @@ MappingQ1Eulerian<dim, EulerVectorType, spacedim>::
 MappingQ1Eulerian (const EulerVectorType  &euler_transform_vectors,
                   const DoFHandler<dim,spacedim> &shiftmap_dof_handler)
                    :
-                  euler_transform_vectors(euler_transform_vectors),
-                  shiftmap_dof_handler(&shiftmap_dof_handler, typeid(*this).name())
+                  euler_transform_vectors(&euler_transform_vectors),
+                  shiftmap_dof_handler(&shiftmap_dof_handler)
 {}
 
 
@@ -48,14 +48,12 @@ compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_
                                   // *before* the mapping object is
                                   // constructed, which is not
                                   // necessarily what we want.
-  Assert (spacedim == shiftmap_dof_handler->get_fe().n_dofs_per_vertex(),
-          ExcWrongNoOfComponents());
-  Assert (shiftmap_dof_handler->get_fe().n_components() == spacedim,
-         ExcWrongNoOfComponents());
 
-  Assert (shiftmap_dof_handler->n_dofs() == euler_transform_vectors.size(),
-          ExcWrongVectorSize(euler_transform_vectors.size(),
-                            shiftmap_dof_handler->n_dofs()));
+//TODO: Only one of these two assertions should be relevant
+  AssertDimension (spacedim, shiftmap_dof_handler->get_fe().n_dofs_per_vertex());
+  AssertDimension(shiftmap_dof_handler->get_fe().n_components(), spacedim);
+
+  AssertDimension (shiftmap_dof_handler->n_dofs(), euler_transform_vectors->size());
 
                                   // cast the
                                   // Triangulation<dim>::cell_iterator
@@ -82,7 +80,7 @@ compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_
                                   // now get the values of the shift
                                   // vectors at the vertices
   Vector<double> mapping_values (shiftmap_dof_handler->get_fe().dofs_per_cell);
-  dof_cell->get_dof_values (euler_transform_vectors, mapping_values);
+  dof_cell->get_dof_values (*euler_transform_vectors, mapping_values);
 
 
   for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
@@ -143,6 +141,7 @@ MappingQ1Eulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
 
 // explicit instantiation
 template class MappingQ1Eulerian<deal_II_dimension, Vector<double> >;
+template class MappingQ1Eulerian<deal_II_dimension, Vector<float> >;
 #ifdef DEAL_II_USE_PETSC
 template class MappingQ1Eulerian<deal_II_dimension, PETScWrappers::Vector>;
 #endif
index 314270e7a61cdcd681a1ab61b8453a70802250b1..f0ded0c90b33a2c48a3236f9dc169fea3f84db5b 100644 (file)
@@ -36,8 +36,8 @@ MappingQEulerian (const unsigned int degree,
                  const DoFHandler<dim,spacedim> &euler_dof_handler)
                :
                MappingQ<dim,spacedim>(degree, true),
-               euler_vector(euler_vector),
-               euler_dof_handler(&euler_dof_handler, typeid(*this).name()),
+               euler_vector(&euler_vector),
+               euler_dof_handler(&euler_dof_handler),
                support_quadrature(degree),
                fe_values(euler_dof_handler.get_fe(),
                          support_quadrature,
@@ -51,7 +51,7 @@ Mapping<dim,spacedim> *
 MappingQEulerian<dim, EulerVectorType, spacedim>::clone () const
 {
   return new MappingQEulerian<dim,EulerVectorType,spacedim>(this->get_degree(),
-                                                  euler_vector,
+                                                  *euler_vector,
                                                   *euler_dof_handler);
 }
 
@@ -111,9 +111,9 @@ compute_mapping_support_points
                                   // with respect to vector size,
 
   const unsigned int n_dofs      = euler_dof_handler->n_dofs();
-  const unsigned int vector_size = euler_vector.size();
+  const unsigned int vector_size = euler_vector->size();
 
-  Assert (n_dofs == vector_size,ExcWrongVectorSize(vector_size,n_dofs));
+  AssertDimension(vector_size,n_dofs);
 
                                   // we then transform our tria iterator
                                   // into a dof iterator so we can
@@ -153,7 +153,7 @@ compute_mapping_support_points
   const unsigned int n_support_pts = support_quadrature.n_quadrature_points;
   const unsigned int n_components  = euler_dof_handler->get_fe().n_components();
 
-  Assert (n_components >= spacedim, ExcWrongNoOfComponents() );
+  Assert (n_components >= spacedim, ExcDimensionMismatch(n_components, spacedim) );
 
   std::vector<Vector<double> > shift_vector(n_support_pts,Vector<double>(n_components));
 
@@ -165,7 +165,7 @@ compute_mapping_support_points
                                   // threads
   Threads::ThreadMutex::ScopedLock lock(fe_values_mutex);
   fe_values.reinit(dof_cell);
-  fe_values.get_function_values(euler_vector,shift_vector);
+  fe_values.get_function_values(*euler_vector, shift_vector);
 
                                   // and finally compute the positions of the
                                   // support points in the deformed

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.