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):
The following picture shows my geometry:
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}} $