@chebfun|qr.m

function [Q, R] = qr(A, econ)
%QR   QR factorization of an array-valued CHEBFUN.
%   [Q, R] = QR(A) or QR(A, 0), where A is a column CHEBFUN with n columns,
%   produces a column CHEBFUN Q with n orthonormal columns and an n x n upper
%   triangular matrix R such that A = Q*R.
%
%   The algorithm used is described in L.N. Trefethen, "Householder
%   triangularization of a quasimatrix", IMA J. Numer. Anal. (30), 887-897
%   (2010).
%
% See also SVD, MRDIVIDE, RANK.

% Copyright 2015 by The University of Oxford and The Chebfun Developers.
% See http://www.chebfun.org/ for Chebfun information.

% Check inputs
if ( (nargin == 2) && (econ ~= 0) )
    error('CHEBFUN:CHEBFUN:qr:twoargs',...
      'Use qr(A) or qr(A, 0) for QR decomposition of an array-valued CHEBFUN.');
end
if ( A(1).isTransposed )
    error('CHEBFUN:CHEBFUN:qr:transpose',...
        'CHEBFUN QR works only for column CHEBFUN objects.')
end
if ( ~all(isfinite(A(1).domain)) )
    error('CHEBFUN:CHEBFUN:qr:infdomain', ...
        'CHEBFUN QR does not support unbounded domains.');
end

numCols = numColumns(A);

% If A has only one column we simply scale it.
if ( numCols == 1 )
    R = sqrt(innerProduct(A, A));
    Q = A./R;
    return
end

% If A is a quasimatrix then try to convert it to an array-valued CHEBFUN and then
% call QR on that. Otherwise, use ABSTRACTQR with a Legendre-Vandermonde matrix
% as the second argument.
if ( numel(A) > 1 )

    [A, isSimple] = quasi2cheb(A);
    if ( isSimple ) % Successfully converted to array-valued CHEBFUN.
        [Q, R] = qr(A, 0);
        return
    end

    % We can't convert the quasimatrix to an array-valued CHEBFUN. Tricky case.
    % TODO: This case is not tested!
    %  DEVELOPER NOTE:
    %   Currently (5th Feb 2016) the only way this is reachable is in the
    %   case of a quasimatrix consisting of SINGFUN or DELTAFUN objects,
    %   neither of which return anything sensible when we attempt to compute
    %   a QR factorization.
    
    % Legendre-Vandermonde matrix converted to also be a quasimatrix:
    L = cheb2quasi(legpoly(0:numCols-1, domain(A), 'norm', 1));
    [Q, R] = abstractQR(A, L, @innerProduct, @normest);     
    return
end

% If A is an array-valued Chebfun, find the number of pieces.
numFuns = numel(A.funs);

% If there is only one piece (no breakpoints), call QR at the fun level.
if ( numFuns == 1 )
    [Q, R] = qr(A.funs{1});
    Q = chebfun({Q});
    return
end

% If there are multiple pieces, perform QR on each piece, stack the
% R matrices from the separate pieces and perform QR on that, yielding
% Qhat and Rhat, then break Qhat back down into sections and fold those
% back into the matrices Q.

% Perform QR on each piece:
R = cell(numFuns, 1);
Q = R;
for k = 1:numFuns
    [Q{k}, R{k}] = qr(A.funs{k});
end

% Compute [Qhat, Rhat] = qr(R):
[Qhat, Rhat] = qr(cell2mat(R));
Rhat = Rhat(1:numCols,:); % Extract required entries.
Qhat = Qhat(:,1:numCols);

% Ensure that the diagonal of Rhat is non-negative:
neg = diag(Rhat) < 0;
Rhat(neg,:) = -1 * Rhat(neg,:);
Qhat(:,neg) = -1 * Qhat(:,neg);

% ... And fold Qhat back in to Q for each segment:
m = cellfun(@(v) size(v, 1), R); % Number of rows per segment.
Qhat = mat2cell(Qhat, m, numCols);
for k = 1:numFuns
    Q{k} = Q{k}*Qhat{k,1};
end

Q = chebfun(Q);
R = Rhat;

end

diff.patch

