Max size of set linear equations to solve? (X=AX+B)

Question

This question has also been asked at Stack Overflow.

This is a very general question regarding the maximum size of a set of linear equations to be solved by today's fastest hardware, in the form:

$$X = AX + B$$

Where we solve for $X$, and

• $A$: $N\times N$ a sparse matrix of floats;

• $B$: $N$-vector of floats.

This becomes $X(I-A) = B$ which is best solved using factorization (and not matrix inversion) as I read here.

Do you know or have a reference to a benchmark or paper which gives some maximum value for $N$ with today's fastest hardware? Assume that I can have that hardware, this is a sort of mental what-if exercise. Most benchmarks I have seen use $N < 10000$. I am thinking about $N>10^7$ or more to be processed within a month.

Please consider not only the computational dimension but also storage for $A$. It can be a problem, e.g., assuming $N = 10^6$, storage would be $4\times 10^{12}$ bytes $\approx$ 4 Terabytes for a totally dense matrix, which is about manageable I guess.

Lastly, can the method to solve the system be parallelized so that I can make the assumption that with parallelization $N$ can become pretty large?

later addition to question:

One possible application of this is: consider N products from factories (e.g. screws, an electric motor, an engine, a car). Most products require parts from other products and are used as parts for other more complex products in addition to some quantities ending in the shops - e.g. an electric motor.

For the $i^{th}$ product we have: $X_i = a_{i,1} X_1 + a_{i,2} X_2 + \cdots + a_{i,N} X_N + B_i$ meaning that the quantity of $X_i$ we aim to produce requires $a_{i,1} X_1$ items from product $X_1$, to $N$, plus $B_i$ which is what ends up in shops for general consumption.

Our aim is to find the solution to this problem given final demands ($B$) and "inter-dependancy" quantities $A$.

I guess $A$ is sparse but this probably depends on applications, e.g. Toyota states that a basic car consists of 30,000 other products (all broken down to screws). This is why I mention sparse but I also go for the worst case of "mild-denseness".

Answer 1

For your applications, dense matrix techniques such as the different factorization techniques will not work. Plus your matrix is sparse anyway so you should be using an iterative solver (e.g. Krylov subspace methods). These easily handle problems of the scale you are interested in and there are several sparse matrix solver libraries implemented for parallel architectures (e.g. Petsc, Trilinos, others).

There are several advantages to using sparse matrix techniques:

1. Storage. Who wants to store a dense matrix? Do you have terabytes of RAM waiting to be used? Most people don't, and you should not want to have to use disk space because that will cripple any sort of algorithm performance. You've already stated that your matrix A is sparse, which generally means that total storage space should only be $O(N)$. So it should only take megabytes to store your matrix and required vectors, not terabytes.

2. Time. Dense matrix factorizations take $O(N^3)$ operations, even for sparse matrices. On a serial processor, a system of equations with $N>10^7$ can easily take on the order of years to compute a factorization. Even when parallelized you still have to deal with storage and communication costs. Sparse methods rely on matrix-vector multiplication, which is an $O(N^2)$ operation for dense matrices but a $O(N)$ operation for sparse matrices. And usually it requires less than $O(N)$ iterations to converge. This will take much less time in general (and can take even less time when a good preconditioner is used).

3. Parallelization. Because iterative methods rely on matrix-vector multiplication and vector-vector addition, they are highly parallelizable (especially on shared memory architectures, but also distributed memory too). Many libraries already have parallel support (again see here) for various architectures.