]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Revise internal functions of MappingQ1.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 27 Jul 2015 18:51:30 +0000 (13:51 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 29 Jul 2015 22:54:25 +0000 (17:54 -0500)
The compute_fill* functions now also take FEValuesData as arguments,
just like the fill_fe_values() functions.

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

index d4b01051b39425920db3cdb2a2e98e8a1726b862..6fc969415f7389d95b2563d7bb77d10b68dada97 100644 (file)
@@ -740,7 +740,10 @@ private:
    *   return value of this function).
    * @param[in] quadrature A reference to the quadrature formula in use
    *   for the current evaluation. This quadrature object is the same
-   *   as the one used when creating the @p internal_data object
+   *   as the one used when creating the @p internal_data object. The
+   *   object is used both to map the location of quadrature points,
+   *   as well as to compute the JxW values for each quadrature
+   *   point (which involves the quadrature weights).
    * @param[in] internal_data A reference to an object previously
    *   created by get_data() and that may be used to store information
    *   the mapping can compute once on the reference cell. See the
@@ -789,7 +792,10 @@ private:
    *   information is requested.
    * @param[in] quadrature A reference to the quadrature formula in use
    *   for the current evaluation. This quadrature object is the same
-   *   as the one used when creating the @p internal_data object
+   *   as the one used when creating the @p internal_data object. The
+   *   object is used both to map the location of quadrature points,
+   *   as well as to compute the JxW values for each quadrature
+   *   point (which involves the quadrature weights).
    * @param[in] internal_data A reference to an object previously
    *   created by get_data() and that may be used to store information
    *   the mapping can compute once on the reference cell. See the
@@ -822,7 +828,10 @@ private:
    *   given cell for which information is requested.
    * @param[in] quadrature A reference to the quadrature formula in use
    *   for the current evaluation. This quadrature object is the same
-   *   as the one used when creating the @p internal_data object
+   *   as the one used when creating the @p internal_data object. The
+   *   object is used both to map the location of quadrature points,
+   *   as well as to compute the JxW values for each quadrature
+   *   point (which involves the quadrature weights).
    * @param[in] internal_data A reference to an object previously
    *   created by get_data() and that may be used to store information
    *   the mapping can compute once on the reference cell. See the
index d3eb9006a9ba073be6018011ed4842e7690840c0..f2b6f0524250f33989725b01b90f779ccc462ac4 100644 (file)
@@ -405,13 +405,8 @@ public:
                           const unsigned int                npts,
                           const DataSetDescriptor           data_set,
                           const std::vector<double>        &weights,
