I have an optimization problem, which is to maximize the following integral over the unit sphere: $$ \max_B \int d\Omega \mathbf{f}^{\dagger}(\theta,\phi) (B^{\dagger} + B) \mathbf{f}(\theta,\phi) $$ subject to: $$ B = Q L^{-1} $$ where $Q$ is an orthogonal matrix, and $L$ is a given lower triangular matrix. The vector $\mathbf{f} = (f_1, \dots, f_n) \in \mathbb{C}^n$ is known. So I need to find the orthogonal matrix $Q$ which maximizes the integral.

The functions $f_i(\theta,\phi)$ can be represented as a finite expansion in spherical harmonics, which means that given a $Q$, the integral can be evaluated exactly using Lebedev quadrature. The only method I can think of to solve this, is a Monte-Carlo approach where I generate random orthogonal matrices $Q$, calculate the integral with Lebedev quadrature, and save the $Q$ which maximizes the integral over some large number of $N$ trials. This seems quite inefficient, and may not yield a good $Q$ in the end anyway.

Can anyone think of a better approach to this optimization problem?

  • $\begingroup$ This is not a quadratic optimization problem. It is linear because $Q$ appears linearly in the objective function. It is only a well-posed problem because the set of orthogonal matrices is bounded. $\endgroup$ Aug 30, 2021 at 20:14
  • $\begingroup$ You might find this answer from the math stack exchange helpful -- it looks like whoever posted that was trying to solve a similar problem. $\endgroup$ Aug 30, 2021 at 21:42
  • $\begingroup$ is B dependent on the angular variables? $\endgroup$ Aug 31, 2021 at 15:02

2 Answers 2


Since the $f_i(\theta,\phi)$ are linear combinations of spherical harmonics, we can write $$ \mathbf{f} = F \mathbf{Y} $$ where $\mathbf{Y}$ is a vector of the orthonormalized spherical harmonics - i.e.: $$ \int d\Omega Y_l^m Y_{l'}^{m'*} = \delta_{ll'} \delta_{mm'} $$ So the integral becomes, $$ \int d\Omega \mathbf{f}^{\dagger} (B^{\dagger} + B) \mathbf{f} = \int d\Omega \mathbf{Y}^{\dagger} F^{\dagger} (B^{\dagger} + B) F \mathbf{Y} $$ Note that since $\mathbf{Y}$ is orthonormal with respect to $\int d\Omega$, we have $$ \int d\Omega \mathbf{Y}^{\dagger} A \mathbf{Y} = Tr(A) $$ for any matrix $A$. This can be seen by explicitely writing out the integrand in components and then doing the trivial integrations using the orthonormality condition. So we have \begin{align} \int d\Omega \mathbf{f}^{\dagger} (B^{\dagger} + B) \mathbf{f} &= Tr(F^{\dagger} (B^{\dagger} + B) F) \\ &= 2 Tr(F^{\dagger} B F) \quad note: (Tr(A+A^{\dagger}) = 2 Tr(A)) \\ &= 2 Tr(Q L^{-1} F F^{\dagger}) \end{align} We want to find the $Q$ which maximizes this trace. By the answer given here, the solution is $$ Q = V U^{\dagger} $$ where $L^{-1} F F^{\dagger} = U S V^{\dagger}$ is the singular value decomposition of the matrix multiplying $Q$ in the trace.


Assuming $B$ is not angle-dependent, you know the coefficients of each $f_i$ a priori, and $\mathbf{f}$ can only be evaluated by summing up the spherical harmonics expansion, simplify your problem by making use of the orthogonality of spherical harmonics instead of numerically integrating.

Let $f_i(\theta,\phi)=\sum_{lm}f_{i,lm}Y_l^m(\theta,\phi)$. Then your cost function becomes, using the summation convention: $$ \begin{align*} g &= \int d\Omega f_{i,lm}^*Y_l^{m*}(\theta,\phi) (B_{ji}^* + B_{ij})f_{j,l'm'}Y_{l'}^{m'}(\theta,\phi)\\ &= f_{i,lm}^*f_{j,l'm'}(B_{ji}^* + B_{ij}) \int d\Omega Y_l^{m*}(\theta,\phi) Y_{l'}^{m'}(\theta,\phi)\\ &= f_{i,lm}^*f_{j,l'm'}(B_{ji}^* + B_{ij}) \delta_{ll',mm'}\\ &= f_{i,lm}^*f_{j,lm}(B_{ji}^* + B_{ij}) \end{align*} $$ Assuming $F_{lm}=(f_{1,lm},\ldots,f_{n,lm})$, this can be written as $$ g=\sum_{lm} F_{lm}^\dagger (B^\dagger+B) F_{lm}. $$ We see that for each scalar element of the sum $F_{lm}^\dagger B^\dagger F_{lm}=(F_{lm}^\dagger B F_{lm})^\dagger$, so we can simplify further to $$ g=2\mathfrak{R}\sum_{lm} F_{lm}^\dagger Q L^{-1} F_{lm}. $$ Since you only have control over the elements of $Q$, you can precompute and store as a vector $x_{lm}=L^{-1}F_{lm}$s. Now, you know the $F_{lm}$ and $x_{lm}$ coefficients, and $Q$ is by definition orthogonal. One way to make the sum large is to make $Q$ span the space of the "largest" $F_{lm},x_{lm}$. You could take the $k \leq n$ largest-norm $F_{lm}$ and $x_{lm}$ vectors, put them in a matrix (horzcat in Matlab), take the SVD, and use its column space matrix as your $Q$.

A similar approach can be easily derived in the case that $\mathbf{f}$ is cheap to evaluate (i.e., with fewer than $\mathcal{O}(L^2)$ operations per quadrature point, where $L$ is the maximal degree of spherical harmonics used), with samples of $\Omega$ taking the place of the $l,m$ indices. I realize that taking the real part in the above expression may change your approach slightly; have not yet thought that through completely.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.