New example

src/03-function-integration.cc

+// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
+// vi: set et ts=4 sw=2 sts=2:
+
+#ifdef HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+
+
+// included standard library headers
+#include <iostream>
+#include <array>
+
+// included dune-common headers
+#include <dune/common/parallel/mpihelper.hh>
+#include <dune/common/exceptions.hh>
+#include <dune/common/fvector.hh>
+
+// included dune-geometry headers
+#include <dune/geometry/quadraturerules.hh>
+
+// included dune-grid headers
+#include <dune/grid/yaspgrid.hh>
+#include <dune/grid/io/file/vtk/vtkwriter.hh>
+
+
+
+
+
+template<class GridView, class F>
+double integrate(const GridView& gridView, const F& f, std::size_t order)
+{
+  static const int dim = GridView::dimension;
+
+  double value = 0.0;
+  for(const auto& e: Dune::elements(gridView))
+  {
+    auto geometryType = e.type();
+
+    const auto& geometry = e.geometry();
+
+    const auto& quad = Dune::template QuadratureRules<double, dim>::rule(geometryType, order);
+
+    for (size_t i=0; i < quad.size(); ++i)
+    {
+      const auto& quadPos = quad[i].position();
+      const auto& weight = quad[i].weight();
+
+      auto globalQuadPos = geometry.global(quadPos);
+
+      // get integration factor
+      const double integrationElement = geometry.integrationElement(quadPos);
+
+      value += f(globalQuadPos) * weight * integrationElement;
+    }
+  }
+  return value;
+}
+
+
+int main(int argc, char** argv)
+{
+  try{
+    // Maybe initialize MPI
+    Dune::MPIHelper& helper = Dune::MPIHelper::instance(argc, argv);
+
+    // Print process rank
+    if(Dune::MPIHelper::isFake)
+      std::cout<< "This is a sequential program." << std::endl;
+    else
+      std::cout<<"I am rank "<<helper.rank()<<" of "<<helper.size()
+        <<" processes!"<<std::endl;
+
+    // fix grid dimension
+    const int dim = 2;
+
+    // use the YaspGrid implementation
+    // YaspGrid = "Yet another structured parallel grid"
+    using Grid = Dune::YaspGrid<dim>;
+
+    // extension of the domain [0,l]
+    Dune::FieldVector<double,dim> l(1);
+
+    // start with 10x10x... elements
+    std::array<int,dim> E;
+    for(auto& e : E)
+      e = 2;
+
+    // create grid
+    Grid grid(l, E);
+
+    // refine the grid globally two times
+    grid.globalRefine(7);
+
+    // create function to integrate
+    auto f = [](const auto& x) {
+      return x.infinity_norm();
+    };
+
+    // get a leaf GridView
+    auto gridView = grid.leafGridView();
+
+    std::cout << integrate(gridView, f, 7) << std::endl;
+
+    using GridView = decltype(gridView);
+    Dune::VTKWriter<GridView> vtkWriter(gridView);
+    vtkWriter.write("01-refined-grid");
+
+    return 0;
+  }
+  catch (Dune::Exception &e){
+    std::cerr << "Dune reported error: " << e << std::endl;
+  }
+  catch (...){
+    std::cerr << "Unknown exception thrown!" << std::endl;
+  }
+}

src/CMakeLists.txt

@@ -6,3 +6,6 @@ target_link_dune_default_libraries("01-basic-application")
 
 add_executable("02-basic-grid-interface" 02-basic-grid-interface.cc)
 target_link_dune_default_libraries("02-basic-grid-interface")
+
+add_executable("03-function-integration" 03-function-integration.cc)
+target_link_dune_default_libraries("03-function-integration")