]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix comments by Wolfgang.
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 10 Aug 2017 23:43:36 +0000 (17:43 -0600)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 10 Aug 2017 23:43:36 +0000 (17:43 -0600)
include/deal.II/meshworker/assemble_flags.h
include/deal.II/meshworker/mesh_loop.h
tests/meshworker/mesh_loop_01.cc

index fa157b17f504a2d8cc801212531604424b5aedee..4fcd6e86e0e514d4f0e23d070a67635b4b89d805 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 1998 - 2016 by the deal.II authors
+// Copyright (C) 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
 //
 // ---------------------------------------------------------------------
 
-#ifndef dealii__assemble_flags_h
-#define dealii__assemble_flags_h
+#ifndef dealii__mesh_worker_assemble_flags_h
+#define dealii__mesh_worker_assemble_flags_h
 
 
 #include <deal.II/base/config.h>
-#include <deal.II/base/table.h>
-#include <deal.II/base/derivative_form.h>
-#include <deal.II/base/point.h>
-#include <deal.II/base/tensor.h>
-
 #include <vector>
 
 
@@ -31,175 +26,173 @@ DEAL_II_NAMESPACE_OPEN
 /*!@addtogroup MeshWorker */
 /*@{*/
 
-/**
- * The enum type given to the MeshLoop function, telling that function which
- * elements need to be assembled.
- *
- * You can select more than one flag by concatenation using the bitwise or
- * operator|(AssembleFlags,AssembleFlags).
- *
- * <h3>Use of these flags flags</h3>
- *
- * More information on the use of this type both in user code as well as
- * internally can be found in the documentation modules on
- * @ref AssembleFlags "The interplay of AssembleFlags, Mapping, and FiniteElement in FEValues"
- * and
- * @ref FE_vs_Mapping_vs_FEValues "How Mapping, FiniteElement, and FEValues work together".
- */
-enum AssembleFlags
+namespace MeshWorker
 {
-  //! No update
-  assemble_default = 0,
-  //! Own cells
-  /**
-   * Assemble on cells.
-   */
-  assemble_own_cells = 0x0001,
-  //! Ghost cells
-  /**
-   * Assemble on ghost cells.
-   */
-  assemble_ghost_cells = 0x0002,
-  //! Own faces once
-  /**
-   * Assemble on own interior faces, visiting each face only once
-   */
-  assemble_own_interior_faces_once = 0x0004,
-  //! Own faces both
-  /**
-   * Assemble on own interior faces, visiting each internal face twice
-   */
-  assemble_own_interior_faces_both = 0x0008,
-  //! Ghost faces once
-  /**
-   * Assemble on faces between a locally owned cell and a ghost cell, making
-   * sure that only one of the processes will assemble these faces (from the
-   * finer side or the process with the lower mpi rank
-   */
-  assemble_ghost_faces_once = 0x0010,
-  //! Ghost faces both
-  /**
-   * Assemble on faces between a locally owned cell and a ghost cell, both
-   * processes will assemble these faces. Note that these faces are never
-   * assembled from both sides on a single process.
-   */
-  assemble_ghost_faces_both = 0x0020,
-  //! Assemble cells before faces
-  /**
-   * Assemble cells integrals before face integrals.
-   */
-  assemble_cells_first = 0x0040,
-  //! Assemble boundary faces
-  /**
-   * Assemble on boundary faces
-   */
-  assemble_boundary_faces = 0x0080,
-  //! Assemble own interior aces
-  /**
-   * Assemble own interior faces, either interior ones or on the boundary.
-   */
-  assemble_own_interior_faces = assemble_own_interior_faces_both | assemble_own_interior_faces_once,
-  //! Assemble own faces
+
   /**
-   * Assemble own faces, either interior ones or on the boundary.
+   * The enum type given to the mesh_loop() function, telling that function which
+   * elements need to be assembled.
+   *
+   * You can select more than one flag by concatenation using the bitwise or
+   * `operator|(AssembleFlags,AssembleFlags)`.
+   *
+   * @author Luca Heltai, 2017.
    */
-  assemble_own_faces = assemble_own_interior_faces | assemble_boundary_faces,
-  //! Assemble ghost faces
+  enum AssembleFlags
+  {
+    //! No update
+    assemble_default = 0,
+    //! Own cells
+    /**
+     * Assemble on cells.
+     */
+    assemble_own_cells = 0x0001,
+    //! Ghost cells
+    /**
+     * Assemble on ghost cells.
+     */
+    assemble_ghost_cells = 0x0002,
+    //! Own faces once
+    /**
+     * Assemble on own interior faces, visiting each face only once
+     */
+    assemble_own_interior_faces_once = 0x0004,
+    //! Own faces both
+    /**
+     * Assemble on own interior faces, visiting each interior face twice
+     */
+    assemble_own_interior_faces_both = 0x0008,
+    //! Ghost faces once
+    /**
+     * Assemble on faces between a locally owned cell and a ghost cell, making
+     * sure that only one of the processes will assemble these faces (from the
+     * finer side or the process with the lower mpi rank)
+     */
+    assemble_ghost_faces_once = 0x0010,
+    //! Ghost faces both
+    /**
+     * Assemble on faces between a locally owned cell and a ghost cell, both
+     * processes will assemble these faces. Note that these faces are never
+     * assembled from both sides on a single process.
+     */
+    assemble_ghost_faces_both = 0x0020,
+    //! Assemble cells before faces
+    /**
+     * Assemble cell integrals before face integrals.
+     */
+    assemble_cells_first = 0x0040,
+    //! Assemble boundary faces
+    /**
+     * Assemble on boundary faces
+     */
+    assemble_boundary_faces = 0x0080,
+    //! Assemble own interior aces
+    /**
+     * Assemble own interior faces, either interior ones or on the boundary.
+     */
+    assemble_own_interior_faces = assemble_own_interior_faces_both | assemble_own_interior_faces_once,
+    //! Assemble own faces
+    /**
+     * Assemble own faces, either interior ones or on the boundary.
+     */
+    assemble_own_faces = assemble_own_interior_faces | assemble_boundary_faces,
+    //! Assemble ghost faces
+    /**
+     * Assemble ghost faces
+     */
+    assemble_ghost_faces = assemble_ghost_faces_both | assemble_ghost_faces_once,
+  };
+
+
   /**
-   * Assemble ghost faces
+   * Output operator which outputs update flags as a set of or'd text values.
+   *
+   * @ref AssembleFlags
    */
-  assemble_ghost_faces = assemble_ghost_faces_both | assemble_ghost_faces_once,
-};
-
-
-/**
- * Output operator which outputs update flags as a set of or'd text values.
- *
- * @ref AssembleFlags
- */
-template <class StreamType>
-inline
-StreamType &operator << (StreamType &s, AssembleFlags u)
-{
-  s << " AssembleFlags";
-  if (u & assemble_own_cells        ) s << "|own_cells"        ;
-  if (u & assemble_own_interior_faces_once   ) s << "|own_faces_once"   ;
-  if (u & assemble_own_interior_faces_both   ) s << "|own_faces_both"   ;
-  if (u & assemble_ghost_cells      ) s << "|ghost_cells"      ;
-  if (u & assemble_ghost_faces_once ) s << "|ghost_faces_once" ;
-  if (u & assemble_ghost_faces_both ) s << "|ghost_faces_both" ;
-  if (u & assemble_boundary_faces   ) s << "|boundary_faces"   ;
-  return s;
-}
+  template <class StreamType>
+  inline
+  StreamType &operator << (StreamType &s, AssembleFlags u)
+  {
+    s << " AssembleFlags";
+    if (u & assemble_own_cells        ) s << "|own_cells"        ;
+    if (u & assemble_own_interior_faces_once   ) s << "|own_faces_once"   ;
+    if (u & assemble_own_interior_faces_both   ) s << "|own_faces_both"   ;
+    if (u & assemble_ghost_cells      ) s << "|ghost_cells"      ;
+    if (u & assemble_ghost_faces_once ) s << "|ghost_faces_once" ;
+    if (u & assemble_ghost_faces_both ) s << "|ghost_faces_both" ;
+    if (u & assemble_boundary_faces   ) s << "|boundary_faces"   ;
+    return s;
+  }
 
 
-/**
- * Global operator which returns an object in which all bits are set which are
- * either set in the first or the second argument. This operator exists since
- * if it did not then the result of the bit-or <tt>operator |</tt> would be an
- * integer which would in turn trigger a compiler warning when we tried to
- * assign it to an object of type AssembleFlags.
- *
- * @ref AssembleFlags
- */
-inline
-AssembleFlags
-operator | (AssembleFlags f1, AssembleFlags f2)
-{
-  return static_cast<AssembleFlags> (
-           static_cast<unsigned int> (f1) |
-           static_cast<unsigned int> (f2));
-}
+  /**
  * Global operator which returns an object in which all bits are set which are
  * either set in the first or the second argument. This operator exists since
  * if it did not then the result of the bit-or <tt>operator |</tt> would be an
  * integer which would in turn trigger a compiler warning when we tried to
  * assign it to an object of type AssembleFlags.
  *
  * @ref AssembleFlags
  */
