When does [the distance from origin to some point in space] vary and when is it fixed? As defined in spherical coordinates, $p =$ the distance from the origin to a point P in space.
In Stewart P1092 16.7.47 (below), $p$ is fixed by the solution to be $a$.
Yet in Stewart P1103 16.9.13, $p$ is a variable. How does one determine whether $p$ varies or is fixed? 
I realise that in 13, the Divergence Theorem effects a triple integral so 3 variables must be integrated. In the interest of spherical coordinates, $\phi \, \& \, \theta$ are used, so $p$ may fit. This is desultory guesswork, so
I'd like to learn a more definitive explanation.

 A: Both problems, 13 and 47, are about the divergence theorem. This theorem reads as follows: Given a flow field ${\bf v}$ in some domain $\Omega\subset{\mathbb R}^3$, and in addition a "body" $B\subset\Omega$ with surface $\partial B$ (oriented outwards) the "integral formula"
$$\int_{\partial B}{\bf v}\cdot d\vec\omega=\int_B{\rm div}({\bf v}){\rm d}({\bf x})\tag{1}$$
is valid. Here the left hand side is a surface integral (involving two variables $u$, $v$, or, e.g., $\phi$, $\theta$, when it has to be computed "the hard way"), and the right hand side is a volume integral (involving three variables $x$, $y$ $z$, or, e.g.,  $r$, $\phi$, $\theta$, when it has to be computed "the hard way").
In problem 47 one has $\Omega={\mathbb R}^3\setminus{\bf 0}$. Therefore we cannot apply the divergence theorem to balls centered at ${\bf 0}$. But we can apply it to spherical shells $B:=\{{\bf x}\>|\>a\leq|{\bf x}|\leq b\}$ with boundary $\partial B=\partial S_b-\partial S_a$, where the minus sign takes care of the fact that the standard orientation of the inner sphere $S_a$ has to inverted when $S_a$ is considered as part of $\partial B$. 
Since it is easily verified that ${\rm div}({\bf v})\equiv0$ in $\Omega$ we obtain
$$\int_{S_b}{\bf v}\cdot d\vec\omega -\int_{S_a}{\bf v}\cdot d\vec\omega=\int_{\partial B}{\bf v}\cdot d\vec\omega=0\ ,$$
which is saying that the integrals $\int_{S_b}{\bf v}\cdot d\vec\omega$ and $\int_{S_a}{\bf v}\cdot d\vec\omega$ over the outward oriented spheres have the same value.
In problem 13 the field  ${\bf v}$ is defined on all of ${\mathbb R}^3$, so that $(1)$ may be applied to the ball $B:=\{{\bf x}\>|\>|{\bf x}|\leq R\}$. You have computed ${\rm div}({\bf v})=12(x^2z+y^2z+z^3)$, whence $(1)$ gives
$$\int_{S_R}{\bf v}\cdot d\vec\omega =12 \int_B (x^2z+y^2z+z^3)\ {\rm d}(x,y,z)\ .$$
Now $B$ is symmetric with respect to ${\bf x}\mapsto -{\bf x}$, and each term in the last integral is odd. It follows that the integral is $0$, by symmetry.
A: It is essentially a matter of convenience: sometimes it looks easier to do the integration over the surface, especially when the field is proportional to $\overrightarrow r$.
Other times it is easier to use Gauss' theorem and switch to $div$ and volume integral. And every time you will have to do some context-related guessing what the author intends you to do, obviously :-).
As you can easily see, $13$ is designed for applying $div$, as by magic a term $2z\rho^2$ appears which makes integrating over a volume easier. Trying to do the sphere surface integral is much more difficult, because one has to calculate the surface integral of $z(4x^4+4y^4+3z^4)$ and substitute $x=Rcos\theta sin\phi$ etcetera, as above.
$47$ is clearly fit for just doing the surface integral, but there's really no problem going the other route here and apply $div$:
$$F_x = \frac{cx}{r^3} = \frac{cx}{(x^2+y^2+z^2)^{3/2}}$$
 etcetera, so
\begin{equation}
\begin{split}
div F &= \frac{\partial F_x}{\partial x} + \frac{\partial F_y}{\partial y} +\frac{\partial F_z}{\partial z}
\\ &= \frac{c}{r^3} - \frac{3c}{2}\frac{x}{r^5}\cdot 2x + \frac{c}{r^3} - \frac{3c}{2}\frac{y}{r^5}\cdot 2y + \frac{c}{r^3} - \frac{3c}{2}\frac{z}{r^5}\cdot 2z
\\ &= \frac{3c}{r^3} - \frac{3c(x^2+y^2+z^2)}{r^5}
\\ &= \frac{3c}{r^3} - \frac{3c}{r^3}
\\ &= 0
\end{split}
\end{equation}
so a volume integral of $div F$ seems to always equal zero - but we have to take into account the value of $div F$ at the origin. As you may know, one of Maxwell's equations states
$$div \overrightarrow E = \frac{\rho}{\epsilon _0}$$
where $\rho$ is the electrical charge density; essentially it says that electrical field lines originate from charges. Here this means that we have to know whether $div F$ has a value in $O$. Obviously it has to: otherwise where are all these field lines coming from? Assuming the $\overrightarrow F$ mentioned is actually just the electrical field, we see that at the origin a (positive) point charge $4\pi \epsilon_0 c$ has to be present (with infinite charge density but obviously related to the point-likeness in such a way as to exactly equal the correct charge), and thus the surface integral of $\overrightarrow F$ over $any$ surface enclosing the origin (not just spheres) equals $4\pi \epsilon_0 c$.
