# Surface integrals in spherical coordinates

If I am given a surface in spherical coordinates $$(r,\theta,\varphi)$$, such that it is parametrised as:

\begin{align} r&=r(\theta,\varphi)\\ \theta&=\theta\\ \varphi&=\varphi \end{align}

What is the area $$S$$ of such surface? Or more specifically, can you show how to get the result: $$S=\int_{0}^{2\pi}\int_{0}^{\pi}\sqrt{r^2+\left(\frac{\partial r}{\partial \theta}\right)^2 + \frac{1}{\sin^2\theta}\left(\frac{\partial r}{\partial \varphi}\right)^2}\;r\sin\theta\;{\rm d}\theta\,{\rm d}\varphi$$

Some definitions that I am using:

$$k$$-surface: Let $$k,N\in\mathbb{N}$$, $$k, $$M\subset \mathbb{R}^N$$ is called a $$k$$-surface, if there exists a non-empty open set $$E\subset \mathbb{R}^k$$ and a map $$\varphi:\mathbb{R}^k\to \mathbb{R}^N$$, such that: (i) $$\varphi(E)=M$$, (ii) $$\varphi\in C^1(E;\mathbb{R}^N)$$, and (iii) the rank of Jacobi matrix of $$\varphi$$ is equal $$k$$ everywhere on $$E$$. The surface is called simple if $$\varphi$$ is also injective on $$E$$ and $$\varphi^{-1}$$ is continuous of $$\varphi(E)$$.

Surface integral of the first kind: Let $$k,N\in\mathbb{N}$$, $$k, $$M\subset \mathbb{R}^N$$ is a simple $$k$$-surface parametrized by the map $$\varphi:\mathbb{R}^k\to \mathbb{R}^N$$, from the open set $$E\subset \mathbb{R}^k$$ and $$f:\mathbb{R}^N\to\mathbb{R}$$ is defined on $$M$$. The surface integral of the first kind is defined by: $$\int_M f\,\mathrm{d}S:=\int_E f(\varphi(t))\sqrt{\det{G(D_\varphi(t))}}\,\mathrm{d}t\,,$$ if the integral on the right exists in the Lebesgue sense and is finite. Here, $$G(A)$$ denotes the Gramm matrix made from columns of $$A$$ and $$D_\varphi$$ is the Jacobi matrix of the map $$\varphi$$. The numeric value of: $$S_k(M):=\int_M f\,\mathrm{d}S\,,$$ is called the $$k$$-dimensional surface area of the $$k$$-surface $$M$$.

Motivation for the question

Now these definitions can be used to calculate e.g. the surface of a unit sphere. When one describes the sphere by the map (omitting one longitudinal line): $$\varphi: (\eta,\psi)\mapsto(\cos\psi\cos\eta,\cos\psi\sin\eta,\sin\psi)$$, where $$(\eta,\psi)\in E=(0,2\pi)\times(-\frac{\pi}{2},\frac{\pi}{2})$$. The Gramm matrix looks like: $$\begin{pmatrix} \cos^2\psi & 0 \\ 0 & 1 \end{pmatrix}$$

one ends up with the (here $$\lambda_2$$ denotes the Lebesgue measure): $$S_2(M) = \int_E 1\sqrt{\det{G(D_\varphi(\eta,\psi))}}\,\mathrm{d}\lambda_2(\eta,\psi)=\int_0^{2\pi}\int_{-\frac{\pi}{2}}^{\frac{\pi}{2}}\cos\psi\,\mathrm{d}\psi\,\mathrm{d}{\eta}=4\pi$$

I did not do anything else but blindly followed the definitions.

The wrong result:

If I do the same approach to the above problem at hand. My map is $$\varphi:\,(\theta,\phi)\mapsto (r(\theta,\phi),\theta,\phi)$$, the Jacobian is then: $$\begin{pmatrix} \frac{\partial r(\theta,\phi)}{\partial \theta} & \frac{\partial r(\theta,\phi)}{\partial \phi} \\ 1 & 0 \\ 0 & 1 \end{pmatrix}$$

The rank of this matrix is 2 as needed for a 2-surface. The Gramm matrix is then: $$G(D_\varphi)= \begin{pmatrix} 1+ \left(\frac{\partial r(\theta,\phi)}{\partial \theta}\right)^2 & \frac{\partial r(\theta,\phi)}{\partial \theta}\frac{\partial r(\theta,\phi)}{\partial \phi} \\ \frac{\partial r(\theta,\phi)}{\partial \theta}\frac{\partial r(\theta,\phi)}{\partial \phi} & 1+ \left(\frac{\partial r(\theta,\phi)}{\partial \phi}\right)^2 \end{pmatrix}$$ with determinant evaluated to: $$\det(G(D_\varphi)) = 1 + \left(\frac{\partial r(\theta,\phi)}{\partial \theta}\right)^2 + \left(\frac{\partial r(\theta,\phi)}{\partial \phi}\right)^2$$

which leads to the surface area:

$$S=\int_{S\subset\mathcal{R}^3}{\rm d}S=\int_{0}^{2\pi}\int_{0}^{\pi}\sqrt{1+\left(\frac{\partial r}{\partial \theta}\right)^2 + \left(\frac{\partial r}{\partial \varphi}\right)^2}\;{\rm d}\theta\,{\rm d}\varphi$$

Which is incorrect. The question is why, and please show the correct way with explanation because as suggested by the example with unit sphere, I did not have to do any transformations and it did work with curvilinear coordinates right away. So the answer provided by Quanto does not address this at all.

Note to the edit:

I have added quite a lot of details since the answer by Quanto, but the root problem is the same - to see by calculation explicitly how one arrives to the correct result and to understand why the my calculation for a sphere works (where i am not using Cartesian coordinates) but fails here.

• The error is very simple. It occurs in the second sentence of “The wrong result” when you write “My map is $\varphi:\,(\theta,\phi)\mapsto (r(\theta,\phi),\theta,\phi)$.” That is not the same map described at the beginning of your post. The map I have quoted is written in Cartesian coordinates and it says that $y=\theta$, $z=\phi$, and $x=r(\theta, \phi)$. That is why the formula you derive is the standard formula for surface area for the graph of a surface in Cartesian coordinates. The map at the beginning of the post, by contrast, is not in Cartesian coordinates. Jan 4, 2021 at 20:27
• Compare this to the map you describe in “Motivation for the question,” which does correctly parametrise the sphere in Cartesian coordinates. When you write “ $\varphi: (\eta,\psi)\mapsto(\cos\psi\cos\eta,\cos\psi\sin\eta,\sin\psi)$,” you are saying $x=\cos\psi\cos\eta$, $y= \cos\psi\sin\eta$, etc. The map in spherical coordinates, by contrast, is just $r=1$. Jan 4, 2021 at 20:31

The map $$\varphi:(\theta,\phi)\rightarrow\Big(f(\theta,\phi),\theta,\phi\Big)$$ is not a parametric representation of the spherical surface $$r=f(\theta,\phi)$$. Rather, $$\varphi$$ is a parametric representation of the surface $$x=f(y,z)$$ since parametric surfaces (and curves) are inherently represented in Cartesian form. What you're doing is equivalent to saying $$\theta \longrightarrow \Big(\theta^2,\theta\Big)$$ is a parametric represention of the polar spiral $$r=\theta^2$$ when, in reality, $$\theta \rightarrow \Big(\theta^2,\theta\Big)$$ is parametric representation of the parabola $$x=y^2$$ while $$\theta \rightarrow \Big(\theta^2 \cos(\theta),\theta^2 \sin(\theta)\Big)$$ is the spiral $$r=\theta^2$$.

Let's refer to your example in which you seek to compute the surface area of a sphere.

In spherical coordinates, the equation of a sphere is $$r=1$$ on the domain $$(\theta,\phi)\in[0,2\pi)\times [0,\pi]$$. You can represent this parametrically as $$(\phi,\theta) \longrightarrow \Big(\sin(\phi)\cos(\theta),\sin(\phi)\sin(\theta),\cos(\phi)\Big)$$ simply by converting from spherical to cartesian coordinates. However, the image of your function $$\varphi:(\phi,\theta)\longrightarrow (1,\phi,\theta)$$ on the same domain is the rectangular patch $$\{1\}\times [0,2\pi)\times [0,\pi]$$ embedded in the vertical plane $$x=1$$ which has area $$2\pi^2$$. This is precisely why $$\int_0^{2\pi}\int_0^{\pi}\sqrt{1+\big(f_{\theta}\big)^2+\big(f_{\phi}\big)^2}d\phi d\theta=2\pi^2 \neq 4\pi$$ whenever $$f(\phi,\theta)=1$$; this integral is calculating the area of the surface $$x=1$$ on $$(y,z)\in [0,2\pi)\times[0,\pi]$$ which is the image of your map $$\varphi$$.

If you want to find the area of the surface given in spherical coordinates by $$r=f(\phi,\theta)$$ defined on domain $$(\phi,\theta)\in \mathcal{U}$$ you will need to express this parametrically as $$\vec{p}(\phi,\theta)=\Big(f(\phi,\theta)\sin(\phi)\cos(\theta),f(\phi,\theta)\sin(\phi)\sin(\theta),f(\phi,\theta)\cos(\phi)\Big)$$ The area of the surface will be $$S=\int \int _{\mathcal{U}}||\vec{p}_\phi \times \vec{p}_\theta||d\phi d\theta$$ If you compute $$||\vec{p}_\phi \times \vec{p}_\theta||$$ you will surely obtain your desired result.

• Upvoted because you use the radial azimuthal polar convention :) Jan 1, 2021 at 20:45

Your surface corresponds to the following idea: $$S:\quad r(\theta,\phi):=R + h(\theta,\phi)\ ,\tag{1}$$ where $$R$$ is the earth radius at sea level, and $$(\theta,\phi)$$ is the height over sea level of the point with geographical coordinates $$(\theta,\phi)$$. You want to know the area of $$S$$. In order to find this area we have to resort to cartesian coordinates, because $$S$$ inherits its area from the euclidean $${\mathbb R}^3$$, and not from the "practical" representation. This means that we have to replace $$(1)$$ by $$S:\quad (\theta,\phi)\mapsto{\bf r}(\theta,\phi)=r(\theta,\phi)\,\bigl(\cos\theta\cos\phi,\cos\theta\sin\phi,\sin\theta\bigr)\ ,\tag{2}$$ where $$(\theta,\phi)\in B:=\bigl[-{\pi\over2},{\pi\over2}\bigr]\times[-\pi,\pi]$$. Note that $$(x,y,z)=\bigl(\cos\theta\cos\phi,\cos\theta\sin\phi,\sin\theta\bigr)$$ is the point on the unit sphere $$S^2$$ having geographical coordinates $$(\theta,\phi)$$.

In calculus 102 it is taught that the area of the surface $$(2)$$ is computed as $${\rm area}(S)=\int_B \bigl|{\bf r}_\theta\times {\bf r}_\phi\bigr|\>{\rm d}(\theta,\phi)=\int_{-\pi/2}^{\pi/2}\int_{-\pi}^\pi \Psi(\theta,\phi)\>d\phi\>d\theta\ ,$$ where $$\Psi(\theta,\phi):=\bigl|{\bf r}_\theta\times {\bf r}_\phi\bigr|_{(\theta,\phi)}$$ has to be computed carefully from $$(2)$$.

In theory one arrives at this expression $$\bigl|{\bf r}_\theta\times {\bf r}_\phi\bigr|_{(\theta,\phi)}$$ by noting that a tiny $$[\theta,\theta+\Delta\theta]\times[\phi,\phi+\Delta\phi]$$ rectangle in the parameter plane is mapped by $$(2)$$ onto a tiny parallelogram of area $$\bigl|{\bf r}_\theta\times {\bf r}_\phi\bigr|_{(\theta,\phi)}\>\Delta\theta\Delta\phi$$.

The standard surface integral, or the first approach, is for the cartesian coordinates, i.e.,

$$S=\int_{S\subset\mathcal{R}^3}{\rm d}S=\int\int \sqrt{1+\left(\frac{\partial z}{\partial x}\right)^2 + \left(\frac{\partial z}{\partial y}\right)^2}\;{\rm d}x\,{\rm d}y$$

It can not be simply recast into one with spherical coordinates.

• If I have a surface $z=\sqrt{x^2+y^2}$ bounded by the interior of a cylinder $x^2+y^2=2x$ I can parametrise the surface $x = 1+r\cos\varphi,\,y=r\sin\varphi,\,z=\sqrt{r^2+2r\cos\varphi+1}$, I get ${\rm d} S=\sqrt{2}r{\rm d}r\,{\rm d}\varphi$ using the first approach quite easily, no additional jacobian or Lame's parameters.. Sep 5, 2019 at 16:27
• I would think so. You should obtain $\sqrt{2}\pi$ as a result. Sep 5, 2019 at 16:43
• yes, but why can't i use this approach for the case in my question? In both problems, you go into curvilinear coordinates. Sep 5, 2019 at 16:44
• There may be confusion. My point is the standard formula works only for Cartesian, not spherical coordinates. Sep 5, 2019 at 16:49
• Mmm, I am well confused. With the determinant of Gram matrix you end up integrating in curvilinear coordinates anyway, at least I think, which happens exactly in the example i gave in the first comment, that is not Cartesian coordinate system... Sep 5, 2019 at 16:53

There is an elementary way to get the areal element, but it's kind of tedious. Along the surface let $$\vec r=\langle x,y,z\rangle=\langle r\sin\theta\cos\phi,r\sin\theta\sin\phi,r\cos\theta\rangle$$ Then \begin{align}d\vec r&=\left\langle\frac{\partial r}{\partial\theta}\sin\theta\cos\phi+r\cos\theta\cos\phi,\frac{\partial r}{\partial\theta}\sin\theta\sin\phi+r\cos\theta\sin\phi,\frac{\partial r}{\partial\theta}\cos\theta-r\sin\theta\right\rangle\,d\theta\\ &\quad+\left\langle\frac{\partial r}{\partial\phi}\sin\theta\cos\phi-r\sin\theta\sin\phi,\frac{\partial r}{\partial\phi}\sin\theta\sin\phi+r\sin\theta\cos\phi,\frac{\partial r}{\partial\phi}\cos\theta\right\rangle\,d\phi\end{align} This writes an infinitesimal displacement along the surface in terms of changes in the independent variables. We find a vector normal to the surface whose magnitude is the area of a parallelogram two of whose sides are the two infinitesimal vectors above: \begin{align}d^2\vec A&=\pm\left\langle\frac{\partial r}{\partial\theta}\sin\theta\cos\phi+r\cos\theta\cos\phi,\frac{\partial r}{\partial\theta}\sin\theta\sin\phi+r\cos\theta\sin\phi,\frac{\partial r}{\partial\theta}\cos\theta-r\sin\theta\right\rangle\,d\theta\\ &\quad\times\left\langle\frac{\partial r}{\partial\phi}\sin\theta\cos\phi-r\sin\theta\sin\phi,\frac{\partial r}{\partial\phi}\sin\theta\sin\phi+r\sin\theta\cos\phi,\frac{\partial r}{\partial\phi}\cos\theta\right\rangle\,d\phi\\ &=\pm\left\langle\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\sin\phi+r\frac{\partial r}{\partial\phi}\cos^2\theta\sin\phi\right.\\ &\quad-\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\sin\phi+r\frac{\partial r}{\partial\phi}\sin^2\theta\sin\phi-r\frac{\partial r}{\partial\theta}\sin\theta\cos\theta\cos\phi+r^2\sin^2\theta\cos\phi,\\ &\quad\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\cos\phi-r\frac{\partial r}{\partial\phi}\sin^2\theta\cos\phi-r\frac{\partial r}{\partial\theta}\sin\theta\cos\theta\sin\phi+r^2\sin^2\theta\sin\phi\\ &\quad-\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\cos\phi-r\frac{\partial r}{\partial\phi}\cos^2\theta\cos\phi,\\ &\quad\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin^2\theta\sin\phi\cos\phi+r\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\sin\phi\cos\phi+r\frac{\partial r}{\partial\theta}\sin^2\theta\cos^2\phi\\ &\quad+r^2\sin\theta\cos\theta\cos^2\phi-\frac{\partial r}{\partial\theta}\frac{\partial r}{\partial\phi}\sin^2\theta\sin\phi\cos\phi-r\frac{\partial r}{\partial\phi}\sin\theta\cos\theta\sin\phi\cos\phi\\ &\quad\left.+r\frac{\partial r}{\partial\theta}\sin^2\theta\sin^2\phi+r^2\sin\theta\cos\theta\sin^2\phi\right\rangle\,d\theta\,d\phi\\ &=\pm\left\langle r\frac{\partial r}{\partial\phi}\sin\phi-r\frac{\partial r}{\partial\theta}\sin\theta\cos\theta\cos\phi+r^2\sin^2\theta\cos\phi,\right.\\ &\quad-r\frac{\partial r}{\partial\phi}\cos\phi-r\frac{\partial r}{\partial\theta}\sin\theta\cos\theta\sin\phi+r^2\sin^2\theta\sin\phi,\\ &\quad\left.r\frac{\partial r}{\partial\theta}\sin^2\theta+r^2\sin\theta\cos\theta\right\rangle\,d\theta\,d\phi\end{align} So now \begin{align}d^2A&=\left\lVert d^2\vec A\right\rVert=\left\{r^2\left(\frac{\partial r}{\partial\phi}\right)^2\sin^2\phi+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^2\theta\cos^2\theta\cos^2\phi+r^4\sin^4\theta\cos^2\phi\right.\\ &\quad+r^2\left(\frac{\partial r}{\partial\phi}\right)^2\cos^2\phi+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^2\theta\cos^2\theta\sin^2\phi+r^4\sin^4\theta\sin^2\phi\\ &\quad-2r^2\left(\frac{\partial r}{\partial\phi}\right)\left(\frac{\partial r}{\partial\theta}\right)\sin\theta\cos\theta\sin\phi\cos\phi+2r^3\left(\frac{\partial r}{\partial\phi}\right)\sin^2\theta\sin\phi\cos\phi\\ &\quad-2r^3\left(\frac{\partial r}{\partial\theta}\right)\sin^3\theta\cos\theta\cos^2\phi+2r^2\left(\frac{\partial r}{\partial\phi}\right)\left(\frac{\partial r}{\partial\theta}\right)\sin\theta\cos\theta\sin\phi\cos\phi\\ &\quad-2r^3\left(\frac{\partial r}{\partial\phi}\right)\sin^2\theta\sin\phi\cos\phi-2r^3\left(\frac{\partial r}{\partial\theta}\right)\sin^3\theta\cos\theta\sin^2\phi\\ &\quad\left.+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^4\theta+2r^3\left(\frac{\partial r}{\partial\theta}\right)\sin^3\theta\cos\theta+r^4\sin^2\theta\cos^2\theta\right\}^{1/2}\,d\theta\,d\phi\\ &=\left\{r^2\left(\frac{\partial r}{\partial\phi}\right)^2+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^2\theta\cos^2\theta+r^4\sin^4\theta-2r^3\left(\frac{\partial r}{\partial\theta}\right)\sin^3\theta\cos\theta\right.\\ &\quad\left.+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^4\theta+2r^3\left(\frac{\partial r}{\partial\theta}\right)\sin^3\theta\cos\theta+r^4\sin^2\theta\cos^2\theta\right\}^{1/2}\,d\theta\,d\phi\\ &=\left\{r^2\left(\frac{\partial r}{\partial\phi}\right)^2+r^2\left(\frac{\partial r}{\partial\theta}\right)^2\sin^2\theta+r^4\sin^2\theta\right\}^{1/2}\,d\theta\,d\phi\\ &=\sqrt{\frac1{\sin^2\theta}\left(\frac{\partial r}{\partial\phi}\right)^2+\left(\frac{\partial r}{\partial\theta}\right)^2+r^2}\,r\sin\theta\,d\theta\,d\phi\end{align} A lot of typing but we have the correct answer if the typo will. The areal element given in the original question was wrong at least because quantities with different units ($$\text{length}^2$$ and no units) were being added under the square root.

EDIT: There is a much more efficient way to organize the arithmetic. Since $$\vec r=\langle r\sin\theta\cos\phi,r\sin\theta\sin\phi,r\cos\theta\rangle$$ we have \begin{align}d\vec r&=\langle\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta\rangle\,dr+\langle r\cos\theta\cos\phi,r\cos\theta\sin\phi,-r\sin\theta\rangle\,d\theta\\ &\quad+\langle-r\sin\theta\sin\phi,r\sin\theta\cos\phi,0\rangle\,d\phi\\ &=\langle\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta\rangle\,dr+\langle \cos\theta\cos\phi,\cos\theta\sin\phi,-\sin\theta\rangle\,rd\theta\\ &\quad+\langle-\sin\phi,\cos\phi,0\rangle\,r\sin\theta d\phi\\ &=\hat rdr+\hat{\theta}rd\theta+\hat{\phi}r\sin\theta d\phi\end{align} Since we are considering $$\theta$$ and $$\phi$$ to be the independent variables, $$dr=\left(\frac{\partial r}{\partial\theta}\right)d\theta+\left(\frac{\partial r}{\partial\phi}\right)d\phi=\frac1r\left(\frac{\partial r}{\partial\theta}\right)rd\theta+\frac1{r\sin\theta}\left(\frac{\partial r}{\partial\phi}\right)r\sin\theta d\phi$$ So along the surface $$d\vec r=\left(\frac1r\left(\frac{\partial r}{\partial\theta}\right)\hat r+\hat{\theta}\right)rd\theta+\left(\frac1{r\sin\theta}\left(\frac{\partial r}{\partial\phi}\right)\hat r+\hat{\phi}\right)r\sin\theta d\phi$$ Then since we can show that $$\hat r\times\hat{\theta}=\hat{\phi}$$, $$\hat{\theta}\times\hat{\phi}=\hat r$$, and $$\hat{\phi}\times\hat r=\hat{\theta}$$, \begin{align}d^2\vec A&=\pm\left(\frac1r\left(\frac{\partial r}{\partial\theta}\right)\hat r+\hat{\theta}\right)rd\theta\times\left(\frac1{r\sin\theta}\left(\frac{\partial r}{\partial\phi}\right)\hat r+\hat{\phi}\right)r\sin\theta d\phi\\ &=\pm\left(-\frac1r\frac{\partial r}{\partial\theta}\hat{\theta}-\frac1{r\sin\theta}\frac{\partial r}{\partial\phi}\hat{\phi}+\hat r\right)r^2\sin\theta d\theta d\phi\end{align} Since we can also show that $$\hat r\cdot\hat r=\hat{\theta}\cdot\hat{\theta}=\hat{\phi}\cdot\hat{\phi}=1$$ and that $$\hat r\cdot\hat{\theta}=\hat{\theta}\cdot\hat{\phi}=\hat{\phi}\cdot\hat r=0$$ we have \begin{align}d^2A&=\left\lVert d^2\vec A\right\lVert\\ &=\sqrt{\left(-\frac1r\frac{\partial r}{\partial\theta}\hat{\theta}-\frac1{r\sin\theta}\frac{\partial r}{\partial\phi}\hat{\phi}+\hat r\right)\cdot\left(-\frac1r\frac{\partial r}{\partial\theta}\hat{\theta}-\frac1{r\sin\theta}\frac{\partial r}{\partial\phi}\hat{\phi}+\hat r\right)}\,r^2\sin\theta d\theta d\phi\\ &=\sqrt{\frac1{r^2}\left(\frac{\partial r}{\partial\theta}\right)^2+\frac1{r^2\sin^2\theta}\left(\frac{\partial r}{\partial\phi}\right)^2+1}\,r^2\sin\theta d\theta d\phi\end{align} Much nicer :)

In the case where the surface can be expressed in both spherical and Cartesian forms the result can be obtained by transforming the Cartesian integral into spherical coordinates. The theorem below gives a proof of this in terms of surface integrals - to obtain the result for surface areas choose an integrand function $$f = 1$$.

The conventions used below are the same as described in this answer to the current question and in this answer to How to prove that the area of a shape is independent of the choice of axes?. The term 'HC' means that the boundary of a set in $$\mathbb{R}^p$$ has 'zero content' [ERA, Definition 24.12, p325].

## Theorem

Let $$S$$ be a surface in $$\mathbb{R}^3 \setminus \{ \mbox{z-axis} \}$$ given in spherical coordinates by $$S = \{ (r(\theta, \varphi), \theta, \varphi) : (\theta, \varphi) \in G \}$$ where $$G \subseteq (0, \pi) \times [0, 2\pi)$$ is open in $$\mathbb{R}^2$$ and $$r : G \rightarrow \mathbb{R}$$ is a positive-valued $$C^1$$ function. Let $$f : S_D \rightarrow \mathbb{R}$$ be a continuous function where $$S_D \subseteq S$$ is defined by $$S_D = \{ (r(\theta, \varphi), \theta, \varphi) : (\theta, \varphi) \in D \}$$, where $$D \subseteq G$$ is compact HC, and let $$I$$ be the following integral : $$I = \int_D f(r(\theta, \varphi), \theta, \varphi) \; \sqrt{r^2 + \left( \frac{\partial r}{\partial \theta} \right)^2 + \frac{1}{\sin^2 \theta}\left( \frac{\partial r}{\partial \varphi} \right)^2} r \sin \theta$$ Let $$CS'$$ denote the associated Cartesian $$xyz$$-coordinate system and suppose that if $$S$$ is expressed in $$CS'$$ coordinates as $$S' \subseteq \mathbb{R}^3$$ then :

(i) the surface $$S'$$ has 1-1 projection onto the $$xy$$-plane of $$CS'$$.

(ii) the tangent plane to $$S'$$ is never vertical within $$CS'$$, ie. never perpendicular to the $$xy$$-plane of $$CS'$$.

Then $$S' = \{ (x, y, s(x, y) : (x, y) \in G' \}$$, where $$G'$$ is open in $$\mathbb{R}^2$$ and $$s : G' \rightarrow \mathbb{R}$$ is $$C^1$$, and if $$f'$$ denotes $$f$$ wrt $$CS'$$ then $$\mathrm{dom} f' = S'_{D'} = \{ (x, y, s(x, y) : (x, y) \in D' \}$$, where $$D' \subseteq G'$$ is compact HC and the surface integral $$I'$$ of $$f'$$ over $$S'_{D'}$$ given by : $$I' = \int_{S'_{D'}} f' = \int_{D'} f'(x, y, s(x, y)) \sqrt{\left(\frac{\partial s}{\partial x}\right)^2 + \left(\frac{\partial s}{\partial y}\right)^2 + 1}$$ is equal to the integral $$I$$.

## Proof

Let $$\alpha$$ be the transform from spherical to Cartesian coordinates, ie. $$\alpha(r, \theta, \varphi) = (r \sin \theta \cos \varphi, r \sin \theta \sin \varphi, r \cos \theta)$$ so $$\alpha : (0, \infty) \times (0, \pi) \times [0, 2\pi) \rightarrow \mathbb{R}^3 \setminus \{ \mbox{z-axis} \}$$ is a bijective $$C^1$$ function with Jacobian $$J_{\alpha} = [ \begin{array}{ccc} \underline{\mathbf{\hat{r}}} & r \, \underline{\mathbf{\hat{\theta}}} & r \sin \theta \, \underline{\mathbf{\hat{\varphi}}} \end{array} ]$$ and $$J_{\alpha}^{-1} = \left[ \begin{array}{c} \underline{\mathbf{\hat{r}}} \\ \underline{\mathbf{\hat{\theta}}} / r \\ \underline{\mathbf{\hat{\varphi}}} / r \sin \theta \end{array} \right]$$.

We firstly show that the tangent plane to $$S$$ at $$P(r, \theta, \varphi)$$ is vertical wrt $$xy$$-plane iff $$\displaystyle r \cos \theta + \frac{\partial r}{\partial \theta} \sin \theta = 0$$. In this answer to the current question it is shown that a normal to $$S$$ is $$\underline{\mathbf{n}} = \nabla k \cdot J_{\alpha}^{-1}$$ where $$k(r, \theta, \varphi) = r - r(\theta, \varphi)$$. We show that the tangent plane at a point $$P$$ on $$S$$ can never be the $$\varphi$$-plane at $$P$$ : if $$\underline{\mathbf{u}} = \cos \varphi \: \underline{\mathbf{i}} + \sin \varphi \: \underline{\mathbf{j}}$$ then $$\underline{\mathbf{n}} \cdot \underline{\mathbf{u}} = \nabla k \cdot J_{\alpha}^{-1} \cdot \underline{\mathbf{u}} = \nabla k \cdot \left[ \begin{array}{c} \underline{\mathbf{\hat{r}}} \cdot \underline{\mathbf{u}} \\ \underline{\mathbf{\hat{\theta}}} \cdot \underline{\mathbf{u}} / r \\ 0 \end{array} \right] = \displaystyle \sin \theta - \frac{1}{r} \frac{\partial r}{\partial \theta} \cos \theta$$ and $$\underline{\mathbf{n}} \cdot \underline{\mathbf{k}} = \nabla k \cdot J_{\alpha}^{-1} \cdot \underline{\mathbf{k}} = \nabla k \cdot \left[ \begin{array}{c} \underline{\mathbf{\hat{r}}} \cdot \underline{\mathbf{k}} \\ \underline{\mathbf{\hat{\theta}}} \cdot \underline{\mathbf{k}} / r \\ 0 \end{array} \right] = \displaystyle \cos \theta + \frac{1}{r} \frac{\partial r}{\partial \theta} \sin \theta$$. But if tangent plane at point $$P$$ on $$S$$ was the $$\varphi$$-plane at $$P$$ then we would have $$0 = (\underline{\mathbf{n}} \cdot \underline{\mathbf{u}})^2 + (\underline{\mathbf{n}} \cdot \underline{\mathbf{k}})^2 = \displaystyle 1 + \frac{1}{r^2} \left(\frac{\partial r}{\partial \theta}\right)^2 > 0$$, a contradiction. Thus in the diagram below, in which the $$\varphi$$-plane is in the plane of the page, the tangent plane at $$P$$ intersects the $$\varphi$$-plane in a line in the page which passes through point $$P$$ and this line is vertical in the diagram (ie. it lies along the line $$PQ$$) iff the tangent plane is vertical wrt $$xy$$-plane. Let $$\theta' = 90^{\circ} - \theta$$, as shown in the diagram. In the $$\varphi$$-plane $$\varphi$$ is fixed and the function $$r = r(\theta, \varphi)$$ traces out a curve in $$S$$ passing through $$P$$ which has a $$\theta\!-\!r$$ gradient at $$P$$ of $$\displaystyle \frac{\partial r}{\partial \theta}$$ and hence a $$\underline{\mathbf{\hat{\theta}}}\!-\!\underline{\mathbf{\hat{r}}}$$ gradient at $$P$$ of $$\displaystyle \frac{1}{r} \frac{\partial r}{\partial \theta}$$, and this is the gradient of the tangent line to the curve at $$P$$ in $$\underline{\mathbf{\hat{\theta}}}\!-\!\underline{\mathbf{\hat{r}}}$$ coordinates. This tangent line lies within the line of intersection of the tangent plane at $$P$$ and the $$\varphi$$-plane and therefore the condition for this line of intersection to lie along $$PQ$$ (ie. the tangent plane is vertical) is that the gradient $$\displaystyle \frac{1}{r} \frac{\partial r}{\partial \theta} = - \tan \theta'$$, ie. $$\displaystyle \frac{1}{r} \frac{\partial r}{\partial \theta} = \frac{-1}{\tan \theta}$$, from which we obtain the relation $$\displaystyle r \cos \theta + \frac{\partial r}{\partial \theta} \sin \theta = 0$$.

By definition we have $$S' = \alpha(S)$$. Define $$G' = \{ (x, y) : \exists \; z \mbox{ such that } (x, y, z) \in S' \}$$. Then by $$S'$$ property (i) we can define $$s(x, y)$$ as the unique such $$z$$, $$\forall \; (x, y) \in G'$$, so $$s : G' \rightarrow \mathbb{R}$$. Then we can write $$S'$$ in two ways as : $$\begin{eqnarray*} S' & = & \{ (x, y, s(x, y)) : (x, y) \in G' \} \\ \mbox{and} \hspace{1em} S' & = & \{ \alpha(r(\theta, \varphi), \theta, \varphi) : (\theta, \varphi) \in G \} \end{eqnarray*}$$ Thus we can define a mapping $$\beta : G \rightarrow G'$$ by $$\beta(\theta, \varphi) =$$ 1st 2 coordinates $$(x, y)$$ of $$\alpha(r(\theta, \varphi), \theta, \varphi)$$. Then we have the following properties for $$\beta, \; \forall \; (\theta, \varphi) \in G$$ :

(a) $$\alpha(r(\theta, \varphi), \theta, \varphi) = (\underline{\mathbf{i}}_2 \cdot \beta(\theta, \varphi), \underline{\mathbf{j}}_2 \cdot \beta(\theta, \varphi), (s \circ \beta)(\theta, \varphi) )$$

(b) $$(s \circ \beta)(\theta, \varphi) = \underline{\mathbf{k}}_3 \cdot \alpha(r(\theta, \varphi), \theta, \varphi)$$

(c) $$\beta(\theta, \varphi) = \left[ \begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 0 \\ \end{array} \right] \alpha(r(\theta, \varphi), \theta, \varphi)$$. Thus $$\beta$$ is a composition of 3 $$C^1$$ functions and hence is itself $$C^1$$. [ERA, 20.9, p235]

(d) $$\beta$$ is a bijection from $$G$$ to $$G'$$ : (i) $$(x, y) \in G' \Rightarrow (x, y, s(x, y)) = \alpha(r(\theta, \varphi), \theta, \varphi)$$ for some $$(\theta, \varphi) \in G$$ and then $$(x, y) = \beta(\theta, \varphi) \therefore \beta$$ is surjective. (ii) $$\beta(\theta_1, \varphi_1) = \beta(\theta_2, \varphi_2) \Rightarrow$$ from (a), $$\alpha(r(\theta_1, \varphi_1), \theta_1, \varphi_1) = \alpha(r(\theta_2, \varphi_2), \theta_2, \varphi_2) \Rightarrow$$ as $$\alpha$$ injective, $$(\theta_1, \varphi_1) = (\theta_2, \varphi_2) \therefore \beta$$ is injective.

(e) From (c) : $$\begin{eqnarray*} \mbox{Jacobian } J_{\beta} & = & \left[ \begin{array}{c@{\hspace{1.2em}}c@{\hspace{1.2em}}c} 1 & 0 & 0 \\[1ex] 0 & 1 & 0 \end{array} \right] \left[ \begin{array}{c@{\hspace{1.2em}}c@{\hspace{1.2em}}c} \underline{\mathbf{\hat{r}}} & r \, \underline{\mathbf{\hat{\theta}}} & r \sin \theta \, \underline{\mathbf{\hat{\varphi}}} \end{array} \right] \left[ \begin{array}{c@{\hspace{1.2em}}c} \frac{\partial r}{\partial \theta} & \frac{\partial r}{\partial \varphi} \\[2ex] 1 & 0 \\[1ex] 0 & 1 \end{array} \right] \\ \Rightarrow J_{\beta} & = & \left[ \begin{array}{c@{\hspace{2em}}c} \cos \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta\right) & \sin \theta\left(\frac{\partial r}{\partial \varphi} \cos \varphi - r \sin \varphi\right) \\[1ex] \sin \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta\right) & \sin \theta \left(\frac{\partial r}{\partial \varphi} \sin \varphi + r \cos \varphi \right) \end{array} \right] \\[2ex] \mbox{and} \hspace{1em} \mathrm{det} \, J_{\beta} & = & \sin \theta \left [ \cos \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta\right)\left(\frac{\partial r}{\partial \varphi} \sin \varphi + r \cos \varphi\right) - \right. \\[2ex] & & \hspace{2em} \left. \sin \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta\right)\left(\frac{\partial r}{\partial \varphi} \cos \varphi - r \sin \varphi\right) \right] \end{eqnarray*}$$ in which the inner terms cancel and the first terms cancel leaving : $$\mathrm{det} \, J_{\beta} = r \sin \theta \, \left(r \cos \theta + \frac{\partial r}{\partial \theta} \, \sin \theta \right)$$

(f) From the above result on tangent planes and property (ii) of $$S'$$ we have $$r \cos \theta + \frac{\partial r}{\partial \theta} \sin \theta \neq 0$$ on $$S$$. Hence from (e), $$\mathrm{det}\,J_{\beta} \neq 0$$ on $$G$$.

(g) $$G' = \beta(G)$$ and (f) imply that $$G'$$ is open in $$\mathbb{R}^2$$. [ERA, 21.9, p256]

(h) Defining $$D' = \beta(D)$$ implies by (f) that $$D'$$ is compact HC. [ERA, Lemma 24.23, p333]

(i) $$\beta^{-1} : G' \rightarrow G$$ is $$C^1$$. [ERA, 21.10, p256]

By definition $$f' = f \circ \alpha^{-1}$$ and $$\mathrm{dom} f' = \{ (x, y, z) : \alpha^{-1}(x, y, z) \in \mathrm{dom} f = S_D \}$$. Thus $$(x, y, z) \in \mathrm{dom} f' \Leftrightarrow (x, y, z) \in \alpha(S_D) \Leftrightarrow (x, y, z) = \alpha(r(\theta, \varphi), \theta, \varphi))$$ for some $$(\theta, \varphi) \in D \Leftrightarrow (x, y) = \beta(\theta, \varphi)$$ and $$z = (s \circ \beta)(\theta, \varphi)$$ for some $$(\theta, \varphi) \in D$$ (from (a) above) $$\Leftrightarrow (x, y) \in D'$$ and $$z = s(x, y)$$ by definition of $$D' \therefore \mathrm{dom} f' = S'_{D'}$$ as required.

To show that $$s : G' \rightarrow \mathbb{R}$$ is $$C^1$$ substitute $$(\theta, \varphi) = \beta^{-1}(x, y)$$ into (b) above, for $$(x, y) \in G'$$, to obtain : $$s(x, y) = \underline{\mathbf{k}}_3 \cdot \alpha(r(\beta^{-1}(x, y)), \; \underline{\mathbf{i}}_2 \cdot \beta^{-1}(x, y), \; \underline{\mathbf{j}}_2 \cdot \beta^{-1}(x, y))$$ $$\therefore$$ since the functions $$\alpha$$, $$r$$, and $$\beta^{-1}$$ are all $$C^1$$ on their respective domains, so $$s$$ as a composition of these is $$C^1$$ on its domain $$G'$$.

Now since $$D' = \beta(D)$$ and $$f' = f \circ \alpha^{-1}$$ we have : $$I' = \int_{\beta(D)} (f \circ \alpha^{-1})(x, y, s(x, y)) \sqrt{\left(\frac{\partial s}{\partial x}\right)^2 + \left(\frac{\partial s}{\partial y}\right)^2 + 1}$$ $$\therefore$$ if we apply the Integral Transform theorem [ERA, 24.26, p335] with the transform function $$\beta$$, so that every $$(x, y) \in \beta(D)$$ within $$I'$$ is substituted with $$\beta(\theta, \varphi)$$ for $$(\theta, \varphi) \in D$$, and that in particular $$(x, y, s(x, y))$$ becomes $$( \underline{\mathbf{i}}_2 \cdot \beta(\theta, \varphi), \underline{\mathbf{j}}_2 \cdot \beta(\theta, \varphi), (s \circ \beta)(\theta, \varphi) )$$, which by (a) above equals $$\alpha(r(\theta, \varphi), \theta, \varphi)$$, then we obtain : $$\begin{eqnarray*} I' & = & {\int_D f(r(\theta, \varphi), \theta, \varphi)\sqrt{\left(\frac{\partial s}{\partial x}(\beta(\theta, \varphi))\right)^2 \cdot \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2 + \left(\frac{\partial s}{\partial y}(\beta(\theta, \varphi))\right)^2 \cdot \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2 + \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2}} \end{eqnarray*}$$ Hence to show $$I' = I$$ it suffices to show that $$\forall\ (\theta, \varphi) \in D$$ : $$\begin{eqnarray} & & \left(\frac{\partial s}{\partial x}(\beta(\theta, \varphi))\right)^2 \cdot \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2 + \left(\frac{\partial s}{\partial y}(\beta(\theta, \varphi))\right)^2 \cdot \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2 + \left(\mathrm{det} \, J_{\beta}(\theta, \varphi)\right)^2 \nonumber \\ & & \hspace{5em} = \left(r^2 + \left(\frac{\partial r}{\partial \theta}\right)^2 + \frac{1}{\sin^2 \theta} \left(\frac{\partial r}{\partial \varphi}\right)^2\right) r^2 \sin^2 \theta \label{suffices} \tag{1} \end{eqnarray}$$ But from (b), applying the Chain Rule [ERA, 20.9, p235] : $$\begin{eqnarray*} & & \left[\begin{array}{c@{\hspace{1.5em}}c} \frac{\partial s}{\partial x}(\beta(\theta, \varphi)) & \frac{\partial s}{\partial y}(\beta(\theta, \varphi)) \end{array} \right] \cdot J_{\beta}(\theta, \varphi) = \mbox{Jacobian of } \underline{\mathbf{k}}_3 \cdot \alpha(r(\theta, \varphi), \theta, \varphi) \\ & & \hspace{2em} = \mbox{3rd row of Jacobian of } \alpha(r(\theta, \varphi), \theta, \varphi) \\ & & \hspace{2em} = \mbox{3rd row of } J_{\alpha} \cdot \left[\begin{array}{c@{\hspace{1.2em}}c} \frac{\partial r}{\partial \theta} & \frac{\partial r}{\partial \varphi} \\ 1 & 0 \\[0.75ex] 0 & 1 \end{array} \right], \hspace{1em} \mbox{again applying Chain Rule} \\ & & \hspace{2em} = \left[ \begin{array}{c@{\hspace{1.5em}}c} \frac{\partial r}{\partial \theta} \cos \theta - r \sin \theta & \frac{\partial r}{\partial \varphi} \cos \theta \end{array} \right] \\ & & \hspace{2em} \mbox{since } J_{\alpha} = \left[ \begin{array}{c@{\hspace{1.5em}}c@{\hspace{1.5em}}c} \sin \theta \cos \varphi & r \cos \theta \cos \varphi & -r \sin \theta \sin \varphi \\[1ex] \sin \theta \sin \varphi & r \cos \theta \sin \varphi & r \sin \theta \cos \varphi \\[1ex] \cos \theta & -r \sin \theta & 0 \end{array} \right] \end{eqnarray*}$$ and thus : $$\begin{eqnarray*} & & \left[ \begin{array}{c@{\hspace{1.5em}}c} \frac{\partial s}{\partial x}(\beta(\theta, \varphi)) \cdot \mathrm{det} \, J_{\beta}(\theta, \varphi) & \frac{\partial s}{\partial y}(\beta(\theta, \varphi)) \cdot \mathrm{det} \, J_{\beta}(\theta, \varphi) \end{array} \right] \\ & & = \left[ \begin{array}{c@{\hspace{1.5em}}c} \frac{\partial r}{\partial \theta} \cos \theta - r \sin \theta & \frac{\partial r}{\partial \varphi} \cos \theta \end{array} \right] \cdot \mathrm{det} \, J_{\beta}(\theta, \varphi) \cdot J_{\beta}(\theta, \varphi)^{-1} \\ & & = \left[ \begin{array}{c@{\hspace{1.5em}}c} \frac{\partial r}{\partial \theta} \cos \theta - r \sin \theta & \frac{\partial r}{\partial \varphi} \cos \theta \end{array} \right] \cdot \left[ \begin{array}{c@{\hspace{1.5em}}c} \sin \theta \left(\frac{\partial r}{\partial \varphi} \sin \varphi + r \cos \varphi \right) & \sin \theta \left(r \sin \varphi - \frac{\partial r}{\partial \varphi} \cos \varphi \right) \\[1ex] -\sin \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta \right) & \cos \varphi \left(\frac{\partial r}{\partial \theta} \sin \theta + r \cos \theta \right) \end{array} \right] \end{eqnarray*}$$ and so : $$\begin{eqnarray*} \frac{\partial s}{\partial x}(\beta(\theta, \varphi)) \cdot \mathrm{det} \, J_{\beta}(\theta, \varphi) & = & r \frac{\partial r}{\partial \theta} \sin \theta \cos \theta \cos \varphi - r \frac{\partial r}{\partial \varphi} \sin \varphi - r^2 \sin^2 \theta \cos \varphi \\[2ex] \mbox{and} \hspace{2em} \frac{\partial s}{\partial y}(\beta(\theta, \varphi)) \cdot \mathrm{det} \, J_{\beta}(\theta, \varphi) & = & r \frac{\partial r}{\partial \theta} \sin \theta \cos \theta \sin \varphi + r \frac{\partial r}{\partial \varphi} \cos \varphi - r^2 \sin^2 \theta \sin \varphi \end{eqnarray*}$$ Hence : $$\begin{eqnarray*} \mbox{LHS (\ref{suffices})} & = & \left( r \frac{\partial r}{\partial \theta} \sin \theta \cos \theta \cos \varphi - r \frac{\partial r}{\partial \varphi} \sin \varphi - r^2 \sin^2 \theta \cos \varphi \right)^2 + \\ & & \left( r \frac{\partial r}{\partial \theta} \sin \theta \cos \theta \sin \varphi + r \frac{\partial r}{\partial \varphi} \cos \varphi - r^2 \sin^2 \theta \sin \varphi \right)^2 + \\ & & \left( r^2 \sin \theta \cos \theta + r \frac{\partial r}{\partial \theta} \sin^2 \theta \right)^2 \end{eqnarray*}$$ ie. a degree 4 polynomial in $$r$$, which we require to equal RHS (\ref{suffices}), ie. $$r^4 \sin^2 \theta + r^2 \left[ \left( \frac{\partial r}{\partial \varphi} \right)^2 + \sin^2 \theta \left( \frac{\partial r}{\partial \theta} \right)^2 \right]$$ Collecting together terms we obtain : $$\begin{eqnarray*} r^4 \mbox{ term} & = & r^4 \sin^4 \theta \cos^2 \varphi + r^4 \sin^4 \theta \sin^2 \varphi + r^4 \sin^2 \theta \cos^2 \theta \\ & = & r^4 \sin^4 \theta + r^4 \sin^2 \theta \cos^2 \theta \\ & = & r^4 \sin^2 \theta, \hspace{1em} \mbox{as required} \\ r^3 \mbox{ term} & = & -\frac{\partial r}{\partial \theta} \sin^3 \theta \cos \theta \cos^2 \varphi + \frac{\partial r}{\partial \varphi} \sin^2 \theta \sin \varphi \cos \varphi + \\ & & -\frac{\partial r}{\partial \theta} \sin^3 \theta \cos \theta \cos^2 \varphi + \frac{\partial r}{\partial \varphi} \sin^2 \theta \sin \varphi \cos \varphi + \\ & & -\frac{\partial r}{\partial \theta} \sin^3 \theta \cos \theta \sin^2 \varphi - \frac{\partial r}{\partial \varphi} \sin^2 \theta \sin \varphi \cos \varphi + \\ & & -\frac{\partial r}{\partial \theta} \sin^3 \theta \cos \theta \sin^2 \varphi - \frac{\partial r}{\partial \varphi} \sin^2 \theta \sin \varphi \cos \varphi + \\ & & \mbox{} 2\frac{\partial r}{\partial \theta} \sin^3 \theta \cos \theta = 0, \hspace{1em} \mbox{as required} \\ r^2 \mbox{ term} & = & \left( \frac{\partial r}{\partial \theta} \right)^2 \sin^2 \theta \cos^2 \theta \cos^2 \varphi - 2\left( \frac{\partial r}{\partial \theta} \right) \left( \frac{\partial r}{\partial \varphi} \right) \sin \theta \cos \theta \sin \varphi \cos \varphi + \\ & & \left( \frac{\partial r}{\partial \varphi} \right)^2 \sin^2 \varphi + \left( \frac{\partial r}{\partial \theta} \right)^2 \sin^2 \theta \cos^2 \theta \sin^2 \varphi + \\ & & 2\left( \frac{\partial r}{\partial \theta} \right) \left( \frac{\partial r}{\partial \varphi} \right) \sin \theta \cos \theta \sin \varphi \cos \varphi + \left( \frac{\partial r}{\partial \varphi} \right)^2 \cos^2 \varphi + \\ & & \left( \frac{\partial r}{\partial \theta} \right)^2 \sin^4 \theta = \left( \frac{\partial r}{\partial \varphi} \right)^2 + \left( \frac{\partial r}{\partial \theta} \right)^2 \sin^2 \theta, \hspace{1em} \mbox{as required} \end{eqnarray*}$$ and the $$r$$ and constant terms are $$0$$, as required.

