]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fill jacobians in FEFaceValues. 465/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 23 Jan 2015 09:38:10 +0000 (10:38 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 23 Jan 2015 09:40:31 +0000 (10:40 +0100)
The previous code did not fill the Jacobian transformations on
FEFaceValues and FESubfaceValues. This is now fixed for the various
mappings.

13 files changed:
doc/news/changes.h
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_cartesian.h
include/deal.II/fe/mapping_q.h
include/deal.II/fe/mapping_q1.h
source/fe/fe_values.cc
source/fe/mapping_cartesian.cc
source/fe/mapping_q.cc
source/fe/mapping_q1.cc
tests/fe/jacobians_face.cc [new file with mode: 0644]
tests/fe/jacobians_face.output [new file with mode: 0644]
tests/fe/jacobians_face_cartesian.cc [new file with mode: 0644]
tests/fe/jacobians_face_cartesian.output [new file with mode: 0644]

index 9a5f4ff10c34b443d1005128ebda427aa195b1ef..47b6ba62e4d0e46eed38d5312e097024783733e5 100644 (file)
@@ -254,6 +254,13 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> Fixed: FEFaceValues and FESubfaceValues did not fill the
+  jacobians and inverse jacobians if requested via the update flags.
+  This is now fixed.
+  <br>
+  (Martin Kronbichler, 2015/01/23)
+  </li>
+
   <li> Fixed: ParameterHandler::read_input() now checks that
   'subsection'/'end' are balanced in the input.
   <br>
index 713d579d96e6792eb5ef01ad808f98de6ceedfac..101fbeabfe754c9d6d000d97e3a84a1e5381aef2 100644 (file)
@@ -679,15 +679,15 @@ private:
    */
   virtual void
   fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                  const Quadrature<dim>                                     &quadrature,
-                  InternalDataBase                                          &internal,
-                  std::vector<Point<spacedim> >                             &quadrature_points,
-                  std::vector<double>                                       &JxW_values,
-                  std::vector<DerivativeForm<1,dim,spacedim>  >       &jacobians,
-                  std::vector<DerivativeForm<2,dim,spacedim>  >       &jacobian_grads,
-                  std::vector<DerivativeForm<1,spacedim,dim>  >       &inverse_jacobians,
-                  std::vector<Point<spacedim> >                             &cell_normal_vectors,
-                  CellSimilarity::Similarity                           &cell_similarity
+                  const Quadrature<dim>                         &quadrature,
+                  InternalDataBase                              &internal,
+                  std::vector<Point<spacedim> >                 &quadrature_points,
+                  std::vector<double>                           &JxW_values,
+                  std::vector<DerivativeForm<1,dim,spacedim>  > &jacobians,
+                  std::vector<DerivativeForm<2,dim,spacedim>  > &jacobian_grads,
+                  std::vector<DerivativeForm<1,spacedim,dim>  > &inverse_jacobians,
+                  std::vector<Point<spacedim> >                 &cell_normal_vectors,
+                  CellSimilarity::Similarity                    &cell_similarity
                  ) const=0;
 
 
@@ -704,27 +704,31 @@ private:
    */
   virtual void
   fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                       const unsigned int                                        face_no,
-                       const Quadrature<dim-1>                                   &quadrature,
-                       InternalDataBase                                          &internal,
-                       std::vector<Point<spacedim> >                             &quadrature_points,
-                       std::vector<double>                                       &JxW_values,
-                       std::vector<Tensor<1,spacedim> >                          &boundary_form,
-                       std::vector<Point<spacedim> >                             &normal_vectors) const = 0;
+                       const unsigned int                            face_no,
+                       const Quadrature<dim-1>                      &quadrature,
+                       InternalDataBase                             &internal,
+                       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 = 0;
 
   /**
    * See above.
    */
   virtual void
   fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                          const unsigned int                        face_no,
-                          const unsigned int                        sub_no,
-                          const Quadrature<dim-1>                  &quadrature,
-                          InternalDataBase                         &internal,
-                          std::vector<Point<spacedim> >        &quadrature_points,
-                          std::vector<double>                      &JxW_values,
-                          std::vector<Tensor<1,spacedim> >     &boundary_form,
-                          std::vector<Point<spacedim> >        &normal_vectors) const = 0;
+                          const unsigned int                face_no,
+                          const unsigned int                sub_no,
+                          const Quadrature<dim-1>          &quadrature,
+                          InternalDataBase                 &internal,
+                          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 = 0;
 
   /**
    * Give class @p FEValues access to the private <tt>get_...data</tt> and
index 03524ab5a987d5f15318e4f0feb2be32ff8a5209..dc7b3baa83b0b9b0fde5705eeb5c730a0d918ef8 100644 (file)
@@ -71,36 +71,40 @@ public:
 
   virtual void
   fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                  const Quadrature<dim>                                     &quadrature,
-                  typename Mapping<dim, spacedim>::InternalDataBase         &mapping_data,
-                  std::vector<Point<spacedim> >                             &quadrature_points,
-                  std::vector<double>                                       &JxW_values,
+                  const Quadrature<dim>                             &quadrature,
+                  typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
+                  std::vector<Point<spacedim> >                     &quadrature_points,
+                  std::vector<double>                               &JxW_values,
                   std::vector<DerivativeForm<1,dim,spacedim> >      &jacobians,
-                  std::vector<DerivativeForm<2,dim,spacedim> >     &jacobian_grads,
-                  std::vector<DerivativeForm<1,spacedim,dim> >   &inverse_jacobians,
+                  std::vector<DerivativeForm<2,dim,spacedim> >      &jacobian_grads,
+                  std::vector<DerivativeForm<1,spacedim,dim> >      &inverse_jacobians,
                   std::vector<Point<spacedim> > &,
-                  CellSimilarity::Similarity                           &cell_similarity) const ;
+                  CellSimilarity::Similarity                        &cell_similarity) const ;
 
 
   virtual void
   fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                       const unsigned int face_no,
-                       const Quadrature<dim-1>& quadrature,
+                       const unsigned int               face_no,
+                       const Quadrature<dim-1>&         quadrature,
                        typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
                        std::vector<Point<dim> >        &quadrature_points,
                        std::vector<double>             &JxW_values,
-                       std::vector<Tensor<1,dim> >        &boundary_form,
-                       std::vector<Point<spacedim> >        &normal_vectors) const ;
+                       std::vector<Tensor<1,dim> >     &boundary_form,
+                       std::vector<Point<spacedim> >   &normal_vectors,
+                       std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                       std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
   virtual void
   fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                          const unsigned int face_no,
-                          const unsigned int sub_no,
-                          const Quadrature<dim-1>& quadrature,
+                          const unsigned int              face_no,
+                          const unsigned int              sub_no,
+                          const Quadrature<dim-1>&        quadrature,
                           typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
                           std::vector<Point<dim> >        &quadrature_points,
                           std::vector<double>             &JxW_values,
-                          std::vector<Tensor<1,dim> >        &boundary_form,
-                          std::vector<Point<spacedim> >        &normal_vectors) const ;
+                          std::vector<Tensor<1,dim> >     &boundary_form,
+                          std::vector<Point<spacedim> >   &normal_vectors,
+                          std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                          std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
 
   virtual void
   transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
index a16cd5d43e2e404deaa1aaf661306988a1e65dea..87a8e628eeab15c897613a2e79acdca9dccc9719 100644 (file)
@@ -195,42 +195,46 @@ protected:
    */
   virtual void
   fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                  const Quadrature<dim>                                     &quadrature,
