How to find electric current density from electric potential in 3d finite element

Stacky
Hello,
The value of electric potential(∅) is known at every node in a 3d finite element mesh. The relation between electric current density(i) and electric potential(∅) is i=k.∇∅, I am writing a code in c, I want to know how to find the gradient of electric potential(∇∅) at every node so as to get current density(i). Please help.

If we subscript the nodes using i,j,k then the bilinear approximation would be,

$$J_x = \frac{\sigma}{2 \Delta} \left[ V_{(i-1),j,k} - V_{(i+1),j,k} \right]$$

and similarly for the y and z components.

Homework Helper
The OP said a finite element mesh. The previous answer looks more appropriate for a regular finite differnce mesh, not an FE mesh with arbitrary geometry.

You can calculate the gradient at each node of each element using the element shape functions. (This will be similar, but simpler, than calculating strains from displacements in a structural FE program).

The problem is that the only output that is consistent with the FE method is gradient averaged over the volume of each element. When several elements meet at a node, you will get a different gradient value at the node for each of the elements.

If you just want to plot the data you can do something fairly simple, for example average of the different nodal values (possibly weighted by the volume of the elements). Or you can use the difference between the element values at each node as a measure of the accuracy of the solution.

If you want to do some mathematical post processing using the gradients, it would be better to do it based on the value within the volume of each element, rather than some more or less arbitrary nodal averaging method.

Ok, it looks like I was over-simplifying it by assuming it was a rectangular grid.

BTW Stacky. The electric field is $E = - \nabla \phi$, so it should be $J = - \sigma (\nabla \phi)$

Stacky
Method suggested by uart is correct, but it suits a regular structured rectangular mesh. Because I am dealing with unstructured mesh, I would prefer to use the method proposed by AlphaZero. Thanks both for their inputs.
Neverthless in strict sense ∇∅ is an elemental solution and we need to use some averaging method to calculate its value at nodes.