+  inline
+  AssembleFlags
+  operator | (AssembleFlags f1, AssembleFlags f2)
+  {
+    return static_cast<AssembleFlags> (
+             static_cast<unsigned int> (f1) |
+             static_cast<unsigned int> (f2));
+  }
 
 
 
 
-/**
- * Global operator which sets the bits from the second argument also in the
- * first one.
- *
- * @ref AssembleFlags
- */
-inline
-AssembleFlags &
-operator |= (AssembleFlags &f1, AssembleFlags f2)
-{
-  f1 = f1 | f2;
-  return f1;
-}
+  /**
  * Global operator which sets the bits from the second argument also in the
  * first one.
  *
  * @ref AssembleFlags
  */
+  inline
+  AssembleFlags &
+  operator |= (AssembleFlags &f1, AssembleFlags f2)
+  {
+    f1 = f1 | f2;
+    return f1;
+  }
 
 
-/**
- * Global operator which returns an object in which all bits are set which are
- * set in the first as well as the second argument. This operator exists since
- * if it did not then the result of the bit-and <tt>operator &</tt> would be
- * an integer which would in turn trigger a compiler warning when we tried to
- * assign it to an object of type AssembleFlags.
- *
- * @ref AssembleFlags
- */
-inline
-AssembleFlags
-operator & (AssembleFlags f1, AssembleFlags f2)
-{
-  return static_cast<AssembleFlags> (
-           static_cast<unsigned int> (f1) &
-           static_cast<unsigned int> (f2));
-}
+  /**
  * Global operator which returns an object in which all bits are set which are
  * set in the first as well as the second argument. This operator exists since
  * if it did not then the result of the bit-and <tt>operator &</tt> would be
  * an integer which would in turn trigger a compiler warning when we tried to
  * assign it to an object of type AssembleFlags.
  *
  * @ref AssembleFlags
  */