-                  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
-                  typename std::vector<Point<spacedim> >                    &quadrature_points,
-                  std::vector<double>                                       &JxW_values,
-                  std::vector<DerivativeForm<1,dim,spacedim> >       &jacobians,
-                  std::vector<DerivativeForm<2,dim,spacedim> >       &jacobian_grads,
-                  std::vector<DerivativeForm<1,spacedim,dim> >      &inverse_jacobians,
-                  std::vector<Point<spacedim> >                             &cell_normal_vectors,
-                  CellSimilarity::Similarity                           &cell_similarity) const ;
+                  const Quadrature<dim>                            &quadrature,
+                  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+                  typename std::vector<Point<spacedim> >           &quadrature_points,
+                  std::vector<double>                              &JxW_values,
+                  std::vector<DerivativeForm<1,dim,spacedim> >     &jacobians,
+                  std::vector<DerivativeForm<2,dim,spacedim> >     &jacobian_grads,
+                  std::vector<DerivativeForm<1,spacedim,dim> >     &inverse_jacobians,
+                  std::vector<Point<spacedim> >                    &cell_normal_vectors,
+                  CellSimilarity::Similarity                       &cell_similarity) const ;
 
   /**
    * Implementation of the interface in Mapping.
    */
   virtual void
   fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                       const unsigned int face_no,
-                       const Quadrature<dim-1>& quadrature,
+                       const unsigned int                           face_no,
+                       const Quadrature<dim-1>&                     quadrature,
                        typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                       typename std::vector<Point<spacedim> >        &quadrature_points,
-                       std::vector<double>             &JxW_values,
-                       typename std::vector<Tensor<1,spacedim> >        &exterior_form,
-                       typename std::vector<Point<spacedim> >        &normal_vectors) const ;
+                       typename std::vector<Point<spacedim> >       &quadrature_points,
+                       std::vector<double>                          &JxW_values,
+                       typename std::vector<Tensor<1,spacedim> >    &exterior_form,
+                       typename std::vector<Point<spacedim> >       &normal_vectors,
+                       std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                       std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
 
   /**
    * Implementation of the interface in Mapping.
    */
   virtual void
   fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                          const unsigned int face_no,
-                          const unsigned int sub_no,
-                          const Quadrature<dim-1>& quadrature,
+                          const unsigned int                           face_no,
+                          const unsigned int                           sub_no,
+                          const Quadrature<dim-1>&                     quadrature,
                           typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                          typename std::vector<Point<spacedim> >        &quadrature_points,
-                          std::vector<double>             &JxW_values,
-                          typename std::vector<Tensor<1,spacedim> >        &exterior_form,
-                          typename std::vector<Point<spacedim> >        &normal_vectors) const ;
+                          typename std::vector<Point<spacedim> >       &quadrature_points,
+                          std::vector<double>                          &JxW_values,
+                          typename std::vector<Tensor<1,spacedim> >    &exterior_form,
+                          typename std::vector<Point<spacedim> >       &normal_vectors,
+                          std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                          std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
 
   /**
    * For <tt>dim=2,3</tt>. Append the support points of all shape functions
index 82d11cf8dc7cf32e75cbbe1965cf06a8b0e5dae9..856864c4071e3f79a5fc371dac42c60643cb8abf 100644 (file)
@@ -306,42 +306,46 @@ public:
    */
   virtual void
   fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                  const Quadrature<dim>                                     &quadrature,
-                  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
-                  typename std::vector<Point<spacedim> >                    &quadrature_points,
-                  std::vector<double>                                       &JxW_values,
-                  std::vector<DerivativeForm<1,dim,spacedim> >        &jacobians,
-                  std::vector<DerivativeForm<2,dim,spacedim> >       &jacobian_grads,
-                  std::vector<DerivativeForm<1,spacedim,dim> >      &inverse_jacobians,
-                  std::vector<Point<spacedim> >                             &cell_normal_vectors,
-                  CellSimilarity::Similarity                           &cell_similarity) const;
+                  const Quadrature<dim>                        &quadrature,
+                  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+                  typename std::vector<Point<spacedim> >       &quadrature_points,
+                  std::vector<double>                          &JxW_values,
+                  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
+                  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
+                  std::vector<Point<spacedim> >                &cell_normal_vectors,
+                  CellSimilarity::Similarity                   &cell_similarity) const;
 
   /**
    * Implementation of the interface in Mapping.
    */
   virtual void
   fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                       const unsigned int                               face_no,
-                       const Quadrature<dim-1>                          &quadrature,
+                       const unsigned int                            face_no,
+                       const Quadrature<dim-1>                      &quadrature,
                        typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                       typename std::vector<Point<spacedim> >                &quadrature_points,
-                       std::vector<double>                              &JxW_values,
-                       typename std::vector<Tensor<1,spacedim> >             &boundary_form,
-                       typename std::vector<Point<spacedim> >           &normal_vectors) const ;
+                       typename std::vector<Point<spacedim> >       &quadrature_points,
+                       std::vector<double>                          &JxW_values,
+                       typename std::vector<Tensor<1,spacedim> >    &boundary_form,
+                       typename std::vector<Point<spacedim> >       &normal_vectors,
+                       std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                       std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
 
   /**
    * Implementation of the interface in Mapping.
    */
   virtual void
   fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                          const unsigned int face_no,
-                          const unsigned int sub_no,
-                          const Quadrature<dim-1>& quadrature,
+                          const unsigned int                           face_no,
+                          const unsigned int                           sub_no,
+                          const Quadrature<dim-1>&                     quadrature,
                           typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                          typename std::vector<Point<spacedim> >        &quadrature_points,
