Taking the Vector Gradient of a Function in $(r)$ where $r = |\mathbf{\vec{r}}|$ I have tried looking for the answer to this online, but probably due to my personal lack of knowledge have been unable to ask the question effectively.
Essentially I would like to understand how to take the vector gradient of the following equation (the equation has been left un-simplified deliberately):
$$u(r_{ij}) = D_0\left[e^{-2\alpha(r_{ij}-r_0)} - 2e^{-\alpha(r_{ij}-r_0)} \right]$$
Where $r_{ij}$ represents distance and is a magnitude and $r_0$, $\alpha$, and $D_0$ are constants. The only problem is I'm not sure what form the vector gradient equation would take when the distance is given as $r_{ij}$. Has anyone seen something similar?
My naive attempt is below using $\nabla_{ij} = \frac{d}{dr_{ij}}$:
$$\boldsymbol\nabla u(r_{ij}) = D_0\left[-2\alpha e^{-2\alpha(r_{ij}-r_0)} + 2\alpha e^{-\alpha(r_{ij}-r_0)} \right]$$
But I assume as its vector quantity a $\vec{r_{ij}}$ must fall out somewhere, I'm just not sure where!
For those who are interested the equation itself is the definition of an atomic pair potential that gives an energy value based on the distance defined by $r_{ij}$. The distance is the magnitude of the cartesian vector joining two atoms and is defined:
$$r_{ij} = |\vec{r_i} - \vec{r_j}|$$
 A: The definition below may not be appropriate, but defining $\mathbf{\vec{r}}$ as:
$$\mathbf{\vec{r}}_{ij} = x \mathbf{\hat{i}} + y \mathbf{\hat{j}} + z \mathbf{\hat{k}}$$
Where $\mathbf{\vec{r}}_{ij}$ defines the vector between atoms $i, j$, and the values $\mathbf{\hat{i}}, \mathbf{\hat{j}}, \mathbf{\hat{k}}$ are unit vectors.
Gives you:
$$r_{ij} = |\mathbf{\vec{r}}_{ij}| = \sqrt{x^2+y^2+z^2} = \left(x^2+y^2+z^2\right)^{1/2}$$
Then defining $\nabla$ as:
$$\nabla = \frac{\partial}{\partial x} \mathbf{\hat{i}} + \frac{\partial}{\partial y} \mathbf{\hat{j}} + \frac{\partial}{\partial z} \mathbf{\hat{k}} $$
Gives you:
$$\nabla u(r_{ij}) = \frac{\partial u(r_{ij})}{\partial x} \mathbf{\hat{i}} + \frac{\partial u(r_{ij})}{\partial y} \mathbf{\hat{j}} + \frac{\partial u(r_{ij})}{\partial z} \mathbf{\hat{k}} $$
This allows you to use the chain rule for each partial derivative, so for each partial derivative we have a variant of what's shown below for the partial with respect to $x$:
$$\frac{\partial u(r_{ij})}{\partial x} = \frac{\partial r_{ij}}{\partial x}\frac{\partial u(r_{ij})}{\partial r_{ij}}$$
Letting $r_{ij} = r$ for notational simplicity we have:
$$\frac{\partial u(r)}{\partial x} = \frac{\partial r}{\partial x}\frac{\partial u(r)}{\partial r}$$
The first part of the right hand side of the equation can be calculated for $x$ as:
$$\begin{eqnarray}\frac{\partial r}{\partial x} = \frac{\partial \left(x^2 + y^2 + z^2\right)^{1/2}}{\partial x} &=& \left(2x\right) \left(\frac{1}{2}\right) \left(x^2 + y^2 + z^2\right)^{-1/2} \\ \\
&=& \frac{x}{\left(x^2 + y^2 + z^2 \right)^{1/2}}\end{eqnarray}$$
As $r = \left(x^2+y^2+z^2\right)^{1/2}$ this simplifies to:
$$\frac{\partial r}{\partial x} = \frac{x}{r}$$
Similar equations for the partial with respect to $y$ and $z$ can be found. This leaves the second part of the right hand side equation:
$$u(r) = D_0\left[ e^{-2\alpha(r-r_0)} -2e^{-\alpha(r-r_0)}\right]$$
Where again we have let $r_{ij} = r$, taking the partial with respect to r we have:
$$\frac{\partial u(r)}{\partial r} = D_0\left[ \left(-2\alpha\right)e^{-2\alpha(r-r_0)} - (-2\alpha)e^{-\alpha(r-r_0)}\right]$$
As $\frac{\partial u(r)}{\partial r}$ will be the same for all partials with respect to $x, y, z$ it can be factored out of the final equation and we have:
$$\nabla u(r) = \frac{\partial u(r)}{\partial r}\left(\frac{\partial r}{\partial x}\mathbf{\hat{i}} + \frac{\partial r}{\partial y}\mathbf{\hat{j}} + \frac{\partial r}{\partial z}\mathbf{\hat{k}}\right)$$
Subbing in the partials with respect to $x, y, z$ along with $\frac{\partial u(r)}{\partial r}$ you get:
$$\begin{eqnarray}\nabla u(r) &=& D_0\left[ \left(-2\alpha\right)e^{-2\alpha(r-r_0)} - (-2\alpha)e^{-\alpha(r-r_0)}\right]\left(\frac{x}{r}\mathbf{\hat{i}} + \frac{y}{r}\mathbf{\hat{j}} + \frac{z}{r}\mathbf{\hat{k}}\right) \\
\\
 &=& \frac{-\left(2 \alpha\right)D_0}{r} \left[ e^{-2\alpha(r-r_0)} - e^{-\alpha(r-r_0)}\right]\left(x\mathbf{\hat{i}} + y\mathbf{\hat{j}} + z\mathbf{\hat{k}}\right)\end{eqnarray}$$

As $\mathbf{\vec{r}}_{ij} = x \mathbf{\hat{i}} + y \mathbf{\hat{j}} + z \mathbf{\hat{k}}$ you have:
$$\nabla u(r) = \frac{-\left(2 \alpha\right)D_0}{r} \left[ e^{-2\alpha(r-r_0)} - e^{-\alpha(r-r_0)}\right]\mathbf{\vec{r}}$$
Which finally simplifies to:
$$\nabla u(r) = -\left(2 \alpha\right)D_0 \left[ e^{-2\alpha(r-r_0)} - e^{-\alpha(r-r_0)}\right]\mathbf{\hat{r}}$$
Where $\mathbf{\hat{r}} = \frac{\mathbf{\vec{r}}}{r}$ and is a unit vector in the direction of $\mathbf{\vec{r}}$, here the notation has also been simplified and $\mathbf{\vec{r}} = \mathbf{\vec{r}}_{ij}$ from the original definition.
