Since Spivak requires absolute convergence, it is not enough to rearrange the series. We need a "blocking" lemma of the following kind:
Lemma. Suppose $\sum a_n$ is a conditionally convergent series. Then for every $\beta\in \mathbb R$ there exists a partition $\mathbb N=\bigcup_{j=1}^\infty B_j$ where each set $B_j$ is finite, and the series $\sum_j \left(\sum_{n\in B_j} a_n\right)$ converges absolutely to $\beta$.
Assume the lemma for now. Our partition should include open sets $U_j=\bigcup_{n\in B_j}(j,j+1)$. These do not cover the integers, though. Cover each integer by a tiny interval $(n-\epsilon,n+\epsilon)$ where $\epsilon$ is small enough so that the interval does not meet the sets $A_n$ and $A_{n-1}$. Note that $f$ is identically zero within such intervals. For each $n$, the set $A_n$ meets exactly one element $U_j$ of our open cover; therefore, the corresponding function $\varphi_j$ from a partition of unity will be identically $1$ on $A_n$. We have $\int f\varphi_j=\sum_{n\in B_j} a_n$ which gives the desired result. $\quad\Box$
As for the lemma, the proof goes like this: first, rearrange the series so that it converges conditionally to $\beta$; then divide the rearranged series $\sum a_{\sigma(k)}$ into blocks $N_i\le \sigma(k)<N_{i+1}$, $i=1,2,3,\dots$ such that the sum within the $i$th block is at most $2^{-i}$ for all $i\ge 2$. (Idea: cut off the tail repeatedly when it gets smaller than $2^{-i}$.)
I know this is an old question, but I thought this explanation might be helpful to some.
By definition (in $\mathbb R^3$):
$$\int_{\partial B(\boldsymbol x,r)}f(\boldsymbol y)dS(\boldsymbol y)=
\int_U f(\boldsymbol y(s,t))\left\|\frac{\partial\boldsymbol y}{\partial s}\times\frac{\partial\boldsymbol y}{\partial t}\right\|dsdt$$
Now, observe that $f(\boldsymbol y)=f(\boldsymbol x+r(\frac{\boldsymbol y-\boldsymbol x}{r}))$, and that if $\boldsymbol y(s,t)$ is a parametrization of $\partial B(\boldsymbol x,r)$ for $(s,t)\in U$, then $\frac{\boldsymbol y(s,t)-\boldsymbol x}{r}$ is a parametrization of $\partial B(\boldsymbol 0,1)$ for $(s,t)\in U$. Finally we observe that
$$\left\|\frac{\partial\boldsymbol y}{\partial s}\times\frac{\partial\boldsymbol y}{\partial t}\right\|=
r^2\left\|\frac{\partial}{\partial s} \left (\frac{\boldsymbol y-\boldsymbol x}{r} \right )\times\frac{\partial }{\partial t} \left (\frac{\boldsymbol y-\boldsymbol x}{r} \right )\right\|$$
So if we let $\boldsymbol z(s,t)=\frac{\boldsymbol y(s,t)-\boldsymbol x}{r}$, then we have
$$\int_U f(\boldsymbol y(s,t))\left\|\frac{\partial\boldsymbol y}{\partial s}\times\frac{\partial\boldsymbol y}{\partial t}\right\|dsdt=
r^2\int_U f(\boldsymbol x +r\boldsymbol z(s,t))\left\|\frac{\partial\boldsymbol z}{\partial s}\times\frac{\partial\boldsymbol z}{\partial t}\right\|dsdt\\=
r^2\int_{\partial B(\boldsymbol 0,1)}f(\boldsymbol x+r\boldsymbol z)dS(\boldsymbol z)$$
Edit by OP
As @user5753974 commented, you can generalize this if you use the fact that in $\mathbb R^n$ $$∫_{∂B(\boldsymbol x,r)}f(\boldsymbol y)dS(\boldsymbol y)=∫_{U}f(\boldsymbol y(\boldsymbol z)) \left \|\det\left (\frac{∂\boldsymbol y}{∂z_1},…,\frac{∂\boldsymbol y}{∂z_{n−1}},\boldsymbol n\right) \right \| d^{n−1}\boldsymbol z,$$ where $\boldsymbol n$ is the normal vector to the surface, and that $\boldsymbol n$ does not change when the surface is scaled and translated.
Best Answer
The following argument in fact uses Fubini's theorem (twice!), but makes simplifying assumptions about the given domain $D$ in the meridian half-plane. I use $(r,z)$ instead of $(y,z)$ as coordinates in this "abstract" plane.
Assume that $D$ is given in the form $$D:=\{(r,z)\ |\ a\leq z\leq b,\ 0<c(z)\leq r\leq d(z)\}\ ,$$ with $a<b$ and $c(\cdot)$, $d(\cdot)$ integrable. Then the three-dimensional body $B$ whose volume we have to compute is given by $$B=\{(x,y,z)\ |\ a\leq z\leq b,\ c(z)\leq \sqrt{x^2+y^2}\leq d(z)\}\ .$$ For given $z\in[a,b]$ let
$$B_z:=\{(x,y)\ |\ c(z)\leq \sqrt{x^2+y^2}\leq d(z)\}$$ denote the intersection of $B$ with the horizontal plane at level $z$. The set $B_z$ is an annulus with inner radius $c(z)$ and outer radius $d(z)$. Therefore its area is given by $${\rm area}(B_z)=\pi\bigl(d^2(z)-c^2(z)\bigr)=2\pi\int_{c(z)}^{d(z)}r \ dr$$ (writing this area as an integral is a trick).
Now comes Fubini: $$\eqalign{{\rm vol}(B)&=\int_B 1\ {\rm d}(x,y,z)=\int_a^b \int_{B_z} 1\ {\rm d}(x,y)\ dz\cr &=\int_a^b {\rm area}(B_z)\ dz =2\pi\int_a^b \int_{c(z)}^{d(z)} r\ dr \ dz \cr &=2\pi \int_D r\ {\rm d}(r,z)\ .\cr}$$