# source:stable/3.9/Ipopt/contrib/MatlabInterface/examples/bayesnet/computeJGHessian.m@1864

Last change on this file since 1864 was 1864, checked in by andreasw, 3 years ago

sync with trunk rev 1862

File size: 1.7 KB
Line
2% This code is published under the Eclipse Public License.
3%
4% Author: Peter Carbonetto
5%         Dept. of Computer Science
6%         University of British Columbia
7%         May 19, 2007
8
9function H = computeJGHessian (qR, qS, sigma, lambda, ...
10                               returnStructureOnly, auxdata)
11  [K C f Rv Rf Sv Sf NS d] = deal(auxdata{:});
12
13  ns  = length(Sv);  % The number of small regions.
14  nqr = length(qR);  % The number of variables associated with the
15                     % marginals defined on the large regions.
16  nqs = length(qS);  % The number of variables associated with the
17                     % marginals defined on the small regions.
18
19  if returnStructureOnly
20    H = speye(nqr + nqs);
21  else
22
23    % Compute the second-order partial derivatives on the subblock HR.
24    HR = spdiag(1 ./ qR);
25
26    % Compute the second-order partial derivatives on the subblock
27    % HS. Repeat for each small region.
28    entries = zeros(1,nqs);
29    is      = 0;
30    for s = 1:ns
31      tablesize = prod(K(Sv{s}));  % The number of entries necessary to
32                                   % store the entire marginal.
33      is          = is(end) + (1:tablesize);
34      entries(is) = (1-d(s)) ./ qS(is);
35    end
36    HS = spdiag(entries);
37
38    % Construct the full Hessian from its blocks.
39    H = sigma * [ HR                spzeros(nqr,nqs);
40                  spzeros(nqs,nqr)  HS             ];
41  end
42
43% ------------------------------------------------------------------
44function A = spdiag (xs)
45  n = length(xs);  % The height and width of the sparse matrix.
46
47  % Produce the sparse, diagonal matrix.
48  A = sparse(1:n,1:n,xs,n,n);
49
Note: See TracBrowser for help on using the repository browser.