%Making active (preconditioned) residuals orthogonal to blockVectorX
if isempty(operatorB)
blockVectorR(:,activeMask) = blockVectorR(:,activeMask) - ...
blockVectorX*(blockVectorX'*blockVectorR(:,activeMask));
else
blockVectorR(:,activeMask) = blockVectorR(:,activeMask) - ...
blockVectorX*(blockVectorBX'*blockVectorR(:,activeMask));
end
https://github.com/lobpcg/blopex/blob/master/blopex_abstract/krylov/lobpcg.c is a somewhat simplified C version, following https://arxiv.org/abs/0705.2626, of my original MATLAB implementation https://github.com/lobpcg/blopex/blob/master/blopex_tools/matlab/lobpcg/lobpcg.m
Comparing my original MATLAB implementation https://github.com/lobpcg/blopex/blob/master/blopex_tools/matlab/lobpcg/lobpcg.m with https://github.com/lobpcg/blopex/blob/master/blopex_abstract/krylov/lobpcg.c code, I see an important step missing: MATLAB lines 601-608
that should have been before line 1026
/*--- orthonormilize blockVectorR */of https://github.com/lobpcg/blopex/blob/master/blopex_abstract/krylov/lobpcg.cOther relevant improvements are implemented in scikit-learn/scikit-learn#12319