Previous Next newton_step_ok.m

ckbs_newton_step Example and Test

Source Code
 
function [ok] = newton_step_ok()
ok = true;
% --------------------------------------------------------
% You can change these parameters
m    = 2;   % number of measurements per time point
n    = 3;   % number of state vector components per time point
N    = 4;   % number of time points
% ---------------------------------------------------------
%  Define the problem
rand('seed', 123);
mu    = .5;
r     = m * N;
p     = n * N;
s     = rand(r, 1) + 1;
y     = rand(p, 1);
u     = rand(r, 1) + 1;
b     = rand(r, 1);
d     = rand(p, 1);
H     = zeros(p, p);
B     = zeros(r, r);
Hdia  = zeros(n, n, N);
Hlow  = zeros(n, n, N);
Bdia  = zeros(m, n, N);
blk_m = m * N + (1 : m);
blk_n = n * N + (1 : n);
for k = N : -1 : 1
        blk_n = blk_n - n;
        blk_m = blk_m - m;
        %
        Bdia(:,:, k)    = rand(m, n);
        tmp             = rand(n, n);
        Hlow(:,:, k)    = tmp;
        Hdia(:,:, k)   = (tmp * tmp') + 4 * eye(n);
        %
        B(blk_m, blk_n) = Bdia(:,:, k);
        H(blk_n, blk_n) = Hdia(:,:, k);
        if k > 1
                H(blk_n, blk_n - n) = Hlow(:,:, k);
        end
        if k < N
                H(blk_n, blk_n + n) = Hlow(:,:, k + 1)';
        end
end
% -------------------------------------------------------------------
[ds, dy, du] = ckbs_newton_step(mu, s, y, u, b, d, Bdia, Hdia, Hlow);
% -------------------------------------------------------------------
F      = [ ...
        s + b + B * y
        H * y + B' * u + d
        s .* u - mu ...
];
dF     = [ ...
        eye(r)      , B           , zeros(r, r)
        zeros(p, r) , H           , B'
        diag(u)     , zeros(r, p) , diag(s) ...
];
delta        = [ ds ; dy ; du ];
ok           = ok & ( max( abs( F + dF * delta ) ) <= 1e-10 );
return
end

Input File: test/newton_step_ok.m