most efficient way to compute nodal strain-rate values

Question

I would like to compute the elements of the strain rate tensor

$$\dot{\epsilon} = \frac{1}{2} \left( \nabla u + \left( \nabla u \right)^T \right)$$

at each nodal value, so with a mesh of $n$ nodes in two dimension, four vectors.