# How to (efficiently) find the set of maxima of an array?

**URL:** <https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423>\
**Category:** Performance\
**Tags:** question, arrays\
**Created:** [December 21, 2021, 8:58am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423 "2021-12-21T08:58:33Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [December 21, 2021, 8:58am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/1 "2021-12-21T08:58:33Z")

</div>

I want to find the set of maxima of an array. More precisely, I want the set of their indices. For example, from

> A = [1 2 3 4 1 2 3 4 1 2 3 4]

I want to obtain something similar to

> 3-element Vector{Int64}:  
> 4  
> 8  
> 12

I was able to get this result by using

> getindex.(findall(A .== maximum(A)),2)
> 
> 3-element Vector{Int64}:  
> 4  
> 8  
> 12

but I would like to have something more efficient. In my code, this instruction needs to be repeated millions of times. My solution is slow and it allocates too much memory.

---

<div class="post-metadata">

**Author:** ![N5N3](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/n5n3/32/17663_2.png) [@N5N3](https://discourse.julialang.org/u/N5N3)\
**Post date:** [December 21, 2021, 9:47am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/2 "2021-12-21T09:47:20Z")

</div>

If you can give a upper bound of output length, and it’s not large (\<= 20 for example), you might try

```julia
findmaxs(A::AbstractArray) = begin
    inds = Vector{eltype(keys(A))}(undef, 20)
    maxval = typemin(eltype(A))
    n = 0
    @inbounds for i in keys(A)
        Ai = A[i]
        Ai < maxval && continue
        if Ai > maxval
            maxval = Ai
            n = 1
            inds[n] = i
        else
            inds[n+=1] = i
        end
    end
    resize!(inds, n)
end

```

For inputs like `A = rand(1:1000, 100, 100)`, this save about 18% time cost on my desktop. (I have to say the space for optimization is quite small). If you accept linear index, replacing `keys` with `eachindex` seems better.

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [December 21, 2021, 9:55am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/3 "2021-12-21T09:55:51Z")

</div>

Using StaticArrays gives your example a boost:

```julia
f1(A) = getindex.(findall(A .== maximum(A)),2)

using StaticArrays
A = SA[1 2 3 4 1 2 3 4 1 2 3 4]
@btime f1($A) # 52 ns (2 allocations: 192 bytes)

A = [1 2 3 4 1 2 3 4 1 2 3 4]
@btime f1($A) # 114 ns (4 allocations: 304 bytes)

```

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [December 21, 2021, 10:58am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/4 "2021-12-21T10:58:23Z")

</div>

Thanks for your help! Your algorithm works quite well. Since I am looping this instruction millions of times, this really makes a difference.

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [December 21, 2021, 11:00am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/5 "2021-12-21T11:00:45Z")

</div>

The problem is I cannot use a static array in my code. The array updates at each iteration, and in each iteration I must search for the set of maxima!

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [December 21, 2021, 11:24am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/6 "2021-12-21T11:24:35Z")

</div>

> [@N5N3](#):
>
> If you can give a upper bound of output length, and it’s not large (\<= 20 for example)

I think this restriction could be lifted as follows:

```julia
function findmaxs2(A::AbstractArray)
    inds = eltype(keys(A))[]
    maxval = typemin(eltype(A))
    n = 0
    push!(inds, CartesianIndex(0,0))
    @inbounds for i in keys(A)
        Ai = A[i]
        Ai < maxval && continue
        if Ai > maxval
            maxval = Ai
            inds[1] = i
            n = 1
        else
            n += 1
            insert!(inds, n, i)
        end
    end
    resize!(inds, n)
end

A = rand(1:1000, 100, 100);
findmaxs(A) == findmaxs2(A) # true

@btime findmaxs($A) # 6.060 μs (1 allocation: 400 bytes)
@btime findmaxs2($A) # 6.180 μs (3 allocations: 880 bytes)

```

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [December 21, 2021, 11:55am UTC](https://discourse.julialang.org/t/how-to-efficiently-find-the-set-of-maxima-of-an-array/73423/7 "2021-12-21T11:55:56Z")

</div>

You may want to consider preallocating the array if that is being called millions of times. And return the number of elements as well, or a view of the preallocated arrays. Just change @N5N3 solution to:

```julia
julia> findmaxs(A::AbstractArray;inds=Vector{eltype(keys(A))}(undef,20)) = begin
           maxval = typemin(eltype(A))
           n = 0
           @inbounds for i in keys(A)
               Ai = A[i]
               Ai < maxval && continue
               if Ai > maxval
                   maxval = Ai
                   n = 1
                   inds[n] = i
               else
                   inds[n+=1] = i
               end
           end
           return @view(inds[1:n])
       end
findmaxs (generic function with 1 method)

julia> A = rand(1:1000, 100, 100);

julia> inds=zeros(eltype(keys(A)),20); # now this can be of any size

julia> @btime findmaxs($A,inds=$inds)
  13.148 μs (0 allocations: 0 bytes)
10-element view(::Vector{CartesianIndex{2}}, 1:10) with eltype CartesianIndex{2}:
 CartesianIndex(21, 7)
 CartesianIndex(98, 23)
 CartesianIndex(63, 37)
 CartesianIndex(24, 49)
 CartesianIndex(54, 60)
 CartesianIndex(55, 62)
 CartesianIndex(19, 75)
 CartesianIndex(52, 76)
 CartesianIndex(63, 95)
 CartesianIndex(21, 99)

```

This is not necessarily faster for _one_ call of the function, but for a loop, it will be.
