zillions of allocations in eig?

[stevengj]

julia> @time inv(rand(1000,1000));
  0.147837 seconds (21 allocations: 23.393 MB, 2.14% gc time)

julia> @time eig(rand(1000,1000));
  4.988959 seconds (9.31 k allocations: 83.650 MB, 0.48% gc time)

The 10,000 allocations required in eig look suspicious, considering that all the work should occur in LAPACK. Is there a type instability in eig somewhere? cc @andreasnoack

[ViralBShah]

The issue is on 0.4 and master. I am pretty sure there was a PR for a type stable eig.

[andreasnoack]

The issue is discussed in #4006. Making eig type stable is possible but problematic since the result the non-symmetric real eigenvalue problem is either real or complex. We could

1. Do the same as sqrt and just fail on complex results. This would probably annoy some people but it would also be expensive as it would then require the real eigenvalue problem to be solved in complex arithmetic whenever there are complex eigenvalues.
2. Always return complex results. Even when the values are actually real. This is probably the only sane solution but it also wastes memory when the result is real.

It should also be noted that the type stability doesn't matter much for non-small problems. Most of the time is spent in LAPACK. Compare

julia> A = rand(1000,1000);

julia> @time eig(A);
  2.258325 seconds (9.33 k allocations: 76.175 MB, 1.38% gc time)

julia> @time LAPACK.geev!('N', 'V', copy(A));
  2.409304 seconds (20 allocations: 15.534 MB, 3.63% gc time)

julia> @time eig(A);
  2.224663 seconds (9.33 k allocations: 76.175 MB, 0.28% gc time)

julia> @time LAPACK.geev!('N', 'V', copy(A));
  2.209915 seconds (20 allocations: 15.534 MB, 0.17% gc time)

julia> @time eig(A);
  2.229527 seconds (9.33 k allocations: 76.175 MB, 0.87% gc time)

julia> @time LAPACK.geev!('N', 'V', copy(A));
  2.277948 seconds (20 allocations: 15.534 MB, 3.28% gc time)

[eschnett]

There could be a special version of eig -- or a Val{Bool} typed argument -- that chooses between the two cases.

[ViralBShah]

I think the reason to do something about this is exactly for the small linear algebra case.

[tkelman]

My #13362 is marginally related (dealing with the number of outputs at the lapack wrapper level) and I've been meaning to get back to that. Half-joking, maybe we could return the imaginary component as a nullable?

[stevengj]

@andreasnoack, even if the result of eig is type-unstable, I don't understand why it needs ten thousand allocations. I'm concerned that there is some type instability inside eig, i.e. there is some loop with a "boxed" variable that is getting heap-allocated over and over again.

[yuyichao]

The loop in eigfact! seems to be allocating new arrays and is also not inlined...

if103:                                            ; preds = %L51
  store %jl_value_t* %103, %jl_value_t** %81, align 8
  store %jl_value_t* inttoptr (i64 140328436789640 to %jl_value_t*), %jl_value_t** %7, align 8
  store %jl_value_t* inttoptr (i64 140328434551472 to %jl_value_t*), %jl_value_t** %82, align 8
  %285 = call %jl_value_t* @jl_box_int64(i64 signext %storemerge148)
  store %jl_value_t* %285, %jl_value_t** %83, align 8
  %286 = call %jl_value_t* @julia__unsafe_getindex_50342(%jl_value_t* inttoptr (i64 140328442456464 to %jl_value_t*), %jl_value_t** %7, i32 4)
  store %jl_value_t* %286, %jl_value_t** %63, align 8
  store %jl_value_t* %103, %jl_value_t** %64, align 8
  br i1 %272, label %L55, label %if106

[andreasnoack]

Yesterday, I looked at the output of @code_warntype and nothing except the
return was flagged but I agree that the number of allocations suggests
something like the problem you describe. I'll take a closer look.

On Saturday, June 4, 2016, Steven G. Johnson [email protected]
wrote:

@andreasnoack https://github.com/andreasnoack, even if the result of
eig is type-unstable, I don't understand why it needs ten thousand
allocations. I'm concerned that there is some type instability inside
eig, i.e. there is some loop with a "boxed" variable that is getting
heap-allocated over and over again.

—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
#16751 (comment),
or mute the thread
https://github.com/notifications/unsubscribe/AAe0qQhk9JWQzqTHCjd5lcfnH1YglP0dks5qIXDjgaJpZM4It6xO
.

[andreasnoack]

Yes. It's creating a ton of temporaries when there are complex values. I'll prepare a PR.