Derivative in dolfinx

Hi, I solved an electrostatic problem and exported the electric potential. The electric potential is smooth and then its numeric derivative was calculated to get the electric field. The problem is that the electric field is not smooth now. So I think if it would be better to directly calculate the derivative in dolfinx. How to calculate derivative of a function u on domain V with regard to z in dolfinx for example?
Figure_1

Please supply a minimal reproducible example that illustrates your issue

Not everything can make a smallest example. I give up. But how to calculate derivative in dolfinx?

Why cannot upload file?

I don’t understand what derivative you want to calculate.
Is it a spatial derivative of a dolfinx.fem.Function?
Is it a derivative of an ufl-Expression with respect to a variable?

If u is a function in say a first order lagrange space, then its derivative in z-direction is a piecewise continuous function on every cell, thus it is no longer continuous.
To get such derivatives, you can use u.dx(i) and interpolate with dolfinx.fem.Expression into a suitable function space (“DG” 0 in the case of u being in P1).

See for instance:
https://jsdokken.com/fenics22-tutorial/heat_eq.html#post-processing-without-projections
and
https://fenics.readthedocs.io/projects/ufl/en/latest/manual/form_language.html?highlight=dx#basic-spatial-derivatives

Dear dokken,
Thanks for your help. The solution u is a second order lagrange function defined as

V = FunctionSpace(domain, ("CG", 3))
u = TrialFunction(V)

Then the electric field Ez along z axis I calculated as:

W = FunctionSpace(domain, ("DG", 2))

expr = Expression(uh.dx(2), W.element.interpolation_points())

w = Function(W)

w.interpolate(expr)

tol = 0.001 # Avoid hitting the outside of the domain
z = np.linspace(0 + tol, 31.6 - tol, 101)
points = np.zeros((3, 101))
points[2] = z
u_values = []
bb_tree = geometry.BoundingBoxTree(domain, domain.topology.dim)

cells = []
points_on_proc = []
cell_candidates = geometry.compute_collisions(bb_tree, points.T)
colliding_cells = geometry.compute_colliding_cells(domain, cell_candidates, points.T)
for i, point in enumerate(points.T):
    if len(colliding_cells.links(i))>0:
        points_on_proc.append(point)
        cells.append(colliding_cells.links(i)[0])

points_on_proc = np.array(points_on_proc, dtype=np.float64)
u_values = w.eval(points_on_proc, cells)

Is it correct?

Yes, this would be one way of doing it. Note that if a point is on a boundary between two cells, there will be some ambiguity as to which value it should take (That is decided by whatever cell is first in colliding_cells for the given point).

Thanks again. The derivative, that is the electric field, is much smooth now.