From 9a984d508b8f80d9610ab80c4d5e4fb8728aabfc Mon Sep 17 00:00:00 2001 From: Pablo Brubeck Date: Mon, 19 Aug 2024 15:10:45 +0100 Subject: [PATCH] Clean up variants --- tsfc/finatinterface.py | 28 ++++++++++++---------------- 1 file changed, 12 insertions(+), 16 deletions(-) diff --git a/tsfc/finatinterface.py b/tsfc/finatinterface.py index 3e832819..d7b4f66c 100644 --- a/tsfc/finatinterface.py +++ b/tsfc/finatinterface.py @@ -143,11 +143,15 @@ def convert_finiteelement(element, **kwargs): return finat.FlattenedDimensions(finat_elem), deps kind = element.variant() - is_interval = element.cell.cellname() == 'interval' if kind is None: kind = 'spectral' # default variant + is_interval = element.cell.cellname() == 'interval' - if element.family() == "Lagrange": + if element.family() in {"Raviart-Thomas", "Nedelec 1st kind H(curl)", + "Brezzi-Douglas-Marini", "Nedelec 2nd kind H(curl)", + "Argyris"}: + lmbda = partial(lmbda, variant=element.variant()) + elif element.family() == "Lagrange": if kind == 'spectral': lmbda = finat.GaussLobattoLegendre elif kind.startswith('integral'): @@ -171,10 +175,6 @@ def convert_finiteelement(element, **kwargs): else: # Let FIAT handle the general case lmbda = partial(finat.Lagrange, variant=kind) - elif element.family() in {"Raviart-Thomas", "Nedelec 1st kind H(curl)", - "Brezzi-Douglas-Marini", "Nedelec 2nd kind H(curl)", - "Argyris"}: - lmbda = partial(lmbda, variant=element.variant()) elif element.family() in ["Discontinuous Lagrange", "Discontinuous Lagrange L2"]: if kind == 'spectral': lmbda = finat.GaussLegendre @@ -186,6 +186,8 @@ def convert_finiteelement(element, **kwargs): lmbda = lambda *args: finat.DiscontinuousElement(finat.FDMLagrange(*args)) elif kind in 'fdm_broken' and is_interval: lmbda = finat.FDMBrokenL2 + elif kind in ['demkowicz', 'fdm']: + lmbda = partial(finat.Legendre, variant=kind) elif kind in ['mgd', 'feec', 'qb', 'mse']: degree = element.degree() shift_axes = kwargs["shift_axes"] @@ -195,16 +197,10 @@ def convert_finiteelement(element, **kwargs): else: # Let FIAT handle the general case lmbda = partial(finat.DiscontinuousLagrange, variant=kind) - elif element.family() == ["DPC", "DPC L2"]: - if element.cell.geometric_dimension() == 2: - element = element.reconstruct(cell=ufl.cell.hypercube(2)) - elif element.cell.geometric_dimension() == 3: - element = element.reconstruct(cell=ufl.cell.hypercube(3)) - elif element.family() == "S": - if element.cell.geometric_dimension() == 2: - element = element.reconstruct(cell=ufl.cell.hypercube(2)) - elif element.cell.geometric_dimension() == 3: - element = element.reconstruct(cell=ufl.cell.hypercube(3)) + elif element.family() == ["DPC", "DPC L2", "S"]: + dim = element.cell.geometric_dimension() + if dim > 1: + element = element.reconstruct(cell=ufl.cell.hypercube(dim)) return lmbda(cell, element.degree()), set()