Out of memory error

I am solving a system of linear algebraic equations with a sparse complex matrix in a loop. At each pass through the loop I form the sparse complex matrix out of real sparse matrices. First pass through the loop is okay, on the second pass through the loop the memory allocation goes crazy. Yet, the types of the matrices apparently have not changed. No clue what is going on at the moment.

With 1.5.1 on Windows 10.

Has anyone seen this?

A bit more information:

ERROR: LoadError: OutOfMemoryError()                                                                                                                                
Stacktrace:                                                                                                                                                         
 [1] Array at .\boot.jl:406 [inlined]                                                                                                                               
 [2] _allocres at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:233 [inlined]                               
 [3] _noshapecheck_map(::SparseArrays.HigherOrderFns.var"#3#4"{typeof(*),SparseArrays.HigherOrderFns.var"#17#20"{Complex{Float64},SparseArrays.HigherOrderFns.var"#21#24"}}, ::SparseArrays.SparseMatrixCSC{Float64,Int64}) at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:164                                                                                                                                                                    
 [4] _shapecheckbc at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:1025 [inlined]                          
 [5] _copy at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:1015 [inlined]                                  
 [6] _copy at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:1020 [inlined]                                  
 [7] copy at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.5\SparseArrays\src\higherorderfns.jl:1011 [inlined]                                   
 [8] materialize at .\broadcast.jl:837 [inlined]                                                                                                                    
 [9] broadcast_preserving_zero_d at .\broadcast.jl:826 [inlined]                                                                                                    
 [10] * at .\arraymath.jl:52 [inlined]    

on this line

U1 .= ((-omega^2+1im*a0)*M + (1+1im*a1)*K)\F;
1 Like

From “Out of memory error”, I take it you eventually run out of memory.

Does it increase on each iteration, or is it already the 2nd iteration that causes problems?

Can you test smaller problems? If the loops complete, is all the memory freed? Or when the function returns? Or do you have to run again for the previous call to free after the second?

I’ve been struggling with excessive memory use as well. I can easily get 100 GB of memory consumed, and GC.gc() that takes half a minute to run without actually freeing anything, despite redefining all returned globals as nothing. Only when I call the offending function again is the memory from the previous call freed, and I still don’t really have a clue about what is going on. But my problem doesn’t have sparse matrices, so it’s probably unrelated.
And instead of OutOfMemoryError(), Linux just kills the Julia process when I run out of memory.

2 Likes

It looks like this:First pass: high water mark around 5.6 GB. Second pass: I hit the wall.

The loop looks like this

    for k in 1:length(frequencies)
        omega = 2*pi*frequencies[k];
        @show typeof(M)
        @show typeof(K)
        #U1 .= (-omega^2*M + 1im*omega*C + K)\F;
        U1 .= ((-omega^2+1im*a0).*M + (1+1im*a1).*K)\F;
        @show typeof(M)
        @show typeof(K)
        push!(frf, U1[sensorndof]...)
        @show frf
        print(".")
    end

The matrices are

typeof(M) = SparseArrays.SparseMatrixCSC{Float64,Int64}                                                                                                             
typeof(K) = SparseArrays.SparseMatrixCSC{Float64,Int64} 

I also checked typeof((-(omega ^ 2) + (1im) * a0) .* M + (1 + (1im) * a1) .* K) = SparseArrays.SparseMatrixCSC{Complex{Float64},Int64} .

What’s F? Might be that the \ forces your matrices to be dense?

1 Like

F is a real vector. BTW, first pass requires only about a gig of memory. Second, way more.

Real as in the abstract type? If I construct these, I get the following:

julia> a = sparse([1.0im 2.0im 3. 4.])
1Ă—4 SparseMatrixCSC{Complex{Float64},Int64} with 4 stored entries:
  [1, 1]  =  0.0+1.0im
  [1, 2]  =  0.0+2.0im
  [1, 3]  =  3.0+0.0im
  [1, 4]  =  4.0+0.0im

julia> F = Real[1. 2. 3. 4.]
1Ă—4 Array{Real,2}:
 1.0  2.0  3.0  4.0

julia> a \ F
4Ă—4 Array{Complex{Float64},2}:
 0.0-1.0im  0.0-2.0im  0.0-3.0im  0.0-4.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im

Note that the result is a dense array and that F has to be a row vector/matrix. The same happens with Float64:

julia> F = [1. 2. 3. 4.]
1Ă—4 Array{Float64,2}:
 1.0  2.0  3.0  4.0

julia> a \ F
4Ă—4 Array{Complex{Float64},2}:
 0.0-1.0im  0.0-2.0im  0.0-3.0im  0.0-4.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im
 0.0+0.0im  0.0+0.0im  0.0+0.0im  0.0+0.0im

Also:

What does @code_warntype say? How big is U1[sensorndof]? Maybe append! is more appropriate than push!.

1 Like

Could it be that you’re multiplying by Inf or NaN ?

That would turn all the zeros into NaN, effectively giving you a dense matrix stored as CSC.

2 Likes

I finally checked the last thing I suspected, the frequencies. Turns out I generated the log space numbers wrong so that the second was Inf. Excellent guess from @Per!

1 Like