From 07838a105471784f2f517248029a47d204a4987f Mon Sep 17 00:00:00 2001
From: Dominic Kempf <dominic.kempf@iwr.uni-heidelberg.de>
Date: Wed, 18 Jan 2017 13:12:14 +0100
Subject: [PATCH] Add testing for full tensor calculations

---
 .../loopy/transformations/collect_rotate.py   |  2 +-
 .../loopy/transformations/vectorview.py       |  1 -
 test/sumfact/poisson/CMakeLists.txt           |  6 ++++
 test/sumfact/poisson/poisson_dg_tensor.mini   | 22 +++++++++++++
 test/sumfact/poisson/poisson_dg_tensor.ufl    | 32 +++++++++++++++++++
 5 files changed, 61 insertions(+), 2 deletions(-)
 create mode 100644 test/sumfact/poisson/poisson_dg_tensor.mini
 create mode 100644 test/sumfact/poisson/poisson_dg_tensor.ufl

diff --git a/python/dune/perftool/loopy/transformations/collect_rotate.py b/python/dune/perftool/loopy/transformations/collect_rotate.py
index d633843a..78070a10 100644
--- a/python/dune/perftool/loopy/transformations/collect_rotate.py
+++ b/python/dune/perftool/loopy/transformations/collect_rotate.py
@@ -188,7 +188,7 @@ def collect_vector_data_rotate(knl):
             else:
                 # Add a vector view to this quantity
                 expr, = quantities[quantity]
-                knl = add_vector_view(knl, quantity)
+                knl = add_vector_view(knl, quantity, flatview=True)
                 replacemap_vec[expr] = prim.Subscript(prim.Variable(get_vector_view_name(quantity)),
                                                       (prim.Variable("vec_index"), prim.Variable(new_iname)),
                                                       )
diff --git a/python/dune/perftool/loopy/transformations/vectorview.py b/python/dune/perftool/loopy/transformations/vectorview.py
index 1450b26e..7fc3ff8f 100644
--- a/python/dune/perftool/loopy/transformations/vectorview.py
+++ b/python/dune/perftool/loopy/transformations/vectorview.py
@@ -59,7 +59,6 @@ def add_vector_view(knl, tmpname, pad_to=None, flatview=False):
         shape = (size, vecsize)
         dim_tags = "c,vec"
     else:
-        assert(temp.shape[-1] == vecsize)
         shape = temp.shape
         # This works around a loopy weirdness (which might as well be a bug)
         # TODO: investigate this!
diff --git a/test/sumfact/poisson/CMakeLists.txt b/test/sumfact/poisson/CMakeLists.txt
index aab2e764..062fc23d 100644
--- a/test/sumfact/poisson/CMakeLists.txt
+++ b/test/sumfact/poisson/CMakeLists.txt
@@ -38,3 +38,9 @@ dune_add_formcompiler_system_test(UFLFILE poisson_dg_2d.ufl
                                   BASENAME sumfact_poisson_fastdg
                                   INIFILE poisson_fastdg.mini
                                   )
+
+# Poisson DG with a full tensor
+dune_add_formcompiler_system_test(UFLFILE poisson_dg_tensor.ufl
+                                  BASENAME sumfact_poisson_dg_tensor
+                                  INIFILE poisson_dg_tensor.mini
+                                  )
diff --git a/test/sumfact/poisson/poisson_dg_tensor.mini b/test/sumfact/poisson/poisson_dg_tensor.mini
new file mode 100644
index 00000000..0ae90c37
--- /dev/null
+++ b/test/sumfact/poisson/poisson_dg_tensor.mini
@@ -0,0 +1,22 @@
+__name = sumfact_poisson_dg_tensor_{__exec_suffix}
+__exec_suffix = {quadvec_suffix}_{gradvec_suffix}
+
+quadvec_suffix = quadvec, nonquadvec | expand quad
+gradvec_suffix = gradvec, nongradvec | expand grad
+
+cells = 8 8 8
+extension = 1. 1. 1.
+
+[wrapper.vtkcompare]
+name = {__name}
+extension = vtu
+
+[formcompiler]
+sumfact = 1
+exact_solution_expression = g
+compare_l2errorsquared = 1e-4
+vectorize_quad = 1, 0 | expand quad
+vectorize_grads = 1, 0 | expand grad
+
+[formcompiler.ufl_variants]
+degree = 2
diff --git a/test/sumfact/poisson/poisson_dg_tensor.ufl b/test/sumfact/poisson/poisson_dg_tensor.ufl
new file mode 100644
index 00000000..ceecebc5
--- /dev/null
+++ b/test/sumfact/poisson/poisson_dg_tensor.ufl
@@ -0,0 +1,32 @@
+cell = hexahedron
+
+x = SpatialCoordinate(cell)
+
+I = Identity(3)
+A = as_matrix([[x[i]*x[j] + I[i,j] for j in range(3)] for i in range(3)])
+g = x[0]**2 + x[1]**2 + x[2]**2
+c = 10.
+f = -6.
+
+V = FiniteElement("DG", cell, 1)
+
+u = TrialFunction(V)
+v = TestFunction(V)
+
+n = FacetNormal(cell)('+')
+
+gamma = 1.0
+theta = -1.0
+
+r = (inner(A*grad(u), grad(v)) + c*u*v)*dx \
+  + inner(n, A*avg(grad(u)))*jump(v)*dS \
+  + gamma*jump(u)*jump(v)*dS \
+  + theta*jump(u)*inner(A*avg(grad(v)), n)*dS \
+  - inner(n, A*grad(u))*v*ds \
+  + gamma*u*v*ds \
+  - theta*u*inner(A*grad(v), n)*ds \
+  - f*v*dx \
+  + theta*g*inner(A*grad(v), n)*ds \
+  - gamma*g*v*ds
+
+forms = [r]
-- 
GitLab