%  This demo illustrates QR with column pivoting for
%  sparsity.   A sparse matrix is generated with
%  columns permuted randomly to give an irregular
%  sparsity pattern to A'A.  The symrcm ordering 
%  is applied to this sparsity pattern to get
%  a permutation P such that  AP  has a good sparsity
%  pattern for obtaining a cholesky factor R that
%  is sparse (in this case banded).   The rows
%  of AP can be processed one at a time via the
%  Row-wise Givens factorization to produce the 
%  QR factorization of AP with the resulting sparsity
%  pattern.
%   
%
%  See George and Heath, "Solution of Sparse Linear
%       Least Squares Problems Using Givens Rotations,"
%       LAA 77, 165-187, (1980).
%
%  D.C. Sorensen
%  22 Sep 03
%
%  Set up Matrix
%
   n = 20;
   a = 0; b = 2*pi;
   h = (b-a)/(n+1);

   x = [a:h:b]';

   m = 100;
   h1 = (b-a)/(m+1);
   t = [a:h1:b]';
   m = length(t);

%
%  set up a sparse A using
%  a B-spline basis for the interval [a,b]
%  partitioned by x
%
   m = length(t); n = length(x);
   A = zeros(m,n);
   for i = 1:m,
      [j, a,b,c,d, ierr] = evalB(t(i),x);
      if(j > 1), A(i,j-1) = a; end
      A(i,j) = b;
      if(j < n), A(i,j+1) = c; end
      if (j < n-1), A(i,j+2) = d; end
   end

%
%  Now give a random permutation to the columns
%  of A to give an irregular sparsity pattern
%
   p = randperm(n);
   A = A(:,p);
%
%  Demonstrate that the sparsity of A'A and the
%  resulting factorization A = QR  are undesirable
%


   spy(A) 
   title('Sparsity of A','fontsize',20)
   pause

   echo on
   M = A'*A;
   spy(M)
   title('Sparsity of A''*A','fontsize',20)
   pause

   [Q,R] = qr(A,0);
   spy(R)
   title('Sparsity of R','fontsize',20)
   pause

%
%  Apply SYMRCM to re-order the columns of A
%  for a more favorable sparsity pattern
%

   p1 = symrcm(M);

   A1 = A(:,p1);
   spy(A1);
   pause


%
%  Show the modified sparsity pattern
%
   
   M1 = A1'*A1;
   spy(M1)
   title('Sparsity of M1 = (AP)''*AP','fontsize',20)
   pause

   [Q1,R1] = qr(A1,0);
   spy(R1)
   title('Sparsity of R1','fontsize',20)

   echo off