]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add a simplified interface for MeshWorker::loop()
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 14 Sep 2012 08:37:12 +0000 (08:37 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 14 Sep 2012 08:37:12 +0000 (08:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@26371 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/meshworker/local_integrator.h [new file with mode: 0644]
deal.II/include/deal.II/meshworker/local_results.h
deal.II/include/deal.II/meshworker/loop.h
deal.II/source/numerics/mesh_worker.cc

diff --git a/deal.II/include/deal.II/meshworker/local_integrator.h b/deal.II/include/deal.II/meshworker/local_integrator.h
new file mode 100644 (file)
index 0000000..4777ff7
--- /dev/null
@@ -0,0 +1,111 @@
+//---------------------------------------------------------------------------
+//    $Id: loop.h 26042 2012-08-20 19:22:01Z bangerth $
+//
+//    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+#ifndef __deal2__mesh_worker_local_integrator_h
+#define __deal2__mesh_worker_local_integrator_h
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/std_cxx1x/function.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace MeshWorker
+{
+  template <int dim, int spacedim, typename number> class DoFInfo;
+  template <int dim, int spacedim> class IntegrationInfo;
+  
+/**
+ * A local integrator object, which can be used to simplify the call
+ * of loop(). Instead of providing the three local integration
+ * functions separately, we bundle them as virtual functions in this
+ * class.
+ *
+ * Additionally, since we cannot have a virtual null function, we
+ * provide flags, which allow us to indicate, whether we want to
+ * integrate on boundary and interior faces. Thes flags are true by
+ * default, but can be modified by applications to speed up the loop.
+ *
+ * @ingroup MeshWorker
+ * @author Guido Kanschat
+ * @date 2012
+ */
+  template <int dim, int spacedim=dim, typename number=double>
+    class LocalIntegrator : public Subscriptor
+  {
+    public:
+                                      /**
+                                       * The constructor setting
+                                       * default values.
+                                       */
+      LocalIntegrator();
+
+                                      /**
+                                       * The empty virtual destructor.
+                                       */
+      ~LocalIntegrator();
+      
+                                      /**
+                                       * Virtual function for
+                                       * integrating on cells.
+                                       */
+      virtual void cell(DoFInfo<dim, spacedim, number>& dinfo,
+                       IntegrationInfo<dim, spacedim>& info) const = 0;
+                                      /**
+                                       * Virtual function for
+                                       * integrating on boundary faces.
+                                       */
+      virtual void boundary(DoFInfo<dim, spacedim, number>& dinfo,
+                           IntegrationInfo<dim, spacedim>& info) const = 0;
+                                      /**
+                                       * Virtual function for
+                                       * integrating on interior faces.
+                                       */
+      virtual void face(DoFInfo<dim, spacedim, number>& dinfo1,
+                       DoFInfo<dim, spacedim, number>& dinfo2,
+                       IntegrationInfo<dim, spacedim>& info1,
+                       IntegrationInfo<dim, spacedim>& info2) const = 0;
+
+                                      /**
+                                       * The flag indicating whether
+                                       * the cell integrator cell()
+                                       * is to be used in the
+                                       * loop. Defaults to
+                                       * <tt>true</tt>.
+                                       */
+      bool use_cell;
+
+                                      /**
+                                       * The flag indicating whether
+                                       * the boundary integrator
+                                       * boundary() is to be
+                                       * used in the loop. Defaults
+                                       * to <tt>true</tt>.
+                                       */
+      bool use_boundary;
+
+                                      /**
+                                       * The flag indicating whether
+                                       * the interior face integrator
+                                       * face() is to be used in the
+                                       * loop. Defaults to
+                                       * <tt>true</tt>.
+                                       */
+      bool use_face;
+      
+  };
+}
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index cb27256847ff8e58f07603b871a80246a7a3c3d4..829f68840c7a7e58c252198976141a2ed9ae2be0 100644 (file)
@@ -31,7 +31,10 @@ template<int,int> class MGDoFHandler;
  * an ubiquitous part of each finite element program.
  *
  * The workhorse of this namespace is the loop() function, which implements a
- * completely generic loop over all mesh cells.
+ * completely generic loop over all mesh cells. Since the calls to
+ * loop() are error-prone due to its generality, for many applications
+ * it is advisable to derive a class from MeshWorker::LocalIntegrator
+ * and use the less general integration_loop() instead.
  *
  * The loop() depends on certain objects handed to it as
  * arguments. These objects are of two types, info objects like
index 1c1f648287f3b120c27078dc33985594d6decf64..c3d2d3224cf4f70d6555e4317aa41f98b4806c79 100644 (file)
 #include <deal.II/base/work_stream.h>
 #include <deal.II/base/template_constraints.h>
 #include <deal.II/grid/tria.h>
+#include <deal.II/meshworker/local_integrator.h>
 #include <deal.II/meshworker/dof_info.h>
 #include <deal.II/meshworker/integration_info.h>
 
+
 #define DEAL_II_MESHWORKER_PARALLEL 1
 
 DEAL_II_NAMESPACE_OPEN
@@ -282,6 +284,9 @@ namespace MeshWorker
   }
 
 /**
+ * @deprecated The simplification in this loop is
+ * insignificant. Therefore, it is recommended to use loop() instead.
+ *
  * Simplified interface for loop() if specialized for integration.
  *
  * @ingroup MeshWorker
@@ -310,6 +315,49 @@ namespace MeshWorker
        assembler,
        cells_first);
   }
+
+
+/**
+ * Simplified interface for loop() if specialized for integration,
+ * using the virtual functions in LocalIntegrator.
+ *
+ * @ingroup MeshWorker
+ * @author Guido Kanschat, 2009
+ */
+  template<int dim, int spacedim, class ITERATOR, class ASSEMBLER>
+  void integration_loop(ITERATOR begin,
+                        typename identity<ITERATOR>::type end,
+                        DoFInfo<dim, spacedim>& dof_info,
+                        IntegrationInfoBox<dim, spacedim>& box,
+                       const LocalIntegrator<dim, spacedim>& integrator,
+                        ASSEMBLER &assembler,
+                        bool cells_first = true)
+  {
+    std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> cell_worker;
+    std_cxx1x::function<void (DoFInfo<dim>&, IntegrationInfo<dim, spacedim>&)> boundary_worker;
+    std_cxx1x::function<void (DoFInfo<dim>&, DoFInfo<dim>&,
+                             IntegrationInfo<dim, spacedim>&,
+                             IntegrationInfo<dim, spacedim>&)> face_worker;
+  if (integrator.use_cell)
+    cell_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::cell, &integrator, std_cxx1x::_1, std_cxx1x::_2);
+  if (integrator.use_boundary)
+    boundary_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::boundary, &integrator, std_cxx1x::_1, std_cxx1x::_2);
+  if (integrator.use_face)
+    face_worker = std_cxx1x::bind(&LocalIntegrator<dim, spacedim>::face, &integrator, std_cxx1x::_1, std_cxx1x::_2, std_cxx1x::_3, std_cxx1x::_4);
+
+  loop<dim, spacedim>
+      (begin, end,
+       dof_info,
+       box,
+       cell_worker,
+       boundary_worker,
+       face_worker,
+       assembler,
+       cells_first);
+  }
+
+
+
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 1ba1bbd28d64c0f8f541055e9911818936ff4cb1..bd92f34bbf4853afa359edea2ee74cc8068e1263 100644 (file)
@@ -11,6 +11,7 @@
 //---------------------------------------------------------------------------
 
 #include <deal.II/meshworker/local_results.h>
+#include <deal.II/meshworker/local_integrator.h>
 #include <deal.II/lac/block_indices.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -51,6 +52,30 @@ namespace MeshWorker
 
   template class LocalResults<float>;
   template class LocalResults<double>;
+
+  template <int dim, int spacedim, typename number>
+  LocalIntegrator<dim, spacedim, number>::LocalIntegrator ()
+                 :
+                 use_cell(true), use_boundary(true), use_face(true)
+  {}
+
+  
+  template <int dim, int spacedim, typename number>
+  LocalIntegrator<dim, spacedim, number>::~LocalIntegrator ()
+  {}
+  
+  template class LocalIntegrator<1,1,float>;
+  template class LocalIntegrator<1,1,double>;
+  template class LocalIntegrator<1,2,float>;
+  template class LocalIntegrator<1,2,double>;
+  template class LocalIntegrator<1,3,float>;
+  template class LocalIntegrator<1,3,double>;
+  template class LocalIntegrator<2,2,float>;
+  template class LocalIntegrator<2,2,double>;
+  template class LocalIntegrator<2,3,float>;
+  template class LocalIntegrator<2,3,double>;
+  template class LocalIntegrator<3,3,float>;
+  template class LocalIntegrator<3,3,double>;
 }
 
 

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.