From ff7b045fc1bc9321b8fa60b8bd088f854f2b644d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Ren=C3=A9=20He=C3=9F?= <rene.hess@iwr.uni-heidelberg.de>
Date: Thu, 25 Aug 2016 12:08:35 +0200
Subject: [PATCH] Compare with exact solution
---
python/dune/perftool/pdelab/driver.py | 11 +++++++++++
test/heatequation/CMakeLists.txt | 10 +++-------
test/heatequation/heatequation.ufl | 3 ++-
test/heatequation/heatequation_explicit_ref.cc | 8 ++++++--
test/heatequation/problem.hh | 6 +++---
5 files changed, 25 insertions(+), 13 deletions(-)
diff --git a/python/dune/perftool/pdelab/driver.py b/python/dune/perftool/pdelab/driver.py
index a3d7e88a..6bed1ab6 100644
--- a/python/dune/perftool/pdelab/driver.py
+++ b/python/dune/perftool/pdelab/driver.py
@@ -1387,6 +1387,17 @@ def time_loop():
elif not linear and not matrix_free:
assert False
+ print_residual()
+ print_matrix()
+ if get_option("exact_solution_expression"):
+ from ufl import MixedElement, VectorElement, TensorElement
+ if isinstance(_driver_data['form'].coefficients()[0].ufl_element(), (MixedElement, VectorElement, TensorElement)):
+ raise NotImplementedError("Comparing to exact solution is olny implemented for scalar elements.")
+
+ if get_option("compare_dofs"):
+ compare_dofs()
+ if get_option("compare_l2errorsquared"):
+ compare_L2_squared()
def generate_driver(formdatas, data):
# The driver module uses a global dictionary for storing necessary data
diff --git a/test/heatequation/CMakeLists.txt b/test/heatequation/CMakeLists.txt
index 3ac71527..543dd2d0 100644
--- a/test/heatequation/CMakeLists.txt
+++ b/test/heatequation/CMakeLists.txt
@@ -1,6 +1,6 @@
add_generated_executable(UFLFILE heatequation.ufl
TARGET heatequation_implicit
- FORM_COMPILER_ARGS
+ FORM_COMPILER_ARGS --exact-solution-expression g --compare-l2errorsquared 1e-7
)
@@ -10,7 +10,7 @@ dune_add_system_test(TARGET heatequation_implicit
add_generated_executable(UFLFILE heatequation.ufl
TARGET heatequation_explicit
- FORM_COMPILER_ARGS --explicit-time-stepping
+ FORM_COMPILER_ARGS --explicit-time-stepping --exact-solution-expression g --compare-l2errorsquared 1e-7
)
@@ -18,14 +18,10 @@ dune_add_system_test(TARGET heatequation_explicit
INIFILE heatequation_explicit.mini
)
+# Hand written reference solutions
add_executable(heatequation_implicit_ref heatequation_implicit_ref.cc)
dune_symlink_to_source_files(FILES heatequation_implicit_ref.ini)
set_target_properties(heatequation_implicit_ref PROPERTIES EXCLUDE_FROM_ALL 1)
add_executable(heatequation_explicit_ref heatequation_explicit_ref.cc)
set_target_properties(heatequation_explicit_ref PROPERTIES EXCLUDE_FROM_ALL 1)
-
-# # Add the reference code with the PDELab localoperator that produced
-# # the reference vtk file
-# add_executable(heatequation_ref heatequation_main.cc)
-# set_target_properties(heatequation_ref PROPERTIES EXCLUDE_FROM_ALL 1)
diff --git a/test/heatequation/heatequation.ufl b/test/heatequation/heatequation.ufl
index 5495df8a..b6d3f65c 100644
--- a/test/heatequation/heatequation.ufl
+++ b/test/heatequation/heatequation.ufl
@@ -1,3 +1,4 @@
+f = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());")
g = Expression("Dune::FieldVector<double,2> c(0.5); c-= x; return std::exp(-1.*c.two_norm2());")
V = FiniteElement("CG", "triangle", 1, dirichlet_expression=g)
@@ -5,6 +6,6 @@ u = TrialFunction(V)
v = TestFunction(V)
mass = (u*v)*dx
-poisson = (inner(grad(u), grad(v)))*dx
+poisson = (inner(grad(u), grad(v)) - f*v)*dx
forms = [mass, poisson]
\ No newline at end of file
diff --git a/test/heatequation/heatequation_explicit_ref.cc b/test/heatequation/heatequation_explicit_ref.cc
index 28470723..9a07f049 100644
--- a/test/heatequation/heatequation_explicit_ref.cc
+++ b/test/heatequation/heatequation_explicit_ref.cc
@@ -53,9 +53,13 @@ public:
//! source term
typename Traits::RangeFieldType
- f (const typename Traits::ElementType& e, const typename Traits::DomainType& x) const
+ f (const typename Traits::ElementType& e, const typename Traits::DomainType& xlocal) const
{
- return 0.0;
+ typename Traits::DomainType x = e.geometry().global(xlocal);
+ Dune::FieldVector<double,2> c(0.5);
+ c-= x;
+ return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());
+ // return 0.0;
}
//! boundary condition type function
diff --git a/test/heatequation/problem.hh b/test/heatequation/problem.hh
index 37c93f23..8c1a13c6 100644
--- a/test/heatequation/problem.hh
+++ b/test/heatequation/problem.hh
@@ -29,9 +29,9 @@ public:
template<typename E, typename X>
Number f (const E& e, const X& local) const
{
- return 0.0;
- // auto x = e.geometry().global(local);
- // Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());
+ // return 0.0;
+ auto x = e.geometry().global(local);
+ Dune::FieldVector<double,2> c(0.5); c-= x; return 4.*(1.-c.two_norm2())*std::exp(-1.*c.two_norm2());
}
//! boundary condition type function (true = Dirichlet)
--
GitLab