]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Instantiate MeshWorker objects for dim != spacedim
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 16 Apr 2017 11:42:49 +0000 (13:42 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 16 Apr 2017 11:48:36 +0000 (13:48 +0200)
doc/news/changes/minor/20170416DanielArndt [new file with mode: 0644]
include/deal.II/meshworker/dof_info.h
include/deal.II/meshworker/integration_info.h
include/deal.II/meshworker/loop.h
include/deal.II/meshworker/vector_selector.h
include/deal.II/meshworker/vector_selector.templates.h
source/meshworker/mesh_worker_info.inst.in
source/meshworker/mesh_worker_vector_selector.inst.in

diff --git a/doc/news/changes/minor/20170416DanielArndt b/doc/news/changes/minor/20170416DanielArndt
new file mode 100644 (file)
index 0000000..f8b22c6
--- /dev/null
@@ -0,0 +1,3 @@
+Fixed: The MeshWorker framework can also be used in case dim != spacedim.
+<br>
+(Daniel Arndt, 2017/04/16)
index d17372380e947a00893614f025ca1cae13736d44..6daf750343eabf7d1995a4b8fba3bf4f2096b511 100644 (file)
@@ -327,7 +327,7 @@ namespace MeshWorker
     const DHFaceIterator &f,
     const unsigned int face_no)
   {
-    face = static_cast<typename Triangulation<dim>::face_iterator> (f);
+    face = static_cast<typename Triangulation<dim, spacedim>::face_iterator> (f);
     face_number = face_no;
     sub_number = numbers::invalid_unsigned_int;
   }
@@ -342,11 +342,11 @@ namespace MeshWorker
     const unsigned int face_no)
   {
     if ((cell.state() != IteratorState::valid)
-        ||  cell != typename Triangulation<dim>::cell_iterator(*c))
+        ||  cell != typename Triangulation<dim, spacedim>::cell_iterator(*c))
       get_indices(c);
     level_cell = c->is_level_cell();
 
-    cell = typename Triangulation<dim>::cell_iterator(*c);
+    cell = typename Triangulation<dim, spacedim>::cell_iterator(*c);
     set_face(f,face_no);
 
     if (block_info)
@@ -364,7 +364,7 @@ namespace MeshWorker
     const unsigned int face_no,
     const unsigned int subface_no)
   {
-    face = static_cast<typename Triangulation<dim>::face_iterator> (f);
+    face = static_cast<typename Triangulation<dim, spacedim>::face_iterator> (f);
     face_number = face_no;
     sub_number = subface_no;
   }
@@ -380,11 +380,11 @@ namespace MeshWorker
     const unsigned int subface_no)
   {
     if (cell.state() != IteratorState::valid
-        || cell != static_cast<typename Triangulation<dim>::cell_iterator> (c))
+        || cell != static_cast<typename Triangulation<dim, spacedim>::cell_iterator> (c))
       get_indices(c);
     level_cell = c->is_level_cell();
 
-    cell = static_cast<typename Triangulation<dim>::cell_iterator> (c);
+    cell = static_cast<typename Triangulation<dim, spacedim>::cell_iterator> (c);
     set_subface(f, face_no, subface_no);
 
     if (block_info)
index c494c1174cfb83e051a4462fa46d24e0f57ac546..e75848a4ed72621e6cabadb57788f90ce0bed003 100644 (file)
@@ -166,7 +166,7 @@ namespace MeshWorker
      * one vector for each component, containing vectors with values for each
      * quadrature point.
      */
-    std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
+    std::vector<std::vector<std::vector<Tensor<1,spacedim> > > > gradients;
 
     /**
      * The vector containing the second derivatives of finite element
@@ -176,7 +176,7 @@ namespace MeshWorker
      * one vector for each component, containing vectors with values for each
      * quadrature point.
      */
