From 21a7a1a5e0d1b0f06fe933526246901868634f20 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Maximilian=20H=C3=B6rl?=
 <maximilian.hoerl@mathematik.uni-stuttgart.de>
Date: Fri, 20 Mar 2020 19:23:44 +0100
Subject: [PATCH] use different problem for inhomogeneousbulkproblem

---
 .../mmdg/problems/inhomogeneousbulkproblem.hh | 28 +++++++++++--------
 1 file changed, 17 insertions(+), 11 deletions(-)

diff --git a/dune/mmdg/problems/inhomogeneousbulkproblem.hh b/dune/mmdg/problems/inhomogeneousbulkproblem.hh
index afa64d2..601d0bb 100644
--- a/dune/mmdg/problems/inhomogeneousbulkproblem.hh
+++ b/dune/mmdg/problems/inhomogeneousbulkproblem.hh
@@ -12,12 +12,13 @@ class InhomogeneousBulkProblem : public DGProblem<Coordinate, Scalar>
     using Matrix = typename Base::Matrix;
 
     //the exact solution at position pos
-    //1d: p(x)     = x,
-    //2d: p(x,y)   = x*y,
-    //3d: p(x,y,z) = x*y*z
-    virtual Scalar exactSolution (const Coordinate& pos) const
+    //1d: p(x)     = x^2,
+    //2d: p(x,y)   = x^2 + y^2,
+    //3d: p(x,y,z) = x^2 + y^2 + z^2
+    Scalar exactSolution (const Coordinate& pos) const
     {
-      return pos * Coordinate(1.0);
+      //return pos * Coordinate(1.0);
+      return pos * pos;
     }
 
     //exact solution is implemented
@@ -27,14 +28,14 @@ class InhomogeneousBulkProblem : public DGProblem<Coordinate, Scalar>
     };
 
     //source term at position pos
-    virtual Scalar q (const Coordinate& pos) const
+    Scalar q (const Coordinate& pos) const
     {
-      return -pos * Coordinate(2.0);
+      return -2.0 * (dim + 3.0 * (pos * pos));
     }
 
     //permeability tensor at position pos
-    virtual Matrix K (const Coordinate& pos) const
-    {
+    Matrix K (const Coordinate& pos) const
+    { //K = diag(1 + x_i^2 : i=1,...,dim)
       Matrix permeability(0.0);
 
       for (int i = 0; i < dim; i++)
@@ -49,14 +50,19 @@ class InhomogeneousBulkProblem : public DGProblem<Coordinate, Scalar>
     //over x * q(x)
     int quadratureOrder_q () const
     {
-      return 2;
+      return 2 * dim + 1;
     }
 
     //returns the recommended quadrature order to compute an integral
     //over x * boundary(x) and K(x) * boundary(x)
     int quadratureOrderBoundary () const
     {
-      return 1 + dim;
+      if constexpr (dim == 1)
+      {
+        return 2;
+      }
+
+      return 4;
     }
 
     //returns the recommended quadrature order to compute an integral
-- 
GitLab