]> https://gitweb.dealii.org/ - dealii.git/commitdiff
compute_fill fixed
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 30 Jan 2004 14:26:09 +0000 (14:26 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 30 Jan 2004 14:26:09 +0000 (14:26 +0000)
git-svn-id: https://svn.dealii.org/trunk@8366 0785d39b-7218-0410-832d-ea1e28bc413d

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

index c6a947f65d033f2881f7ac9f14331a1d4f0f80b5..245404d652d2a20d44632ca17236e149b90a2fbc 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003 by the deal.II authors
+//    Copyright (C) 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,8 @@
  * Mapping of an axis-parallel cell.
  *
  * This class maps the unit cell to a grid cell with surfaces parallel
- * to the coordinate lines/planes. It is specifically developed for
+ * to the coordinate lines/planes. The mapping is therefore a scaling
+ * along the coordinate directions. It is specifically developed for
  * cartesian meshes. Apply this mapping to a general mesh to get
  * strange results.
  *
@@ -37,37 +38,21 @@ template <int dim>
 class MappingCartesian : public Mapping<dim>
 {
   public:
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual
     typename Mapping<dim>::InternalDataBase *
     get_data (const UpdateFlags,
              const Quadrature<dim>& quadrature) const;
 
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual
     typename Mapping<dim>::InternalDataBase *
     get_face_data (const UpdateFlags flags,
                   const Quadrature<dim-1>& quadrature) const;
 
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual
     typename Mapping<dim>::InternalDataBase *
     get_subface_data (const UpdateFlags flags,
                      const Quadrature<dim-1>& quadrature) const;
 
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     fill_fe_values (const typename DoFHandler<dim>::cell_iterator &cell,
                    const Quadrature<dim>& quadrature,
@@ -75,10 +60,6 @@ class MappingCartesian : public Mapping<dim>
                    std::vector<Point<dim> >        &quadrature_points,
                    std::vector<double>             &JxW_values) const ;
 
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     fill_fe_face_values (const typename DoFHandler<dim>::cell_iterator &cell,
                         const unsigned int face_no,
@@ -88,11 +69,6 @@ class MappingCartesian : public Mapping<dim>
                         std::vector<double>             &JxW_values,
                         std::vector<Tensor<1,dim> >        &boundary_form,
                         std::vector<Point<dim> >        &normal_vectors) const ;
-
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     fill_fe_subface_values (const typename DoFHandler<dim>::cell_iterator &cell,
                            const unsigned int face_no,
@@ -104,53 +80,30 @@ class MappingCartesian : public Mapping<dim>
                            std::vector<Tensor<1,dim> >        &boundary_form,
                            std::vector<Point<dim> >        &normal_vectors) const ;
 
-
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     transform_covariant (Tensor<1,dim>          *begin,
                         Tensor<1,dim>          *end,
                         const Tensor<1,dim>    *src,
                         const typename Mapping<dim>::InternalDataBase &internal) const;
 
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     transform_covariant (Tensor<2,dim>          *begin,
                         Tensor<2,dim>          *end,
                         const Tensor<2,dim>    *src,
                         const typename Mapping<dim>::InternalDataBase &internal) const;
     
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     transform_contravariant (Tensor<1,dim>          *begin,
                             Tensor<1,dim>          *end,
                             const Tensor<1,dim>    *src,
                             const typename Mapping<dim>::InternalDataBase &internal) const;
     
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
     virtual void
     transform_contravariant (Tensor<2,dim>          *begin,
                             Tensor<2,dim>          *end,
                             const Tensor<2,dim>    *src,
                             const typename Mapping<dim>::InternalDataBase &internal) const;
 
-                                    /**
-                                     * Transforms the point @p{p} on
-                                     * the unit cell to the point
-                                     * @p{p_real} on the real cell
-                                     * @p{cell} and returns @p{p_real}.
-                                     */
     virtual Point<dim>
     transform_unit_to_real_cell (
       const typename Triangulation<dim>::cell_iterator &cell,
@@ -171,16 +124,7 @@ class MappingCartesian : public Mapping<dim>
       const typename Triangulation<dim>::cell_iterator &cell,
       const Point<dim>                                 &p) const;
     
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
-    virtual UpdateFlags update_once (const UpdateFlags) const;
-    
-                                    /**
-                                     * Implementation of the interface in
-                                     * @ref{Mapping}.
-                                     */
+    virtual UpdateFlags update_once (const UpdateFlags) const;    
     virtual UpdateFlags update_each (const UpdateFlags) const;
     
                                     /**
@@ -195,7 +139,7 @@ class MappingCartesian : public Mapping<dim>
   protected:
                                     /** 
                                      * Storage for internal data of
-                                     * d-linear transformation.
+                                     * the scaling.
                                      */
     class InternalData : public Mapping<dim>::InternalDataBase
     {
@@ -216,22 +160,22 @@ class MappingCartesian : public Mapping<dim>
                                         /**
                                          * Length of the cell in
                                          * different coordinate
-                                         * directions, @p{h_x},
-                                         * @p{h_y}, @p{h_z}.
+                                         * directions, <i>h<sub>x</sub></i>,
+                                         * <i>h<sub>y</sub></i>, <i>h<sub>z</sub></i>.
                                          */
        Tensor<1,dim> length;
 
                                         /**
                                          * Vector of all quadrature
                                          * points. Especially, all
-                                         * points of all faces.
+                                         * points on all faces.
                                          */
        std::vector<Point<dim> > quadrature_points;
     };
     
                                     /**
                                      * Do the computation for the
-                                     * @p{fill_*} functions.
+                                     * <tt>fill_*</tt> functions.
                                      */
     void compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
                       const unsigned int face_no,
@@ -253,6 +197,7 @@ class MappingCartesian : public Mapping<dim>
 
 /* -------------- declaration of explicit specializations ------------- */
 
+/// @if NoDoc
 
 template <> void MappingCartesian<1>::fill_fe_face_values (
   const DoFHandler<1>::cell_iterator &,
@@ -275,6 +220,6 @@ template <> void MappingCartesian<1>::fill_fe_subface_values (
   std::vector<Tensor<1,1> >&,
   std::vector<Point<1> >&) const;
 
-  
+/// @endif  
 
 #endif
index c750a1657729600d53e78cdfa0627507deab6eda..434ca0fdca30995da027b1f080f7e44e160e995f 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+//    Copyright (C) 1998, 1999, 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
@@ -136,26 +136,6 @@ MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterat
 {
   const UpdateFlags update_flags(data.current_update_flags());
 
-  const unsigned int npts = quadrature_points.size ();
-
-  const typename QProjector<dim>::DataSetDescriptor offset
-    = (face_no == invalid_face_number
-       ?
-       QProjector<dim>::DataSetDescriptor::cell()
-       :
-       (sub_no == invalid_face_number
-        ?
-                                        // called from FEFaceValues
-        QProjector<dim>::DataSetDescriptor::face (face_no,
-                                                  cell->face_orientation(face_no),
-                                                  quadrature_points.size())
-        :
-                                        // called from FESubfaceValues
-        QProjector<dim>::DataSetDescriptor::sub_face (face_no, sub_no,
-                                                      cell->face_orientation(face_no),
-                                                      quadrature_points.size())
-        ));
-  
                                    // some more sanity checks
   if (face_no != invalid_face_number)
     {
@@ -213,9 +193,28 @@ MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterat
                                   // each direction
   if (update_flags & update_q_points)
     {
-      Assert (quadrature_points.size() == npts,
-             ExcDimensionMismatch(quadrature_points.size(), npts));
-      for (unsigned int i=0; i<npts; ++i)
+      const typename QProjector<dim>::DataSetDescriptor offset
+       = (face_no == invalid_face_number
+          ?
+          QProjector<dim>::DataSetDescriptor::cell()
+          :
+          (sub_no == invalid_face_number
+           ?
+                                            // called from FEFaceValues
+           QProjector<dim>::DataSetDescriptor::face (face_no,
+                                                     cell->face_orientation(face_no),
+                                                     quadrature_points.size())
+           :
+                                            // called from FESubfaceValues
+           QProjector<dim>::DataSetDescriptor::sub_face (face_no, sub_no,
+                                                         cell->face_orientation(face_no),
+                                                         quadrature_points.size())
+          ));
+
+//This assertion was nonsense, since we let npts = quadrature_points.size() above
+//      Assert (quadrature_points.size() == npts,
+//           ExcDimensionMismatch(quadrature_points.size(), npts));
+      for (unsigned int i=0; i<quadrature_points.size(); ++i)
        {
          quadrature_points[i] = start;
          for (unsigned int d=0; d<dim; ++d)
@@ -234,8 +233,11 @@ MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterat
                                   // value
   if (update_flags & update_normal_vectors)
     {
-      Assert (normal_vectors.size() == npts,
-             ExcDimensionMismatch(normal_vectors.size(), npts));
+// We would like to have an assertion like this, but if
+// we do not compute quadrature points, npts is zero.
+      
+//      Assert (normal_vectors.size() == npts,
+//           ExcDimensionMismatch(normal_vectors.size(), npts));
       Assert (face_no < GeometryInfo<dim>::faces_per_cell,
               ExcInternalError());
       

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.