]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new filtered iterator LocallyOwnedLevelCell
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 24 Jun 2013 13:17:37 +0000 (13:17 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 24 Jun 2013 13:17:37 +0000 (13:17 +0000)
git-svn-id: https://svn.dealii.org/trunk@29871 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-50/step-50.cc
deal.II/include/deal.II/grid/filtered_iterator.h
deal.II/include/deal.II/grid/tria_accessor.h
deal.II/include/deal.II/grid/tria_accessor.templates.h

index c3ffb9c06456585d432ca173ce074680ca88eec2..937e93846b75366e25a81df0f33db8da39c165a8 100644 (file)
@@ -695,51 +695,35 @@ namespace Step50
                                        local_dof_indices,
                                        mg_matrices[cell->level()]);
 
-          // The next step is again slightly more
-          // obscure (but explained in the @ref
-          // mg_paper): We need the remainder of
-          // the operator that we just copied
-          // into the <code>mg_matrices</code>
-          // object, namely the part on the
-          // interface between cells at the
-          // current level and cells one level
-          // coarser. This matrix exists in two
-          // directions: for interior DoFs (index
-          // $i$) of the current level to those
-          // sitting on the interface (index
-          // $j$), and the other way around. Of
-          // course, since we have a symmetric
-          // operator, one of these matrices is
-          // the transpose of the other.
+          // The next step is again slightly more obscure (but
+          // explained in the @ref mg_paper): We need the remainder of
+          // the operator that we just copied into the
+          // <code>mg_matrices</code> object, namely the part on the
+          // interface between cells at the current level and cells
+          // one level coarser. This matrix exists in two directions:
+          // for interior DoFs (index $i$) of the current level to
+          // those sitting on the interface (index $j$), and the other
+          // way around. Of course, since we have a symmetric
+          // operator, one of these matrices is the transpose of the
+          // other.
           //
-          // The way we assemble these matrices
-          // is as follows: since the are formed
-          // from parts of the local
-          // contributions, we first delete all
-          // those parts of the local
-          // contributions that we are not
-          // interested in, namely all those
-          // elements of the local matrix for
-          // which not $i$ is an interface DoF
-          // and $j$ is not. The result is one of
-          // the two matrices that we are
-          // interested in, and we then copy it
-          // into the
-          // <code>mg_interface_matrices</code>
-          // object. The
-          // <code>boundary_interface_constraints</code>
-          // object at the same time makes sure
-          // that we delete contributions from
-          // all degrees of freedom that are not
-          // only on the interface but also on
-          // the external boundary of the domain.
+          // The way we assemble these matrices is as follows: since
+          // the are formed from parts of the local contributions, we
+          // first delete all those parts of the local contributions
+          // that we are not interested in, namely all those elements
+          // of the local matrix for which not $i$ is an interface DoF
+          // and $j$ is not. The result is one of the two matrices
+          // that we are interested in, and we then copy it into the
+          // <code>mg_interface_matrices</code> object. The
+          // <code>boundary_interface_constraints</code> object at the
+          // same time makes sure that we delete contributions from
+          // all degrees of freedom that are not only on the interface
+          // but also on the external boundary of the domain.
           //
