Solve $\left( \sum_i J_i J_i^T \right) x = \sum_i J_i \cdot \alpha_i$ efficiently

101 Views Asked by At

Is it possible to solve for $x$

$$ \left( \sum_i^N J_i J_i^T \right) x = \sum_i^N J_i \alpha_i$$ where $\alpha_i \in \mathbb{R}$ is a scalar and $J_i \in \mathbb{R}^6$ is a vector

in such a way that I do not have to do a summation over all $\mathbb{R}^{6 \times 6}$ matrices ($J_i J_i^T \in \mathbb{R}^{6 \times 6}$).


What is the reason?

  • $N \approx 10^5$ is quite large.
  • On a GPU this involves 2 sum reductions.
  • The left-hand side involves a $\mathbb{R}^{6 \times 6}$ sum reduction - the right-hand side a $ \mathbb{R}^{6}$ sum reduction: $$ \left( \sum_i^N \mathbb{R}^{6 \times 6}_i \right) x = \sum_i^N \mathbb{R}_i^{6} \mathbb{R}_i$$
  • Reducing $\mathbb{R}^{6 \times 6}$ matrices $N$ times drains the GPU bandwidth.

What is my goal?

  • Doing a sum reduction with a lower dimensionality. For example $$ x = \sum_i^N \mathbb{R}^{6}$$

Important notes

  • The GPU is memory bound but not compute bound.
  • Consider $\left( \sum_i^N J_i J_i^T \right)$ as positive semi-definite.

Thanks

3

There are 3 best solutions below

1
On

If you can find a vector $X$ that is perpendicular to all but one of the $J_i$'s, then you can left multiply by $X^T$ to get $$ \sum_{i} X^T J_i J_i^T x = \sum_{i} X^T J_i g, $$ which will simplify to $$ X^T J_k J_k^T x = X^T J_k g $$ where we assume that $X^T J_k \neq 0$. This simplifies to $$ J_k^T x = g $$ which, as in your previous question, has infinitely many solutions.

1
On

$\sum_i J_i J_i^T$ is a positive semidefinite matrix, positive definite if the inverse exists. The Cholesky decomposition is then then both efficient and numerically stable.

3
On

I don't understand your rationale exactly, but what about conjugate gradient? This will solve your problem exactly in exactly $6$ steps (and $6$ gradient evaluations/ matrix $A=\sum J_iJ_i^T$ vector multiplies). You don't need to actually compute $A$ itself, you only need to compute $Ax$, for vectors, which can be done by computing the individual $J_i\times(J_i\cdot x)$ terms in the sum.