I am attempting to calculate the functional derivative of a functional $E[\rho] = \int G(\rho(\mathbf{r}),\nabla\rho(\mathbf{r}),\mathbf{r})d\mathbf{r},$ where $G(\rho(\mathbf{r}),\nabla\rho(\mathbf{r}),\mathbf{r})=\rho(\mathbf{r})^{4/3}\left(\alpha-\frac{(\nabla\rho(\mathbf{r})\cdot\nabla\rho(\mathbf{r}))^{3/4}}{137 \rho(\mathbf{r})^{2}}\right),$ and $\alpha$ is a constant. This is for use in a computational chemistry code.
To find the functional derivative I think I should use the Euler-Lagrange equation, $\frac{\delta G}{\delta \rho}=\frac{\partial G}{\partial \rho} - \nabla\cdot\frac{\partial G}{\partial \nabla \rho}, $ as given on the Wikipedia article on functional derivatives.
What I am struggling with is the second term in the E-L equation. Firstly, I am not sure how to approach the partial derivative with respect to $\nabla\rho$. So far, I have use the chain rule to obtain $ \frac{\partial G}{\partial \nabla\rho}=-\frac{3}{4\times 137 \rho^{2/3}}\frac{1}{(\nabla\rho(\mathbf{r})\cdot\nabla\rho(\mathbf{r}))^{1/4}}\left(\frac{\partial}{\partial \nabla\rho}(\nabla\rho(\mathbf{r})\cdot\nabla\rho(\mathbf{r}))\right), $ but I am not sure how to proceed with the differentiation of the dot product. Furthermore, it appears from the E-L equation that I must then find the divergence of this partial derivative. I think that the result of $\frac{\partial G}{\partial \nabla\rho}$ will be a scalar function, so am not sure how the divergence can be applied here.
I would appreciate some advice on how to tackle the partial derivative and subsequent divergence. Perhaps I am missing something, or there is a flaw in my reasoning.