Calculate flux from solution with Fenicsx

[dokken]

@bragostin

There is only one thing I still struggle to understand, is why I should not use CG1 in Expression(k*grad(T), V.element.interpolation_points())considering that the problem variable T was defined with

V = dolfinx.fem.FunctionSpace(mesh, ('CG', 1))
T = ufl.TrialFunction(V)

When you differentiate a piecewise linear continuous function, you end you with cellwise constant basis functions. See for instance Introduction to the Finite Element Method (FEM) — Simula Summer School in Computational Physiology