firedrakeproject / firedrake

Firedrake is an automated system for the portable solution of partial differential equations using the finite element method (FEM)
https://firedrakeproject.org
Other
517 stars 160 forks source link

Problem with FunctionSpaces with similar FiniteElements but different variants #1610

Open tommbendall opened 4 years ago

tommbendall commented 4 years ago

If two instances of FunctionSpace are made with FiniteElements that are similar but have a different variant, I think the two FunctionSpaces aren't completely distinguished. For instance, they will have the same set of boundary nodes (even if these should be different).

Example:

from firedrake import *
mesh = RectangleMesh(4, 4, 4, 4, quadrilateral=True)
cell = mesh.ufl_cell().cellname()
DG_equi_elt = FiniteElement("DG", cell, 1, variant="equispaced")
DG_spec_elt = FiniteElement("DG", cell, 1, variant="spectral")

DG_spec = FunctionSpace(mesh, DG_spec_elt)
print(DG_spec.boundary_nodes("on_boundary", "geometric"))

DG_equi = FunctionSpace(mesh, DG_equi_elt)
print(DG_equi.boundary_nodes("on_boundary", "geometric"))

produces the output

[ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 20 21 22 23 24 25 26 27
 36 37 38 39 40 41 42 43 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63]
[ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 20 21 22 23 24 25 26 27
 36 37 38 39 40 41 42 43 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63]

but we expect these to differ as the equispaced elements should have basis functions which are zero on the boundaries. If the order of building function spaces is swapped, i.e.

DG_equi = FunctionSpace(mesh, DG_equi_elt)
print(DG_equi.boundary_nodes("on_boundary", "geometric"))

DG_spec = FunctionSpace(mesh, DG_spec_elt)
print(DG_spec.boundary_nodes("on_boundary", "geometric"))

then the output is

[ 0  1  2  4  6  8 10 12 14 20 22 24 26 27 36 38 39 41 43 50 51 53 55 58
 59 61 62 63]
[ 0  1  2  4  6  8 10 12 14 20 22 24 26 27 36 38 39 41 43 50 51 53 55 58
 59 61 62 63]
wence- commented 4 years ago

Thanks, this is an issue in how we stash "shared" data on function spaces. The introduction of variants is not safe since the variant is not taken into account in the cache key. Perhaps we should cache on the finat element instead in this case...

wence- commented 1 month ago

Unassigning myself, but just noting that this is still an issue (although geometric boundary nodes have been removed).