Reduce allocs for beamformers

[microcassidy]

Hello,
Thanks for making such an easy-to-use library.

I have made some changes to reduce the allocations for the Capon, MUSIC, and Bartlett beamformer functions.
Perhaps it might be useful to you too.

All the best,

M

[mchitre]

Thanks for the contribution. Will review changes shortly.

[microcassidy]

Another thing that I noticed though was that Precompile tools runs the test workload for the package. The time to first execution is longer because test is a dependency of both the library and the test. If you remove the test decency from the source, test is only called on running test, so it ends up being quicker for development and compilation on both sides. The only thing that's slower is running test

[microcassidy]

https://github.com/org-arl/SignalAnalysis.jl/actions/runs/11425445079/job/31864162141#step:5:152

[mchitre]

I added the tests in precompile a while back to force pre-compilation caching of commonly used methods. I do recognize that it slows down the installation of the package, and in case of running tests, runs them twice. I intend to fix that by changing the running of tests during package init to just loading the precompile statements instead.

[mchitre]

https://github.com/org-arl/SignalAnalysis.jl/actions/runs/11425445079/job/31864162141#step:5:152

Will fix shortly.

[mchitre]

Looks great. Made minor suggestions to replace allocation with zeros followed by filling array with direct computation of array.

[microcassidy]

Thanks. I'll take another look when I've got a chance.

…

On Tue, Oct 22, 2024, 17:48 Mandar Chitre ***@***.***> wrote: ***@***.**** requested changes on this pull request. Looks great. Made minor suggestions to replace allocation with zeros followed by filling array with direct computation of array. ------------------------------ In src/array.jl <#46 (comment)> : > + RS = R * sv + pseudospectrum = zeros(Float64, size(sv, 2)) + pseudospectrum .= abs.(dot.(eachcol(sv), eachcol(RS))) Would this work? ⬇️ Suggested change - RS = R * sv - pseudospectrum = zeros(Float64, size(sv, 2)) - pseudospectrum .= abs.(dot.(eachcol(sv), eachcol(RS))) + abs.(dot.(eachcol(sv), eachcol(R * sv))) Verification and benchmark: julia> sv = randn(64,256); julia> R = let A = randn(64,64); (A + A')/2; end julia> abs.(dot.(eachcol(sv), eachcol(R * sv))) ≈ [abs.(sv[:,j]' * R * sv[:,j]) for j ∈ 1:size(sv,2)]true julia> julia> @Btime abs.(dot.(eachcol($sv), eachcol($R * $sv))); 12.584 μs (6 allocations: 130.14 KiB) julia> @Btime [abs.($sv[:,j]' * $R * $sv[:,j]) for j ∈ 1:size($sv,2)]; 107.792 μs (1539 allocations: 434.06 KiB) ------------------------------ In src/array.jl <#46 (comment)> : > + pseudospectrum = zeros(Float64, size(sv, 2)) + RiS = inv(Hermitian(R)) * sv + pseudospectrum .= 1 ./ abs.(dot.(eachcol(sv), eachcol(RiS))) ⬇️ Suggested change - pseudospectrum = zeros(Float64, size(sv, 2)) - RiS = inv(Hermitian(R)) * sv - pseudospectrum .= 1 ./ abs.(dot.(eachcol(sv), eachcol(RiS))) + 1 ./ abs.(dot.(eachcol(sv), eachcol(inv(Hermitian(R)) * sv))) Verification and benchmark: julia> 1 ./ abs.(dot.(eachcol(sv), eachcol(inv(Hermitian(R)) * sv))) ≈ [1 ./ abs.(sv[:,j]' * inv(R) * sv[:,j]) for j ∈ 1:size(sv,2)]true julia> @Btime 1 ./ abs.(dot.(eachcol($sv), eachcol(inv($R) * $sv))); 59.042 μs (15 allocations: 194.88 KiB) julia> @Btime [1 ./ abs.($sv[:,j]' * inv($R) * $sv[:,j]) for j ∈ 1:size($sv,2)]; 12.500 ms (3843 allocations: 16.61 MiB) ------------------------------ In src/array.jl <#46 (comment)> : > + pseudospectrum = zeros(Float64, size(sv, 2)) + _, Q = eigen(Hermitian(R), 1:(size(R, 1)-opt.nsignals)) + QS = Q' * sv + pseudospectrum .= 1 ./ real.(dot.(eachcol(QS), eachcol(QS))) ⬇️ Suggested change - pseudospectrum = zeros(Float64, size(sv, 2)) - _, Q = eigen(Hermitian(R), 1:(size(R, 1)-opt.nsignals)) - QS = Q' * sv - pseudospectrum .= 1 ./ real.(dot.(eachcol(QS), eachcol(QS))) + _, Q = eigen(Hermitian(R), 1:(size(R, 1)-opt.nsignals)) + QS = Q' * sv + 1 ./ real.(dot.(eachcol(QS), eachcol(QS))) Verification: julia> F = eigen(R); ndx = sortperm(F.values; rev=true); Q = F.vectors[:, ndx[1+1:end]]; Q = Q * Q'; x = [1 ./ abs.(sv[:,j]' * Q * sv[:,j]) for j ∈ 1:size(sv,2)]; julia> _, Q = eigen(Hermitian(R), 1:(size(R, 1)-1)); QS = Q' * sv; y = 1 ./ real.(dot.(eachcol(QS), eachcol(QS))); julia> x ≈ ytrue — Reply to this email directly, view it on GitHub <#46 (review)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABZZTWONULCLLDL7ZZYEOLTZ4XYL5AVCNFSM6AAAAABQILXQNOVHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDGOBUGA2DCOBTGA> . You are receiving this because you authored the thread.Message ID: ***@***.***>

[mchitre]

Almost there. One minor formatting change requested to be reverted.