diff --git a/python/dune/perftool/options.py b/python/dune/perftool/options.py
index 1aeb7f272d52b7c23ac361b9df132fa1381d36d6..e4b51c3a15f191f86691b1b0816e1a3057f23073 100644
--- a/python/dune/perftool/options.py
+++ b/python/dune/perftool/options.py
@@ -14,9 +14,9 @@ def get_form_compiler_arguments():
parser.add_argument('--version', action='version', version='%(prog)s 0.1')
parser.add_argument("--numerical-jacobian", action="store_true", help="use numerical jacobians (only makes sense, if uflpdelab for some reason fails to generate analytic jacobians)")
parser.add_argument(
- "--dirichlet-is-exact-solution",
- action="store_true",
- help="test if solution is correct by providing the exact solution as dirichlet boundary condition")
+ "--exact-solution-expression",
+ nargs = 2,
+ help="requires two arguments: the name of the expression corresponding to solution in the ufl file and the maximal allowed error (between difference of the dof vectors of num. solution and interpolated exact solution")
parser.add_argument("--interactive", action="store_true", help="whether the optimization process should be guided interactively (also useful for debugging)")
# These are the options that I deemed necessary in uflpdelab
diff --git a/python/dune/perftool/pdelab/driver.py b/python/dune/perftool/pdelab/driver.py
index a20571f850473d9bf20f012cad1b6be4a7e5c7c2..f7a27ff09e204568ade79e48b19ac0907ce5d0d0 100644
--- a/python/dune/perftool/pdelab/driver.py
+++ b/python/dune/perftool/pdelab/driver.py
@@ -706,6 +706,7 @@ def name_boundary_lambda(boundary, name):
def define_boundary_function(boundary, name):
gv = name_leafview()
lambdaname = name_boundary_lambda(boundary, name)
+ include_file('dune/pdelab/function/callableadapter.hh', filetag='driver')
return "auto {} = Dune::PDELab::makeGridFunctionFromCallable({}, {});".format(name,
gv,
lambdaname,
@@ -744,6 +745,20 @@ def name_boundary_function(expr):
return name
+@symbol
+def name_solution_function(expr):
+ # Get expression representing solution
+ from dune.perftool.generation import get_global_context_value
+ solution_expression_name = get_option("exact_solution_expression")[0]
+ solution_expression = get_global_context_value("data").object_by_name[solution_expression_name]
+
+ # Create function just like it is done for boundary_function
+ name = "solution_function"
+ define_boundary_function(solution_expression, name)
+
+ return name
+
+
@preamble
def interpolate_vector(name):
define_vector(name)
@@ -756,6 +771,18 @@ def interpolate_vector(name):
)
+@preamble
+def interpolate_solution_expression(name):
+ define_vector(name)
+ element = _form.coefficients()[0].ufl_element()
+ sol = name_solution_function(element)
+ gfs = name_gfs(element)
+ return "Dune::PDELab::interpolate({}, {}, {});".format(sol,
+ gfs,
+ name,
+ )
+
+
def maybe_interpolate_vector(name):
element = _form.coefficients()[0].ufl_element()
if has_constraints(element):
@@ -769,11 +796,13 @@ def name_vector():
maybe_interpolate_vector("x")
return "x"
+
@symbol
-def name_solution_from_dirichlet_vector():
- interpolate_vector("solution")
+def name_vector_from_solution_expression():
+ interpolate_solution_expression("solution")
return "solution"
+
@preamble
def typedef_linearsolver(name):
include_file("dune/pdelab/backend/istl.hh", filetag="driver")
@@ -895,9 +924,9 @@ def name_vtkfile():
@preamble
-def exact_solution_from_dirichlet():
+def exact_solution():
v = name_vector()
- solution = name_solution_from_dirichlet_vector();
+ solution = name_vector_from_solution_expression()
return ["using Dune::PDELab::Backend::native;",
"double maxerror(0.0);",
@@ -905,13 +934,12 @@ def exact_solution_from_dirichlet():
" if (std::abs(native({})[i]-native({})[i]) > maxerror)".format(v, solution),
" maxerror = std::abs(native({})[i]-native({})[i]);".format(v, solution),
"",
- "std::cout << \"Maxerror: \" << maxerror << std::endl;",
+ "std::cout << std::endl << \"Maxerror: \" << maxerror << std::endl;",
"",
- "if (maxerror>1e-14)",
+ "if (maxerror>{})".format(get_option("exact_solution_expression")[1]),
" throw;"]
-
@preamble
def print_residual():
ini = name_initree()
@@ -991,8 +1019,9 @@ def vtkoutput():
dune_solve()
print_residual()
print_matrix()
- if get_option("dirichlet_is_exact_solution"):
- exact_solution_from_dirichlet()
+
+ if get_option("exact_solution_expression"):
+ exact_solution()
return ["Dune::PDELab::addSolutionToVTKWriter({}, {}, {}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vec, predicate),
"{}.write({}, Dune::VTK::ascii);".format(vtkwriter, vtkfile)]
diff --git a/test/nonlinear/CMakeLists.txt b/test/nonlinear/CMakeLists.txt
index ff8e162462606d16ab8d5952b1f70675f48a7b1b..5ebedaa84ecbcecfa73d5c17ed4d24597aa7038c 100644
--- a/test/nonlinear/CMakeLists.txt
+++ b/test/nonlinear/CMakeLists.txt
@@ -9,6 +9,15 @@ dune_add_system_test(TARGET nonlinear_symdiff
SCRIPT dune_vtkcompare.py
)
+add_generated_executable(UFLFILE nonlinear_dg.ufl
+ TARGET nonlinear_dg
+ FORM_COMPILER_ARGS --exact-solution-expression g 1e-2
+ )
+
+dune_add_system_test(TARGET nonlinear_dg
+ INIFILE nonlinear_dg.mini
+ )
+
# Add the reference code with the PDELab localoperator that produced
# the reference vtk file
add_executable(nonlinear_ref reference_main.cc)
diff --git a/test/nonlinear/nonlinear_dg.ufl b/test/nonlinear/nonlinear_dg.ufl
new file mode 100644
index 0000000000000000000000000000000000000000..bf606f9cdb7b3726b457688db4921c4aacfa99ab
--- /dev/null
+++ b/test/nonlinear/nonlinear_dg.ufl
@@ -0,0 +1,26 @@
+f = Expression("return -2.0*x.size();")
+g = Expression("return x.two_norm2();", on_intersection=True)
+
+cell = triangle
+V = FiniteElement("DG", cell, 1)
+
+u = TrialFunction(V)
+v = TestFunction(V)
+
+n = FacetNormal(cell)('+')
+
+gamma = 1.0
+theta = 1.0
+
+r = inner(grad(u), grad(v))*dx \
+ + inner(n, avg(grad(u)))*jump(v)*dS \
+ + gamma*jump(u)*jump(v)*dS \
+ + theta*jump(u)*inner(avg(grad(v)), n)*dS \
+ - inner(n, grad(u))*v*ds \
+ + gamma*u*v*ds \
+ - theta*u*inner(grad(v), n)*ds \
+ - f*v*dx \
+ + theta*g*inner(grad(v), n)*ds \
+ - gamma*g*v*ds
+
+forms = [r]