+  inline
+  AssembleFlags
+  operator & (AssembleFlags f1, AssembleFlags f2)
+  {
+    return static_cast<AssembleFlags> (
+             static_cast<unsigned int> (f1) &
+             static_cast<unsigned int> (f2));
+  }
 
 
-/**
- * Global operator which clears all the bits in the first argument if they are
- * not also set in the second argument.
- *
- * @ref AssembleFlags
- */
-inline
-AssembleFlags &
-operator &= (AssembleFlags &f1, AssembleFlags f2)
-{
-  f1 = f1 & f2;
-  return f1;
+  /**
+   * Global operator which clears all the bits in the first argument if they are
+   * not also set in the second argument.
+   *
+   * @ref AssembleFlags
+   */
+  inline
+  AssembleFlags &
+  operator &= (AssembleFlags &f1, AssembleFlags f2)
+  {
+    f1 = f1 & f2;
+    return f1;
+  }
 }
 
 /*@}*/
index 7631af58d541d507cd7b4a4ed7877e9c16d04782..145e436f78dacf57f1b71953073424eda4601efc 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2006 - 2016 by the deal.II authors
+// Copyright (C) 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -14,8 +14,8 @@
 // ---------------------------------------------------------------------
 
 
-#ifndef dealii__mesh_loop_h
-#define dealii__mesh_loop_h
+#ifndef dealii__mesh_worker_mesh_loop_h
+#define dealii__mesh_worker_mesh_loop_h
 
 #include <deal.II/base/config.h>
 #include <deal.II/base/work_stream.h>