--- untitled 8
+++ (clipboard)
@@ -13,14 +13,6 @@
 % Copyright 2015 by The University of Oxford and The Chebfun Developers.
 % See http://www.chebfun.org/ for Chebfun information.
 
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-% Developer note:
-%  If A contains only a single FUN, then FUN/QR is used directly. If A has
-%  multiple pieces but each of these are simple CHEBTECH objects, then
-%  QRSIMPLE() is called. This violates OOP principles, but is _much_ more
-%  efficient.
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
 % Check inputs
 if ( (nargin == 2) && (econ ~= 0) )
     error('CHEBFUN:CHEBFUN:qr:twoargs',...
@@ -37,74 +29,78 @@
 
 numCols = numColumns(A);
 
+% If A has only one column we simply scale it.
 if ( numCols == 1 )
-    % If f has only one column we simply scale it.
     R = sqrt(innerProduct(A, A));
     Q = A./R;
-    
-elseif ( numel(A) > 1 )
-    % Quasimatrix case:
+    return
+end
+
+% If A is a quasimatrix then try to convert it to an array-valued CHEBFUN and then
+% call QR on that. Otherwise, use ABSTRACTQR with a Legendre-Vandermonde matrix
+% as the second argument.
+if ( numel(A) > 1 )
 
-    % Try to convert to an array-valued CHEBFUN:
     [A, isSimple] = quasi2cheb(A);
-    
     if ( isSimple ) % Successfully converted to array-valued CHEBFUN.
-        
         [Q, R] = qr(A, 0);
-        
-    else            % We can't convert to array valued. Tricky case.
-        
-        % TODO: This case is not tested!
-        %  DEVELOPER NOTE:
-        %   Currently (5th Feb 2016) the only way this is reachable is in the
-        %   case of a quasimatrix consisting of SINGFUN or DELTAFUN objects,
-        %   neither of which return anything sensible when we attempt to compute
-        %   a QR factorization.
-
-        % Legendre-Vandermonde matrix:
-        L = legpoly(0:numCols-1, domain(A), 'norm', 1);
-        % Convert so that L is also quasimatrix:
-        L = cheb2quasi(L);
-        % Call abstract QR:
-        [Q, R] = abstractQR(A, L, @innerProduct, @normest);
-        
+        return
     end
-    
-elseif ( numel(A.funs) == 1 )
-    % No breakpoints = easy case.
-    
-    % Call QR at the FUN level:
-    [Q, R] = qr(A.funs{1});
-    Q = chebfun({Q});
 
-else
+    % We can't convert the quasimatrix to an array-valued CHEBFUN. Tricky case.
+    % TODO: This case is not tested!
+    %  DEVELOPER NOTE:
+    %   Currently (5th Feb 2016) the only way this is reachable is in the
+    %   case of a quasimatrix consisting of SINGFUN or DELTAFUN objects,
+    %   neither of which return anything sensible when we attempt to compute
+    %   a QR factorization.
+    
+    % Legendre-Vandermonde matrix converted to also be a quasimatrix:
+    L = cheb2quasi(legpoly(0:numCols-1, domain(A), 'norm', 1));
+    [Q, R] = abstractQR(A, L, @innerProduct, @normest);     
+    return
+end
 
-    numFuns = numel(A.funs);
+% If A is an array-valued Chebfun, find the number of pieces.
+numFuns = numel(A.funs);
 
-    % Perform QR on each piece:
-    R = cell(numFuns, 1);
-    Q = R;
-    for k = 1:numFuns
-        [Q{k}, R{k}] = qr(A.funs{k});
-    end
+% If there is only one piece (no breakpoints), call QR at the fun level.
+if ( numFuns == 1 )
+    [Q, R] = qr(A.funs{1});
+    Q = chebfun({Q});
+    return
+end
 
-    % Compute [Qhat, Rhat] = qr(Q):
-    m = cellfun(@(v) size(v, 1), R);
-    [Qhat, Rhat] = qr(cell2mat(R));
-    Rhat = Rhat(1:numCols,:); % Extract require entries.
-    Qhat = mat2cell(Qhat(:,1:numCols), m, numCols);
-
-    % Ensure diagonal has positive sign ...
-    s = sign(diag(Rhat)); s(~s) = 1; 
-    S = spdiags(s, 0, numCols, numCols); 
-    R = S*Rhat;
-    % ... and fold Qhat back in to Q:
-    for k = 1:numFuns
-        Q{k} = Q{k}*(Qhat{k,1}*S);
-    end
-    
-    Q = chebfun(Q);
+% If there are multiple pieces, perform QR on each piece, stack the
+% R matrices from the separate pieces and perform QR on that, yielding
+% Qhat and Rhat, then break Qhat back down into sections and fold those
+% back into the matrices Q.
+
+% Perform QR on each piece:
+R = cell(numFuns, 1);
+Q = R;
+for k = 1:numFuns
+    [Q{k}, R{k}] = qr(A.funs{k});
+end
 
+% Compute [Qhat, Rhat] = qr(R):
+[Qhat, Rhat] = qr(cell2mat(R));
+Rhat = Rhat(1:numCols,:); % Extract required entries.
+Qhat = Qhat(:,1:numCols);
+
+% Ensure that the diagonal of Rhat is non-negative:
+neg = diag(Rhat) < 0;
+Rhat(neg,:) = -1 * Rhat(neg,:);
+Qhat(:,neg) = -1 * Qhat(:,neg);
+
+% ... And fold Qhat back in to Q for each segment:
+m = cellfun(@(v) size(v, 1), R); % Number of rows per segment.
+Qhat = mat2cell(Qhat, m, numCols);
+for k = 1:numFuns
+    Q{k} = Q{k}*Qhat{k,1};
 end
 
+Q = chebfun(Q);
+R = Rhat;
+
 end

Hello. I'm sorry. I'm weak in coding. Can you submit the code for the clenshaw method to my email?
Thank you very much
my email:
milad.mahboub71@gmail.com