I am working on a problem where I had to find the following expression: $$ \ell = \operatorname{Tr}(P'HP)$$
I already modified my model formulation using cholesky decomposition for PSD matrices and came up with equivalent expression with $L$ instead of $P.$
$$ \ell\,' = \operatorname{Tr}(L'HL)$$ Where,
$L$ is lower triangular and $H$ is symmetric. How can I compute this fast? I do not want to compute the full product and then take trace since they are huge matrices.
I tried to decompose $H$ into a lower and an upper triangular matrix, but it doesn't help.
I think that the complexity of the calculation is $\Theta(n^3)$.
Let $L=[L_{i,j}],H=[H_{i,j}]$. Thus $tr(L^THL)=\sum_{i,j,k}L_{j,i}H_{j,k}L_{k,i}=\sum_{j,k}<L_j,L_k>H_{j,k}$ where $<L_j,L_k>$ is the scalar product of $L_j$ and $L_k$ -the rows of $L$ of indices $j,k$-. Since the $(H_{j,k})$ do not depend on $L$, the searched complexity is the one of the calculation of the $(<L_j,L_k>)_{j\leq k}$. I do not see how to do better than $\Theta(n^3)$; indeed, a priori, there are no relations linking the $(<L_j,L_k>)_{j\leq k}$.
EDIT. About the user1551's post. Each calculation of $x^TL^THLx$ has complexity $O(n^2)$. Let $n=200$; a sample of $10^5$ tests gives an error $\approx 1.5$%; a sample of $200$ tests does not give any significant digit. Unfortunately, it seems that the number of tests is necessarily much greater than $O(n)$.