Fix `enstrophy` computation for `DGMulti`

The original implementation of the `enstrophy` analysis callback introduced in [PR1239](https://github.com/trixi-framework/Trixi.jl/commit/3f48f64d07ab90cc42c324a430e3b372162d91e4#diff-af5fa59f044a6b62471fd6e4189b4e136c339901623a229a05d5155049044fbc) does not work and was not tested.

It is removed in #2140 .

A starting point for a fix would be (note that there are already couple changes to the original version)

```julia
function integrate(func::typeof(enstrophy), u,
mesh::DGMultiMesh,
equations, equations_parabolic::CompressibleNavierStokesDiffusion3D,
dg::DGMulti,
cache, cache_parabolic; normalize=true)
gradients_x, gradients_y, gradients_z = cache_parabolic.gradients
# allocate local storage for gradients.
# TODO: can we avoid allocating here?
local_gradient_quadrature_values = ntuple(_ -> copy(cache_parabolic.local_u_values_threaded), 3)
integral = zero(eltype(eltype(u)))
for e in eachelement(mesh, dg)
u_quadrature_values = cache_parabolic.local_u_values_threaded[Threads.threadid()]
gradient_x_quadrature_values = local_gradient_quadrature_values[1][Threads.threadid()]
gradient_y_quadrature_values = local_gradient_quadrature_values[2][Threads.threadid()]
gradient_z_quadrature_values = local_gradient_quadrature_values[3][Threads.threadid()]
# interpolate to quadrature on each element
apply_to_each_field(mul_by(dg.basis.Vq), u_quadrature_values, view(u, :, e))
apply_to_each_field(mul_by(dg.basis.Vq), gradient_x_quadrature_values, view(gradients_x, :, e))
apply_to_each_field(mul_by(dg.basis.Vq), gradient_y_quadrature_values, view(gradients_y, :, e))
apply_to_each_field(mul_by(dg.basis.Vq), gradient_z_quadrature_values, view(gradients_z, :, e))
# integrate over the element
for i in eachindex(u_quadrature_values)
gradients_i = SVector(gradient_x_quadrature_values[i],
gradient_y_quadrature_values[i],
gradient_z_quadrature_values[i])
integral += mesh.md.wJq[i, e] * func(u_quadrature_values[i], gradients_i, equations_parabolic)
end
end
return integral
end
```

Start with the enstrophy analysis callback and the integrate path for DGMulti, using the removed implementation in PR1239 and the changes referenced in #2140 as context. Check how DGMultiMesh, CompressibleNavierStokesDiffusion3D, and the parabolic gradients are handled. Done means enstrophy computes correctly for DGMulti and the previously untested behavior is covered by a suitable check.

Written by the indexing model from the issue text.