Divergence in curvilinear coordinates In a curvilinear system with coordinates: $u_1, u_1, u_3$.
I want to understand why:
$$\nabla\cdot \vec{V} \neq  \frac{1}{h_1} \frac{\partial V_1}{\partial u_1} +\frac{1}{h_2} \frac{\partial V_2}{\partial u_2} + \frac{1}{h_3}\frac{\partial V_3}{\partial u_3}. \tag1$$
My reasoning to explain why eq.(1) is not true is the following: 
Taking the first term of the dot product:
$$\frac{\hat{e_1} }{h_1}\frac{\partial}{\partial u_1} (V_1 \hat{e_1}) =  \frac{\hat{e_1}}{h_1} \big( \hat{e_1} \frac{\partial V_1}{\partial u_1} + V_1 \frac{\partial \hat{e_1}}{\partial u_1} \big)$$
This seems rare to me since the definition of the dot product is:
$$A \cdot B = A_1B_1 + A_2B_2 + A_3B_3 $$
And using that definition the divergence should be eq.(1) (with an equality).
I would appreciate any references that explain this, since all I can find is the formula for the divergence in curvilinear systems but without a deduction. 
 A: Defining a "del" vector $\nabla=\mathbf{e}_i\frac{\partial}{\partial x^i}$ (summation on repeated indices assumed) is sort-of a mistake. This isn't a vector (field). It is just a nice notational device, because in this case, in cartesian coordinates, we have $$ \text{grad}(f)=\mathbf{e}_i\frac{\partial f}{\partial x^ i},$$ which looks like $\nabla f$, where $\nabla$ is the "vector" defined above and we have $$ \text{div}(\mathbf{A})=\frac{\partial A^ i}{\partial x^i}, $$ which sort of looks like $\nabla\cdot\mathbf{A}$, etc. But this is just a useful way to remember these formulas.

