Hi
You may write it as
def tau(u, mu):
n = u.geometric_dimension()
return mu*(grad(u)+grad(u).T-2./n*div(u)*Identity(n))
The 2/3 part in your form only works for the deviatoric tensor in 3D, not in 2D. A more general form can also make use of the UFL functions sym and dev, like this:
def tau(u, mu):
return mu*dev(2*sym(grad(u)))
where
sym(grad(u)) = 0.5(grad(u)+grad(u).T)
and
dev(A) = A - 1/n trace(A) Identity(n)