1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
|
% STK_KREQ_QR [STK internal]
%
% CALL: kreq = stk_kreq_qr (model, xi, xt)
%
% constructs an stk_kreq_qr object.
% Copyright Notice
%
% Copyright (C) 2016, 2017 CentraleSupelec
% Copyright (C) 2013, 2014 SUPELEC
%
% Author: Julien Bect <julien.bect@centralesupelec.fr>
% Copying Permission Statement
%
% This file is part of
%
% STK: a Small (Matlab/Octave) Toolbox for Kriging
% (https://github.com/stk-kriging/stk/)
%
% STK is free software: you can redistribute it and/or modify it under
% the terms of the GNU General Public License as published by the Free
% Software Foundation, either version 3 of the License, or (at your
% option) any later version.
%
% STK is distributed in the hope that it will be useful, but WITHOUT
% ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
% or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
% License for more details.
%
% You should have received a copy of the GNU General Public License
% along with STK. If not, see <http://www.gnu.org/licenses/>.
function kreq = stk_kreq_qr (model, xi, xt)
if nargin == 0
kreq = struct ('n', 0, 'r', 0, 'P_scaling', [], ...
'LS_Q', [], 'LS_R', [], 'RS', [], 'lambda_mu', []);
else
[Kii, Pi] = stk_make_matcov (model, xi);
[n, r] = size (Pi);
% heuristics: scale Pi to (try to) avoid conditioning issues later
P_scaling = compute_P_scaling (Kii, Pi);
Pi = bsxfun (@times, Pi, P_scaling);
% kriging matrix (left-hand side of the kriging equation)
LS = [[Kii, Pi]; [Pi', zeros(r)]];
% orthogonal-triangular decomposition
[Q, R] = qr (LS);
% Check for duplicated rows if R is close to singular
u = abs (diag (R));
if (~ stk_isnoisy (model)) && (min (u) < eps * max (u))
stk_assert_no_duplicates (xi);
end
kreq = struct ('n', n, 'r', r, 'P_scaling', P_scaling, ...
'LS_Q', Q, 'LS_R', R, 'RS', [], 'lambda_mu', []);
end
kreq = class (kreq, 'stk_kreq_qr');
% prepare the right-hand side of the kriging equation
if nargin > 2
[Kti, Pt] = stk_make_matcov (model, xt, xi);
kreq = stk_set_righthandside (kreq, Kti, Pt);
end
end % function
%!test stk_test_class ('stk_kreq_qr')
|
