function a = r8bto_dif2 ( m, l )

%*****************************************************************************80
%
%% R8BTO_DIF2 sets up an R8BTO second difference matrix.
%
%  Discussion:
%
%    To get the second difference matrix, it is assumed that M is 1%
%
%    The "indicator matrix" simply has a value like I*10+J at every
%    entry of a dense matrix, or at every entry of a compressed storage
%    matrix for which storage is allocated. 
%
%    The R8BTO storage format is for a block Toeplitz matrix. The matrix
%    can be regarded as an L by L array of blocks, each of size M by M.
%    The full matrix has order N = M * L.  The L by L matrix is Toeplitz,
%    that is, along its diagonal, the blocks repeat.
%
%    Storage for the matrix consists of the L blocks of the first row,
%    followed by the L-1 blocks of the first column (skipping the first row).
%    These items are stored in the natural way in an (M,M,2*L-1) array.
%
%  Example:
%
%    M = 2, L = 3
%
%    1 2 | 3 4 | 5 6
%    5 5 | 6 6 | 7 7
%    ----+-----+-----
%    7 8 | 1 2 | 3 4
%    8 8 | 5 5 | 6 6
%    ----+-----+-----
%    9 0 | 7 8 | 1 2
%    9 9 | 8 8 | 5 5
%
%    X = (/ 1, 2, 3, 4, 5, 6 /)
%
%    B = (/ 91, 134, 73, 125, 97, 129 /)
%
%  Licensing:
%
%    This code is distributed under the MIT license. 
%
%  Modified:
%
%    04 July 2016
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M, the order of the blocks of the matrix A.
%    For the second different matrix, M should be 1.
%
%    Input, integer L, the number of blocks in a row or column.
%
%    Output, real A(M,M,2*L-1), the R8BTO matrix.
%
  a = zeros ( m, m, 2 * l - 1 );
%
%  Blocks 1 to L form the first row.
%
  j = 0;

  for k = 1 : l

    if ( k == 1 )
      value = 2.0;
    elseif ( k == 2 )
      value = -1.0;
    else
      value = 0.0;
    end

    for j2 = 1 : m
      j = j + 1;
      for i = 1 : m
        a(i,j2,k) = value;
      end
    end
  end
%
%  Blocks L+1 through 2*L-1 form the remainder of the first column.
%
  i = m;

  for k = l + 1 : 2 * l - 1

    if ( k == l + 1 )
      value = -1.0;
    else
      value = 0.0;
    end

    for i2 = 1 : m
      i = i + 1;
      for j = 1 : m
        a(i2,j,k) = value;
      end
    end

  end

  return
end