Assembling an identity matrix in Mixed spaces

FEniCS version: 2019.1.0 installed via Conda

Consider the following MWE on mixed spaces:

from dolfin import *

mesh = UnitSquareMesh(1,1)
P1 = FiniteElement('CG', triangle, 1)
P0 = FiniteElement('DG', triangle, 0)
element = MixedElement([P1, P0])
V = FunctionSpace(mesh, element)

v1, v2 = TestFunctions(V)
u1, u2 = TrialFunctions(V)

# Identity matrix of size dim(P0) x dim(P0) assembled at a particular location in the global matrix
cell = Cell(mesh, 0)
Eq = assemble((1/cell.volume())*u2*v2*dx)

print(Eq.array())

Output:

[[0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 1. 0. 0.]
 [0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 1.]]

The above MWE successfully assembles an Identity matrix of dimension dim(\mathbb{P}_0) \times dim(\mathbb{P}_0) at the correct global positions.

Question: Is there a way to extend it to \mathbb{P}_1 finite element case?

One approach is discussed here https://fenicsproject.org/qa/5372/assemble-an-identity-matrix/,

from dolfin import *

mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, 'CG', 1)
u = TrialFunction(V)
v = TestFunction(V)

form = inner(u, v)*dx
I = assemble(form)
I.zero()
I.set_diagonal(interpolate(Constant(1), V).vector())

However I can not add the obtained matrix I to a variational form. Is there a way to do this before assembling it?

Hi, a diagonal mass matrix can be obtained with P1 element by changing the integration rule of the measure to a “vertex” scheme, see here:
https://comet-fenics.readthedocs.io/en/latest/demo/tips_and_tricks/mass_lumping.html#Second-method

However, please note that you will not necessarily end up with an identity matrix, depending on the number of elements sharing a common dof.