0
$\begingroup$

I am looking for a computationally cheap way to compute $x$ such that $$(L L^T + \mu^2 I)x = y$$ where $L \in \mathbb{R}^{n \times n}$ is a lower triangular definite positive matrix (with some very small eigenvalues), $y \in \mathbb{R}^n$ and $\mu \in \mathbb{R}$ are known. If it is necessary, I can assume that $$\mu \ll 1$$ But $\mu^2$ is larger than the smallest eigenvalue of $LL^T$.

Basically, I would like to make the most of my knowledge of the Cholesky decomposition $L L^T$. Eventually, I hope to be able to compute $x$ in $\mathcal{O}(n^2)$. Approximate approaches are also welcomed.

I have seen here that this does not seem to be doable in a more general situation, but I hope the smallness of $\mu$ may help...

Any idea, reference or warning? Thanks for your help.

$\endgroup$

1 Answer 1

1
$\begingroup$

Let $A=LL^T,\lambda=\mu^2,f:X\rightarrow X^{-1}$. That follows is an approximate approach in $O(n^2)$ that is valid only if $\lambda$ is small with respect to $\inf(spectrum(A))$. $Df_A(H)=-A^{-1}HA^{-1},D^2f_A(H,K)=A^{-1}KA^{-1}HA^{-1}+A^{-1}HA^{-1}KA^{-1}$. Thus, according to the Taylor formula, $(A+\lambda I)^{-1}\approx A^{-1}-\lambda A^{-2}+\lambda^2 A^{-3}$. Do not calculate $A^{-2},A^{-3}$ but solve $LL^Tx_0=y,LL^TLL^Tx_1=y,LL^TLL^TLL^Tx_2=y$, that is $LL^Tx_0=y,LL^Tx_1=x_0,LL^Tx_2=x_1$.

$\endgroup$
4
  • $\begingroup$ Unfortunately, the role of $\lambda = \mu^2$ is to regularize the problem: $\lambda$ is larger (resp. smaller) than the smallest (resp. largest) eigenvalue of $A$. So I don't think we can apply this Taylor expansion... Thanks for your help, I edited the question accordingly. $\endgroup$ Commented Jan 8, 2014 at 11:46
  • $\begingroup$ Hi Mathieu, to regularize the problem, addind a scalar matrix, is a bad method. It is much better to do as follows: use this robust Cholesky factorization: math.berkeley.edu/~cinnawu/hss.pdf $\endgroup$
    – loup blanc
    Commented Jan 8, 2014 at 17:21
  • $\begingroup$ Unfortunately (for you), according your comment, your question is an absolute non-sense. Indeed, if your problem is ill-conditioned, then you know the matrix $L$ with a very poor precision ; then what is the interest to regularize the problem in these conditions ? $\endgroup$
    – loup blanc
    Commented Jan 8, 2014 at 22:25
  • $\begingroup$ I don't think we mean the same thing by regularization. I am doing this in the framework of an adaptive algo where I can compute exactly $L$ in $\mathcal{O}(n^2)$. The problem is that it has eigenvalues very close to $0$. Therefore, when I want to invert it I get extremely large values for $x$. I am willing to loose a bit of precision in this inversion if it gives me smaller values for $x$. The regularization needs to be done, not for the Cholesky decomposition, but for the inversion. Thanks for your time. $\endgroup$ Commented Jan 9, 2014 at 9:56

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .