Mercurial > repos > public > sbplib

comparison +sbp/D2VariableCompatible.m @ 1217:a3d9567d9004 feature/d2_compatible

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Improve efficiency of construction of D2(b) by precomputing as much as possible
author Martin Almquist <malmquist@stanford.edu>
date Wed, 13 Nov 2019 14:54:22 -0800
parents 2924b3a9b921
children 60c875c18de3
comparison
equal deleted inserted replaced
1198:2924b3a9b921 1217:a3d9567d9004
46 otherwise 46 otherwise
47 error('Invalid operator order %d.',order); 47 error('Invalid operator order %d.',order);
48 end 48 end
49 obj.borrowing.H11 = obj.H(1,1)/obj.h; % First element in H/h, 49 obj.borrowing.H11 = obj.H(1,1)/obj.h; % First element in H/h,
50 obj.borrowing.M.d1 = obj.H(1,1)/obj.h; % First element in H/h is borrowing also for M 50 obj.borrowing.M.d1 = obj.H(1,1)/obj.h; % First element in H/h is borrowing also for M
51 obj.borrowing.R.delta_D = inf; % Because delta_D is zero, one can borrow infinitely much. 51 obj.borrowing.R.delta_D = inf;
52 % This sets penalties of the form 1/borrowing to 0, which is
53 % the desired behaviour.
54 obj.m = m; 52 obj.m = m;
55 obj.M = []; 53 obj.M = [];
54
56 55
57 D1 = obj.D1; 56 D1 = obj.D1;
58 e_r = obj.e_r; 57 e_r = obj.e_r;
59 e_l = obj.e_l; 58 e_l = obj.e_l;
60 59
61 % D2 = Hinv * (-M + br*er*d1r^T - bl*el*d1l^T); 60 % D2 = Hinv * (-M + br*er*d1r^T - bl*el*d1l^T);
62 % Replace d1' by e'*D1 in D2. 61 % Replace d1' by e'*D1 in D2.
63 D2_compatible = @(b) D2(b) - obj.HI*(b(m)*e_r*d1_r' - b(m)*e_r*e_r'*D1) ... 62 correction_l = obj.HI*(e_l*d1_l' - e_l*e_l'*D1);
64 + obj.HI*(b(1)*e_l*d1_l' - b(1)*e_l*e_l'*D1); 63 correction_r = - obj.HI*(e_r*d1_r' - e_r*e_r'*D1);
64
65 D2_compatible = @(b) D2(b) + b(1)*correction_l + b(m)*correction_r;
65 66
66 obj.D2 = D2_compatible; 67 obj.D2 = D2_compatible;
67 obj.d1_l = (e_l'*D1)'; 68 obj.d1_l = (e_l'*D1)';
68 obj.d1_r = (e_r'*D1)'; 69 obj.d1_r = (e_r'*D1)';
69 70
71 function str = string(obj) 72 function str = string(obj)
72 str = [class(obj) '_' num2str(obj.order)]; 73 str = [class(obj) '_' num2str(obj.order)];
73 end 74 end
74 end 75 end
75 76
77
76 end 78 end
77 79
78 80
79 81
80 82