Parallelization of SLE solution with sparse matrices

Hello,

in the framework of a research project I implemented a 2D CFD-Code in Matlab. The implicit formulation of the linearized Navier-Stokes-Equations ends up in a system of linear equations represented by a sparse unsymmetric matrix on the left hand side. For its solution at every iteration step the built-in function mldivide (A\b) has been used which basically uses the UMFPACK-Library. I want to start a quite large calculation campaign on a linux cluster (https://www.lrz.de/services/compute/linux-cluster/overview/) with manifold geometry and parameter variations. As calculation time is an important issue here I wanted to profit from the cluster architecture, i.e. the possibility to execute parallelized code. Unfortunately UMFPACK is a sequential solver and the parallel computing toolbox does not offer the possibility to treat the sparse matrix solution in a parallel manner. As the number of rows/columns of the square matrix is about 200,000 the treatment as full matrix does not appear reasonable in terms of memory usage (sparse 13MB, full 320GB). So I have been looking for other libraries (GNU) that could help me to speed up the calculations, but without striking success: The library MUMPS offers a Matlab interface but only runs in single-thread mode when combined with Matlab. The library PARDISO was the next candidate to be examined. Unfortunately it requires glibc version 2.14 on linux where the available cluster provides only version 2.11. The third potential option is the library WSMP from IBM. But there is no Matlab interface available for it. Trying to implement an interface by myself exceeds my current programming abilities and therefore would lead to exorbitant time expense.

So here is my question if anybody has another idea to solve the sparse matrix (200,000 x 200,000) in a parallelized manner in order to benefit from cluster hardware?

Thank you very much in advance for any tips!

Korbi

[Matt J]

Does A change over the iterations, or just b? If A has to be recomputed in each iteration, are you sure that A\b operation is the bottleneck? The recomputation of A might take similar time or greater.

[Korbi]

Hi, thanks for you interest. Yes, the matrix A changes over the iterations. A time analysis revealed that more than 90% of computation time is spent on the solution of the SLEs consisting of one big system (200k x 200k) and three smaller ones (60k x 60k). These systems correspond to nonlinear equations which are physically coupled, but have to be linearized and lead to the aforementioned SLEs. The calculation and assembly of the respective matrices A performs quite fast as the Matlab internal parallel treatment seems to be very efficient with vectorized code.

Answer 1

[Matt J]

the parallel computing toolbox does not offer the possibility to treat the sparse matrix solution in a parallel manner

Maybe the mldivide operation itself cannot be parallelized, but since you have multiple iterations to do, why wouldn't it be helpful to parallelize the iterations, like in the following

    c=zeros(4,10);

    parfor i=1:10

        A=sparse(rand(4));
        b=sparse(rand(4,1));
        c(:,i)=A\b;

    end

[Matt J]

It might also be worth trying to solve each system iteratively, using PCG for example. You could have many PCG instances running in parallel, using parfor as above.

Which solver is more time-efficient (iterative or analytical) can depend a lot on the specific sparsity pattern.

[Korbi]

Thanks for your suggestion. The matrix of the subsequent SLE unfortunately depends on the previous solution x (c in your code). So, A(i+1) = function(x(i)). Therefore a parallelization with parfor at this point is not possible.

[Rylan]

An interesting question.