Non-Unit Bidiagonal Matrices for Factorization of Vandermonde Matrices

A non-unit bidiagonal matrix and its inverse with simple structures are introduced. These matrices can be constructed easily using the entries of a given non-zero vector without any computations among the entries. The matrix transforms the given vector to a column of the identity matrix. The given vector can be computed back without any round off error using the inverse matrix. Since a Vandermonde matrix can also be constructed using given n quantities, it is established in this paper that Vandermonde matrices can be factorized in a simple way by applying these bidiagonal matrices. Also it is demonstrated that factors of Vandermonde and inverse Vandermonde matrices can be conveniently presented using the matrices introduced here.


Introduction
Bjorck and Pereyra in 1970 used in their classical work [1] unit bidiagonal matrices with constant off-diagonal entries and diagonal matrices for the LU representation of the inverse of Vandermonde matrices. A non-unit bidiagonal matrix with row-wise constant entries having opposite signs is also used for representing the factors. A recent extension of this approach is adapted to Vandermonde like matrices in Nicholas Higham's book [2]. Regarding the stability of confluent Vandermonde systems, weak stability and weakly stable algorithm concepts are presented in [2]. Weakly stable algorithms solve the dual of non-confluent or confluent Vandermonde or Vandermonde like systems with good accuracy in floating point arithmetic, when there will be not much subtractive cancellations in the inverse Vandermonde UL representation. The desirable criterion for making a minimal amount of subtractive cancellation is that those individual factors of U and L have alternating sign pattern for A=(a ij );-1 (i+j) a ij ≥ 0. The lower triangular components of L are bidiagonal matrices with row-wise constant entries and alternate sign patterns. Higham reports that these components will maintain alternating sign pattern if the points are distinct and arranged in increasing order. Note that Higham does not consider the properties of these matrices outside this stability and accuracy domain and later extended their use for deriving stable factors for Vandermonde like matrices [2] in line with the Bjorck and Pereyra factorization of Vandermonde matrices. It can also be noted that in these two works, the linear transformation that maps a given vector to a column of the identity matrix is not at all considered and inverse of this transformation is not utilized for the factorization of Vandermonde System matrices. From a totally different background and new perception we are going to introduce the lower bidiagonal version of these matrices and present several interesting features with such matrices. An interesting quoting from Gasca and Pena [3] is as follows. " At our knowledge, the uniqueness of different factorizations, which is a consequence of the uniqueness of elimination process, is a novelty in this type of results". It is applicable in this context of introducing the proposed bidiagonal matrices for matrix factorization.
The organization of the paper is as follows. First we will introduce the bidiagonal matrix, its inverse and basic features which make it an ideal choice for factorization of matrices. After that we will discuss the factorization of Vandermonde Matrices, representation of factors and solution of general Vandermonde Linear Systems. This is followed by computational cost of the approach, numerical illustration of solving the Vandermonde linear systems by virtue of the factorization procedure and concluding remarks.

The bidiagonal matrix and its features.
Let a non-zero vector x=[x 1 ,x 2 …x n ] T ; x i ≠0, i=1,2…n be given. Consider the lower bidiagonal matrix and its inverse defined as below. (2.1) Typical examples for the case n=3 is as below.
If we look at the columns in (2.2), these are the given vector itself and its projection to the subspaces of dimension k=n-1, n-2,…,1. These columns constitute a basis and hence can represent any given vector in a unique way. Since the first column itself is the very same vector, the linear combination can be only the entries from 1 e . This forms the elementary theory behind the factorization. Clearly 1 B( x )x e and B( x ) e x 11 . If in the given vector x, x ;k , , , j k 1 2 are zeros and x ; k j , ,n k 0 1then the first j rows in Bxcan be set identical to that of the identity matrix and then 1 j B( x )x e and B( x ) e x j 11 . In general Bxhas to be appropriately tuned with the rows and columns of the identity matrix so that mapping of x to a column of the identity matrix is possible. In any case, the mapping will be to another vector y whose entries will be consisting of ±1 and zeros. Accept that as discussed in J.H. Wilkinson [4], a negligibly small error, say t e 2 , where the computer has t digits mantissa, is bound to occur. Still the mapping and inverse mapping will be always without any round off errors because of the structure of the matrix (2.2) and the presence of unity. See Plamen Koev [5] that relative accuracy will be affected when floating point subtractions are involved as cancellation of significant digits during subtraction of intermediate quantities. This is applicable to the proposed factorization also. But the intermediate quantity is exactly maintained in the inverse of the operator matrix. This structural feature can thus contribute to preserve the relative accuracy. Recall the remarks of Higham [2] about the association of these matrices with his definition of weak stability in maintaining accuracy and stability. These operator matrices can be tuned appropriately with the columns of the identity matrix in presence of zero elements of a column instead of row or column exchanges. The matrices in (2.1) and (2.2) are the results of applying a sequence of column or row operations in corresponding diagonal matrices and these can be illustrated as follows.

