Open in the MATLAB editor: f08ck_example
function f08ck_example fprintf('f08ck example results\n\n'); % Find x that minimizes norm(c-Ax) subject to Bx = d . m = int64(6); n = int64(4); p = int64(2); a = [-0.57, -1.28, -0.39, 0.25; -1.93, 1.08, -0.31, -2.14; 2.30, 0.24, 0.40, -0.35; -1.93, 0.64, -0.66, 0.08; 0.15, 0.30, 0.15, -2.13; -0.02, 1.03, -1.43, 0.50]; b = [1, 0, -1, 0; 0, 1, 0, -1]; c = [-1.50; -2.14; 1.23; -0.54; -1.68; 0.82]; d = [0; 0]; % Compute the generalized RQ factorization of (B,A) as % A = ZRQ, B = TQ [TQ, taub, ZR, taua, info] = ... f08zf(b, a); % Set Qx = y. The problem reduces to: % minimize (Ry - Z^Tc) subject to Ty = d % Update c = Z^T*c -> minimize (Ry-c) [cup, info] = f08ag( ... 'Left','Transpose',ZR,taua,c); % Solve Ty = d for last p elements T12 = triu(TQ(1:p,n-p+1:n)); [y2, info] = f07te( ... 'Upper', 'No transpose', 'Non-unit', T12, d); % (from Ry-c) R11*y1 + R12*y2 = c1 --> R11*y1 = c1 - R12*y2 % Update c1 c1 = cup(1:n-p) - ZR(1:n-p,n-p+1:n)*y2; % Solve R11*y1 = c1 for y1 R11 = triu(ZR(1:n-p,1:n-p)); [y1, info] = f07te( ... 'Upper', 'No transpose', 'Non-unit', R11, c1); % Contruct y and backtransform for x = Q^Ty y = [y1;y2]; [~, x, info] = f08ck( ... 'Left', 'Transpose', TQ, taub, y); fprintf('Constrained least squares solution\n'); disp(x); fprintf('Residuals computed directly\n'); res = a*x - c; disp(res); fprintf('Residual norm\n'); disp(norm(res));
f08ck example results Constrained least squares solution 0.4890 0.9975 0.4890 0.9975 Residuals computed directly 0.0030 -0.0129 -0.0193 -0.0084 0.0012 -0.0029 Residual norm 0.0251