From 5fc9a781b0524f012ae9f24264cd695c315cf9dc Mon Sep 17 00:00:00 2001 From: Pablo Brubeck Date: Mon, 28 Oct 2024 18:25:02 +0000 Subject: [PATCH] test Hu-Zhang --- tests/regression/test_stress_els.py | 31 ++++++++++++++++++----------- 1 file changed, 19 insertions(+), 12 deletions(-) diff --git a/tests/regression/test_stress_els.py b/tests/regression/test_stress_els.py index 16ca97bc58..3b524652f8 100755 --- a/tests/regression/test_stress_els.py +++ b/tests/regression/test_stress_els.py @@ -7,7 +7,7 @@ convergence_orders = lambda x: -np.log2(relative_magnitudes(x)) -@pytest.fixture(scope='module', params=["conforming", "nonconforming", "macro"]) +@pytest.fixture(params=["conforming", "nonconforming", "macro", "high-order"]) def stress_element(request): if request.param == "conforming": return FiniteElement("AWc", triangle, 3) @@ -15,11 +15,13 @@ def stress_element(request): return FiniteElement("AWnc", triangle, 2) elif request.param == "macro": return FiniteElement("JM", triangle, 1) + elif request.param == "high-order": + return FiniteElement("HZ", triangle, 3) else: raise ValueError("Unknown family") -@pytest.fixture(scope='module') +@pytest.fixture def mesh_hierarchy(request): N_base = 2 mesh = UnitSquareMesh(N_base, N_base) @@ -64,7 +66,10 @@ def epsilon(u): l2_u = [] l2_sigma = [] l2_div_sigma = [] - displacement_element = VectorElement("DG", cell=mesh.ufl_cell(), degree=1, variant="integral") + sdegree = stress_element.degree() + sfamily = stress_element.family() + vdegree = sdegree-1 if sfamily == "Hu-Zhang" else 1 + displacement_element = VectorElement("DG", cell=mesh.ufl_cell(), degree=vdegree, variant="integral") element = MixedElement([stress_element, displacement_element]) for msh in mesh_hierarchy[1:]: x, y = SpatialCoordinate(msh) @@ -126,21 +131,23 @@ def epsilon(u): l2_sigma.append(error_sigma) l2_div_sigma.append(error_div_sigma) + sdegree = V[0].ufl_element().degree() + vdegree = V[1].ufl_element().degree() if stress_element.family().startswith("Conforming"): - assert min(convergence_orders(l2_u)) > 1.9 - assert min(convergence_orders(l2_sigma)) > 2.9 - assert min(convergence_orders(l2_div_sigma)) > 1.9 + rates = (vdegree+1, sdegree, sdegree-1) elif stress_element.family().startswith("Nonconforming"): - assert min(convergence_orders(l2_u)) > 1.9 - assert min(convergence_orders(l2_sigma)) > 1 - assert min(convergence_orders(l2_div_sigma)) > 1.9 + rates = (vdegree+1, sdegree-1, sdegree) elif stress_element.family().startswith("Johnson-Mercier"): - assert min(convergence_orders(l2_u)) > 1.9 - assert min(convergence_orders(l2_sigma)) > 1 - assert min(convergence_orders(l2_div_sigma)) > 0.9 + rates = (vdegree+1, sdegree, sdegree) + elif stress_element.family().startswith("Hu-Zhang"): + rates = (vdegree+1, sdegree+1, sdegree) else: raise ValueError("Don't know what the convergence should be") + assert min(convergence_orders(l2_u)) > rates[0] * 0.9 + assert min(convergence_orders(l2_sigma)) > rates[1] * 0.9 + assert min(convergence_orders(l2_div_sigma)) > rates[2] * 0.9 + def test_mass_conditioning(stress_element, mesh_hierarchy): mass_cond = []