-          // The last part to remember is how to
-          // get the other matrix. Since it is
-          // only the transpose, we will later
-          // (in the <code>solve()</code>
-          // function) be able to just pass the
-          // transpose matrix where necessary.
+          // The last part to remember is how to get the other
+          // matrix. Since it is only the transpose, we will later (in
+          // the <code>solve()</code> function) be able to just pass
+          // the transpose matrix where necessary.
           for (unsigned int i=0; i<dofs_per_cell; ++i)
             for (unsigned int j=0; j<dofs_per_cell; ++j)
               if ( !(interface_dofs[cell->level()][local_dof_indices[i]]==true &&
@@ -763,51 +747,40 @@ namespace Step50
 
   // @sect4{LaplaceProblem::solve}
 
-  // This is the other function that is
-  // significantly different in support of the
-  // multigrid solver (or, in fact, the
-  // preconditioner for which we use the
-  // multigrid method).
+  // This is the other function that is significantly different in
+  // support of the multigrid solver (or, in fact, the preconditioner
+  // for which we use the multigrid method).
   //
-  // Let us start out by setting up two of the
-  // components of multilevel methods: transfer
-  // operators between levels, and a solver on
-  // the coarsest level. In finite element
-  // methods, the transfer operators are
-  // derived from the finite element function
-  // spaces involved and can often be computed
-  // in a generic way independent of the
-  // problem under consideration. In that case,
-  // we can use the MGTransferPrebuilt class
-  // that, given the constraints on the global
-  // level and an MGDoFHandler object computes
-  // the matrices corresponding to these
-  // transfer operators.
+  // Let us start out by setting up two of the components of
+  // multilevel methods: transfer operators between levels, and a
+  // solver on the coarsest level. In finite element methods, the
+  // transfer operators are derived from the finite element function
+  // spaces involved and can often be computed in a generic way
+  // independent of the problem under consideration. In that case, we
+  // can use the MGTransferPrebuilt class that, given the constraints
+  // on the global level and an MGDoFHandler object computes the
+  // matrices corresponding to these transfer operators.
   //
-  // The second part of the following lines
-  // deals with the coarse grid solver. Since
-  // our coarse grid is very coarse indeed, we
-  // decide for a direct solver (a Householder
-  // decomposition of the coarsest level
-  // matrix), even if its implementation is not
-  // particularly sophisticated. If our coarse
-  // mesh had many more cells than the five we
-  // have here, something better suited would
-  // obviously be necessary here.
+  // The second part of the following lines deals with the coarse grid
+  // solver. Since our coarse grid is very coarse indeed, we decide
+  // for a direct solver (a Householder decomposition of the coarsest
+  // level matrix), even if its implementation is not particularly
+  // sophisticated. If our coarse mesh had many more cells than the
+  // five we have here, something better suited would obviously be
+  // necessary here.
   template <int dim>
   void LaplaceProblem<dim>::solve ()
   {
 
-    // Create the object that deals with the transfer
-    // between different refinement levels. We need to
-    // pass it the hanging node constraints.
+    // Create the object that deals with the transfer between
+    // different refinement levels. We need to pass it the hanging
+    // node constraints.
     MGTransferPrebuilt<vector_t> mg_transfer(hanging_node_constraints, mg_constrained_dofs);
-    // Now the prolongation matrix has to be built.
-    // This matrix needs to take the boundary values on
-    // each level into account and needs to know about
-    // the indices at the refinement egdes. The
-    // <code>MGConstraints</code> knows about that so
-    // pass it as an argument.
+    // Now the prolongation matrix has to be built.  This matrix needs
+    // to take the boundary values on each level into account and
+    // needs to know about the indices at the refinement egdes. The
+    // <code>MGConstraints</code> knows about that so pass it as an
+    // argument.
     mg_transfer.build_matrices(mg_dof_handler);
 
     matrix_t & coarse_matrix = mg_matrices[0];
@@ -822,38 +795,27 @@ namespace Step50
         coarse_matrix,
         id);
 
-    // The next component of a multilevel
-    // solver or preconditioner is that we need
-    // a smoother on each level. A common
-    // choice for this is to use the
-    // application of a relaxation method (such
-    // as the SOR, Jacobi or Richardson method). The
-    // MGSmootherPrecondition class provides
-    // support for this kind of
-    // smoother. Here, we opt for the
-    // application of a single SOR
-    // iteration. To this end, we define an
-    // appropriate <code>typedef</code> and
-    // then setup a smoother object.
+    // The next component of a multilevel solver or preconditioner is
+    // that we need a smoother on each level. A common choice for this
+    // is to use the application of a relaxation method (such as the
+    // SOR, Jacobi or Richardson method). The MGSmootherPrecondition
+    // class provides support for this kind of smoother. Here, we opt
+    // for the application of a single SOR iteration. To this end, we
+    // define an appropriate <code>typedef</code> and then setup a
+    // smoother object.
     //
-    // The last step is to initialize the
-    // smoother object with our level matrices
-    // and to set some smoothing parameters.
-    // The <code>initialize()</code> function
-    // can optionally take additional arguments
-    // that will be passed to the smoother
-    // object on each level. In the current
-    // case for the SOR smoother, this could,
-    // for example, include a relaxation
-    // parameter. However, we here leave these
-    // at their default values. The call to
-    // <code>set_steps()</code> indicates that
-    // we will use two pre- and two
-    // post-smoothing steps on each level; to
-    // use a variable number of smoother steps
-    // on different levels, more options can be
-    // set in the constructor call to the
-    // <code>mg_smoother</code> object.
+    // The last step is to initialize the smoother object with our
+    // level matrices and to set some smoothing parameters.  The
+    // <code>initialize()</code> function can optionally take
+    // additional arguments that will be passed to the smoother object
+    // on each level. In the current case for the SOR smoother, this
+    // could, for example, include a relaxation parameter. However, we
+    // here leave these at their default values. The call to
+    // <code>set_steps()</code> indicates that we will use two pre-
+    // and two post-smoothing steps on each level; to use a variable
+    // number of smoother steps on different levels, more options can
+    // be set in the constructor call to the <code>mg_smoother</code>
+    // object.
     //
     // The last step results from the fact that
     // we use the SOR method as a smoother -
index 2117169cb5e68fb2e8fcb44112c920f822887809..e321f635c450acd6c3afcd84c2ef77b1b263705e 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2011, 2012 by the deal.II authors
+//    Copyright (C) 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2011, 2012, 2013 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -22,17 +22,16 @@ DEAL_II_NAMESPACE_OPEN
 
 /**
  * In this namespace a number of classes is declared that may be used
- * as filters in the FilteredIterator class. The filters either
- * check for binary information (for example, the IteratorFilters::Active filter
- * class checks whether the object pointed to is active), or for
- * valued information by comparison with prescribed values (for
- * example, the LevelEqualTo filter class checks whether the
- * level of the object pointed to by the iterator under consideration
- * is equal to a value that was given to the filter upon construction.
+ * as filters in the FilteredIterator class. The filters either check
+ * for binary information (for example, the IteratorFilters::Active
+ * filter class checks whether the object pointed to is active), or
+ * for valued information by comparison with prescribed values (for
+ * example, the LevelEqualTo filter class checks whether the level of
+ * the object pointed to by the iterator under consideration is equal
+ * to a value that was given to the filter upon construction.
  *
  * For examples of use of these classes as well as requirements on
- * filters see the general description of the FilteredIterator
- * class.
+ * filters see the general description of the FilteredIterator class.
  *
  * @ingroup Iterators
  * @author Wolfgang Bangerth, 2002
@@ -40,10 +39,8 @@ DEAL_II_NAMESPACE_OPEN
 namespace IteratorFilters
 {
   /**
-   * Filter that evaluates to true if
-   * either the iterator points to an
-   * active object or an iterator
-   * past the end.
+   * Filter that evaluates to true if either the iterator points to an
+   * active object or an iterator past the end.
    *
    * @ingroup Iterators
    */
@@ -51,21 +48,17 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Evaluate the iterator and
-     * return true if the object is
-     * active or past the end.
+     * Evaluate the iterator and return true if the object is active
+     * or past the end.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
   };
 
   /**
-   * Filter that evaluates to true if
-   * either the iterator points to an
-   * object for which the user flag
-   * is set or an iterator past the
-   * end. See @ref GlossUserFlags
-  * for information about user flags.
+   * Filter that evaluates to true if either the iterator points to an
+   * object for which the user flag is set or an iterator past the
+   * end. See @ref GlossUserFlags for information about user flags.
    *
    * @ingroup Iterators
    */
@@ -73,10 +66,8 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Evaluate the iterator and
-     * return true if the object
-     * has a set user flag or past
-     * the end.
+     * Evaluate the iterator and return true if the object has a set
+     * user flag or past the end.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
@@ -84,12 +75,9 @@ namespace IteratorFilters
 
 
   /**
-   * Filter that evaluates to true if
-   * either the iterator points to an
-   * object for which the user flag
-   * is not set or an iterator past
-   * the end. Inverse filter to the
-   * previous class.
+   * Filter that evaluates to true if either the iterator points to an
+   * object for which the user flag is not set or an iterator past the
+   * end. Inverse filter to the previous class.
    *
    * @ingroup Iterators
    */
@@ -97,10 +85,8 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Evaluate the iterator and
-     * return true if the object
-     * has an unset user flag or
-     * past the end.
+     * Evaluate the iterator and return true if the object has an
+     * unset user flag or past the end.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
@@ -108,12 +94,9 @@ namespace IteratorFilters
 
 
   /**
-   * Filter for iterators that
-   * evaluates to true if either the
-   * iterator is past the end or the
-   * level of the object pointed to
-   * is equal to a value given to the
-   * constructor.
+   * Filter for iterators that evaluates to true if either the
+   * iterator is past the end or the level of the object pointed to is
+   * equal to a value given to the constructor.
    *
    * @ingroup Iterators
    */
@@ -121,27 +104,22 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Constructor. Store the level
-     * which iterators shall have
-     * to be evaluated to true.
+     * Constructor. Store the level which iterators shall have to be
+     * evaluated to true.
      */
     LevelEqualTo (const unsigned int level);
 
     /**
-     * Evaluation operator. Returns
-     * true if either the level of
-     * the object pointed to is
-     * equal to the stored value or
-     * the iterator is past the
-     * end.
+     * Evaluation operator. Returns true if either the level of the
+     * object pointed to is equal to the stored value or the iterator
+     * is past the end.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
 
   protected:
     /**
-     * Stored value to compare the
-     * level with.
+     * Stored value to compare the level with.
      */
     const unsigned int level;
   };
