From 2888162c748acb2d5d8a12f52897f6d17a07b446 Mon Sep 17 00:00:00 2001 From: Pablo Brubeck Date: Fri, 25 Oct 2024 12:12:44 +0100 Subject: [PATCH] add continuity test --- tests/regression/test_tensor_elements.py | 39 ++++++++++++++++++++++++ 1 file changed, 39 insertions(+) create mode 100644 tests/regression/test_tensor_elements.py diff --git a/tests/regression/test_tensor_elements.py b/tests/regression/test_tensor_elements.py new file mode 100644 index 0000000000..a0bb5352b6 --- /dev/null +++ b/tests/regression/test_tensor_elements.py @@ -0,0 +1,39 @@ +from firedrake import * +import pytest + + +rg = RandomGenerator(PCG64(seed=1)) + + +@pytest.fixture(params=("square", "cube")) +def mesh(request): + if request.param == "square": + return UnitSquareMesh(2, 2) + elif request.param == "cube": + return UnitCubeMesh(2, 2, 2) + else: + raise ValueError("Unrecognized mesh type") + + +@pytest.mark.parametrize("degree", (0, 1, 2)) +@pytest.mark.parametrize("family", ("Regge", "HHJ")) +def test_tensor_continuity(mesh, family, degree): + V = FunctionSpace(mesh, family, degree) + u = rg.beta(V, 1.0, 2.0) + + dim = mesh.topological_dimension() + n = FacetNormal(mesh) + + space = V.ufl_element().sobolev_space + if space == HDivDiv: + utrace = dot(n, dot(u, n)) + elif space == HEin: + if dim == 2: + t = dot(as_matrix([[0, -1], [1, 0]]), n) + utrace = dot(t, dot(u, t)) + else: + t = as_matrix([[0, n[2], -n[1]], [-n[2], 0, n[0]], [n[1], -n[0], 0]]) + utrace = dot(t, dot(u, t.T)) + + ujump = utrace('+') - utrace('-') + assert assemble(inner(ujump, ujump)*dS) < 1E-10