Now let's note several things. We can define a matrix whose $i,j$-th element is $\partial A^ i/\partial x^ j$. Let's call this matrix $\nabla\otimes\mathbf{A}$. We clearly have $$ \text{div}(\mathbf{A})=\text{Tr}(\nabla\otimes\mathbf{A})=\sum_i\frac{\partial A^ i}{\partial x^ i}. $$ (I'm gonna suspend the automatic summation from now on.)
If we have an orthogonal curvilinear coordinate system $(u^1,u^2,u^3)$, we can define the following quantities:

*

*Coordinate basis vectors $$ \mathbf{g}_i=\frac{\partial \mathbf{r}}{\partial u^ i}, $$


*"Reciprocal" coordinate basis vectors $$ \mathbf{g}^i=\nabla u^ i, $$


*An orthonormal frame $$ \hat{\mathbf{e}}_i=\frac{1}{h_i}\mathbf{g}_i, $$ where $h_i=\sqrt{\mathbf{g}_i\cdot\mathbf{g}_i}$.
These quantities have the property that $$ \mathbf{g}^ i\cdot\mathbf{g}_j=\delta^ i_j, $$ since ($x^ i$ are cartesian coordinates and $\mathbf{e}_i$ are cartesian basis vectors) we have $$ \mathbf{g}^i\cdot\mathbf{g}_j=\nabla u^i\cdot\frac{\partial\mathbf{r}}{\partial u^ j}=\sum_k\frac{\partial u^ i}{\partial x^k}\mathbf{e}_k\cdot\sum_l\frac{\partial x^l}{\partial u^ j}\mathbf{e}_l=\sum_{kl}\frac{\partial u^ i}{\partial x^k}\frac{\partial x^l}{\partial u_j}\delta_{kl}=\sum_k \frac{\partial u^ i}{\partial x^k}\frac{\partial x^k}{\partial u_j}=\delta^i_j.$$
We further introduce a matrix $g_{ij}$ given as $g_{ij}=\mathbf{g}_i\cdot\mathbf{g}_j$. Because of orthogonality, we have $g_{ij}=h_i^2\delta_{ij}$. We also introduce $g^{ij}=\mathbf{g}^i\cdot\mathbf{g}^j$, feel free to check that $g^{ij}=\frac{1}{h_i^2}\delta_{ij}$, thus $g_{ij}$ and $g^{ij}$ are inverse matrices.
If $\mathbf{A}$ is a vector field, we associate three kinds of components with $\mathbf{A}$. We can have

*

*$\mathbf{A}=\sum_i A^ i\mathbf{g}_i$, and the $A^ i$ are called contravariant components,


*$\mathbf{A}=\sum_i A_i\mathbf{g}^i$, and the $A_i$ are called covariant components and


*$\mathbf{A}=\sum_i \hat{A}_i\hat{\mathbf{e}}_i$ and the $\hat{A}_i$ are called physical components (here whether you put the index upwards or downwards doesn't matter!).
It is clear that we have $A^ i=\mathbf{g}^i\cdot\mathbf{A}$ and $A_i=\mathbf{g}_i\cdot\mathbf{A}$, and what is less clear is that differentiation with respect to the coordinates $\partial/\partial u^ i$ produce covariant components, eg. we have for any scalar function $f$, $\nabla f=\sum_i\mathbf{g}^i\frac{\partial f}{\partial u^ i}$ (check the definition of $\mathbf{g}^i$!).
We define Christoffel symbols of the second kind as $$ \Gamma^ k_{ij}=\frac{1}{2}\sum_l g^{kl}\left(\frac{\partial g_{jl}}{\partial u^ i}+\frac{\partial g_{il}}{\partial u^ j}-\frac{\partial g_{ij}}{\partial u^ l}\right). $$
Check that $\Gamma^k_{ij}=\mathbf{g}^k\cdot\frac{\partial \mathbf{g}_i}{\partial u^ j}$ (I really don't want to derive this identity right now).
In curvilinear coordinates, the basis vectors also depend on positions, so every time you differentiate a vector field, you need to make sure to take the variation of the basis vectors also into account, so we calculate the divergence as $$ \text{div}(\mathbf{A})=\sum_i\mathbf{g}^i\cdot\frac{\partial\mathbf{A}}{\partial u^i}=\sum_{ij}\mathbf{g}^i\left(\frac{\partial A^j}{\partial u^ i}\mathbf{g}_j+A^j\frac{\partial \mathbf{g}_j}{\partial u^ i}\right)=\sum_{ij}\left(\frac{\partial A^ i}{\partial u^i}+\Gamma^i_{ji}A^j\right). $$
We now calculate $$\Gamma^i_{ji}=\sum\frac{1}{2}g^{il}\left(\frac{\partial g_{jl}}{\partial u^i}+\frac{\partial g_{il}}{\partial u^j}-\frac{\partial g_{ij}}{\partial u^l}\right)=\sum\frac{1}{2}g^{il}\frac{\partial g_{il}}{\partial u^j},$$ which, in matrix notation ($(g_{ij})=g$) can be written as $$ \Gamma^i_{ji}=\frac{1}{2}\text{Tr}\left(g^{-1}\frac{\partial g}{\partial u^j}\right), $$ which, using Jacobi's formula, can be written as $$ \sum_i\Gamma^ i_{ji}=\frac{1}{2}\frac{1}{\det g}\frac{\partial}{\partial u^ j}\det g=\frac{1}{2}\frac{\partial}{\partial u^ j}\ln\det g=\frac{\partial}{\partial u^j}\ln\sqrt{\det g}=\frac{1}{\sqrt{\det g}}\frac{\partial}{\partial u^j}\sqrt{\det g}. $$
Plugging this back into the divergence formula gives $$\text{div}(\mathbf{A})=\sum_{ij}\left(\frac{\partial A^ i}{\partial u^i}+\Gamma^i_{ji}A^j\right)=\sum_{ij}\left(\frac{\partial A^i}{\partial u^i}+\frac{1}{\sqrt{\det g}}\frac{\partial}{\partial u^j}\sqrt{\det g}A^j\right)= \\ =\sum_{i}\frac{1}{\sqrt{\det g}}\frac{\partial}{\partial u^i}(\sqrt{\det g}A^i). $$
But since $g=\text{diag}(h_1^2,h_2^2,h_3^2)$
, we have $\det g=(h_1h_2h_3)^2$, so $$ \text{div}(\mathbf{A})=\frac{1}{h_1h_2h_3}\left\{\frac{\partial}{\partial u^1}(h_1h_2h_3A^1)+\frac{\partial}{\partial u^2}(h_1h_2h_3A^2)+\frac{\partial}{\partial u^3}(h_1h_2h_3A^3)\right\}= \\ =\frac{1}{h_1h_2h_3}\left\{\frac{\partial}{\partial u^1}(h_2h_3\hat{A}_1)+\frac{\partial}{\partial u^2}(h_1h_3\hat{A}_2)+\frac{\partial}{\partial u^3}(h_1h_2\hat{A}_3)\right\}, $$
which is the formula for divergence in terms of the physical components.
A: Any good book on vector or tensor analysis should provide you with a first principles derivation for some specific curvilinear coordinate systems (spherical and cylindrical systems for example).  For the more general case, you need to understand that your definition for the dot product of two vectors is valid only for Cartesian coordinates.  Why is that?  Because in general the dot product is defined as the multiplicative operation that converts two vector quantities into a scalar quantity.  Likewise the divergence is defined as the first order differential operator that yields a scalar quantity when operating on a vector.  In tensor analysis this operator is called the covariant derivative and for curvilinear coordinate systems the usual partial derivatives are augmented by terms called Christoffel symbols.  These extra terms render the divergence a scalar quantity.  This answer should provide you with keywords that will take you to a complete understanding. 
A: By definition $\nabla \cdot \vec{E}=\lim_{\Delta V\to 0 }\frac{\int \vec{E}\cdot \hat{n}dA}{\Delta V} =\lim_{\Delta V\to 0} \frac{\int E_in_ih_jh_kdx_jdx_k +\int E_jn_jh_ih_kdx_idx_k+  \int E_kn_kh_jh_idx_jdx_i}{\int h_ih_jh_kdx_idx_jdx_k}$
$$= \frac{1}{h_ih_jh_k}(\frac{\partial(E_ih_jh_k) }{\partial x_i }+\frac{\partial(E_jh_ih_k) }{\partial x_j }+ \frac{\partial(E_kh_jh_i) }{\partial x_k } )$$
So for example in cylindrical coordinates:
Cancel all the $dx_i$'s in the numerator with the denominator,term by term and differentiate with the respect to the remaing differential. Treat where applicable the $h_i$'s as constants.
$$\nabla \cdot \vec{E}=\frac{1}{r}\frac{\partial (rE_r)}{\partial r}+\frac{1}{r}\frac{\partial E_\theta }{\partial  \theta}+\frac{\partial E_z}{\partial z}$$
The  scaling factors can make things complicated. Sometimes they don't effect the differentiation and so cancel upon division.