@@ -149,15 +127,10 @@ namespace IteratorFilters
 
 
   /**
-   * Filter for iterators that
-   * evaluates to true if either the
-   * iterator is past the end or the
-   * subdomain id of the object
-   * pointed to is equal to a value
-   * given to the constructor,
-   * assuming that the iterator
-   * allows querying for a subdomain
-   * id).
+   * Filter for iterators that evaluates to true if either the
+   * iterator is past the end or the subdomain id of the object
+   * pointed to is equal to a value given to the constructor, assuming
+   * that the iterator allows querying for a subdomain id).
    *
    * @ingroup Iterators
    */
@@ -165,28 +138,22 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Constructor. Store the
-     * subdomain which iterators
-     * shall have to be evaluated
-     * to true.
+     * Constructor. Store the subdomain which iterators shall have to
+     * be evaluated to true.
      */
     SubdomainEqualTo (const types::subdomain_id subdomain_id);
 
     /**
-     * Evaluation operator. Returns
-     * true if either the subdomain
-     * of the object pointed to is
-     * equal to the stored value or
-     * the iterator is past the
-     * end.
+     * Evaluation operator. Returns true if either the subdomain of
+     * the object pointed to is equal to the stored value or the
+     * iterator is past the end.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
 
   protected:
     /**
-     * Stored value to compare the
-     * subdomain with.
+     * Stored value to compare the subdomain with.
      */
     const types::subdomain_id subdomain_id;
   };
@@ -194,14 +161,12 @@ namespace IteratorFilters
 
 
   /**
-   * Filter for iterators that evaluates to
-   * true if a cell is owned by the current
-   * processor, i.e., if it is a
-   * @ref GlossLocallyOwnedCell "locally owned cell".
+   * Filter for iterators that evaluates to true if a cell is owned by
+   * the current processor, i.e., if it is a @ref
+   * GlossLocallyOwnedCell "locally owned cell".
    *
-   * This class is used in step-32, in
-   * connection with the methods of the @ref
-   * distributed module.
+   * This class is used in step-32, in connection with the methods of
+   * the @ref distributed module.
    *
    * @ingroup Iterators
    */
@@ -209,8 +174,26 @@ namespace IteratorFilters
   {
   public:
     /**
-     * Evaluation operator. Returns true if
-     * the cell is locally owned.
+     * Evaluation operator. Returns true if the cell is locally owned.
+     */
+    template <class Iterator>
+    bool operator () (const Iterator &i) const;
+  };
+
+
+
+  /**
+   * Filter for iterators that evaluates to true if th level subdomain
+   * id of a cell is equal to the current processor id.
+   *
+   * @ingroup Iterators
+   */
+  class LocallyOwnedLevelCell
+  {
+  public:
+    /**
+     * Evaluation operator. Returns true if the level subdomain id of
+     * the cell is equal to the current processor id.
      */
     template <class Iterator>
     bool operator () (const Iterator &i) const;
@@ -1142,6 +1125,17 @@ namespace IteratorFilters
   {
     return (i->is_locally_owned());
   }
+
+
+// ---------------- IteratorFilters::LocallyOwnedLevelCell ---------
+
+  template <class Iterator>
+  inline
+  bool
+  LocallyOwnedLevelCell::operator () (const Iterator &i) const
+  {
+    return (i->is_locally_owned_on_level());
+  }
 }
 
 
index 750e4738ff80a1583759ad296bbda962a83186b5..fea050957357e11be55ee2855636ba73a42ef9bf 100644 (file)
@@ -2890,6 +2890,12 @@ public:
    */
   bool is_locally_owned () const;
 
+  /**
+   * Return true if either the Triangulation is not distributed or if
+   * level_subdomain_id() is equal to the id of the current processor.
+   */
+  bool is_locally_owned_on_level () const;
+
   /**
    * Return whether this cell
    * exists in the global mesh but
index b712d5f072fec7a95776024873a4142b801983f4..970817688088049b843e74d05a0c589f2583145b 100644 (file)
@@ -2925,6 +2925,25 @@ CellAccessor<dim,spacedim>::is_locally_owned () const
 }
 
 
+template <int dim, int spacedim>
+inline
+bool
+CellAccessor<dim,spacedim>::is_locally_owned_on_level () const
+{
+#ifndef DEAL_II_WITH_P4EST
+  return true;
+#else
+  const parallel::distributed::Triangulation<dim,spacedim> *pdt
+    = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(this->tria);
+
+  if (pdt == 0)
+    return true;
+  else
+    return (this->level_subdomain_id() == pdt->locally_owned_subdomain());
+#endif
+}
+
+
 template <int dim, int spacedim>
 inline
 bool

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.