1

Let $y \in \mathbb{R}^n$ be a given point, and $\alpha \in \mathbb{R}^n_+$ and $r > 0$ be given parameters. What is an efficient way to project $y$ onto the set $$S := \left\lbrace x \in \mathbb{R}^n \: : \: x \geq 0, \sum_i x_i = 1, \sum_i \alpha_i x^2_i \leq r \right\rbrace$$ that is the intersection of a (scaled) Euclidean ball and a simplex (assuming $S \neq \emptyset$)?

I can think of three plausible approaches:

  1. Solve the projection problem directly using a convex QCQP solver (such as CPLEX/Gurobi)
  2. Use the method of alternating projections by projecting alternately on the simplex and the scaled Euclidean ball, but I am not sure that this will yield the projection even in the limit
  3. Following the approach outlined here, we can write the Lagrangian as $$L(x,\mu,\lambda) := \frac{1}{2} \left\lVert x - y \right\rVert^2 + \mu \left( \sum_i x_i - 1 \right) + \lambda \left(\sum_i \alpha_i x^2_i - r \right)$$ and the dual function as $$g(x,\mu,\lambda) := \min_{x \geq 0} L(x,\mu,\lambda).$$ The solution to the dual problem is then $$x^*_i = \left(\frac{y_i - \mu}{1+2\alpha_i\lambda}\right)_+, \quad \forall i \in \{1,\cdots,n\}.$$ The KKT conditions yield $$\sum_i \left(\frac{y_i - \mu}{1+2\alpha_i\lambda}\right)_+ = 1, \quad \lambda \geq 0, \quad \lambda \left(\sum_i \alpha_i (x^*_i)^2 - r \right) = 0, \quad \sum_i \alpha_i (x^*_i)^2 \leq r.$$ We can first check if $\lambda = 0$ satisfies the KKT conditions; if not, we could use bisection over $\lambda$ (with a large upper bound for $\lambda$), which fixes the value of $\lambda$ at a positive value and tries to solve the nonlinear system of equations $$\sum_i \left(\frac{y_i - \mu}{1+2\alpha_i\lambda}\right)_+ = 1, \quad \sum_i \alpha_i \left(\frac{y_i - \mu}{1+2\alpha_i\lambda}\right)^2 - r = 0,$$ for $\mu$ using bisection or the approach suggested here.

Is there a faster approach that can exploit the problem structure? I wish to solve this projection problem on the order of a million times within my algorithm, so speed is critical.

EDIT: The dimension $n$ is around $1000$. For $n = 1000$, Gurobi takes $\approx 0.1$ seconds on my laptop for each projection step. I'm hoping to reduce the time at least by an order of magnitude.

ProAmateur
  • 1,828
  • How large is $n$? – littleO Apr 18 '20 at 23:33
  • @littleO $n$ is around $1000$ – ProAmateur Apr 18 '20 at 23:34
  • Thanks. How fast are you hoping to be able to solve a million instances of this problem? – littleO Apr 18 '20 at 23:35
  • @littleO Gurobi takes $\approx 0.1$ seconds for each projection step. I'm hoping to reduce the time at least by an order of magnitude. – ProAmateur Apr 18 '20 at 23:44
  • Another thought is that Boyd gave some talks a few years ago about solving convex optimization problems in milliseconds or microseconds using custom interior point methods. I think one of his students from that time posted some code. It doesn't seem easy, but you could look into that and see if it works in your case. You can Google for Boyd's work on "real time embedded convex optimization". – littleO Apr 19 '20 at 00:14
  • @littleO Thanks for the comments! I will look into the methods you mentioned – ProAmateur Apr 19 '20 at 00:17
  • A few thoughts. You could try using the Douglas-Rachford method to minimize $I_A(x) + I_B(x) + \frac12 | x - y |^2$, where $A$ is the probability simplex and $B$ is the scaled Euclidean ball. Here $I_A$ is the convex indicator function of the set $A$. It might be easier to take $A = { x ,\mid \sum_i x_i = 1}$ and $B = { x \mid x\geq 0, \sum_i \alpha_i x_i^2 \leq \alpha}$. It's possible that projecting onto $A$ and $B$ will be easier that way. – littleO Apr 19 '20 at 00:19
  • 2
    I think @littleO is referring to CVXGEN http://www.cvxgen.com/docs/index.html .That does not handle QCQP (or SOCP). – Mark L. Stone Apr 19 '20 at 14:00

0 Answers0