I am trying to calculate the element-wise mass balance (i.e., \int_{\Omega_element} div(u) d\Omega) over every cell and project the output to DG0 space. Has this been done or is there a way (e.g., with par_loop() ) to implement this in firedrake?
This is how I achieved it in fenics:
-------------------------------------------------------------------
cell_domains = CellFunction('size_t', mesh)
cell_domains.set_all(0)
dx = Measure('dx')[cell_domains]
area = [cell.volume()   for cell in cells(mesh)]
for i in range(mesh.num_cells()):
    Error.vector()[i]  = (1./area[i]) * assemble(div(solution)*dx(i))
-------------------------------------------------------------------
Regards,
Sarraf