Y Saukh - C- factorization method for sparse matrices - страница 1

1  2 

CR- factorization method for sparse matrices

Sergiy Ye. Saukh National Academy of Sciences of Ukraine saukh@svitonline.com


Sergiy Ye. Saukh. CR- factorization method for sparse matrices. The

column-row (CR-) factorization method is offered. The CR-factorization method differs from the LU- factorization method by property of adaptive to the placement of pivoting entries. Suggested method allows to give up implementation of actions by transposition of columns and rows in the process of matrix factorization and therefore to accelerate the solving of the large-scale systems of linear algebraic equalizations.


Triangular decomposition (or lu - factorization) method is one of the most frequently used variant of Gaussian elimination method. In this case matrix A is presented as the product of lower triangular matrix L and upper triangular matrix U. The necessary conditions to realize lu - factorization are regularity of A matrix and nonzero entry ukk of matrix U on the k - th step of factorization.

The triangular decomposition of regular matrix A can always be achieved by permuting of rows and columns so that conditions ukk ф 0 are held true. For that

purpose permutation matrices P and Q are used such that P A Q = L U. In the conditions of calculations with finite number of digits the Gaussian elimination method realization features require such selection of pivoting entries aX] that the

growth of absolute values of entries of L and U matrices is staying limited during the factorization process [3-8]. In the case of the sparse matrix A the permutations influence on the level of filling by nonzero entries of L and U matrices are also taken into account. This level may raise significantly in comparison with initial value of the level of filling by nonzero entries of matrix A. So the selection of the pivoting entries is based on the combination of demands to support of numerical stability of the Gaussian elimination and to minimize the filling levels of the L and U matrices. In accordance with these demands the different strategies of pivoting entries selection are realized in modern software [3-6].

For storage nonzero entries of sparse matrices are used special formats of data so called row-wise or linked list representation formats. Using these formats significantly complicate procedures of the permutation of rows and columns [3­6]. Although the permutation operations don't need really displacements of the

values of entries in the computer memory, but they require searching and changing index values of the rows and columns. In general the sub arrays of permuting column and row indexes don't coincide by lengths and it complicates the realization of index displacements. Our researches show the commensurability of the volumes of operations which are necessary for execution of the permutations and the lu - factorization. We suggest the new method of column-row factorization or CR - factorization briefly, which is distinguished on principle from the method of lu - factorization by the adaptation to the selection of pivoting entries. It permits to refuse from the row and column permutations in the process of matrix factorization and therefore to accelerate the process of calculation of matrix factors.

CR-factorization method

The basic relation of the method is the expression:

in which the initial n x n matrix A is presented by the sum of n x n matrix A ^(i, j)) and the product of nxi column-vector Cj and lxn row-vector Ri.

If the entry aij of matrix A doesn't equal to zero then from the relation (l) can determine all the entries of the column Cj and row Ri accurate to a certain factor, under the assumption made that the entries of row i and column j of matrix A f (i, j)) equal to zero. Indeterminacy "accurate to a certain factor" in determination of entry column Cj and row Ri is conditioned by impossibility to find two required values cX] and rX] from one equation cX] rX] = aX].

For the elimination of indeterminacy in ci j and ri j values determination it ought to accept one of the eventual assumptions, for example,

Notice that the assumptions (2) and (3) are analogous to the assumptions which are accepted in point of diagonal entries of one of the L and U factor matrices in lu - factorization method. The assumption (4) also takes place in the special variant of lu - factorization method for the symmetrical positive determinate matrices and are known as the Cholesky elimination method ( llt - factorization) or the square root method [7, 8].

The relation (l) determines the partial one-step column-row factorization of A matrix relative to the aX] pivoting entry. Here the column C] and row Rx

are multipliers or factors, and A ^(i,j)) matrix is a residual matrix. Its entries,

A = C j R x + A i( (i, j))


r . = l

(2) (3) (4)

except those which are in i - th row and j - th column, form an active (n - l)x(n -1) submatrix.

Analyzing the A f(i, j)) matrix structure in the relation (1), we notice that the obtained entries of column Cj and row Ri naturally disposes in the matrix A f (i, j)) on the place of zero entries its i - th row and j - th column. This superposition of the entries of column Cj and row Ri with the entries of matrix A f (i, j)) leads to the superposed matrix forming as

A^r((x, j)) = Ai((x, j))©(cj ,Rx) (5)

where the symbol © is the superposition operation.

Thus the relations (1) - (5) determine the first step of a row-column factorization of matrix A relatively ai j pivoting entry without the rows and

columns permutation. The following analogous step of factorization is realized relative to the pivoting entry, which is chosen in the active part of matrix A f(i, j)). It coincides with the analogous part of matrix A^R{(i, j)).

If the pivoting entry is apq -cpirjq ф 0 which is situated on the intersection

of p - row and q - column of matrix A f (i, j)) then taking into account the identity

(1), we can write

A i((i, j)) = CqRp + A 2((i, j), (p, q)) (6)