BMSA Volume 12
Then we have The matrices It is evident that the matrices in (2.1) and (2.2) are derived from the corresponding diagonal matrices by elementary column operations and whenever a diagonal element is zero it is equivalent to the cancellation of the column operations with the particular diagonal element. Thus the column is reverted to the corresponding column of the identity matrix in (2.3) and (2.4).

Proposition 2.2
From proposition (2.1) and from (2.3)  This is an interesting lower triangular matrix and this construction (2.8) is possible only when the entries of x are distinct and non-zero. A typical 4 x 4 matrix of (2.8) is as below. (2.9) The matrix in (2.9) has an eigen system where in the general case eigen vector corresponding to x1 is [x 1 , x 2 …x n ] T , eigen vector corresponding to x2 is [0 x 2 …x n ] T and so on and that corresponding to xn is [0 0…x n ] T .The diagonal entries will constitute the terms (xj -xj+1) of each column of the matrix (2.9) and the fractional terms will be determined by the entries of its eigen vectors. For example, ( -xj+ ) will be the terms corresponding to its n th power whereas thefractional terms will not be changing. Thus merely by looking at the matrix, one can easily derive the eigen system. The attraction is that its inverse and any power can be easily arrived at without any computations. In the open interval (0,1) , this system attains the minimum and maximum when the off-diagonal entries are uniformly approaching zero. This matrix corresponds to all strictly monotonic decreasing and increasing sequences in the interval (0,1) and correspondence among such sequence matrices are realized by similarity transformation using appropriate diagonal matrices.Recall the remarks by Higham [2] that entries of the inverse Vandermonde lower triangular components should be distinct and in ascending order. This is a pointer to the association of the eigen vectors of matrices (2.9) to such special matrices which are basically generated out of given n distinct quantities.

Proposition 2.3
Given a non-zero n-vector x=[x 1 ,x 2 …x n ] T ; x i ≠0, i=1,2…n then 2 n-1 bidiagonal matrices can be constructed with absolute values of the entries same as that of type (2.1), all of which will map x to e 1 .
Proof: Let B be a lower bidiagonal matrix and consider the equation In (2.10) let α 1 and α 2 be two adjacent entries of a row of B. Assume that α 2 is a diagonal element and α 1 is the corresponding sub-diagonal element in B. For the first row in B, there is only one unique choice as α 1 =; α 2 =1/x 1. For the rest of the rows, assigning one of these unknowns a value, the other can be obtained. So for the remaining 2( n-1) entries, there are infinitely many choices. Here the choice with respect to the diagonal element is a 2 =1/x k ;k=1,2…n. Then the off diagonal elements will be obviously a 2 =-1/x k-1 ;k=2,3…n. Accordingly with this choice we have settled for the matrix (2.1). But a 1 =1/x k -1 ; a 2 =-1/x k also will satisfy equation (2.10). Hence with respect to each of the n-1 rows, the entries can be filled in 2 ways and thus the result follows.
One has the freedom to select a bidiagonal matrix of choice. Then there is a chance that resulting reduction process will be handicapped with the problem of inverting the bidiagonal matrix at every step in addition to disturbing the structural property of the given matrix. For example, P.V Sankar and A.K Sen [7] report that the factorization algorithm proposed by them has the problem of inverting triangular matrices at every step. In the case of the proposed scheme, the operator matrix and its inverse can be easily constructed. These constructions do not call for any additional computations among the entries as in Neville or Gauss. With respect to the number of iterative steps to eliminate column elements, the proposed matrix completes it simultaneously in a single step. Obviously the operator transforms the given vector to a column of the most stable identity matrix and in a stable way. In short, from the infinite set of bidiagonal matrices of (2.10), an ideal bidiagonal matrix for factorization of a given matrix is presented. The detailed results for factorizing a given matrix using these non-unit bi-diagonal matrices are discussed in Nair [6].
The operator B(x) has a close association with Vandermonde matrices. Factors of a Vandermonde matrix can be conveniently presented using B(x)-1 and B(x) matrices. Gaussian Elimination (GE) [5] is applied to a general Vandermonde matrix and as a convenient way to present the resultant factors, B(x) matrix, where x=[1 1...1]T has been used by H. Van de Vel [7], but without referring to the general operator with arbitrary non-zero x , its properties and their applications to the Vandermonde System. Let x0, x1, x2,…., xn є R, be non-zero and distinct. Then a non-singular Vandermonde (n+1)X(n+1) matrix V can be defined as below.

