rank-one update with normalization

cholupdate_norm.m

% H = A'*A;
% D = diag(1 ./ sqrt(diag(H)));
% L = chol(D*H*D)';
% cholupd_norm(L, A_new_row, H)

%% http://math.stackexchange.com/a/1345569/61331
function [Lx] = cholupdate_norm(L, x, H)
    n = 1 ./ sqrt(diag(H));
    D = diag(n);
    nx = 1 ./ sqrt(diag(H + x*x'));
    Dx = diag(nx);
    Lx = cholupdate(Dx*inv(D)*L, Dx*x);
end

%% https://en.wikipedia.org/wiki/Cholesky_decomposition#Rank-one_update
function [L] = cholupdate(L,x)
    p = length(x);
    for k=1:p
        r = sqrt(L(k,k)^2 + x(k)^2);
        c = r / L(k, k);
        s = x(k) / L(k, k);
        L(k, k) = r;
        L(k+1:p,k) = (L(k+1:p,k) + s*x(k+1:p)) / c;
        x(k+1:p) = c*x(k+1:p) - s*L(k+1:p,k);
    end
end