where n x l vector-column Cq and lx n vector-row R p contain the zero entries ciq = 0 and rpj = 0 and n x n matrix A2((i,j), (p,q)) has two zero rows and two zero columns with indexes i, p and j, q correspondingly. We receive after the substitution (6) to (l):

A = CjRx + CqRp + A 2{(i, j), (p, q)) (7)

or in the combined form

A C2R((i, j), (p, q)) = A 2{(i, j), (p,q))®(Cj ,R, ,Cq ,R p ). (8)

Using the expressions (l) and (6), which determine partial one- or two-steps column-row factorization, we receive the full formula of CR - factorization of regular n x n matrix A for the case when the pivoting entry co-ordinates (i, j) are chosen dynamically in the process of factorization and form certain ordered set P:

A 2 C j  R x . (9)

(x,j )e P

Here the summing is realized in the entry sequence in set P. Notice that in contrast to (l) and (7) the right part of this expression don't contain the summand type matrix A n(p) of entries those n rows and n columns in the

intersection of which the pivoting entries are chosen. These entries are replaced by zero in the process of forming A n(p) for n steps. Just these places vector­

columns Cj and vector-rows Ri are situated at the forming of the coincidence matrix Acr{p) . So A n(p) = 0.

Obviously the expression (9) can be presented in the matrix form:

A = C(P).R(P), (l0)

where the matrices c(p) and r(p) are aggregated from the vector-columns Cj

and vector-rows  Ri   accordingly.  The situation order of these vectors

corresponds to the order of their indexes in p set. Therefore the situation orders of vector-columns Cj and vector-rows Ri in the coincidence matrix

ACR{P) = CCR ©RCR (ll)

correspond to their j and i indexes so in fact in the matrix ACR(p) are superposed two matrices Ccr(p) and Rcr(p) instead of the matrices c(p) and


In general case the matrices c(p) and r(p) don't take triangular form.

Only when the pivoting entries are chosen sequentially in the co-ordinate p((l,l),(2,2),...,(n,n)) in the process CR - factorization, then forming vector-

columns Cj and vector-rows Ri form the lower triangular and upper triangular

matrices c{p) = L and r(p) = U, which are situated in the coincidence matrix

ACR(p) so that CCR(p) = c(p) = L and rcr(p) = r(p) = U . The basic structural

peculiarity of general matrices c(p) and r(p) is such as it in principle can be

converted to the triangular form by the way of rows and columns permutations.

The results of the experimental researches

It is evident that the sample-comparison operations over the indexes into rows and columns lists of nonzero matrix entries, which are executed in the process of rows and columns permutation, cannot be directly correlated with the arithmetical operations executed under the nonzero entries in the factorization process. So the number of sample-comparison operations should be estimated experimentally by the time expenses for the fulfillment of those operations and to determine their specific weight in the general time expenses on the matrices factorization. The description of the experiments conditions and its results are represented below.

For calculations the computer Intel P4 (Chipset Intel 865 PE, FSB 800 MHz, CPU 3.0GHz with HT, Dual Channel Memory lGB: 2 x 5l2MB, DDR400) running under operational system Microsoft Windows XP are used. The program code was written by the author on the C++ language in the Microsoft Visual Studio.Net and was taken as a base for the experimental researches. Sequentially three algorithms are realized for solution of linear

algebraic equation A X = B with a vector of a right part B = A ■ 1, where 1 is unit vector. For CR - factorization of matrix A the first algorithm realizes the method with the improved generalized Markovits's strategy of the pivoting entries choice [4]. The second algorithm realizes the method of CR - factorization of the same matrix A, using the set P of pivoting entries co-ordinates on their choice order, which was formed by the first algorithm. The third algorithm differs from the second one only the lu - factorization method realization. Notice that in the both realization of CR - factorization method for matrix A the assumption (3) is used.

The implementation of a written program code allows to estimate the summary time expenses for the search of X vector-solution by the given vector B and to differentiate the time expenses on CR - and lu - factorization.

In all our experiments the equality fulfillment nnz(cCR © RCR) = nnz(l © u) is supported, this is the coincidence of the numbers of nonzero entries in the c and R factor matrices and in the L and u factor matrices which are received by the methods of CR - and lu - factorization.

For the supporting of correct experimental researches the influence of control procedures of nonzero matrix entries located in working arrays was slacked maximally. First of all the sizes of arrays of nonzero entries and their rows and columns indexes were chosen equal 45 ■ l06. This is so large that not to need to initiate the special procedure "cleaning the refuse" in them [4, 5]. In addition the improving extended Markovits's strategy is used in all tests, where the search of pivoting entries is limited by the set Sp of p rows of active

submatrix with the minimal numbers of nonzero entries. The pivoting is chosen the entry aXj with the least Markovits's value but that the absolute value of

which not more than in т times less than maximal entry  on absolute value

a max ^ Sp .

Table l. The test matrices, their dimensions n, the numbers of nonzero entries in initial matrices nnz(a) and in the coincidence factor matrices nnz(c © r),

mistakes є and the test systems solution time TCR .




nnz (c © r)



ASIC 680k






ASIC 680ks





1  2 

Похожие статьи

Y Saukh - C- factorization method for sparse matrices