SPDIAGS

 

  SPDIAGS Extract and create sparse band and diagonal matrices.

 

 

  SPDIAGS, which generalizes the builtin function "diag", deals with

  three matrices, in various combinations, as both input and output.

  A is an m-by-n matrix, usually (but not necessarily) sparse,

  with its nonzero, or specified, elements located on p diagonals.

  B is a min(m,n)-by-p matrix, usually (but not necessarily) full,

  whose columns are the diagonals of A.

  d is a vector of length p whose integer components specify

  the diagonals in A.

 

 

  Roughly, A, B and d are related by

      for k = 1:p

          B(:,k) = diag(A,d(k))

      end

 

 

  Four different operations, distinguished by the number of

  input arguments, are possible with SPDIAGS:

 

 

  Extract all nonzero diagonals:

      [B,d] = spdiags(A);

  Extract specified diagonals:

      B = spdiags(A,d);

  Replace specified diagonals:

      A = spdiags(B,d,A);

  Create a sparse matrix from its diagonals:

      A = spdiags(B,d,m,n);

 

 

  The precise relationship among A, B and d is:

     if m >= n

        for k = 1:p

           for j = max(1,1+d(k)):min(n,m+d(k))

              B(j,k) = A(j-d(k),j);

           end

        end  

     if m < n

        for k = 1:p

           for i = max(1,1-d(k)):min(m,n-d(k))

              B(i,k) = A(i,i+d(k));

           end

        end  

     end

  Some elements of B, corresponding to positions "outside" of A,

  are not defined by these loops.  They are not referenced when

  B is input and are set to zero when B is output.

 

 

  For example, this generates a sparse tridiagonal representation

  of the classic second difference operator on n points.

 

 

      e = ones(n,1);

      A = spdiags([e -2*e e], -1:1, n, n)

 

 

  This turns it into Wilkinson's test matrix (see WILKINSON5TV0YJ(n)).

 

 

      A = spdiags(abs(-(n-1)/2:(n-1)/2)',0,A)

 

 

  Finally, this recovers the 3 diagonals.

 

 

      B = spdiags(A)

 

 

  The second example is not square.

 

 

      A = [ 11    0   13    0

             0   22    0   24

             0    0   33    0

            41    0    0   44

             0   52    0    0

             0    0   63    0

             0    0    0   74]

 

 

  has m = 7, n = 4 and p = 3.

 

 

  The statement [B,d] = spdiags(A)  produces d = [-3 0 2]'  and

       

      B = [ 41   11    0

            52   22    0

            63   33   13

            74   44   24 ]

 

 

  Conversely, with the above B and d, the expression

  spdiags(B,d,7,4) reproduces the original A.

 

 

  See also DIAG16RMFU2.