-    std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
+    std::vector<std::vector<std::vector<Tensor<2,spacedim> > > > hessians;
 
     /**
      * Reinitialize internal data structures for use on a cell.
@@ -666,13 +666,13 @@ namespace MeshWorker
         if (info.sub_number != numbers::invalid_unsigned_int)
           {
             // This is a subface
-            FESubfaceValues<dim> &fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
+            FESubfaceValues<dim,spacedim> &fe = dynamic_cast<FESubfaceValues<dim,spacedim>&> (febase);
             fe.reinit(info.cell, info.face_number, info.sub_number);
           }
         else if (info.face_number != numbers::invalid_unsigned_int)
           {
             // This is a face
-            FEFaceValues<dim> &fe = dynamic_cast<FEFaceValues<dim>&> (febase);
+            FEFaceValues<dim,spacedim> &fe = dynamic_cast<FEFaceValues<dim,spacedim>&> (febase);
             fe.reinit(info.cell, info.face_number);
           }
         else
index d062db5dd81a6bcbcacb9a898824b0e1e5762227..7e025aae406198736ff608220763cd66c231331d 100644 (file)
@@ -438,9 +438,9 @@ namespace MeshWorker
                         ASSEMBLER &assembler,
                         const LoopControl &lctrl = LoopControl())
   {
-    std::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> cell_worker;
-    std::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> boundary_worker;
-    std::function<void (DoFInfo<dim> &, DoFInfo<dim> &,
+    std::function<void (DoFInfo<dim, spacedim>&, IntegrationInfo<dim, spacedim>&)> cell_worker;
+    std::function<void (DoFInfo<dim, spacedim>&, IntegrationInfo<dim, spacedim>&)> boundary_worker;
+    std::function<void (DoFInfo<dim, spacedim>&, DoFInfo<dim, spacedim> &,
                         IntegrationInfo<dim, spacedim> &,
                         IntegrationInfo<dim, spacedim> &)> face_worker;
     if (integrator.use_cell)
index bee6022bf53271f5828ce1493007a5956a6b81b4..1917b9d5f1ae2af07355b2436e7731cd8ff32b79 100644 (file)
@@ -241,8 +241,8 @@ namespace MeshWorker
      */
     virtual void fill(
       std::vector<std::vector<std::vector<Number> > > &values,
-      std::vector<std::vector<std::vector<Tensor<1,dim,Number> > > > &gradients,
-      std::vector<std::vector<std::vector<Tensor<2,dim,Number> > > > &hessians,
+      std::vector<std::vector<std::vector<Tensor<1,spacedim,Number> > > > &gradients,
+      std::vector<std::vector<std::vector<Tensor<2,spacedim,Number> > > > &hessians,
       const FEValuesBase<dim,spacedim> &fe,
       const std::vector<types::global_dof_index> &index,
       const unsigned int component,
@@ -258,8 +258,8 @@ namespace MeshWorker
      */
     virtual void mg_fill(
       std::vector<std::vector<std::vector<Number> > > &values,
-      std::vector<std::vector<std::vector<Tensor<1,dim,Number> > > > &gradients,
-      std::vector<std::vector<std::vector<Tensor<2,dim,Number> > > > &hessians,
+      std::vector<std::vector<std::vector<Tensor<1,spacedim,Number> > > > &gradients,
+      std::vector<std::vector<std::vector<Tensor<2,spacedim,Number> > > > &hessians,
       const FEValuesBase<dim,spacedim> &fe,
       const unsigned int level,
       const std::vector<types::global_dof_index> &index,
@@ -315,8 +315,8 @@ namespace MeshWorker
 
     virtual void fill(
       std::vector<std::vector<std::vector<typename VectorType::value_type> > > &values,
-      std::vector<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > &gradients,
-      std::vector<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > &hessians,
+      std::vector<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > &gradients,
+      std::vector<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > &hessians,
       const FEValuesBase<dim,spacedim> &fe,
       const std::vector<types::global_dof_index> &index,
       const unsigned int component,
@@ -326,8 +326,8 @@ namespace MeshWorker
 
     virtual void mg_fill(
       std::vector<std::vector<std::vector<typename VectorType::value_type> > > &values,
-      std::vector<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > &gradients,
-      std::vector<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > &hessians,
+      std::vector<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > &gradients,
+      std::vector<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > &hessians,
       const FEValuesBase<dim,spacedim> &fe,
       const unsigned int level,
       const std::vector<types::global_dof_index> &index,
index 32a393b58b2fd508d0d226bd42afc9602a85ee96..24c73a890044aa676238c44a95dc6978dd7c31a9 100644 (file)
@@ -55,8 +55,8 @@ namespace MeshWorker
   void
   VectorDataBase<dim, spacedim, Number>::fill(
     std::vector<std::vector<std::vector<Number> > > &,
-    std::vector<std::vector<std::vector<Tensor<1,dim,Number> > > > &,
-    std::vector<std::vector<std::vector<Tensor<2,dim,Number> > > > &,
+    std::vector<std::vector<std::vector<Tensor<1,spacedim,Number> > > > &,
+    std::vector<std::vector<std::vector<Tensor<2,spacedim,Number> > > > &,
     const FEValuesBase<dim,spacedim> &,
     const std::vector<types::global_dof_index> &,
     const unsigned int,
@@ -72,8 +72,8 @@ namespace MeshWorker
   void
   VectorDataBase<dim, spacedim, Number>::mg_fill(
     std::vector<std::vector<std::vector<Number> > > &,
-    std::vector<std::vector<std::vector<Tensor<1,dim,Number> > > > &,
-    std::vector<std::vector<std::vector<Tensor<2,dim,Number> > > > &,
+    std::vector<std::vector<std::vector<Tensor<1,spacedim,Number> > > > &,
+    std::vector<std::vector<std::vector<Tensor<2,spacedim,Number> > > > &,
     const FEValuesBase<dim,spacedim> &,
     const unsigned int,
     const std::vector<types::global_dof_index> &,
@@ -123,8 +123,8 @@ namespace MeshWorker
   void
   VectorData<VectorType, dim, spacedim>::fill(
     std::vector<std::vector<std::vector<typename VectorType::value_type> > > &values,
-    std::vector<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > &gradients,
-    std::vector<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > &hessians,
+    std::vector<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > &gradients,
+    std::vector<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > &hessians,
     const FEValuesBase<dim,spacedim> &fe,
     const std::vector<types::global_dof_index> &index,
     const unsigned int component,
@@ -147,14 +147,14 @@ namespace MeshWorker
     for (unsigned int i=0; i<this->n_gradients(); ++i)
       {
         const VectorType *src = data.read_ptr<VectorType>(this->gradient_index(i));
-        VectorSlice<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > dst(gradients[i], component, n_comp);
+        VectorSlice<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > dst(gradients[i], component, n_comp);
         fe.get_function_gradients(*src, make_slice(index, start, size), dst, true);
       }
 
     for (unsigned int i=0; i<this->n_hessians(); ++i)
       {
         const VectorType *src = data.read_ptr<VectorType>(this->hessian_index(i));
-        VectorSlice<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > dst(hessians[i], component, n_comp);
+        VectorSlice<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > dst(hessians[i], component, n_comp);
         fe.get_function_hessians(*src, make_slice(index, start, size), dst, true);
       }
   }
@@ -207,8 +207,8 @@ namespace MeshWorker
   void
   VectorData<VectorType, dim, spacedim>::mg_fill
   (std::vector<std::vector<std::vector<typename VectorType::value_type> > >                &values,
-   std::vector<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > &gradients,
-   std::vector<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > &hessians,
+   std::vector<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > &gradients,
+   std::vector<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > &hessians,
    const FEValuesBase<dim,spacedim>           &fe,
    const unsigned int                         level,
    const std::vector<types::global_dof_index> &index,
@@ -233,14 +233,14 @@ namespace MeshWorker
     for (unsigned int i=0; i<this->n_gradients(); ++i)
       {
         const MGLevelObject<VectorType> *src = data.read_ptr<MGLevelObject<VectorType> >(this->value_index(i));
-        VectorSlice<std::vector<std::vector<Tensor<1,dim,typename VectorType::value_type> > > > dst(gradients[i], component, n_comp);
+        VectorSlice<std::vector<std::vector<Tensor<1,spacedim,typename VectorType::value_type> > > > dst(gradients[i], component, n_comp);
         fe.get_function_gradients((*src)[level], make_slice(index, start, size), dst, true);
       }
 
     for (unsigned int i=0; i<this->n_hessians(); ++i)
       {
         const MGLevelObject<VectorType> *src = data.read_ptr<MGLevelObject<VectorType> >(this->value_index(i));
-        VectorSlice<std::vector<std::vector<Tensor<2,dim,typename VectorType::value_type> > > > dst(hessians[i], component, n_comp);
+        VectorSlice<std::vector<std::vector<Tensor<2,spacedim,typename VectorType::value_type> > > > dst(hessians[i], component, n_comp);
         fe.get_function_hessians((*src)[level], make_slice(index, start, size), dst, true);
       }
   }
index 20f3cadd8e1013e2835673285ec3e1276ae1abe1..8d83c58d2a0d7bd0fa97c369c32fc2acc10dd652 100644 (file)
 // ---------------------------------------------------------------------
 
 
-for (deal_II_dimension : DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
 {
     namespace MeshWorker
     \{
-    template class IntegrationInfo<deal_II_dimension, deal_II_dimension>;
-    template class IntegrationInfoBox<deal_II_dimension, deal_II_dimension>;
+#if deal_II_dimension <= deal_II_space_dimension
+    template class IntegrationInfo<deal_II_dimension, deal_II_space_dimension>;
+    template class IntegrationInfoBox<deal_II_dimension, deal_II_space_dimension>;
 
-    template class DoFInfo<deal_II_dimension,deal_II_dimension,float>;
+    template class DoFInfo<deal_II_dimension, deal_II_space_dimension,float>;
     template class DoFInfoBox<deal_II_dimension,
-                              DoFInfo<deal_II_dimension,deal_II_dimension,float> >;
+                              DoFInfo<deal_II_dimension, deal_II_space_dimension,float> >;
 
-    template void IntegrationInfo<deal_II_dimension>::fill_local_data(
-        const DoFInfo<deal_II_dimension, deal_II_dimension, float>&, bool);
+    template void IntegrationInfo<deal_II_dimension, deal_II_space_dimension>::fill_local_data(
+        const DoFInfo<deal_II_dimension, deal_II_space_dimension, float>&, bool);
 
-    template class DoFInfo<deal_II_dimension,deal_II_dimension,double>;
+    template class DoFInfo<deal_II_dimension, deal_II_space_dimension,double>;
     template class DoFInfoBox<deal_II_dimension,
-                              DoFInfo<deal_II_dimension,deal_II_dimension,double> >;
+                              DoFInfo<deal_II_dimension, deal_II_space_dimension,double> >;
 
-    template void IntegrationInfo<deal_II_dimension>::fill_local_data(
-        const DoFInfo<deal_II_dimension, deal_II_dimension, double>&, bool);
-
-//   template void IntegrationInfo<deal_II_dimension>
-//   ::initialize<FEValues<deal_II_dimension> >(
-//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-//     const Quadrature<FEValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-//   template void IntegrationInfo<deal_II_dimension>
-//   ::initialize<FEFaceValues<deal_II_dimension> >(
-//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-//     const Quadrature<FEFaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-//   template void IntegrationInfo<deal_II_dimension>
-//   ::initialize<FESubfaceValues<deal_II_dimension> >(
-//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-//     const Quadrature<FESubfaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
+    template void IntegrationInfo<deal_II_dimension, deal_II_space_dimension>::fill_local_data(
+        const DoFInfo<deal_II_dimension, deal_II_space_dimension, double>&, bool);
+#endif
     \}
 }
index 0cbd1e740184d52c2c8b576d31e1ff976c85f845..e6490ba12847b94d5e1549d2c8973eef69e169c1 100644 (file)
 
 
 
-for (deal_II_dimension : DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
 {
     namespace MeshWorker
     \{
-    template class VectorDataBase<deal_II_dimension>;
+#if deal_II_dimension <= deal_II_space_dimension
+    template class VectorDataBase<deal_II_dimension, deal_II_space_dimension>;
+#endif
     \}
 }
 
-for (VECTOR : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS)
+for (VECTOR : SERIAL_VECTORS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
 {
     namespace MeshWorker
     \{
-    template class VectorData<VECTOR,deal_II_dimension>;
-    template class VectorData<const VECTOR,deal_II_dimension>;
-    template class MGVectorData<VECTOR,deal_II_dimension>;
-    template class MGVectorData<const VECTOR,deal_II_dimension>;
+#if deal_II_dimension <= deal_II_space_dimension
+    template class VectorData<VECTOR, deal_II_dimension, deal_II_space_dimension>;
+    template class VectorData<const VECTOR, deal_II_dimension, deal_II_space_dimension>;
+    template class MGVectorData<VECTOR, deal_II_dimension, deal_II_space_dimension>;
+    template class MGVectorData<const VECTOR, deal_II_dimension, deal_II_space_dimension>;
+#endif
     \}
 }

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.