I was exploring quadrature rules and manual integration in dolfinx while encountered something that seems like a problem. When manually assembling any matrix, it is necessary to know which basis value computed in some point (quadrature for example, but not necessary) corresponds to which cell dof. From documentation it looks like function T_apply (dolfinx.fem — DOLFINx 0.10.0.0 documentation) should do what I’m looking for and rearrange basis functions values to match dof ordering. However, in my MWE I do not notice any difference after calling T_apply:
import basix
import dolfinx
import numpy as np
from mpi4py import MPI
# Define mesh
mesh = dolfinx.mesh.create_unit_square(MPI.COMM_SELF, 2, 2, cell_type=dolfinx.mesh.CellType.triangle)
V = dolfinx.fem.functionspace(mesh, ("Lagrange", 1))
# Define quadrature rule
order=2
quad_points, weights = basix.make_quadrature(cell=V.mesh.basix_cell(),degree=order)
# Tabuate basis function values
basis_values = V.element.basix_element.tabulate(0, quad_points)[0, :, :, 0]
# Compute permutation data
mesh.topology.create_entity_permutations()
cell_permutations = mesh.topology.get_cell_permutation_info()
print(cell_permutations)
# Apply permutations to get basis function values in real cells
basis_values_extended = np.tile(basis_values.flatten(),len(cell_permutations))
print("before \n")
print(basis_values_extended.reshape((len(cell_permutations),*basis_values.shape)))
V.element.T_apply(basis_values_extended,cell_permutations,mesh.topology.dim)
print("after \n")
print(basis_values_extended.reshape((len(cell_permutations),*basis_values.shape)))
I would really appreciate any explanation of T_apply behavior and errors that I’ve possibly made. Note that using mesh.topology.dim**2 as last argument for T_apply introduces no changes at all which seems odd. Same result is true for unit cube meshed with tetrahedrons.
Regarding specs, here are versions of related packages that I use:
Lagrange elements are not such an element, as the dofmap is permitted at the creation of the function space to account for a consistent ordering.
What you seem like you would like is the element-dof-layout: dolfinx.cpp.fem — DOLFINx 0.10.0.0 documentation
which tells you how each dof in V.cell_dofs(i) relates to the element entities. You can either look at the closure (for instance all dofs related to a facet, the edges of a facet and the vertices), or just those associated with such an entity.
Do I understand correctly, that if (for 2D case) calling V.dofmap.dof_layout.entity_closure_dofs(2,0) yieds [0,1,2], then no matter how each cell is oriented, dof indices returned by V.dofmap.cell_dofs(i) are arranged in consistency with reference element? For example, if cell dofs are [0,3,2], than tabulating basis function values in any point inside this cell will return triplet where first value corresponds to dof 0, second to dof 3 and last one for dof 2?
Instead of calling closure dofs, I would call entity_dofs(2,0), then (2,1), then (2,2) to get which dof is associated with the 1st,2nd and 3rd vertex (as you avoid assumptions on the order in entity closure dofs).
With respect to the order, if the above returns (0),(1),(2) then the vertices of the reference element aligns with the tabulated basis functions.
If cell dofs(i) align as this, then you have that the first vertex of cell i has number 0 (global number, local to process), the second has 3, etc.
I suppose, you meant entity_dofs(0,0), entity_dofs(0,1) and entity_dofs(0,2)? Yes, they return [0],[1] and [2]. However, I was a bit confused by this function as none of the higher dimension entities are indexed. Calling entity_dofs(1,0) as well as entity_dofs(2,0) will result in []. Same I see for tetrahedrons, the function can accept expected number of entities (up to 4 faces for example), but will always return empty list for any dim except 0. What is the idea behind this function then?
Yes, sorry. I’m on my phone atm, so typing sometime is a bit too fast.
I guess you are «only» considering first order Lagrange functions. DOLFINx supports many other elements of various degrees, where dogs are associated with entities. You can for instance see this for a third order element at:
