Integrating a function over part of a spherical surface

definite integralsintegration

Let's consider the following expressions:

$$
\begin{aligned}
C_{p}(\sigma, \phi) &= 2 \left(\lambda \cos{\left(\sigma \right)} + \nu \sin{\left(\phi \right)} \sin{\left(\sigma \right)} – \omega \sin{\left(\sigma \right)} \cos{\left(\phi \right)}\right)^{2} \\
f(\sigma, \phi) &= C_{p} \sin^{2}{\left(\sigma \right)} \cos{\left(\phi \right)} \\
&= 2 \left(\lambda \cos{\left(\sigma \right)} + \nu \sin{\left(\phi \right)} \sin{\left(\sigma \right)} – \omega \sin{\left(\sigma \right)} \cos{\left(\phi \right)}\right)^{2} \sin^{2}{\left(\sigma \right)} \cos{\left(\phi \right)}
\end{aligned}
$$

where $\lambda, \nu, \omega$ are constants.

I'd like to integrate (analytically) the function $f(\sigma, \phi)$ over some part of a spherical surface. The spherical coordinate system is shown in the following picture (in the following, consider a sphere of unit radius):

Spherical Coordinate System

The following picture shows my geometry:

enter image description here

Here:

  • There is a "clipping" plane that is cutting in half the sphere, with some angle with respect to the x-axis. The intersection between this plane and the sphere is where $C_{p} = 0$.
  • The green arrow indicates the normal vector to the plane. It is placed in that position for rendering purposes, otherwise the clipping plane would have hidden it.
  • P1 is the first point of contact between the clipping plane and the sphere and P2 is the last point of contact. The values $\sigma_{1}$ and $\sigma_{2}$ are known.
  • The blue region indicates a spherical cap, where the limits of integration are $[0, 2\pi]$ in the $\phi$ direction, and $[0, \sigma_{1}]$ in the $\sigma$ direction. Here, I've already computed the analytical solution of $\int_{0}^{\sigma_{1}} \int_{0}^{2\pi} f(\sigma, \phi) d\phi d\sigma$.
  • The yellow-ish region represents the surface where I'm unable to perform the integration.
  • The two arcs (one magenta, the other orange) connecting P1 and P2 represents the two solutions to $C_{p}=0$. Moving from $\sigma_{1}$ to $\sigma_{2}$, the limiting $\phi$ coordinates varies according to:
    $$
    \begin{aligned}
    \phi_{i} &= 2 \tan^{-1}{\left(\frac{-\nu \sin{\sigma} + \sqrt{(\nu^{2} + \omega^{2}) \sin^{2}{\sigma} – \lambda^{2} \cos^{2}{\sigma}}}{\lambda \cos{\sigma} + \omega \sin{\sigma}}\right)} \\
    \phi_{f} &= 2 \tan^{-1}{\left(\frac{-\nu \sin{\sigma} – \sqrt{(\nu^{2} + \omega^{2}) \sin^{2}{\sigma} – \lambda^{2} \cos^{2}{\sigma}}}{\lambda \cos{\sigma} + \omega \sin{\sigma}}\right)}
    \end{aligned}
    $$

Numerical integration over the yellow-ish region produces the correct result. However, it is slow.

Is it possible to compute an analytical solution over this region? How would I set up such integration? Can you provide the analytical solution so that I can verify if I understand the process correctly (I know it will be a very long expression…)?

Best Answer

I'll illustrate the method that can be used with a slightly different function from the one given in the question. For the purpose of this demonstration, I'm going to consider the function

$ f(\sigma, \phi) = C_p \sin(\sigma) \cos(\phi) $

Now, on the unit sphere, with the parameterization given and illustrated in the diagrams, we have

$x = \cos( \sigma)$

$y = \sin (\sigma )\sin (\phi)$

$z = - \sin( \sigma) \cos( \phi)$

Define the unit vector $u = [x, y, z]^T $, then

$C_p = 2 (v^T u)^2$ where

$v = [\lambda, \nu, \omega]^T$

Hence, my function is

$ f = 2 (v^T u)^2 (- z ) $

The surface integral is

$I =\displaystyle \iint_S f dS$

where $dS = \sin(\sigma) d\sigma d\phi$,

Since the region for the surface integral is the hemi-sphere that is pointed to by the vector $v$, then this suggests that we define a rotation of coordinates $r = R r'$ , where $R$ is a $3\times3$ rotation matrix whose columns point in the direction of the new axes. We will take the third column (the $z'$ axis) to be the unit vector along the vector $v$, while the first and second columns are arbitrary, as long as they together with the third column form a rotation matrix.

Written explicitly, the matrix $R = [u_1, u_2, u_3] $ where $u_1, u_2, u_3$ are the columns of matrix $R$, and we have $u_3 = \dfrac{v}{\|v\|} $

Since the $r = R r'$, then $r = x' u_1 + y' u_2 + z' u_3$. Where $r$ is the position vector in the original frame and $r'$ is the position vector in the new frame.

Now,

$v^T u = v^T R u' = (R^T v)^T u' = [0, 0, \|v\|] u' = \|v\| u'_z $

In terms of the new frame, points on the unit sphere are defined in terms of two angle $\theta, \psi$ as follows

$u' = \begin{bmatrix} \sin(\theta) \cos(\psi) \\ \sin(\theta) \sin(\psi) \\ \cos(\theta) \end{bmatrix} $

Therefore, $v^T u = \|v\| \cos(\theta) $

Recalling that the function we want to integrate is

$ f( \sigma, \phi ) = - 2 (v^T u)^2 z$

And from the expansion of $r$ in terms of $r'$, we get

$z = x' u_{13} + y' u_{23} + z' u_{33} = c_1 \sin( \theta) cos( \psi) + c_2 \sin (\theta) \sin(\psi) + c_3 \cos(\theta) $

Where $c_1 = u_{13} , c_2 = u_{23} , c_3 = u_{33} $

Hence, our function to be integrated over the hemisphere surface is

$ f(\theta, \psi) = -2 \|v\|^2 \cos^2(\theta) (c_1 \sin( \theta) \cos( \psi) + c_2 \sin (\theta) \sin(\psi) + c_3 \cos(\theta)) $

And the integral is

$ I = \displaystyle \large -\int_0^{2\pi} \int_0^\frac{\pi}{2} f(\theta, \psi) \sin(\theta) \text{d}\theta \text{d} \phi $

Integrating with respect to $\psi$ first, the terms involving $\cos(\psi)$ and $ \sin(\psi) $ cancel out, while the remaining term is multiplied by $2 \pi$

Hence, the integral is now

$ I = \displaystyle \large -(2\pi) \int_0^\frac{\pi}{2} 2 \|v\|^2 c_3 \cos^3(\theta) \sin(\theta) \text{d}\theta $

And this evaluates to

$ I = \large -\pi c_3 \| v \|^2$

And since $ v = [ \lambda, \nu, \omega ] $ , the above is simply

$ I = - \large \pi \omega ( \lambda^2 + \nu^2 + \omega^2 )^{\frac{1}{2}} $