@@ -42,51 +42,44 @@ namespace MeshWorker
    * @ingroup MeshWorker
    * @author Luca Heltai, 2017
    */
-  template <class Iterator,
-            class CellWorker, class Copyer,
-            class ScratchData, class CopyData,
-            class BoundaryWorker, class FaceWorker>
-  void mesh_loop(const Iterator &begin,
-                 const typename identity<Iterator>::type &end,
-
-//                 const std::function<void (const Iterator &, ScratchData &, CopyData &)> &cell_worker,
-//                 const std::function<void (const CopyData &)> &copyer,
-                 CellWorker cell_worker,
-                 Copyer copyer,
+  template <class CellIteratorType,
+            class ScratchData, class CopyData>
+  void mesh_loop(const CellIteratorType &begin,
+                 const typename identity<CellIteratorType>::type &end,
+
+                 const typename identity<std::function<void (const CellIteratorType &, ScratchData &, CopyData &)>>::type &cell_worker,
+                 const typename identity<std::function<void (const CopyData &)>>::type &copyer,
 
                  ScratchData &scratch_data,
                  CopyData &copy_data,
 
                  const AssembleFlags flags = assemble_own_cells,
 
-//                 const std::function<void (const Iterator &, const unsigned int&, ScratchData &, CopyData &)> &boundary_worker=
-//                  std::function<void (const Iterator &, const unsigned int&, ScratchData &, CopyData &)>(),
-
-//                 const std::function<void (const Iterator &, const unsigned int &, const unsigned int &,
-//                                           const Iterator &, const unsigned int &, const unsigned int &,
-//                                           ScratchData &, CopyData &)> &face_worker=
-//                  std::function<void (const Iterator &, unsigned int, unsigned int,
-//                                      const Iterator &, unsigned int, unsigned int,
-//                                      ScratchData &, CopyData &)>(),
+                 const typename identity<std::function<void (const CellIteratorType &, const unsigned int&, ScratchData &, CopyData &)>>::type &boundary_worker=
+                  std::function<void (const CellIteratorType &, const unsigned int&, ScratchData &, CopyData &)>(),
 
-                 BoundaryWorker boundary_worker = BoundaryWorker(),
-                 FaceWorker face_worker = FaceWorker(),
+                 const typename identity<std::function<void (const CellIteratorType &, const unsigned int &, const unsigned int &,
+                                           const CellIteratorType &, const unsigned int &, const unsigned int &,
+                                           ScratchData &, CopyData &)>>::type &face_worker=
+                  std::function<void (const CellIteratorType &, const unsigned int &, const unsigned int &,
+                                      const CellIteratorType &, const unsigned int &, const unsigned int &,
+                                      ScratchData &, CopyData &)>(),
 
                  const unsigned int   queue_length = 2*MultithreadInfo::n_threads(),
                  const unsigned int   chunk_size = 8)
   {
-    auto cell_action = [&] (const Iterator &cell, ScratchData &scratch, CopyData &copy)
+    auto cell_action = [&] (const CellIteratorType &cell, ScratchData &scratch, CopyData &copy)
     {
       const bool ignore_subdomain = (cell->get_triangulation().locally_owned_subdomain()
                                      == numbers::invalid_subdomain_id);
 
-      types::subdomain_id csid = /*(cell->is_level_cell())
+      types::subdomain_id current_subdomain_id = /*(cell->is_level_cell())
                                  ? cell->level_subdomain_id()
                                  : */cell->subdomain_id();
 
-      const bool own_cell = ignore_subdomain || (csid == cell->get_triangulation().locally_owned_subdomain());
+      const bool own_cell = ignore_subdomain || (current_subdomain_id == cell->get_triangulation().locally_owned_subdomain());
 
-      if ((!ignore_subdomain) && (csid == numbers::artificial_subdomain_id))
+      if ((!ignore_subdomain) && (current_subdomain_id == numbers::artificial_subdomain_id))
         return;
 
       if ( (flags & (assemble_cells_first)) &&
@@ -95,9 +88,9 @@ namespace MeshWorker
         cell_worker(cell, scratch, copy);
 
       if (flags & assemble_own_faces)
-        for (unsigned int face_no=0; face_no < GeometryInfo<Iterator::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
+        for (unsigned int face_no=0; face_no < GeometryInfo<CellIteratorType::AccessorType::Container::dimension>::faces_per_cell; ++face_no)
           {
-            typename Iterator::AccessorType::Container::face_iterator face = cell->face(face_no);
+            typename CellIteratorType::AccessorType::Container::face_iterator face = cell->face(face_no);
             if (cell->at_boundary(face_no) && !cell->has_periodic_neighbor(face_no))
               {
                 // only integrate boundary faces of own cells
@@ -109,17 +102,17 @@ namespace MeshWorker
             else if (flags & assemble_own_interior_faces)
               {
                 // Interior face
-                TriaIterator<typename Iterator::AccessorType> neighbor = cell->neighbor_or_periodic_neighbor(face_no);
+                TriaIterator<typename CellIteratorType::AccessorType> neighbor = cell->neighbor_or_periodic_neighbor(face_no);
 
-                types::subdomain_id neighbid = numbers::artificial_subdomain_id;
+                types::subdomain_id neighbor_subdomain_id = numbers::artificial_subdomain_id;
 //                if (neighbor->is_level_cell())
 //                  neighbid = neighbor->level_subdomain_id();
                 //subdomain id is only valid for active cells
                 /*else */if (neighbor->active())
-                  neighbid = neighbor->subdomain_id();
+                  neighbor_subdomain_id = neighbor->subdomain_id();
 
                 const bool own_neighbor = ignore_subdomain ||
-                                          (neighbid == cell->get_triangulation().locally_owned_subdomain());
+                                          (neighbor_subdomain_id == cell->get_triangulation().locally_owned_subdomain());
 
                 // skip all faces between two ghost cells
                 if (!own_cell && !own_neighbor)
@@ -152,7 +145,7 @@ namespace MeshWorker
                       = periodic_neighbor?
                         cell->periodic_neighbor_of_coarser_periodic_neighbor(face_no):
                         cell->neighbor_of_coarser_neighbor(face_no);
-                    const typename Iterator::AccessorType::Container::face_iterator nface
+                    const typename CellIteratorType::AccessorType::Container::face_iterator nface
                       = neighbor->face(neighbor_face_no.first);
 
                     face_worker(cell, face_no, numbers::invalid_unsigned_int,
@@ -171,9 +164,7 @@ namespace MeshWorker
                     // If iterator is active and neighbor is refined, skip
                     // internal face.
                     if (internal::is_active_iterator(cell) && neighbor->has_children())
-                      {
-                        continue;
-                      }
+                      continue;
 
                     // Now neighbor is on same level, double-check this:
                     Assert(cell->level()==neighbor->level(), ExcInternalError());
@@ -197,7 +188,7 @@ namespace MeshWorker
                     // face.
                     if (own_cell && !own_neighbor
                         && (flags & assemble_ghost_faces_once)
-                        && (neighbid < csid))
+                        && (neighbor_subdomain_id < current_subdomain_id))
                       continue;
 
                     const unsigned int neighbor_face_no = periodic_neighbor?
index d4da345483b427df858646d3a744a6e3e30a66ed..dd1e2270770b3499d43428bef0d8549b8703a072 100644 (file)
@@ -22,6 +22,7 @@
 
 #include <fstream>
 
+using namespace MeshWorker;
 
 int main()
 {

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.