-                          std::vector<double>             &JxW_values,
-                          typename std::vector<Tensor<1,spacedim> >        &boundary_form,
-                          typename std::vector<Point<spacedim> >        &normal_vectors) const ;
+                          typename std::vector<Point<spacedim> >       &quadrature_points,
+                          std::vector<double>                          &JxW_values,
+                          typename std::vector<Tensor<1,spacedim> >    &boundary_form,
+                          typename std::vector<Point<spacedim> >       &normal_vectors,
+                          std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                          std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const ;
 
   /**
    * Compute shape values and/or derivatives.
@@ -388,16 +392,18 @@ public:
    * Do the computation for the <tt>fill_*</tt> functions.
    */
   void compute_fill_face (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-                          const unsigned int      face_no,
-                          const unsigned int      subface_no,
-                          const unsigned int      npts,
-                          const DataSetDescriptor data_set,
-                          const std::vector<double>   &weights,
-                          InternalData           &mapping_data,
+                          const unsigned int                face_no,
+                          const unsigned int                subface_no,
+                          const unsigned int                npts,
+                          const DataSetDescriptor           data_set,
+                          const std::vector<double>        &weights,
+                          InternalData                     &mapping_data,
                           std::vector<Point<spacedim> >    &quadrature_points,
-                          std::vector<double>         &JxW_values,
+                          std::vector<double>              &JxW_values,
                           std::vector<Tensor<1,spacedim> > &boundary_form,
-                          std::vector<Point<spacedim> > &normal_vectors) const;
+                          std::vector<Point<spacedim> >    &normal_vectors,
+                          std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                          std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const;
 
   /**
    * Compute shape values and/or derivatives.
index 9aec4f71900ff4d9843e6630d26c74a0bdd63fd9..d133464e5efa2b85617851c33cd9394139140ac9 100644 (file)
@@ -3702,19 +3702,22 @@ void FEFaceValues<dim,spacedim>::reinit (const typename Triangulation<dim,spaced
 template <int dim, int spacedim>
 void FEFaceValues<dim,spacedim>::do_reinit (const unsigned int face_no)
 {
-  // first of all, set the
-  // present_face_index (if
-  // available)
+  // first of all, set the present_face_index (if available)
   const typename Triangulation<dim,spacedim>::cell_iterator cell=*this->present_cell;
   this->present_face_index=cell->face_index(face_no);
 
+  Assert(!(this->update_flags & update_jacobian_grads),
+         ExcNotImplemented());
+
   this->get_mapping().fill_fe_face_values(*this->present_cell, face_no,
                                           this->quadrature,
                                           *this->mapping_data,
                                           this->quadrature_points,
                                           this->JxW_values,
                                           this->boundary_forms,
-                                          this->normal_vectors);
+                                          this->normal_vectors,
+                                          this->jacobians,
+                                          this->inverse_jacobians);
 
   this->get_fe().fill_fe_face_values(this->get_mapping(),
                                      *this->present_cell, face_no,
@@ -3947,8 +3950,10 @@ void FESubfaceValues<dim,spacedim>::do_reinit (const unsigned int face_no,
       this->present_face_index=subface_index;
     }
 
-  // now ask the mapping and the finite element
-  // to do the actual work
+  Assert(!(this->update_flags & update_jacobian_grads),
+         ExcNotImplemented());
+
+  // now ask the mapping and the finite element to do the actual work
   this->get_mapping().fill_fe_subface_values(*this->present_cell,
                                              face_no, subface_no,
                                              this->quadrature,
@@ -3956,7 +3961,9 @@ void FESubfaceValues<dim,spacedim>::do_reinit (const unsigned int face_no,
                                              this->quadrature_points,
                                              this->JxW_values,
                                              this->boundary_forms,
-                                             this->normal_vectors);
+                                             this->normal_vectors,
+                                             this->jacobians,
+                                             this->inverse_jacobians);
 
   this->get_fe().fill_fe_subface_values(this->get_mapping(),
                                         *this->present_cell,
index c8fb5b03ac4eb02a322304ed6f201998b3d3cee5..02c7a0814dd24d89444053dc1204422cee749671 100644 (file)
@@ -390,88 +390,19 @@ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
 
 
 
-// template <>
-// void
-// MappingCartesian<1,1>::fill_fe_face_values (
-//   const typename Triangulation<1,1>::cell_iterator &,
-//   const unsigned int,
-//   const Quadrature<0>&,
-//   typename Mapping<1,1>::InternalDataBase&,
-//   std::vector<Point<1> >&,
-//   std::vector<double>&,
-//   std::vector<Tensor<1,1> >&,
-//   std::vector<Point<1> >&) const
-// {
-//   Assert(false, ExcNotImplemented());
-// }
-
-
-// template <>
-// void
-// MappingCartesian<1,2>::fill_fe_face_values (
-//   const typename Triangulation<1,2>::cell_iterator &,
-//   const unsigned int,
-//   const Quadrature<0>&,
-//   typename Mapping<1,2>::InternalDataBase&,
-//   std::vector<Point<1> >&,
-//   std::vector<double>&,
-//   std::vector<Tensor<1,1> >&,
-//   std::vector<Point<2> >&) const
-// {
-//   Assert(false, ExcNotImplemented());
-// }
-
-
-
-// template <>
-// void
-// MappingCartesian<1,1>::fill_fe_subface_values (
-//   const typename Triangulation<1,1>::cell_iterator &,
-//   const unsigned int,
-//   const unsigned int,
-//   const Quadrature<0>&,
-//   typename Mapping<1,1>::InternalDataBase&,
-//   std::vector<Point<1> >&,
-//   std::vector<double>&,
-//   std::vector<Tensor<1,1> >&,
-//   std::vector<Point<1> >&) const
-// {
-//   Assert(false, ExcNotImplemented());
-// }
-
-
-
-// template <>
-// void
-// MappingCartesian<1,2>::fill_fe_subface_values (
-//   const typename Triangulation<1,2>::cell_iterator &,
-//   const unsigned int,
-//   const unsigned int,
-//   const Quadrature<0>&,
-//   typename Mapping<1,2>::InternalDataBase&,
-//   std::vector<Point<1> >&,
-//   std::vector<double>&,
-//   std::vector<Tensor<1,1> >&,
-//   std::vector<Point<2> >&) const
-// {
-//   Assert(false, ExcNotImplemented());
-// }
-
-
-
-// Implementation for dim != 1
-
 template<int dim, int spacedim>
 void
 MappingCartesian<dim, spacedim>::fill_fe_face_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const unsigned int            face_no,
-  const Quadrature<dim-1>      &q,
+  const unsigned int             face_no,
+  const Quadrature<dim-1>       &q,
   typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  std::vector<Point<dim> >     &quadrature_points,
-  std::vector<double>          &JxW_values,
-  std::vector<Tensor<1,dim> >  &boundary_forms,
-  std::vector<Point<spacedim> >     &normal_vectors) const
+  std::vector<Point<dim> >      &quadrature_points,
+  std::vector<double>           &JxW_values,
+  std::vector<Tensor<1,dim> >   &boundary_forms,
+  std::vector<Point<spacedim> > &normal_vectors,
+  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   // convert data object to internal
   // data for this class. fails with
@@ -511,6 +442,22 @@ MappingCartesian<dim, spacedim>::fill_fe_face_values (
         J *= data.length[d];
       data.volume_element = J;
     }
+
+  if (data.current_update_flags() & update_jacobians)
+    for (unsigned int i=0; i<jacobians.size(); ++i)
+      {
+        jacobians[i] = DerivativeForm<1,dim,spacedim>();
+        for (unsigned int d=0; d<dim; ++d)
+          jacobians[i][d][d] = data.length[d];
+      }
+
+  if (data.current_update_flags() & update_inverse_jacobians)
+    for (unsigned int i=0; i<inverse_jacobians.size(); ++i)
+      {
+        inverse_jacobians[i] = DerivativeForm<1,dim,spacedim>();
+        for (unsigned int d=0; d<dim; ++d)
+          inverse_jacobians[i][d][d] = 1./data.length[d];
+      }
 }
 
 
@@ -519,14 +466,16 @@ template<int dim, int spacedim>
 void
 MappingCartesian<dim, spacedim>::fill_fe_subface_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const unsigned int       face_no,
-  const unsigned int       sub_no,
-  const Quadrature<dim-1> &q,
+  const unsigned int             face_no,
+  const unsigned int             sub_no,
+  const Quadrature<dim-1>       &q,
   typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
-  std::vector<Point<dim> >     &quadrature_points,
-  std::vector<double>          &JxW_values,
-  std::vector<Tensor<1,dim> >  &boundary_forms,
-  std::vector<Point<spacedim> >     &normal_vectors) const
+  std::vector<Point<dim> >      &quadrature_points,
+  std::vector<double>           &JxW_values,
+  std::vector<Tensor<1,dim> >   &boundary_forms,
+  std::vector<Point<spacedim> > &normal_vectors,
+  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   // convert data object to internal
   // data for this class. fails with
@@ -580,6 +529,22 @@ MappingCartesian<dim, spacedim>::fill_fe_subface_values (
         J *= data.length[d];
       data.volume_element = J;
     }
+
+  if (data.current_update_flags() & update_jacobians)
+    for (unsigned int i=0; i<jacobians.size(); ++i)
+      {
+        jacobians[i] = DerivativeForm<1,dim,spacedim>();
+        for (unsigned int d=0; d<dim; ++d)
+          jacobians[i][d][d] = data.length[d];
+      }
+
+  if (data.current_update_flags() & update_inverse_jacobians)
+    for (unsigned int i=0; i<inverse_jacobians.size(); ++i)
+      {
+        inverse_jacobians[i] = DerivativeForm<1,spacedim,dim>();
+        for (unsigned int d=0; d<dim; ++d)
+          inverse_jacobians[i][d][d] = 1./data.length[d];
+      }
 }
 
 
index 1bd07590d29944ba4a071e811e0367e7a815ab6b..43bc3f7814f39ed689147375dc094394bea97478 100644 (file)
@@ -348,13 +348,15 @@ template<int dim, int spacedim>
 void
 MappingQ<dim,spacedim>::fill_fe_face_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const unsigned int       face_no,
-  const Quadrature<dim-1> &q,
+  const unsigned int                face_no,
+  const Quadrature<dim-1>          &q,
   typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-  std::vector<Point<spacedim> >     &quadrature_points,
-  std::vector<double>          &JxW_values,
-  std::vector<Tensor<1,spacedim> >  &exterior_forms,
-  std::vector<Point<spacedim> >     &normal_vectors) const
+  std::vector<Point<spacedim> >    &quadrature_points,
+  std::vector<double>              &JxW_values,
+  std::vector<Tensor<1,spacedim> > &exterior_forms,
+  std::vector<Point<spacedim> >    &normal_vectors,
+  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   // convert data object to internal data for this class. fails with an
   // exception if that is not possible
@@ -391,7 +393,8 @@ MappingQ<dim,spacedim>::fill_fe_face_values (
                            q.get_weights(),
                            *p_data,
                            quadrature_points, JxW_values,
-                           exterior_forms, normal_vectors);
+                           exterior_forms, normal_vectors, jacobians,
+                           inverse_jacobians);
 }
 
 
@@ -405,7 +408,9 @@ MappingQ<dim,spacedim>::fill_fe_subface_values (const typename Triangulation<dim
                                                 std::vector<Point<spacedim> >     &quadrature_points,
                                                 std::vector<double>          &JxW_values,
                                                 std::vector<Tensor<1,spacedim> >  &exterior_forms,
-                                                std::vector<Point<spacedim> >     &normal_vectors) const
+                                                std::vector<Point<spacedim> >     &normal_vectors,
+                                                std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                                                std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   // convert data object to internal data for this class. fails with an
   // exception if that is not possible
@@ -443,7 +448,8 @@ MappingQ<dim,spacedim>::fill_fe_subface_values (const typename Triangulation<dim
                            q.get_weights(),
                            *p_data,
                            quadrature_points, JxW_values,
-                           exterior_forms, normal_vectors);
+                           exterior_forms, normal_vectors, jacobians,
+                           inverse_jacobians);
 }
 
 
index d7f90a67a17a12ade8b83d83e9d1e17a7b3e3150..399bdaa2fb21b8768607f5a03c47a910271c5b8a 100644 (file)
@@ -732,15 +732,15 @@ template<int dim, int spacedim>
 void
 MappingQ1<dim,spacedim>::fill_fe_values (
   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
-  const Quadrature<dim>                                     &q,
-  typename Mapping<dim,spacedim>::InternalDataBase          &mapping_data,
-  std::vector<Point<spacedim> >                             &quadrature_points,
-  std::vector<double>                                       &JxW_values,
-  std::vector<DerivativeForm<1,dim,spacedim>   >      &jacobians,
-  std::vector<DerivativeForm<2,dim,spacedim>    >     &jacobian_grads,
-  std::vector<DerivativeForm<1,spacedim,dim>   >      &inverse_jacobians,
-  std::vector<Point<spacedim> >                             &normal_vectors,
-  CellSimilarity::Similarity                                &cell_similarity) const
+  const Quadrature<dim>                        &q,
+  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+  std::vector<Point<spacedim> >                &quadrature_points,
+  std::vector<double>                          &JxW_values,
+  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
+  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
+  std::vector<Point<spacedim> >                &normal_vectors,
+  CellSimilarity::Similarity                   &cell_similarity) const
 {
   // ensure that the following static_cast is really correct:
   Assert (dynamic_cast<InternalData *>(&mapping_data) != 0,
@@ -924,7 +924,9 @@ namespace internal
                        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<Point<spacedim> >    &normal_vectors,
+                       std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                       std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians)
     {
       const UpdateFlags update_flags(data.current_update_flags());
 
@@ -1039,6 +1041,14 @@ namespace internal
                 if (update_flags & update_normal_vectors)
                   normal_vectors[i] = boundary_forms[i] / boundary_forms[i].norm();
               }
+
+          if (update_flags & update_jacobians)
+            for (unsigned int point=0; point<n_q_points; ++point)
+              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();
         }
     }
   }
@@ -1057,17 +1067,20 @@ MappingQ1<dim,spacedim>::compute_fill_face (
   const DataSetDescriptor          data_set,
   const std::vector<double>        &weights,
   InternalData                     &data,
-  std::vector<Point<spacedim> >         &quadrature_points,
+  std::vector<Point<spacedim> >    &quadrature_points,
   std::vector<double>              &JxW_values,
   std::vector<Tensor<1,spacedim> > &boundary_forms,
-  std::vector<Point<spacedim> >    &normal_vectors) const
+  std::vector<Point<spacedim> >    &normal_vectors,
+  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   compute_fill (cell, n_q_points, data_set, CellSimilarity::none,
                 data, quadrature_points);
   internal::compute_fill_face (*this,
                                cell, face_no, subface_no, n_q_points,
                                weights, data,
-                               JxW_values, boundary_forms, normal_vectors);
+                               JxW_values, boundary_forms, normal_vectors,
+                               jacobians, inverse_jacobians);
 }
 
 
@@ -1078,10 +1091,12 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
                      const unsigned int                               face_no,
                      const Quadrature<dim-1>                          &q,
                      typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
-                     std::vector<Point<spacedim> >                         &quadrature_points,
+                     std::vector<Point<spacedim> >                    &quadrature_points,
                      std::vector<double>                              &JxW_values,
-                     std::vector<Tensor<1,spacedim> >                      &boundary_forms,
-                     std::vector<Point<spacedim> >                    &normal_vectors) const
+                     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
 {
   // ensure that the following cast
   // is really correct:
@@ -1103,7 +1118,9 @@ fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &
                      quadrature_points,
                      JxW_values,
                      boundary_forms,
-                     normal_vectors);
+                     normal_vectors,
+                     jacobians,
+                     inverse_jacobians);
 }
 
 
@@ -1117,9 +1134,11 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
                         const Quadrature<dim-1> &q,
                         typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
                         std::vector<Point<spacedim> >     &quadrature_points,
-                        std::vector<double>          &JxW_values,
+                        std::vector<double>               &JxW_values,
                         std::vector<Tensor<1,spacedim> >  &boundary_forms,
-                        std::vector<Point<spacedim> >     &normal_vectors) const
+                        std::vector<Point<spacedim> >     &normal_vectors,
+                        std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
+                        std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians) const
 {
   // ensure that the following cast
   // is really correct:
@@ -1142,7 +1161,9 @@ fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterato
                      quadrature_points,
                      JxW_values,
                      boundary_forms,
-                     normal_vectors);
+                     normal_vectors,
+                     jacobians,
+                     inverse_jacobians);
 }
 
 
diff --git a/tests/fe/jacobians_face.cc b/tests/fe/jacobians_face.cc
new file mode 100644 (file)
index 0000000..62ea558
--- /dev/null
@@ -0,0 +1,123 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2014 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Show the Jacobians and inverse Jacobians of FEFaceValues and
+// FESubfaceValues on a hyperball mesh with one quadrature point
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_q.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+
+template<int dim>
+void test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_ball (tria);
+  static const HyperBallBoundary<dim> boundary;
+  tria.set_boundary (0, boundary);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  MappingQ<dim> mapping(3);
+  FE_Nothing<dim> dummy;
+  // choose a point that is not right in the middle of the cell so that the
+  // Jacobian contains many nonzero entries
+  Point<dim-1> quad_p;
+  for (int d=0; d<dim-1; ++d)
+    quad_p(d) = 0.42 + 0.11 * d;
+  Quadrature<dim-1> quad(quad_p);
+  FEFaceValues<dim> fe_val (mapping, dummy, quad,
+                            update_jacobians | update_inverse_jacobians);
+  FESubfaceValues<dim> fe_sub_val (mapping, dummy, quad,
+                                   update_jacobians | update_inverse_jacobians);
+  deallog << dim << "d Jacobians:" << std::endl;
+  typename Triangulation<dim>::active_cell_iterator
+  cell = tria.begin_active(), endc = tria.end();
+  for ( ; cell != endc; ++cell)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+        fe_val.reinit (cell, f);
+
+        for (unsigned int d=0; d<dim; ++d)
+          for (unsigned int e=0; e<dim; ++e)
+            deallog << fe_val.jacobian(0)[d][e] << " ";
+        deallog << std::endl;
+
+        // Also check the Jacobian with FESubfaceValues
+        if (cell->at_boundary(f) == false &&
+            cell->neighbor(f)->level() < cell->level())
+          {
+            fe_sub_val.reinit(cell->neighbor(f), cell->neighbor_face_no(f),
+                              cell->neighbor_of_coarser_neighbor(f).second);
+
+            for (unsigned int d=0; d<dim; ++d)
+              for (unsigned int e=0; e<dim; ++e)
+                deallog << fe_sub_val.jacobian(0)[d][e] << " ";
+            deallog << std::endl;
+          }
+      }
+  deallog << std::endl;
+
+  deallog << dim << "d inverse Jacobians:" << std::endl;
+  cell = tria.begin_active();
+  endc = tria.end();
+  for ( ; cell != endc; ++cell)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+        fe_val.reinit (cell, f);
+
+        for (unsigned int d=0; d<dim; ++d)
+          for (unsigned int e=0; e<dim; ++e)
+            deallog << fe_val.inverse_jacobian(0)[d][e] << " ";
+        deallog << std::endl;
+
+        // Also check the inverse Jacobian with FESubfaceValues
+        if (cell->at_boundary(f) == false &&
+            cell->neighbor(f)->level() < cell->level())
+          {
+            fe_sub_val.reinit(cell->neighbor(f), cell->neighbor_face_no(f),
+                              cell->neighbor_of_coarser_neighbor(f).second);
+
+            for (unsigned int d=0; d<dim; ++d)
+              for (unsigned int e=0; e<dim; ++e)
+                deallog << fe_sub_val.inverse_jacobian(0)[d][e] << " ";
+            deallog << std::endl;
+          }
+
+    }
+  deallog << std::endl;
+}
+
+
+int
+main()
+{
+  std::ofstream logfile ("output");
+  deallog << std::setprecision(4) << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test<2>();
+  test<3>();
+
+  return 0;
+}
diff --git a/tests/fe/jacobians_face.output b/tests/fe/jacobians_face.output
new file mode 100644 (file)
index 0000000..7c4bebc
--- /dev/null
@@ -0,0 +1,301 @@
+
+DEAL::2d Jacobians:
+DEAL::1.2396 -0.1855 0.1241 1.5579 
+DEAL::0.5337 0 0.0751 0.5858 
+DEAL::0.4142 -0.2889 0.4142 1.0601 
+DEAL::0.4142 0.2889 -0.4142 1.0601 
+DEAL::0.5858 0 0 0.5858 
+DEAL::0.5858 0 0 0.5858 
+DEAL::0.5858 0 0 0.5858 
+DEAL::0.5858 0 0 0.5858 
+DEAL::0.2889 -0.4142 1.0601 0.4142 
+DEAL::-0.2889 -0.4142 1.0601 -0.4142 
+DEAL::0.1855 -1.2396 1.5579 0.1241 
+DEAL::0 -0.5337 0.5858 0.0751 
+DEAL::0.1241 1.5579 -1.2396 0.1855 
+DEAL::0.0751 0.5858 -0.5337 0 
+DEAL::0.4142 1.0601 -0.4142 0.2889 
+DEAL::-0.4142 1.0601 -0.4142 -0.2889 
+DEAL::0.5904 0.2071 -0.3010 0.2071 
+DEAL::0.4142 -0.6003 0.4142 1.1510 
+DEAL::0.6484 0 -0.1611 0.3536 
+DEAL::0.7056 0.2097 -0.3449 0.4656 
+DEAL::0.5000 0.1332 -0.1464 0.2971 
+DEAL::0.6484 0 0.1611 0.3536 
+DEAL::0.5904 -0.2071 0.3010 0.2071 
+DEAL::0.6003 -0.4142 1.1510 0.4142 
+DEAL::0.7428 -0.1629 0.2551 0.4947 
+DEAL::0.5000 -0.0979 0.1464 0.3214 
+DEAL::0.4130 0.2071 -0.0849 0.2071 
+DEAL::0.4142 -0.1053 0.4142 0.8181 
+DEAL::0.4130 0 -0.0849 0.3536 
+DEAL::0.5000 0.1201 -0.1464 0.2686 
+DEAL::0.2929 0.1201 0 0.2686 
+DEAL::0.5858 0 0 0.5858 
+DEAL::0.4130 0 0.0849 0.3536 
+DEAL::0.4130 -0.2071 0.0849 0.2071 
+DEAL::0.1053 -0.4142 0.8181 0.4142 
+DEAL::0.5000 -0.0870 0.1464 0.2920 
+DEAL::0.2929 -0.0870 0 0.2920 
+DEAL::0.5858 0 0 0.5858 
+DEAL::
+DEAL::2d inverse Jacobians:
+DEAL::0.7972 0.0949 -0.0635 0.6343 
+DEAL::1.8738 0 -0.2403 1.7071 
+DEAL::1.8971 0.5171 -0.7413 0.7413 
+DEAL::1.8971 -0.5171 0.7413 0.7413 
+DEAL::1.7071 0 0 1.7071 
+DEAL::1.7071 0 0 1.7071 
+DEAL::1.7071 0 0 1.7071 
+DEAL::1.7071 0 0 1.7071 
+DEAL::0.7413 0.7413 -1.8971 0.5171 
+DEAL::-0.7413 0.7413 -1.8971 -0.5171 
+DEAL::0.0635 0.6343 -0.7972 0.0949 
+DEAL::0.2403 1.7071 -1.8738 0 
+DEAL::0.0949 -0.7972 0.6343 0.0635 
+DEAL::0 -1.8738 1.7071 0.2403 
+DEAL::0.5171 -1.8971 0.7413 0.7413 
+DEAL::-0.5171 -1.8971 0.7413 -0.7413 
+DEAL::1.1218 -1.1218 1.6305 3.1980 
+DEAL::1.5867 0.8275 -0.5710 0.5710 
+DEAL::1.5424 0 0.7030 2.8284 
+DEAL::1.1615 -0.5233 0.8604 1.7603 
+DEAL::1.7679 -0.7923 0.8715 2.9754 
+DEAL::1.5424 0 -0.7030 2.8284 
+DEAL::1.1218 1.1218 -1.6305 3.1980 
+DEAL::0.5710 0.5710 -1.5867 0.8275 
+DEAL::1.2096 0.3982 -0.6237 1.8162 
+DEAL::1.8362 0.5594 -0.8367 2.8565 
+DEAL::2.0082 -2.0082 0.8236 4.0048 
+DEAL::2.1389 0.2753 -1.0830 1.0830 
+DEAL::2.4212 0 0.5817 2.8284 
+DEAL::1.7684 -0.7908 0.9641 3.2917 
+DEAL::3.4142 -1.5268 0 3.7228 
+DEAL::1.7071 0 0 1.7071 
+DEAL::2.4212 0 -0.5817 2.8284 
+DEAL::2.0082 2.0082 -0.8236 4.0048 
+DEAL::1.0830 1.0830 -2.1389 0.2753 
+DEAL::1.8395 0.5479 -0.9224 3.1494 
+DEAL::3.4142 1.0169 0 3.4241 
+DEAL::1.7071 0 0 1.7071 
+DEAL::
+DEAL::3d Jacobians:
+DEAL::0.8135 -0.0139 0.4404 -0.0123 0.7662 0.0586 -0.1278 -0.0139 0.4404 
+DEAL::0.8135 0.0139 -0.4404 0.0123 0.7662 0.0586 0.1278 -0.0139 0.4404 
+DEAL::0.8487 0.0077 -0.0231 0.0077 0.9111 0.4803 0.0077 -0.2042 0.4803 
+DEAL::0.8487 -0.0077 -0.0231 -0.0077 0.9111 -0.4803 0.0077 0.2042 0.4803 
+DEAL::1.5498 0.0063 0.1758 0.0093 1.6994 -0.0880 -0.1757 0.0651 1.6471 
+DEAL::0.4226 0 0.0764 0 0.4226 -0.0323 0 0 0.5391 
+DEAL::0.2042 -0.4803 -0.0077 0.9111 0.4803 0.0077 0.0077 -0.0231 0.8487 
+DEAL::-0.2042 -0.4803 -0.0077 0.9111 -0.4803 -0.0077 -0.0077 -0.0231 0.8487 
+DEAL::-0.0651 -1.6471 0.1757 1.6994 -0.0880 0.0093 0.0063 0.1758 1.5498 
+DEAL::0 -0.5391 0 0.4226 0.0864 0 0 -0.0286 0.4226 
+DEAL::0.0139 -0.4404 0.1278 0.7662 0.0586 -0.0123 -0.0139 0.4404 0.8135 
+DEAL::0.0139 -0.4404 -0.1278 0.7662 0.0586 0.0123 0.0139 -0.4404 0.8135 
+DEAL::0.9111 0.4803 0.0077 0.0077 -0.0231 0.8487 0.2042 -0.4803 -0.0077 
+DEAL::0.8204 -0.4404 0.0139 0.0134 0.0586 0.7662 -0.1278 -0.4404 0.0139 
+DEAL::1.6994 -0.0880 0.0093 0.0063 0.1758 1.5498 -0.0651 -1.6471 0.1757 
+DEAL::0.4226 0.0864 0 0 -0.0286 0.4226 0 -0.5391 0 
+DEAL::0.7662 0.0586 -0.0123 -0.0139 0.4404 0.8135 0.0139 -0.4404 0.1278 
+DEAL::0.7662 0.0586 0.0123 0.0139 -0.4404 0.8135 0.0139 -0.4404 -0.1278 
+DEAL::1.6471 -0.1757 0.0651 0.1758 1.5498 0.0063 -0.0880 0.0093 1.6994 
+DEAL::0.5391 0 0 0.0764 0.4226 0 -0.0323 0 0.4226 
+DEAL::0.4404 -0.1278 -0.0139 0.4404 0.8135 -0.0139 0.0586 -0.0123 0.7662 
+DEAL::0.4404 0.1278 -0.0139 -0.4404 0.8135 0.0139 0.0586 0.0123 0.7662 
+DEAL::0.4404 -0.0139 -0.1278 0.0586 0.7662 -0.0134 0.4404 -0.0139 0.8204 
+DEAL::0.4803 0.0077 0.2042 -0.0231 0.8487 -0.0077 -0.4803 -0.0077 0.9111 
+DEAL::0.8204 0.4404 -0.0139 -0.1278 0.4404 -0.0139 -0.0134 0.0586 0.7662 
+DEAL::0.9111 -0.4803 -0.0077 0.2042 0.4803 0.0077 -0.0077 -0.0231 0.8487 
+DEAL::1.6994 -0.0880 0.0093 0.0651 1.6471 -0.1757 0.0063 0.1758 1.5498 
+DEAL::0.4226 -0.0323 0 0 0.5391 0 0 0.0764 0.4226 
+DEAL::0.8487 -0.0231 0.0075 0.0077 0.4803 -0.2042 0.0077 0.4803 0.9092 
+DEAL::0.7662 0.0586 0.0123 -0.0139 0.4404 0.1278 0.0139 -0.4404 0.8135 
+DEAL::0.1758 1.5498 0.0063 -1.6471 0.1757 -0.0651 -0.0880 0.0093 1.6994 
+DEAL::-0.0286 0.4226 0 -0.5391 0 0 0.0864 0 0.4226 
+DEAL::0.4404 0.8135 -0.0139 -0.4404 0.1278 0.0139 0.0586 -0.0123 0.7662 
+DEAL::-0.4803 0.9092 -0.0077 -0.4803 -0.2042 -0.0077 -0.0231 -0.0075 0.8487 
+DEAL::-0.0231 0.8487 0.0077 -0.4803 -0.0077 0.2042 0.4803 0.0077 0.9111 
+DEAL::0.0586 0.7662 0.0134 -0.4404 0.0139 -0.1278 -0.4404 0.0139 0.8204 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4796 0 0 0.2569 0.4226 0 0.2223 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0.2223 0 0 0.4796 0 0 0.2569 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0 0.2569 0 0.4226 0.2223 0 0 0.4796 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0 -0.4796 0 0.4226 0.2761 0 0 0.2067 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 -0.2478 0 0 0.4737 0 0 0.2543 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0 -0.1895 0 0.4226 0.2320 0 0 0.4968 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4968 0 0 -0.1895 0.4226 0 0.2320 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2067 0.4226 0 -0.4796 0 0 0.2761 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0 0.2543 0 0.4226 -0.2478 0 0 0.4737 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0 -0.4968 0 0.4226 -0.2068 0 0 0.2127 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::-0.2316 0.4226 0 -0.4737 0 0 0.2720 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0 -0.1875 0 0.4226 -0.2582 0 0 0.4902 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4737 0 0 0.2543 0.4226 0 -0.2478 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0.2320 0 0 0.4968 0 0 -0.1895 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0.2761 0 0 0.2067 0.4226 0 -0.4796 0 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0 -0.4737 0 0.4226 0.2720 0 0 -0.2316 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 -0.2582 0 0 0.4902 0 0 -0.1875 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 -0.2068 0 0 0.2127 0.4226 0 -0.4968 0 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4902 0 0 -0.1875 0.4226 0 -0.2582 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2127 0.4226 0 -0.4968 0 0 -0.2068 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 0.2720 0 0 -0.2316 0.4226 0 -0.4737 0 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0 -0.4902 0 0.4226 -0.2035 0 0 -0.2381 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::-0.2381 0.4226 0 -0.4902 0 0 -0.2035 0 0.4226 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.2113 0 0 0 0.2113 0 0 0 0.2113 
+DEAL::0.4226 -0.2035 0 0 -0.2381 0.4226 0 -0.4902 0 
+DEAL::
+DEAL::3d inverse Jacobians:
+DEAL::1.0623 0 -1.0623 -0.0066 1.3019 -0.1668 0.3081 0.0411 1.9570 
+DEAL::1.0623 0 1.0623 0.0066 1.3019 -0.1668 -0.3081 0.0411 1.9570 
+DEAL::1.1777 0.0023 0.0544 0 0.8966 -0.8966 -0.0189 0.3813 1.7000 
+DEAL::1.1777 -0.0023 0.0544 0 0.8966 0.8966 -0.0189 -0.3813 1.7000 
+DEAL::0.6375 0.0002 -0.0680 0.0000 0.5872 0.0314 0.0680 -0.0232 0.5986 
+DEAL::2.3660 0 -0.3353 0 2.3660 0.1420 0 0 1.8549 
+DEAL::0.8966 0.8966 0 -1.7000 0.3813 -0.0189 -0.0544 0.0023 1.1777 
+DEAL::-0.8966 0.8966 0 -1.7000 -0.3813 -0.0189 -0.0544 -0.0023 1.1777 
+DEAL::-0.0314 0.5872 0.0000 -0.5986 -0.0232 0.0680 0.0680 0.0002 0.6375 
+DEAL::0.3790 2.3660 0 -1.8549 0 0 -0.1256 0 2.3660 
+DEAL::0.1668 1.3019 -0.0066 -1.9570 0.0411 0.3081 1.0623 0 1.0623 
+DEAL::0.1668 1.3019 0.0066 -1.9570 0.0411 -0.3081 -1.0623 0 1.0623 
+DEAL::0.8966 0 0.8966 0.3813 -0.0189 -1.7000 0.0023 1.1777 -0.0544 
+DEAL::1.0545 0 -1.0545 -0.3059 0.0411 -1.9592 0.0050 1.3019 0.1684 
+DEAL::0.5872 0.0000 -0.0314 -0.0232 0.0680 -0.5986 0.0002 0.6375 0.0680 
+DEAL::2.3660 0 0.3790 0 0 -1.8549 0 2.3660 -0.1256 
+DEAL::1.3019 -0.0066 0.1668 0.0411 0.3081 -1.9570 0 1.0623 1.0623 
+DEAL::1.3019 0.0066 0.1668 0.0411 -0.3081 -1.9570 0 1.0623 -1.0623 
+DEAL::0.5986 0.0680 -0.0232 -0.0680 0.6375 0.0002 0.0314 0.0000 0.5872 
+DEAL::1.8549 0 0 -0.3353 2.3660 0 0.1420 0 2.3660 
+DEAL::1.9570 0.3081 0.0411 -1.0623 1.0623 0 -0.1668 -0.0066 1.3019 
+DEAL::1.9570 -0.3081 0.0411 1.0623 1.0623 0 -0.1668 0.0066 1.3019 
+DEAL::1.9592 0.0411 0.3059 -0.1684 1.3019 -0.0050 -1.0545 0 1.0545 
+DEAL::1.7000 -0.0189 -0.3813 0.0544 1.1777 -0.0023 0.8966 0 0.8966 
+DEAL::1.0545 -1.0545 0 0.3059 1.9592 0.0411 -0.0050 -0.1684 1.3019 
+DEAL::0.8966 0.8966 0 -0.3813 1.7000 -0.0189 -0.0023 0.0544 1.1777 
+DEAL::0.5872 0.0314 0.0000 -0.0232 0.5986 0.0680 0.0002 -0.0680 0.6375 
+DEAL::2.3660 0.1420 0 0 1.8549 0 0 -0.3353 2.3660 
+DEAL::1.1777 0.0543 0.0024 -0.0189 1.6994 0.3819 0 -0.8981 0.8981 
+DEAL::1.3019 -0.1668 0.0066 0.0411 1.9570 -0.3081 0 1.0623 1.0623 
+DEAL::0.0680 -0.5986 -0.0232 0.6375 0.0680 0.0002 0.0000 -0.0314 0.5872 
+DEAL::0 -1.8549 0 2.3660 -0.1256 0 0 0.3790 2.3660 
+DEAL::0.3081 -1.9570 0.0411 1.0623 1.0623 0 -0.0066 0.1668 1.3019 
+DEAL::-0.3819 -1.6994 -0.0189 0.8981 -0.8981 0 -0.0024 -0.0543 1.1777 
+DEAL::-0.0189 -1.7000 0.3813 1.1777 -0.0544 0.0023 0 0.8966 0.8966 
+DEAL::0.0411 -1.9592 -0.3059 1.3019 0.1684 0.0050 0 -1.0545 1.0545 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.0849 0 0 -1.2674 2.3660 0 -1.0966 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 -1.0966 0 0 2.0849 0 0 -1.2674 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 -1.2674 0 2.3660 -1.0966 0 0 2.0849 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::1.3621 2.3660 0 -2.0849 0 0 1.0199 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 1.2378 0 0 2.1111 0 0 -1.2703 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 0.9027 0 2.3660 -1.1048 0 0 2.0129 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.0129 0 0 0.9027 2.3660 0 -1.1048 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::0 -2.0849 0 2.3660 1.0199 0 0 1.3621 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 -1.2703 0 2.3660 1.2378 0 0 2.1111 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::-0.9849 2.3660 0 -2.0129 0 0 1.0129 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::0 -2.1111 0 2.3660 -1.1569 0 0 1.3588 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 0.9052 0 2.3660 1.2464 0 0 2.0401 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.1111 0 0 -1.2703 2.3660 0 1.2378 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 -1.1048 0 0 2.0129 0 0 0.9027 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 1.3621 0 0 -2.0849 0 2.3660 1.0199 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::1.3588 2.3660 0 -2.1111 0 0 -1.1569 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 1.2464 0 0 2.0401 0 0 0.9052 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 -0.9849 0 0 -2.0129 0 2.3660 1.0129 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.0401 0 0 0.9052 2.3660 0 1.2464 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::0 -2.0129 0 2.3660 1.0129 0 0 -0.9849 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 1.3588 0 0 -2.1111 0 2.3660 -1.1569 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::-0.9820 2.3660 0 -2.0401 0 0 -1.1494 0 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::0 -2.0401 0 2.3660 -1.1494 0 0 -0.9820 2.3660 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::4.7321 0 0 0 4.7321 0 0 0 4.7321 
+DEAL::2.3660 0 -0.9820 0 0 -2.0401 0 2.3660 -1.1494 
+DEAL::
diff --git a/tests/fe/jacobians_face_cartesian.cc b/tests/fe/jacobians_face_cartesian.cc
new file mode 100644 (file)
index 0000000..e2ba285
--- /dev/null
@@ -0,0 +1,122 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2014 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Show the Jacobians and inverse Jacobians of FEFaceValues and
+// FESubfaceValues on a Cartesian mesh with MappingCartesian with one
+// quadrature point
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_cartesian.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+
+template<int dim>
+void test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube (tria);
+  tria.refine_global(1);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  MappingCartesian<dim> mapping;
+  FE_Nothing<dim> dummy;
+  // choose a point that is not right in the middle of the cell so that the
+  // Jacobian contains many nonzero entries
+  Point<dim-1> quad_p;
+  for (int d=0; d<dim-1; ++d)
+    quad_p(d) = 0.42 + 0.11 * d;
+  Quadrature<dim-1> quad(quad_p);
+  FEFaceValues<dim> fe_val (mapping, dummy, quad,
+                            update_jacobians | update_inverse_jacobians);
+  FESubfaceValues<dim> fe_sub_val (mapping, dummy, quad,
+                                   update_jacobians | update_inverse_jacobians);
+  deallog << dim << "d Jacobians:" << std::endl;
+  typename Triangulation<dim>::active_cell_iterator
+  cell = tria.begin_active(), endc = tria.end();
+  for ( ; cell != endc; ++cell)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+        fe_val.reinit (cell, f);
+
+        for (unsigned int d=0; d<dim; ++d)
+          for (unsigned int e=0; e<dim; ++e)
+            deallog << fe_val.jacobian(0)[d][e] << " ";
+        deallog << std::endl;
+
+        // Also check the Jacobian with FESubfaceValues
+        if (cell->at_boundary(f) == false &&
+            cell->neighbor(f)->level() < cell->level())
+          {
+            fe_sub_val.reinit(cell->neighbor(f), cell->neighbor_face_no(f),
+                              cell->neighbor_of_coarser_neighbor(f).second);
+
+            for (unsigned int d=0; d<dim; ++d)
+              for (unsigned int e=0; e<dim; ++e)
+                deallog << fe_sub_val.jacobian(0)[d][e] << " ";
+            deallog << std::endl;
+          }
+      }
+  deallog << std::endl;
+
+  deallog << dim << "d inverse Jacobians:" << std::endl;
+  cell = tria.begin_active();
+  endc = tria.end();
+  for ( ; cell != endc; ++cell)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+        fe_val.reinit (cell, f);
+
+        for (unsigned int d=0; d<dim; ++d)
+          for (unsigned int e=0; e<dim; ++e)
+            deallog << fe_val.inverse_jacobian(0)[d][e] << " ";
+        deallog << std::endl;
+
+        // Also check the inverse Jacobian with FESubfaceValues
+        if (cell->at_boundary(f) == false &&
+            cell->neighbor(f)->level() < cell->level())
+          {
+            fe_sub_val.reinit(cell->neighbor(f), cell->neighbor_face_no(f),
+                              cell->neighbor_of_coarser_neighbor(f).second);
+
+            for (unsigned int d=0; d<dim; ++d)
+              for (unsigned int e=0; e<dim; ++e)
+                deallog << fe_sub_val.inverse_jacobian(0)[d][e] << " ";
+            deallog << std::endl;
+          }
+
+    }
+  deallog << std::endl;
+}
+
+
+int
+main()
+{
+  std::ofstream logfile ("output");
+  deallog << std::setprecision(4) << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test<2>();
+  test<3>();
+
+  return 0;
+}
diff --git a/tests/fe/jacobians_face_cartesian.output b/tests/fe/jacobians_face_cartesian.output
new file mode 100644 (file)
index 0000000..7df785f
--- /dev/null
@@ -0,0 +1,277 @@
+
+DEAL::2d Jacobians:
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.5000 0 0 0.5000 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.2500 0 0 0.2500 
+DEAL::0.5000 0 0 0.5000 
+DEAL::
+DEAL::2d inverse Jacobians:
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::4.0000 0 0 4.0000 
+DEAL::2.0000 0 0 2.0000 
+DEAL::
+DEAL::3d Jacobians:
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.2500 0 0 0 0.2500 0 0 0 0.2500 
+DEAL::0.5000 0 0 0 0.5000 0 0 0 0.5000 
+DEAL::
+DEAL::3d inverse Jacobians:
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::4.0000 0 0 0 4.0000 0 0 0 4.0000 
+DEAL::2.0000 0 0 0 2.0000 0 0 0 2.0000 
+DEAL::

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.