I want to derive $\nabla \times \nabla \times \mathbf{V}$ in cylindrical coordinates. My variables are $(u,v,w)$ in $(z,r, \theta)$ (axissymetric, radial and azimuthal directions). I computed first stage of this operator: \begin{equation} \nabla \times (\nabla \times V) = \nabla \times \begin{bmatrix} \partial_r w - \frac{1}{r} (\partial_{\theta} v -w) \\ \frac{1}{r} \partial_{\theta} u - \partial_z w \\ \partial_z v - \partial_r u \end{bmatrix} \end{equation} Additional terms are due to covariant derivative in cylindrical coordinates: \begin{equation} A_{j,k}=\frac{1}{\sqrt{g_{kk}}} \frac{\partial A_j}{\partial x_k} - \Gamma^i_{jk} A_i \end{equation} Which will translate to: $\partial_{\theta} w = \frac{1}{r} \partial_{\theta} w + \frac{v}{r}$ and $\partial_{\theta} v = \frac{1}{r} \partial_{\theta} v - \frac{w}{r}$
Up to now I believe everything is OK. (Strong believe) But how we deal with covarant derivates when we compute second derivative?: $\partial_{\theta}(\partial_z w) = \frac{1}{r} \partial_{z\theta} w + \frac{1}{r} \partial_z v$
Here I took the derivative in z before derivative in r. But for radial direction it is not so obvious: \begin{equation} \partial_{\theta}(\partial_r w) = \partial_r(\partial_{\theta} w) = \partial_r (\frac{1}{r} \partial_{\theta} w + \frac{v}{r}) = -\frac{1}{r^2} \partial_{\theta} w + \frac{1}{r} \partial_{r\theta} w - \frac{v}{r^2} + \frac{1}{r} \partial_r v \end{equation} compared to \begin{equation} \partial_{\theta}(\partial_r w) = \frac{1}{r} \partial_{r\theta} w + \frac{1}{r} \partial_r v \end{equation}
Which approach is correct ?