  function P1Grade()
% Checks out LowRankLS
% Generate a small example...


% Solve min || FG'x - b || for xLS
clc
% A small problem to check correctness..
m = 10; n = 6; r = 3;
F = randn(m,r);
G = randn(n,r);
b = randn(m,1);
format long
xLS_YOU = LowRankLS(F,G,b);
normxLS_YOU = norm(xLS_YOU);
xLS_CVL = LowRankLSCVL1(F,G,b);
normxLS_CVL = norm(xLS_CVL);
TheTwoSolutions = [xLS_YOU xLS_CVL]
TheTwoNorms = [normxLS_YOU normxLS_CVL]
% A big problem to make sure there is no superfluous computation
m = 100000; n = 100000; r = 10;
F = randn(m,r);
G = randn(n,r);
b = randn(m,1);
tic
x = LowRankLS(F,G,b);
t_YOU = toc;
tic
x = LowRankLSCVL(F,G,b);
t_CVL = toc;
fprintf('\ntYOU/tCVL = %10.3f\n',t_YOU/t_CVL)



  function xLS = LowRankLSSVD(F,G,b)
% F is m-by-r and rank(F) = r
% G is n-by-r and rank(G) = r
% b is m-by-1
% xLS is the minimum norm solution to min ||(FG')x - b ||

[m,r] = size(F);
M = F*G';
[U,S,V] = svd(M,0);
xLS = V(:,1:r)*((U(:,1:r)'*b)./diag(S(1:r,1:r)));

  function xLS = LowRankLSCVL1(F,G,b)
% F is m-by-r and rank(F) = r
% G is n-by-r and rank(G) = r
% b is m-by-1
% xLS is the minimum norm solution to min ||(FG')x - b ||

% Obtain the LS solution to min ||Fy - b||:
[Q,R] = qr(F,0);
y = R\(Q'*b);
% Now get the minimum norm solution of the underdetermined system
% Gx = y. Use QR of G'...
[Q,R] = qr(G,0);
xLS = Q*(R'\y);

  function xLS = LowRankLSCVL(F,G,b)
% F is m-by-r and rank(F) = r
% G is n-by-r and rank(G) = r
% b is m-by-1
% xLS is the minimum norm solution to min ||(FG')x - b ||


% Suppose F = QF*RF and G = QG*RG are QR factorizations. It follows that

%          || F*G'*x - b || = || (QF*RF)*(QG*RG)'*x - b ||
%                           = || (RF*RG')*xtilde - btilde ||

% where xtilde = QG'*x and btilde = QF'*b
%
% If RF1 = RF(1:r,1:r), RG1 = RG(1:r,1:r), y = xtilde(1:r), and z = xtilde(r+1:n),
% then

%   || (RF*RG')*xtilde - btilde ||^2 = 
%                   || [RF1*RG1' 0; 0 0][y;z] - btilde||^2
%               || RF1*RG1'*y - btilde(1:r)||^2 + ||btilde(r+1:m)||^2
%
% To minimize this y must satisfy 
%
%            RF1*RG1'y  = btilde(1:r) = QF(:,1:r)'*b
%
% and z can be anything. Since x = QG*[y;z] and QG is orthogonal, it 
% follows that we must set z = 0 to obtain the minimum norm solution. Thus,
%
%        xLS = QG(:,1:r)*y = QG(:,1:r)*(RG1'\(RF1\(QF(:,1:r)'*b)))
%
% We can live with the thin qr factorizations of F and G...

[QF1,RF1] = qr(F,0);
[QG1,RG1] = qr(G,0);
xLS = QG1*(RG1'\(RF1\(QF1'*b)));

