# OutOfMemoryError() while solving a linear system

I am trying to solve a linear system with the direct solver (backslash operator). My matrix is sparse,square and unsymmetric. The system is defined as `Ax=b`, where A is defined to be a SparseMatrixCSC with dimensions = 803800 x 803800 `nnz(A_matrix)=808613995`. the `b` vector is dense and defined to be like that. My A matrix is very close to be a banded matrix, but it is not, due to some elements on the corners, so I believe that the solver is LU factorizing the matrix and then solving the system. I am running my code on a 16 GB RAM laptop. My question is how to know the size limit of my linear system which I can solve on this machine given the sparsity pattern of my matrix?

Think of the memory consumption of such a system: you would need

``````julia> (808613995 * 16) / (2^30)
12.049287483096123
``````

gigabytes to represent the nonzero entries and the row indexes (`Float64`, and `Int`).

Maybe someone can suggest something, but I am not sure that doing this in 16GB is feasible.

Consider even with iterative methods with lower memory requirements, eg

but even these would need more RAM I assume.

2 Likes

Thanks !

This is true. Unless you have a good preconditioner, IterativeSolver will take ages and also massive memory.

Normally, a Krylov method only allocates few n-vectors, where n is the dimension of your matrix. With n = 803800, a n-vector only needs 6.13Mo. So iterative methods could be relevant in your case.

I can recommend you to test DQGMRES method, you can define the number of n-vectors that you want to use during the Arnoldi process. It’s a kind of limited memory variant of GMRES.

``````(x, stats) = dqgmres(A, b, memory=400, itmax=10000, verbose=true)
``````

You can adjust parameters (set verbose to false and/or reduce/increase memory).

Do you have your linear system available somewhere or a way to generate it that I can give a try too ?

1 Like

Thanks for your suggestion! I tried it but it is not actually working with A being sparse. I tried to convert A back to the dense form using `Matrix(A)` (which I know is not a good idea), and I got an OutOfMemoryError.

Have you tried with Pardiso.jl and MKL?

1 Like

@AhmedAlreweny, method should work with A being sparse. What is the error that is returned?

``````using SparseArrays
A = spzeros(100,100)
for i = 1 : 100
A[i,i] = 100
end
b = ones(100)
dqgmres(A, b)

([0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01  …  0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01],
Simple stats
solved: true
inconsistent: false
residuals:  [ 1.0e+01  5.4e-16 ]
Aresiduals: [ ]
status: solution good enough given atol and rtol
)
``````

Otherwise, you can use MUMPS that has been developed for linear systems of high dimensions.

Thanks again! This is the error I got

dqgmres(A_matrix, b_vector)
ERROR: MethodError: no method matching dqgmres(::SparseMatrixCSC{Float64,Int64}, ::Array{Float64,2})
Closest candidates are:
dqgmres(::AbstractArray{T<:AbstractFloat,2}, ::SparseVector{T<:AbstractFloat,Ti} where Ti<:Integer, ::Any…; kwargs…) where T<:AbstractFloat

Your right-hand side b is not a vector but a matrix `Array{Float64,2}` with one column.
You can convert it to a vector with `b = b[:]`.