--- /dev/null
+//---------------------------------------------------------------------------
+// $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
* 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
#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
}
/**
+ * @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
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
//---------------------------------------------------------------------------
#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
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>;
}