1function L = cholupdown (Lold, sigma, w)
2%CHOLUPDOWN Cholesky update/downdate (Bischof, Pan, and Tang method)
3% Example:
4%   L = cholupdown (Lold, sigma, w)
5% See also: cs_demo
6
7% Copyright 2006-2012, Timothy A. Davis, http://www.suitesparse.com
8
9beta = 1 ;
10n = size (Lold,1) ;
11L = Lold ;
12% x = weros (n,1) ;
13worig = w ;
14
15for k = 1:n
16
17    alpha = w(k) / L(k,k) ;
18    beta_new = sqrt (beta^2 + sigma*alpha^2) ;
19    gamma = alpha / (beta_new * beta) ;
20
21    if (sigma < 0)
22
23        % downdate
24        bratio = beta_new / beta ;
25        w (k+1:n) = w (k+1:n) - alpha * L (k+1:n,k) ;
26        L (k,k) = bratio * L (k,k) ;
27        L (k+1:n,k) = bratio * L (k+1:n,k) - gamma*w(k+1:n) ;
28
29    else
30
31        % update
32        bratio = beta / beta_new ;
33
34%       wold = w (k+1:n) ;
35%       w (k+1:n) = w (k+1:n) - alpha * L (k+1:n,k) ;
36%       L (k    ,k) = bratio * L (k    ,k) + gamma*w(k) ;
37%       L (k+1:n,k) = bratio * L (k+1:n,k) + gamma*wold ;
38
39        L (k,k) = bratio * L (k,k) + gamma*w(k) ;
40        for i = k+1:n
41
42            wold = w (i) ;
43            w (i) = w (i) - alpha * L (i,k) ;
44            L (i,k) = bratio * L (i,k) + gamma*wold ;
45
46        end
47
48    end
49
50    w (k) = alpha ;
51
52    beta = beta_new ;
53
54end
55
56norm (w-(Lold\worig))
57