[cig-commits] r1396 - in trunk/aspect: include/aspect/mesh_refinement source/mesh_refinement

[bangerth at dealii.org]

Author: bangerth
Date: 2012-11-29 08:08:50 -0700 (Thu, 29 Nov 2012)
New Revision: 1396

Added:
   trunk/aspect/include/aspect/mesh_refinement/composition.h
   trunk/aspect/source/mesh_refinement/composition.cc
Modified:
   trunk/aspect/include/aspect/mesh_refinement/interface.h
Log:
Implement the first mesh refinement criterion.

Copied: trunk/aspect/include/aspect/mesh_refinement/composition.h (from rev 1394, trunk/aspect/include/aspect/gravity_model/radial.h)
===================================================================
--- trunk/aspect/include/aspect/mesh_refinement/composition.h	                        (rev 0)
+++ trunk/aspect/include/aspect/mesh_refinement/composition.h	2012-11-29 15:08:50 UTC (rev 1396)
@@ -0,0 +1,61 @@
+/*
+  Copyright (C) 2011, 2012 by the authors of the ASPECT code.
+
+  This file is part of ASPECT.
+
+  ASPECT is free software; you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation; either version 2, or (at your option)
+  any later version.
+
+  ASPECT is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with ASPECT; see the file doc/COPYING.  If not see
+  <http://www.gnu.org/licenses/>.
+*/
+/*  $Id$  */
+
+
+#ifndef __aspect__mesh_refinement_composition_h
+#define __aspect__mesh_refinement_composition_h
+
+#include <aspect/mesh_refinement/interface.h>
+#include <aspect/simulator.h>
+
+namespace aspect
+{
+  namespace MeshRefinement
+  {
+
+    /**
+     * A class that implements a mesh refinement criterion based on
+     * the compositional fields (if available).
+     *
+     * @ingroup MeshRefinement
+     */
+    template <int dim>
+    class Composition : public Interface<dim>,
+                        public SimulatorAccess<dim>
+    {
+      public:
+      /**
+       * Execute this mesh refinement criterion.
+       *
+       * @param[out] error_indicators A vector that for every active
+       * cell of the current mesh
+       * (which may be a partition of a distributed mesh) provides an error
+       * indicator. This vector will already have the correct size when the
+       * function is called.
+       */
+      virtual
+      void
+      execute (Vector<float> &error_indicators) const;
+    };
+  }
+}
+
+#endif

Modified: trunk/aspect/include/aspect/mesh_refinement/interface.h
===================================================================
--- trunk/aspect/include/aspect/mesh_refinement/interface.h	2012-11-29 14:21:19 UTC (rev 1395)
+++ trunk/aspect/include/aspect/mesh_refinement/interface.h	2012-11-29 15:08:50 UTC (rev 1396)
@@ -82,7 +82,7 @@
          */
         virtual
         void
-        execute (Vector<float> &error_indicators) = 0;
+        execute (Vector<float> &error_indicators) const = 0;
 
         /**
          * Declare the parameters this class takes through input files.

Copied: trunk/aspect/source/mesh_refinement/composition.cc (from rev 1394, trunk/aspect/source/gravity_model/radial.cc)
===================================================================
--- trunk/aspect/source/mesh_refinement/composition.cc	                        (rev 0)
+++ trunk/aspect/source/mesh_refinement/composition.cc	2012-11-29 15:08:50 UTC (rev 1396)
@@ -0,0 +1,76 @@
+/*
+  Copyright (C) 2011, 2012 by the authors of the ASPECT code.
+
+  This file is part of ASPECT.
+
+  ASPECT is free software; you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation; either version 2, or (at your option)
+  any later version.
+
+  ASPECT is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with ASPECT; see the file doc/COPYING.  If not see
+  <http://www.gnu.org/licenses/>.
+*/
+/*  $Id$  */
+
+
+#include <aspect/mesh_refinement/composition.h>
+
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/numerics/error_estimator.h>
+
+namespace aspect
+{
+  namespace MeshRefinement
+  {
+    template <int dim>
+    void
+    Composition<dim>::execute(Vector<float> &indicators) const
+    {
+      AssertThrow (this->n_compositional_fields() >= 1,
+                   ExcMessage ("This refinement criterion can not be used when no "
+                       "compositional fields are active!"));
+      indicators = 0;
+
+      for (unsigned int c=0; c<this->n_compositional_fields(); ++c)
+        {
+          Vector<float> this_indicator (indicators.size());
+
+          std::vector<bool> composition_component (dim+2+this->n_compositional_fields(), false);
+          composition_component[dim+2+c] = true;
+          KellyErrorEstimator<dim>::estimate (this->get_dof_handler(),
+//TODO: Replace the 2 by something reasonable, adjusted to the polynomial degree
+                                              QGauss<dim-1>(2),
+                                              typename FunctionMap<dim>::type(),
+                                              this->get_solution(),
+                                              this_indicator,
+                                              composition_component,
+                                              0,
+                                              0,
+                                              this->get_triangulation().locally_owned_subdomain());
+          indicators += this_indicator;
+        }
+    }
+  }
+}
+
+// explicit instantiations
+namespace aspect
+{
+  namespace MeshRefinement
+  {
+    ASPECT_REGISTER_MESH_REFINEMENT_CRITERION(Composition,
+                                  "composition",
+                                  "A mesh refinement criterion that takes computes "
+                                  "refinement indicators from the compositional fields. "
+                                  "If there is more than one compositional field, then "
+                                  "it simply takes the sum of the indicators computed "
+                                  "from each of the compositional field.")
+  }
+}