www.gusucode.com > SAE RBM 程序MATLAB源码代码实现的一个关于sae的例子 > minFunc/mchol.m
function [l,d,perm] = mchol(A,mu) % Compute a modified LDL factorization of A % (MEX ME!) if nargin < 2 mu = 1e-12; end n = size(A,1); l = eye(n); d = zeros(n,1); perm = 1:n; for i = 1:n c(i,i) = A(i,i); end % Compute modification parameters gamma = max(abs(diag(A))); xi = max(max(abs(setdiag(A,0)))); delta = mu*max(gamma+xi,1); if n > 1 beta = sqrt(max([gamma xi/sqrt(n^2-1) mu])); else beta = sqrt(max([gamma mu])); end for j = 1:n % Find q that results in Best Permutation with j [maxVal maxPos] = max(abs(diag(c(j:end,j:end)))); q = maxPos+j-1; % Permute d,c,l,a d([j q]) = d([q j]); perm([j q]) = perm([q j]); c([j q],:) = c([q j],:); c(:,[j q]) = c(:,[q j]); l([j q],:) = l([q j],:); l(:,[j q]) = l(:,[q j]); A([j q],:) = A([q j],:); A(:,[j q]) = A(:,[q j]); for s = 1:j-1 l(j,s) = c(j,s)/d(s); end for i = j+1:n c(i,j) = A(i,j) - sum(l(j,1:j-1).*c(i,1:j-1)); end theta = 0; if j < n && j > 1 theta = max(abs(c(j+1:n,j))); end d(j) = max([abs(c(j,j)) (theta/beta)^2 delta]); if j < n for i = j+1:n c(i,i) = c(i,i) - (c(i,j)^2)/d(j); end end end