References

[SPH] Spherical coordinate system, https://en.wikipedia.org/wiki/Spherical_coordinate_system

[ERA] Robert G. Bartle (1964), The Elements of Real Analysis, John Wiley & Sons.

The result can be derived from the theorem below on surface integrals, by setting the integrand function $$f$$ to unity.

I have used the following conventions below :

1. The physics convention $$(r, \theta, \varphi)$$ for spherical coordinates, with $$\theta$$ the polar angle and $$\varphi$$ the azimuthal angle [SPH].

2. The local spherical coordinate system basis vectors $$\underline{\mathbf{\hat{r}}}, \, \underline{\mathbf{\hat{\theta}}}, \, \underline{\mathbf{\hat{\varphi}}}$$ given by $$\underline{\mathbf{\hat{r}}} = \sin \theta \cos \varphi \: \underline{\mathbf{i}} + \sin \theta \sin \varphi \: \underline{\mathbf{j}} + \cos \theta \: \underline{\mathbf{k}}$$, $$\underline{\mathbf{\hat{\theta}}} = \cos \theta \cos \varphi \: \underline{\mathbf{i}} + \cos \theta \sin \varphi \: \underline{\mathbf{j}} - \sin \theta \: \underline{\mathbf{k}}$$, and $$\underline{\mathbf{\hat{\varphi}}} = -\sin \varphi \: \underline{\mathbf{i}} + \cos \varphi \: \underline{\mathbf{j}}$$.

3. The definition of the Riemann integral on $$\mathbb{R}^p$$ in [ERA, p320], which is defined natively on $$\mathbb{R}^p$$ and written $$\int_D f$$. Below we use the notation of form $$\int_D f(x, y, z)$$ to indicate the notation used for points in $$D$$, but without the $$dx\,dy\,dz$$ to avoid confusion with iterated integrals which are not used below.

4. $$C^1$$ applied to a function means it has continuous first partial derivatives.

5. The term 'HC' applied to a set in $$\mathbb{R}^p$$ means the set 'has content' as defined in [ERA, Definition 24.12, 325]. This means the boundary of $$D$$ has 'zero content', as intuitively would a line segment in $$\mathbb{R}^2$$, or a plane in $$\mathbb{R}^3$$. Many regions encountered in practice have this property.

6. A product partition of a compact interval in $$\mathbb{R}^p$$ is a partition formed from individual partitions on its $$p$$ component intervals.

Further relevant results from [ERA] are described in more detail at the start of this answer to How to prove that the area of a shape is independent of the choice of axes?.

## Theorem

Let $$S$$ be a surface in $$\mathbb{R}^3 \setminus \{ \mbox{z-axis} \}$$ given in spherical coordinates by $$S = \{ (r(\theta, \varphi), \theta, \varphi) : (\theta, \varphi) \in G \}$$ where $$G \subseteq (0, \pi) \times [0, 2\pi)$$ is open in $$\mathbb{R}^2$$ and $$r : G \rightarrow \mathbb{R}$$ is a positive-valued $$C^1$$ function. Then the surface integral of continuous function $$f : S_D \rightarrow \mathbb{R}$$ over $$S_D = \{ (r(\theta, \varphi), \theta, \varphi) : (\theta, \varphi) \in D \}$$, where $$D \subseteq G$$ is compact HC, is given by : $$\int_{S_D} f = \int_D f(r(\theta, \varphi), \theta, \varphi) \sqrt{r^2 + \left(\frac{\partial r}{\partial \theta}\right)^2 + \frac{1}{\sin^2 \theta} \left(\frac{\partial r}{\partial \varphi}\right)^2 } \; r \sin \theta$$

## Proof

We first of all show that at any point on $$S$$ the tangent plane is well-defined with a non-zero normal vector $$\underline{\mathbf{n}}$$ which is never perpendicular to the vector $$\underline{\mathbf{\hat{r}}}$$. Define $$C^1$$ function $$k : \mathbb{R}^+ \times G \rightarrow \mathbb{R}$$ by $$k(r, \theta, \varphi) = r - r(\theta, \varphi)$$, so that at all points of $$S$$ $$k$$ is zero. Take a point $$P$$ on $$S$$ and suppose $$\underline{\mathbf{r}} : [a, b] \rightarrow \mathbb{R}^3$$ is a $$C^1$$ curve in Cartesian coordinates on $$S$$ passing through $$P$$, $$\underline{\mathbf{r}}(t_o) = P$$, say. Define function $$l(t) = k(\alpha^{-1}(\underline{\mathbf{r}}(t)))$$, for $$t \in [a, b]$$, where $$\alpha$$ is the transform from spherical to Cartesian coordinates. Then $$l$$ is the zero function. But by the Chain Rule [ERA, 20.9, p235], restricting $$\alpha$$ to $$\mathbb{R}^3$$ minus the $$\varphi = 0$$ half plane (ie. to $$(r, \theta, \varphi) \in (0, \infty) \times (0, \pi) \times (0, 2\pi)$$) so that both $$\alpha$$ and $$\alpha^{-1}$$ are $$C^1$$ functions (note $$\alpha^{-1}$$ is not $$C^1$$ on $$\mathbb{R}^3 \setminus \{ \mbox{z-axis} \}$$ since it is discontinuous where the azimuthal angle $$\varphi$$ wraps around), we have : $$\begin{eqnarray*} l'(t) & = & \nabla k(\alpha^{-1} (\underline{\mathbf{r}}(t))) \cdot J_{\alpha^{-1}}(\underline{\mathbf{r}}(t)) \cdot \underline{\mathbf{r}}'(t), \hspace{1em} \forall t \in [a, b], \\ & & J \mbox{ denoting the Jacobian matrix} \\ \mbox{ie.} \hspace{2em} 0 & = & \nabla k(\alpha^{-1} (\underline{\mathbf{r}}(t))) \cdot J_{\alpha}(\alpha^{-1}(\underline{\mathbf{r}}(t)))^{-1} \cdot \underline{\mathbf{r}}'(t), \hspace{1em} \forall t \in [a, b], \\ & & \mbox{using [ERA, 21.10, p256]} \\ \mbox{and so in particular} \hspace{2em} 0 & = & (\nabla k \cdot J_{\alpha}^{-1}) \cdot \underline{\mathbf{r}}'(t_0), \end{eqnarray*}$$ $$\nabla k$$ and $$J_{\alpha}^{-1}$$ being evaluated at the point $$\alpha^{-1}(\underline{\mathbf{r}}(t_0)) = P(r, \theta, \varphi)$$. Since this applies to all curves $$\underline{\mathbf{r}}$$ in $$S$$ through $$P$$, $$\underline{\mathbf{n}} = \nabla k \cdot J_\alpha^{-1}$$ is then a normal to the tangent plane to $$S$$ at $$P$$, provided $$\underline{\mathbf{n}} \neq \underline{\mathbf{0}}$$. But $$J_{\alpha}^{-1} = \left[ \begin{array}{c} \underline{\mathbf{\hat{r}}} \\ \underline{\mathbf{\hat{\theta}}} / r \\ \underline{\mathbf{\hat{\varphi}}} / r \sin \theta \end{array} \right]$$ (valid everywhere away from $$z$$-axis), and we have $$\displaystyle{\nabla k = \left(1, -\frac{\partial r}{\partial \theta}, -\frac{\partial r}{\partial \varphi} \right)}$$, so $$\underline{\mathbf{n}} \cdot \underline{\mathbf{\hat{r}}} = 1$$. Thus $$\underline{\mathbf{n}} \neq \underline{\mathbf{0}}$$ as required. And we also have the required $$\underline{\mathbf{n}} \cdot \underline{\mathbf{\hat{r}}} \neq 0$$. We had to omit a $$\varphi = 0$$ half plane above, however the above conclusion still remains valid for the whole of $$S$$ as we can repeat the argument omitting a different $$\varphi$$ half plane, $$\varphi = \lambda$$ say, and temporarily changing the range of $$\varphi$$ to $$(\lambda, \lambda + 2\pi)$$.

Let $$I \subseteq \mathbb{R}^2$$ be any compact interval containing $$D$$, as in the ERA-$$\mathbb{R}^p$$ integral definition [ERA, p320] . For each $$n$$ let $$P_n$$ be a product partition of $$I$$, such that $$\| P_n \| \rightarrow 0$$ as $$n \rightarrow \infty$$ (for example define $$P_n$$ by dividing up both component intervals of $$I$$ into $$n$$ equal parts). The subintervals in $$P_n$$ then have form $$J = [\theta, \theta + \Delta \theta] \times [\varphi, \varphi + \Delta \varphi]$$ which defines a spherical segment, $$C$$ say, and the domain $$S_D$$ of the surface integral is then partitioned by these spherical segments. For each $$J \in P_n$$ with $$J \cap D \neq \emptyset$$ choose a tag point $$(\theta_J, \varphi_J) \in J \cap D$$, and for each $$J \in P_n$$ with $$J \cap D = \emptyset$$ choose any tag point $$(\theta_J, \varphi_J) \in J$$. To ensure that for every $$J$$ with $$J \cap D \neq \emptyset$$ we have $$J \subseteq G$$ we can consider $$n$$ sufficiently large that $$\| P_n \| < d$$, where $$d > 0$$ is the minimum distance between the disjoint closed sets $$D$$ and $$\mathbb{R}^3 \setminus G$$ ($$D$$ being bounded) - this then ensures that $$J \subseteq (0, \pi) \times [0, 2\pi)$$ so that $$J$$ as a spherical segment is guaranteed to be a 4-sided figure as shown in the diagram.

When $$J \cap D \neq \emptyset$$ the point $$\underline{\mathbf{v_J}} = (r(\theta_J, \varphi_J), \theta_J, \varphi_J)$$ of $$S_D$$ lies within the spherical segment $$C$$ as shown in the figure. The points $$P, Q, R, U$$ in the figure all lie on a common sphere center $$O$$ which intersects $$C$$ in the spherical element $$PQRU$$. When $$n$$ is large, $$\Delta \theta$$ and $$\Delta \varphi$$ are small and in a reasonable vicinity of $$\underline{\mathbf{v_J}}$$ the spherical segment column $$C$$ can be approximated by a rectangular column $$C_R$$ which is situated on the local basis vectors $$\underline{\mathbf{\hat{r}}}, \, \underline{\mathbf{\hat{\theta}}}, \, \underline{\mathbf{\hat{\varphi}}}$$ (these 'free vectors' being evaluated at $$\underline{\mathbf{v_J}}$$), with $$\underline{\mathbf{\hat{r}}}$$ along the length of the column, $$\underline{\mathbf{\hat{\theta}}}$$ along the sides $$PU$$ and $$QR$$ of the column, and $$\underline{\mathbf{\hat{\varphi}}}$$ along the sides $$PQ$$ and $$UR$$ of the column. Since $$C_R$$ is in the vicinity of $$\underline{\mathbf{v_J}}$$, the width of $$C_R$$ is approximately $$r \, \Delta \theta$$ along the $$\underline{\mathbf{\hat{\theta}}}$$ direction, and approximately $$r \sin \theta \, \Delta \varphi$$ along the $$\underline{\mathbf{\hat{\varphi}}}$$ direction ($$r$$ and $$\theta$$ being evaluated at $$\underline{\mathbf{v_J}}$$).

From the result above the tangent plane $$T_J$$ to $$S$$ at $$\underline{\mathbf{v_J}}$$ cannot be 'vertical' wrt the column $$C_R$$, ie. its normal cannot be perpendicular to $$\underline{\mathbf{\hat{r}}}$$ and so it must intersect all 4 sides of rectangular column $$C_R$$ in a planar element $$\Delta T_J$$ which by Euclid XI-16 must be a parallelogram. Let $$\Delta T_J$$ be contained by the vectors $$\underline{\mathbf{a}}$$ and $$\underline{\mathbf{b}}$$, where $$\underline{\mathbf{a}}$$ lies in the two $$\underline{\mathbf{\hat{\theta}}}$$ faces of $$C_R$$ and has positive $$\underline{\mathbf{\hat{\theta}}}$$ component of $$r \, \Delta \theta$$, and $$\underline{\mathbf{b}}$$ lies in the two $$\underline{\mathbf{\hat{\varphi}}}$$ faces of $$C_R$$ and has positive $$\underline{\mathbf{\hat{\varphi}}}$$ component of $$r \sin \theta \, \Delta \varphi$$.

Thus we have an approximation to the total surface $$S_D$$ by the collection of all the non-overlapping parallelograms $$\Delta T_J$$ ranging over all $$J \in P_n$$ with $$J \cap D \neq \emptyset$$. When $$J$$ is near the boundary of $$D$$ we may not have $$J \subseteq D$$ and then the corresponding $$\Delta T_J$$ will protrude outwith the $$\theta\!-\!\varphi$$ boundary of $$S_D$$. However, intuitively we see that in the limit as $$n \rightarrow \infty$$ the error due to this becomes vanishingly small. Thus we can take the $$n$$th approximation $$I_n$$ to the surface integral $$\int_{S_D} f$$ to be : $$I_n = \sum_{\substack{J \in P_n, \\ J \cap D \neq \emptyset }} f(r(\theta_J, \varphi_J), \theta_J, \varphi_J) \; \Delta T_J$$ where $$\Delta T_J = |\underline{\mathbf{a}} \times \underline{\mathbf{b}}|$$ denotes the area of the parallelogram element corresponding to the partition subinterval $$J \subseteq I$$.

Consider now a plane $$\Gamma$$ through $$\underline{\mathbf{v_J}}$$ which is parallel to the two $$\underline{\mathbf{\hat{\theta}}}$$ faces of $$C_R$$. On $$\Gamma$$, $$\varphi$$ must be constant with value $$\varphi_J$$. Thus on $$\Gamma \cap S$$ we will have $$r = r(\theta, \varphi_J)$$ as a function of $$\theta$$ only. Thus a curve is traced out on $$\Gamma$$ by $$S$$, of $$r$$ versus $$\theta$$, with $$\theta$$ varying about $$\theta_J$$. In $$\Gamma$$ near $$\underline{\mathbf{v_J}}$$ we have the local 'horizontal' coordinate direction $$\underline{\mathbf{\hat{\theta}}}$$ and the local 'vertical' coordinate direction $$\underline{\mathbf{\hat{r}}}$$ (these vectors being evaluated at $$\underline{\mathbf{v_J}}$$). And we also note since $$G$$ is open in $$\mathbb{R}^2$$ that $$\theta_J$$ is a cluster point of $$G$$ in the direction of $$\theta$$ at $$\underline{\mathbf{v_J}}$$. Thus $$\displaystyle{\frac{\partial r}{\partial \theta}}$$ gives the $$\theta\!-\!r$$ gradient of this curve in $$\Gamma$$ at $$\underline{\mathbf{v_J}}$$. However to get the gradient wrt the local coordinate vectors $$\underline{\mathbf{\hat{\theta}}}$$ and $$\underline{\mathbf{\hat{r}}}$$ in $$\Gamma$$ we need to consider the distances along $$\underline{\mathbf{\hat{\theta}}}$$ and $$\underline{\mathbf{\hat{r}}}$$. Thus we must multiply the $$\partial \theta$$ by $$r$$, therefore the $$\underline{\mathbf{\hat{\theta}}}\!-\!\underline{\mathbf{\hat{r}}}$$ gradient of the curve at $$\underline{\mathbf{v_J}}$$ is $$\displaystyle{\frac{1}{r} \, \frac{\partial r}{\partial \theta}}$$ so a tangent vector to the curve at $$\underline{\mathbf{v_J}}$$ is : $$\underline{\mathbf{\hat{\theta}}} + \frac{1}{r} \: \frac{\partial r}{\partial \theta} \;\; \underline{\mathbf{\hat{r}}}$$ and this is then a vector in $$\Gamma \cap T_J$$, as $$T_J$$ by definition contains all tangents to curves in $$S$$ at the point $$\underline{\mathbf{v_J}}$$. But if we now apply Euclid XI-16 to the parallel planes which are $$\Gamma$$ and the two $$\underline{\mathbf{\hat{\theta}}}$$ faces of $$C_R$$ and the plane $$T_J$$ that cuts across them we see that the above tangent vector must be parallel to the vector $$\underline{\mathbf{a}}$$, by the definition of $$\underline{\mathbf{a}}$$. But the $$\underline{\mathbf{\hat{\theta}}}$$ component of $$\underline{\mathbf{a}}$$ is known to be $$r \, \Delta \theta$$, thus : $$\underline{\mathbf{a}} = r \Delta \theta \; \underline{\mathbf{\hat{\theta}}} + \left(\frac{\partial r}{\partial \theta}\right) \, \Delta \theta \; \underline{\mathbf{\hat{r}}}$$

Consider now a plane $$\Delta$$ through $$\underline{\mathbf{v_J}}$$ which is parallel to the two $$\underline{\mathbf{\hat{\varphi}}}$$ faces of $$C_R$$. On $$\Delta$$, $$\theta$$ must be constant with value $$\theta_J$$. Thus on $$\Delta \cap S$$ we will have $$r = r(\theta_J, \varphi)$$ as a function of $$\varphi$$ only. Thus a curve is traced out on $$\Delta$$ by $$S$$, of $$r$$ versus $$\varphi$$, with $$\varphi$$ varying about $$\varphi_J$$. In $$\Delta$$ near $$\underline{\mathbf{v_J}}$$ we have the local 'horizontal' coordinate direction $$\underline{\mathbf{\hat{\varphi}}}$$ and the local 'vertical' coordinate direction $$\underline{\mathbf{\hat{r}}}$$ (these vectors being evaluated at $$\underline{\mathbf{v_J}}$$). And we also note since $$G$$ is open in $$\mathbb{R}^2$$ that $$\varphi_J$$ is a cluster point of $$G$$ in the direction of $$\varphi$$ at $$\underline{\mathbf{v_J}}$$. Thus $$\displaystyle{\frac{\partial r}{\partial \varphi}}$$ gives the $$\varphi\!-\!r$$ gradient of this curve in $$\Delta$$ at $$\underline{\mathbf{v_J}}$$. However to get the gradient wrt local coordinate vectors $$\underline{\mathbf{\hat{\varphi}}}$$ and $$\underline{\mathbf{\hat{r}}}$$ in $$\Delta$$ we need to consider the distances along $$\underline{\mathbf{\hat{\varphi}}}$$ and $$\underline{\mathbf{\hat{r}}}$$. Thus we must multiply the $$\partial \varphi$$ by $$r \sin \theta$$, therefore the $$\underline{\mathbf{\hat{\varphi}}}\!-\!\underline{\mathbf{\hat{r}}}$$ gradient of the curve at $$\underline{\mathbf{v_J}}$$ is $$\displaystyle{\frac{1}{r \sin \theta} \, \frac{\partial r}{\partial \varphi}}$$ (this is where we use the condition $$\theta \in (0, \pi)$$), so a tangent vector to the curve at $$\underline{\mathbf{v_J}}$$ is: $$\underline{\mathbf{\hat{\varphi}}} + \frac{1}{r \sin \theta} \: \frac{\partial r}{\partial \varphi} \;\; \underline{\mathbf{\hat{r}}}$$ and this is then a vector in $$\Delta \cap T_J$$, as $$T_J$$ by definition contains all tangents to curves in $$S$$ at the point $$\underline{\mathbf{v_J}}$$. But if we now apply Euclid XI-16 to the parallel planes which are $$\Delta$$ and the two $$\underline{\mathbf{\hat{\varphi}}}$$ faces of $$C_R$$ and the plane $$T_J$$ that cuts across them we see that the above tangent vector must be parallel to the vector $$\underline{\mathbf{b}}$$, by the definition of $$\underline{\mathbf{b}}$$. But the $$\underline{\mathbf{\hat{\varphi}}}$$ component of $$\underline{\mathbf{b}}$$ is known to be $$r \sin \theta \, \Delta \varphi$$, thus : $$\underline{\mathbf{b}} = r \sin \theta \, \Delta \varphi \; \underline{\mathbf{\hat{\varphi}}} + \left(\frac{\partial r}{\partial \varphi}\right) \, \Delta \varphi \; \underline{\mathbf{\hat{r}}}$$

Now since $$\underline{\mathbf{\hat{r}}} \times \underline{\mathbf{\hat{\theta}}} = \underline{\mathbf{\hat{\varphi}}}$$, $$\underline{\mathbf{\hat{\theta}}} \times \underline{\mathbf{\hat{\varphi}}} = \underline{\mathbf{\hat{r}}}$$, and $$\underline{\mathbf{\hat{\varphi}}} \times \underline{\mathbf{\hat{r}}} = \underline{\mathbf{\hat{\theta}}}$$, we have : $$\underline{\mathbf{a}} \times \underline{\mathbf{b}} = \left( r^2 \sin \theta \; \underline{\mathbf{\hat{r}}} - r \sin \theta \, \frac{\partial r}{\partial \theta} \; \underline{\mathbf{\hat{\theta}}} - r \frac{\partial r}{\partial \varphi} \; \underline{\mathbf{\hat{\varphi}}} \right) \Delta \theta \, \Delta \varphi$$ Thus $$\Delta T_J = | \underline{\mathbf{a}} \times \underline{\mathbf{b}} | = \sqrt{r^2 + \left( \frac{\partial r}{\partial \theta} \right)^2 + \frac{1}{\sin^2 \theta} \left(\frac{\partial r}{\partial \varphi}\right)^2 } \; \; r \sin \theta \; \Delta \theta \, \Delta \varphi$$ where all the terms are evaluated at the point $$\underline{\mathbf{v_J}}$$. In the terminology of [ERA, p316], $$\Delta \theta \, \Delta \varphi$$ is the content $$A(J)$$ of subinterval $$J$$, so we can write $$\Delta T_J = g(\theta_J, \varphi_J) A(J)$$ where $$g : G \rightarrow \mathbb{R}$$ is continuous and we then have : $$I_n = \sum_{\substack{J \in P_n, \\ J \cap D \neq \emptyset }} f(r(\theta_J, \varphi_J), \theta_J, \varphi_J) \, g(\theta_J, \varphi_J) \, A(J) = \sum_{\substack{J \in P_n, \\ J \cap D \neq \emptyset }} h(\theta_J, \varphi_J) \, A(J)$$ where $$h : D \rightarrow \mathbb{R}$$ is continuous.

As in the ERA-$$\mathbb{R}^p$$ integral definition [ERA, p320] we now extend $$h$$ by zero to the bounding box $$I \supseteq D$$. Then for any $$J \in P_n$$ with $$J \cap D = \emptyset$$ our tag $$(\theta_J, \varphi_J) \in J$$ will satisfy $$h(\theta_J, \varphi_J) = 0$$, thus we have : $$I_n = \sum_{J \in P_n} h(\theta_J, \varphi_j) \, A(J)$$ ie. $$I_n$$ is a Riemann sum $$S(P_n; h)$$ for $$\int_D h$$, [ERA, p320]. But since $$\| P_n \| \rightarrow 0$$ we have $$S(P_n; h) \rightarrow \int_D h$$ (that result follows from the equivalence of partition refinement based convergence and partition norm based convergence for the ERA-$$\mathbb{R}^p$$ integral). Thus $$\int_{S_D} f = \lim I_n = \int_D h$$ as required.

Note that $$(I_n)$$ converges to the same value whatever choices we make for (i) bounding box $$I \supseteq D$$, (ii) partitions $$P_n$$ of $$I$$ satisfying $$\| P_n \| \rightarrow 0$$, and (iii) choice of tags $$(\theta_J, \varphi_J)$$ for the subintervals $$J \in P_n$$, which represents quite a considerable variation in the approximating sequences $$(I_n)$$ to $$\int_{S_D} f$$. Furthermore the existence of the integral $$\int_D h$$ is assured from the 2nd integrability theorem [ERA, 24.11, p324] as $$h$$ is continuous and $$D$$ is compact HC.

For a surface expressible in both spherical and Cartesian coordinates it is possible to obtain the above spherical formula for the surface integral from the corresponding Cartesian formula by transforming the integral [ERA, 24.26, p335]. This provides confirmation that the additional approximations that are made in the spherical case (ie. approximating the spherical segment column $$C$$ by the rectangular column $$C_R$$) still result in the correct value in the limit. If $$\gamma$$ is the angle between the surface normal $$\underline{\mathbf{n}}$$ and $$\underline{\mathbf{\hat{r}}}$$ then from above we have $$\cos \gamma = 1 / |\underline{\mathbf{n}}|$$, with $$\underline{\mathbf{n}}$$ continuous for $$(\theta, \varphi)$$ in the compact set $$D$$ and thus $$\cos \gamma$$ attains a maximum value on $$D$$ so $$\gamma$$ has a fixed limit below $$90^{\circ}$$, ie. there is a fixed limit to the 'steepness' of the tangent plane wrt $$\underline{\mathbf{\hat{r}}}$$, which is a main cause of the error in approximating $$C$$ by $$C_R$$. Then as $$\Delta \theta$$ and $$\Delta \varphi$$ get less and less that error effect tends to zero.

References

[SPH] Spherical coordinate system, https://en.wikipedia.org/wiki/Spherical_coordinate_system [ERA] Robert G. Bartle (1964), The Elements of Real Analysis, John Wiley & Sons.