increases again. An illustration of this is shown in Figure 10.14. This plot suggests that

there is an optimal value for which is far from the smallest admissible one. For small ,

8 w

the diagonal dominance of is weak and, as a result, the computed IC factorization

§ 8 w

a `%§ a `

is a poor approximation to the matrix . In other words, is close to § §

a ` §E

the original matrix , but the IC(0) factorization is far from . For large , the oppo-

§

a `

site is true. The matrix has a large deviation from , but its IC(0) factorization

§ § a`

may be quite good. Therefore, the general shape of the curve shown in the ¬gure is not too

surprising.

˜ Ug

To implement the algorithm, the matrix need not be formed explicitly. All

§

that is required is to be able to access one row of at a time. This row can be computed, §

˜S

used, and then discarded. In the following, the -th row of is denoted by . The £ S

algorithm is row-oriented and all vectors denote row vectors. It is adapted from the ILU(0)

˜

¨

factorization of a sparse matrix, i.e., Algorithm 10.4, but it actually computes the

˜

§

¦

factorization instead of an or factorization. The main difference with Algorithm

10.4 is that the loop in line 7 is now restricted to because of symmetry. If only the S

˜ `3 ¥ ¦

elements are stored row-wise, then the rows of which are needed in this loop are

˜

´

¦

not directly available. Denote the -th row of by . These rows are accessible by

¥

adding a column data structure for the matrix which is updated dynamically. A linked

list data structure can be used for this purpose. With this in mind, the IC(0) algorithm will

have the following structure.

m5• &1 ¨0 ¤ £ ˜ 75§"0§ ¡w¤h¨ n¦ ¤ £¢

¶¢ ¦ ¢

8˜ H ˜ ¤ ¨°¡ $ H&’ ”8 ” &’ ¤ ©¦§

©8 $ $

¨ ¡¨ ¨¦ ¤

§ © §¥ c#&) $

©1

¢

¡

£ ff

1. Initial step: Set ,

f¢

s® D¬BBp"% H ª

ff G

¬ ¬%

2. For Do: „

%

3. Obtain all the nonzero inner products

« § S £ §& SS ¨i% D¬BD% H E% G n¥ %„a S &£ % £ ` S w

¬ ¬

4. , and

G

¦ ( ¢ S ¡ i ¥ „¬BD¬B% a `

5. Set

˜ a ` G ¨ ´ % ¬ G

6. For and if Do:

a 7`

7. Extract row

¢ S w S

8. Compute

¦ a ¥ 0 ` ¬ ¬ G m¥

%„B D¬B%

9. For and if Do:

a 7` %

´ S ¡¨ S S

10. Compute

11. EndDo

12. EndDo

SS SS

13. Set ,

S¢

G

14. EndDo

´

Note that initially the row in the algorithm is de¬ned as the ¬rst row of . All vectors

f

in the algorithm are row vectors.

The step represented by lines 3 and 4, which computes the inner products of row

number with all previous rows, needs particular attention. If the inner products

˜£ ˜ ¬¬¬ ˜

S £ S £ % DBB% S £ « £

£ %S

f fd

are computed separately, the total cost of the incomplete factorization would be of the

«®

order of steps and the algorithm would be of little practical value. However, most of

these inner products are equal to zero because of sparsity. This indicates that it may be

possible to compute only those nonzero inner products at a much lower cost. Indeed, if is ¨

S

R

¨ ra G R7`

®¯ ¨

the column of the inner products , then is the product of the rectangular

˜S ¤£¨ % D¬BB% ˜ £ ¨

G

¬¬f

matrix whose rows are by the vector , i.e., £

b

S S

'd

f d

f ¦G) ¤ v0¥ £

‘)

¨ ¬ S £ 'd S

f

This is a sparse matrix-by-sparse vector product which was discussed in Section 10.5. It

is best performed as a linear combination of the columns of which are sparse. The S

d

f

only dif¬culty with this implementation is that it requires both the row data structure of

and of its transpose. A standard way to handle this problem is by building a linked-list data

structure for the transpose. There is a similar problem for accessing the transpose of ,

as mentioned earlier. Therefore, two linked lists are needed: one for the matrix and the

other for the matrix. These linked lists avoid the storage of an additional real array for

the matrices involved and simplify the process of updating the matrix when new rows

are obtained. It is important to note that these linked lists are used only in the preprocessing

phase and are discarded once the incomplete factorization terminates.

’ ”8 ” ’ ¤ §

$ © } ”vG

’ ˜8’ ’ ”p¥ &)0©

8 1 ¨