-                          const InternalData                     &mapping_data,
-                          std::vector<Point<spacedim> >    &quadrature_points,
-                          std::vector<double>              &JxW_values,
-                          std::vector<Tensor<1,spacedim> > &boundary_form,
-                          std::vector<Point<spacedim> >    &normal_vectors,
-                          std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
-                          std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const;
+                          const InternalData                     &internal_data,
+                          FEValuesData<dim,spacedim>                                &output_data) const;
 
   /**
    * Compute shape values and/or derivatives.
index f3bcb5394f1a021e0bf2d8df8f20b6457efabdc2..034c0540c83c1813054cd7b7c6c646a105679e9f 100644 (file)
@@ -343,12 +343,7 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
                                  n_q_points),
                            quadrature.get_weights(),
                            *p_data,
-                           output_data.quadrature_points,
-                           output_data.JxW_values,
-                           output_data.boundary_forms,
-                           output_data.normal_vectors,
-                           output_data.jacobians,
-                           output_data.inverse_jacobians);
+                           output_data);
 }
 
 
@@ -398,12 +393,7 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
                                     cell->subface_case(face_no)),
                            quadrature.get_weights(),
                            *p_data,
-                           output_data.quadrature_points,
-                           output_data.JxW_values,
-                           output_data.boundary_forms,
-                           output_data.normal_vectors,
-                           output_data.jacobians,
-                           output_data.inverse_jacobians);
+                           output_data);
 }
 
 
index c4cb6be7fb1d13566ee6f7b3b5efc658343baf14..b5bc80917d55d5ced5d39c1e3f975adda03373ed 100644 (file)
@@ -924,21 +924,17 @@ namespace internal
                        const unsigned int               n_q_points,
                        const std::vector<double>        &weights,
                        const typename dealii::MappingQ1<dim,spacedim>::InternalData &data,
-                       std::vector<double>              &JxW_values,
-                       std::vector<Tensor<1,spacedim> > &boundary_forms,
-                       std::vector<Point<spacedim> >    &normal_vectors,
-                       std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
-                       std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians)
+                       FEValuesData<dim,spacedim>                                   &output_data)
     {
       const UpdateFlags update_flags(data.current_update_flags());
 
       if (update_flags & update_boundary_forms)
         {
-          AssertDimension (boundary_forms.size(), n_q_points);
+          AssertDimension (output_data.boundary_forms.size(), n_q_points);
           if (update_flags & update_normal_vectors)
-            AssertDimension (normal_vectors.size(), n_q_points);
+            AssertDimension (output_data.normal_vectors.size(), n_q_points);
           if (update_flags & update_JxW_values)
-            AssertDimension (JxW_values.size(), n_q_points);
+            AssertDimension (output_data.JxW_values.size(), n_q_points);
 
           // map the unit tangentials to the real cell. checking for d!=dim-1
           // eliminates compiler warnings regarding unsigned int expressions <
@@ -970,14 +966,14 @@ namespace internal
                     // fields (because it has only dim-1 components), but we
                     // can still compute the boundary form by simply
                     // looking at the number of the face
-                    boundary_forms[i][0] = (face_no == 0 ?
-                                            -1 : +1);
+                    output_data.boundary_forms[i][0] = (face_no == 0 ?
+                                                        -1 : +1);
                     break;
                   case 2:
-                    cross_product (boundary_forms[i], data.aux[0][i]);
+                    cross_product (output_data.boundary_forms[i], data.aux[0][i]);
                     break;
                   case 3:
-                    cross_product (boundary_forms[i], data.aux[0][i], data.aux[1][i]);
+                    cross_product (output_data.boundary_forms[i], data.aux[0][i], data.aux[1][i]);
                     break;
                   default:
                     Assert(false, ExcNotImplemented());
@@ -998,9 +994,9 @@ namespace internal
                   if (dim==1)
                     {
                       // J is a tangent vector
-                      boundary_forms[point] = data.contravariant[point].transpose()[0];
-                      boundary_forms[point] /=
-                        (face_no == 0 ? -1. : +1.) * boundary_forms[point].norm();
+                      output_data.boundary_forms[point] = data.contravariant[point].transpose()[0];
+                      output_data.boundary_forms[point] /=
+                        (face_no == 0 ? -1. : +1.) * output_data.boundary_forms[point].norm();
 
                     }
 
@@ -1014,7 +1010,7 @@ namespace internal
 
                       // then compute the face normal from the face tangent
                       // and the cell normal:
-                      cross_product (boundary_forms[point],
+                      cross_product (output_data.boundary_forms[point],
                                      data.aux[0][point], cell_normal);
 
                     }
@@ -1026,31 +1022,32 @@ namespace internal
 
           if (update_flags & (update_normal_vectors
                               | update_JxW_values))
-            for (unsigned int i=0; i<boundary_forms.size(); ++i)
+            for (unsigned int i=0; i<output_data.boundary_forms.size(); ++i)
               {
                 if (update_flags & update_JxW_values)
                   {
-                    JxW_values[i] = boundary_forms[i].norm() * weights[i];
+                    output_data.JxW_values[i] = output_data.boundary_forms[i].norm() * weights[i];
 
                     if (subface_no!=numbers::invalid_unsigned_int)
                       {
                         const double area_ratio=GeometryInfo<dim>::subface_ratio(
                                                   cell->subface_case(face_no), subface_no);
-                        JxW_values[i] *= area_ratio;
+                        output_data.JxW_values[i] *= area_ratio;
                       }
                   }
 
                 if (update_flags & update_normal_vectors)
-                  normal_vectors[i] = Point<spacedim>(boundary_forms[i] / boundary_forms[i].norm());
+                  output_data.normal_vectors[i] = Point<spacedim>(output_data.boundary_forms[i] /
+                                                                  output_data.boundary_forms[i].norm());
               }
 
           if (update_flags & update_jacobians)
             for (unsigned int point=0; point<n_q_points; ++point)
-              jacobians[point] = data.contravariant[point];
+              output_data.jacobians[point] = data.contravariant[point];
 
           if (update_flags & update_inverse_jacobians)
             for (unsigned int point=0; point<n_q_points; ++point)
-              inverse_jacobians[point] = data.covariant[point].transpose();
+              output_data.inverse_jacobians[point] = data.covariant[point].transpose();
         }
     }
   }
@@ -1061,28 +1058,23 @@ namespace internal
 
 template<int dim, int spacedim>
 void
-MappingQ1<dim,spacedim>::compute_fill_face (
-  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const unsigned int               face_no,
-  const unsigned int               subface_no,
-  const unsigned int               n_q_points,
-  const DataSetDescriptor          data_set,
-  const std::vector<double>        &weights,
-  const InternalData                     &data,
-  std::vector<Point<spacedim> >    &quadrature_points,
-  std::vector<double>              &JxW_values,
-  std::vector<Tensor<1,spacedim> > &boundary_forms,
-  std::vector<Point<spacedim> >    &normal_vectors,
-  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
-  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
+MappingQ1<dim,spacedim>::
+compute_fill_face (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                   const unsigned int               face_no,
+                   const unsigned int               subface_no,
+                   const unsigned int               n_q_points,
+                   const DataSetDescriptor          data_set,
+                   const std::vector<double>        &weights,
+                   const InternalData                     &internal_data,
+                   FEValuesData<dim,spacedim>                                &output_data) const
 {
   compute_fill (cell, n_q_points, data_set, CellSimilarity::none,
-                data, quadrature_points);
+                internal_data,
+                output_data.quadrature_points);
   internal::compute_fill_face (*this,
                                cell, face_no, subface_no, n_q_points,
-                               weights, data,
-                               JxW_values, boundary_forms, normal_vectors,
-                               jacobians, inverse_jacobians);
+                               weights, internal_data,
+                               output_data);
 }
 
 
@@ -1112,12 +1104,7 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
                                               n_q_points),
                      quadrature.get_weights(),
                      data,
-                     output_data.quadrature_points,
-                     output_data.JxW_values,
-                     output_data.boundary_forms,
-                     output_data.normal_vectors,
-                     output_data.jacobians,
-                     output_data.inverse_jacobians);
+                     output_data);
 }
 
 
@@ -1150,12 +1137,7 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
                                                  cell->subface_case(face_no)),
                      quadrature.get_weights(),
                      data,
-                     output_data.quadrature_points,
-                     output_data.JxW_values,
-                     output_data.boundary_forms,
-                     output_data.normal_vectors,
-                     output_data.jacobians,
-                     output_data.inverse_jacobians);
+                     output_data);
 }
 
 

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.