BMSA Volume 12
Since Vandermonde matrices are generated out of n+1 given quantities, its factors, determinant, inverse and solution to Vandermonde linear systems etc. can be expressed as a function of these given set of quantities using the matrices in (2.1) and (2.2) in a convenient way. These are the targets to be achieved this paper.
Following are the Vandermonde linear systems to be solved in this paper.
To solve (2.11) and (2.12), the first attempt is to express the factors of the Vandermonde matrix V and its inverse V -1

Factors of Vandermonde Matrix V
A 4X4 generalized Vandermonde matrix shall be factorized which will give the required ideas regarding the general factors of these systems. So consider STEP-1 Observe that factors T1 and T1-1 of order n are the basic matrices which generate the matrices B(x) and its inverse from the corresponding diagonal matrices as described in (2.5).

STEP-2
Here observe that operators T 2 and are constructed out of forward differences and are of the form x i -x j ; n-1 ≥ i ≥ j ≥ 0.

STEP-3
Here also observe that operators T 3 and are constructed out of forward differences and are of the form x ixj ; n-1 ≥ i ≥ j ≥ 0.

But it is the individual factors of L that can be easily and nicely formulated as below.
A factor Ti of VT can be expressed as (3.1) where and Ei is a nXn matrix whose upper left i rows and columns are identical with the identity matrix and remaining lower right n-i rows and columns are zeros and We know, the order and structure of V are fully defined by the given number of distinct quantities x 0 ,x 1 ,x 2 ,….,x n-1 є R. Hence all the factorization results can be extended to higher order by analogy or mathematical induction. So we can take (3.1) as a general expression for the lower factors of V T . Equation (3.1) clearly establishes the close relationship between the operator B(x) and the Vandermonde systems. ii) The first supra off-diagonal elements u i(i+1) ; i=1,2,…,n-1 are given by ; Polynomials of degree 1 iii) Second supra off-diagonal elements u i(i+2) ; i=1,2,…,n-2 are given by polynomial expressions of order-2 of the form + x i-2 x i-1 + ; i=2,3,n-1. iv) In general, i th supra off-diagonal elements will be polynomial of order i and n-i such elements will constitute this diagonal. Based on these observations, consider the factors of UT for Vandermonde matrix of order 4x4 derived above.
All these factors, U i ; i=1,2,3,4 are unit lower triangular matrices and fourth factor U 4 =I, the identity matrix. The factor U i is defined by x i-1 , (i-1) th element of the given vector [x 0 x 1 …. x n-1 ] T ; i=1,2,….,n, defining the Vandermonde matrix. The first i-1 columns of U i are identical to corresponding column of the identity matrix, as in the case of lower factors, which is expected also. For j=i to n, j th column of U i is identical with the first column of a Vandermonde matrix of dimension n-j+1. So in general appropriate columns of Vandermonde matrix itself and that of the identity matrix are generating U T . Hence we can represent the upper factor as (3.4) Here Vk(xk) is a unit lower triangular matrix of dimension n-k+1, each of its columns are identical with that of the first column of a Vandermonde matrix of dimension n-i+1 with the single variable xi -1 ;i=k,k-1,…1. I k-1 is a (k-1)X(k-1) matrix whose columns are identical with that of the identity matrix.

8
BMSA Volume 12 (3.5) Now to complete the discussion on factors of V and V -1 , the inverse of U has to be derived for dealing with the inverse of the Vandermonde matrix. It is very easy to obtain this matrix as below. V k (1) is the nXn lower triangular matrix which is the base matrix generating B(x) of equation 2.1.

So inverse of Vk(1) is the bi-diagonal matrix
By analogy, we get (3.6) So it is demonstrated that B(x) and B(x) -1 are matrices fundamentally associated with the Vandermonde systems. Now the factorization has come out with the details to spell out how the Vandermonde matrix and its inverse are composed of and the structure of the associated components. So we are in a position to solve the Vandermonde Systems (2.11) and (2.12). Before that let it be pointed out that Vandermonde determinant can be easily derived from the factorization as below.

Vandermonde Determinant.
From the L factor, the determinant of the Vandermonde matrix will be given by That is product of all the terms (x i -x j ) with i >j ; i,j =1,2,….,n. W.Gautschi [9] has treated the norm estimates of Vandermonde inverse matrices and optimal conditioning of these matrices [10] based on the convenient expression (3.7) of Vandermonde determinant.