# Matrix multiplication but sum of min() instead of sum of product

**URL:** <https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174>\
**Category:** General Usage\
**Tags:** statistics, network-analysis\
**Created:** [September 14, 2021, 9:01pm UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174 "2021-09-14T21:01:28Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![sardinecan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sardinecan/32/29141_2.png) [@sardinecan](https://discourse.julialang.org/u/sardinecan)\
**Post date:** [September 14, 2021, 9:01pm UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/1 "2021-09-14T21:01:28Z")

</div>

Hi everyone,

I’m new to Julia and I’m using it to do some network analysis on historical data. As I have valued bimodal data for my graphs, I need to “normalize" those values.

Basically, I need to do some “matrix multiplications” but instead of the sum of product, I need the sum of min()

Here is a small example :

```julia
# a Matrix m
m = [0 1 2
     3 4 5]

# and the transpose of m
mt = [0 3
      1 4
      2 5]

```

Instead of the sum of product `m * mt`

```julia
result = [5 14
         14 50]
# (0*0 + 1*1 + 2*2) (0*3 + 1*4 + 2*5)
# (3*0 + 4*1 + 5*2) (3*3 + 4*4 + 5*5)

```

I need the sum of the min()

```julia
result = [3 3
          3 12
# (min(0,0) + min(1, 1) + min(2, 2)) (min(0, 3) + min(1, 4) + min(2,5))
# (min(3,0) + min(4, 1) + min(5, 2)) (min(3, 3) + min(4, 4) + min(5, 5))

```

I tried something like

```julia
[min(m[i], mt[j]) for i in eachrow(m), j in eachcol(mt)]

```

but it doesn’t work and I think I’m missing something, certainly with the syntax…

Thank you in advance for your help, Best,  
Josselin.

---

<div class="post-metadata">

**Author:** ![vkv](https://avatars.discourse-cdn.com/v4/letter/v/d9b06d/32.png) [@vkv](https://discourse.julialang.org/u/vkv)\
**Post date:** [September 14, 2021, 9:19pm UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/2 "2021-09-14T21:19:16Z")

</div>

Try this

```julia
result = zeros(size(m,1), size(mt,2))

for i in 1:size(m, 1)
    for j in 1:size(mt, 2)
        for k in 1:size(mt, 1)
           result[i,j] += min(m[i,k], mt[k,j])
        end
    end
end
```

---

<div class="post-metadata">

**Author:** ![viraltux](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/viraltux/32/15236_2.png) [@viraltux](https://discourse.julialang.org/u/viraltux)\
**Post date:** [September 14, 2021, 10:42pm UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/3 "2021-09-14T22:42:09Z")

</div>

one-liners! welcome @sardinecan!

```julia
[sum(mapslices(minimum,[i j],dims=2)) for i in eachrow(m),j in eachcol(mt)]

2×2 Matrix{Int64}:
 3 3
 3 12

```

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [September 15, 2021, 10:55am UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/4 "2021-09-15T10:55:44Z")

</div>

This sounds similar to  
[https://github.com/TensorBFS/TropicalGEMM.jl](https://github.com/TensorBFS/TropicalGEMM.jl)  
although not the same. LoopVectorization or Tullio could probably speed this use case up significantly.

---

<div class="post-metadata">

**Author:** ![sardinecan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sardinecan/32/29141_2.png) [@sardinecan](https://discourse.julialang.org/u/sardinecan)\
**Post date:** [September 15, 2021, 11:31am UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/5 "2021-09-15T11:31:03Z")

</div>

@vkv, @viraltux it works perfectly ! Thank you very much !!

---

<div class="post-metadata">

**Author:** ![sardinecan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sardinecan/32/29141_2.png) [@sardinecan](https://discourse.julialang.org/u/sardinecan)\
**Post date:** [September 15, 2021, 11:32am UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/6 "2021-09-15T11:32:24Z")

</div>

@baggepinnen thank you for the package, I found the Tensor pkg for python but I missed this one !!

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [September 15, 2021, 11:32am UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/7 "2021-09-15T11:32:48Z")

</div>

You can take the code from [https://github.com/JuliaSIMD/LoopVectorization.jl/blob/master/docs/src/examples/matrix\_multiplication.md](https://github.com/JuliaSIMD/LoopVectorization.jl/blob/master/docs/src/examples/matrix_multiplication.md)

```julia
function A_mul_B!(C, A, B)
    @turbo for n ∈ indices((C,B), 2), m ∈ indices((C,A), 1)
        Cmn = zero(eltype(C))
        for k ∈ indices((A,B), (2,1))
            Cmn += A[m,k] * B[k,n]
        end
        C[m,n] = Cmn
    end
end

```

Just replace `A[m,k] * B[k,n]` with `min(A[m,k], B[k,n])`. It’s about as fast as a BLAS call. If you replace `@turbo` with `@tturbo`, you might get even further speedup.

**Edit:** I see that you are operating on really tiny matrices, so LoopVectorization may be overkill. Perhaps try StaticArrays, if you need performance.

---

<div class="post-metadata">

**Author:** ![learned\_fool](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/learned_fool/32/21934_2.png) [@learned\_fool](https://discourse.julialang.org/u/learned_fool)\
**Post date:** [September 15, 2021, 11:36am UTC](https://discourse.julialang.org/t/matrix-multiplication-but-sum-of-min-instead-of-sum-of-product/68174/8 "2021-09-15T11:36:21Z")

</div>

Here is a comparison of the results from LV and a normal implementation; LV was about 30x faster on my system. the threaded version was 100x better, but then it’s subjective.

my `versioninfo()` being

```julia
julia> versioninfo()
Julia Version 1.6.2
Commit 1b93d53fc4 (2021-07-14 15:36 UTC)
Platform Info:
  OS: Linux (x86_64-pc-linux-gnu)
  CPU: Intel(R) Core(TM) i7-9750H CPU @ 2.60GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-11.0.1 (ORCJIT, skylake)
Environment:
  JULIA_NUM_THREADS = 12

```

```julia
julia> m = rand(100,1000);

julia> mt = rand(1000,100);

julia> result = zeros(100, 100);

julia> using BenchmarkTools, LoopVectorization
    
julia> function min_loops_lv!(m, mt, result)
       
          @turbo for i in 1:size(m ,1)
               for j in 1:size(mt, 2)
                   for k in 1:size(mt, 1)
                       result[i, j] += min(m[i,k], mt[k, j])
                   end
               end
           end
       end
min_loops_lv! (generic function with 1 method)

julia> function min_loops_tlv!(m, mt, result)
       
          @tturbo for i in 1:size(m ,1)
               for j in 1:size(mt, 2)
                   for k in 1:size(mt, 1)
                       result[i, j] += min(m[i,k], mt[k, j])
                   end
               end
           end
       end
min_loops_tlv! (generic function with 1 method)

julia> function min_loops_normal!(m, mt, result)
       
           for i in 1:size(m ,1)
               for j in 1:size(mt, 2)
                   for k in 1:size(mt, 1)
                       result[i, j] += min(m[i,k], mt[k, j])
                   end
               end
           end
       end
min_loops_normal! (generic function with 1 method)

julia> @benchmark min_loops_normal!($m, $mt, $result)
BenchmarkTools.Trial: 230 samples with 1 evaluation.
 Range (min … max): 21.075 ms … 23.152 ms ┊ GC (min … max): 0.00% … 0.00%
 Time (median): 21.735 ms ┊ GC (median): 0.00%
 Time (mean ± σ): 21.784 ms ± 331.447 μs ┊ GC (mean ± σ): 0.00% ± 0.00%

              ▃ ▄ ▇▄ ▆▄▆▃█▄▄▂▄▅                                 
  ▃▁▁▁▃▁▆▃▅▅▇▇█▆█▅█████████████▅▆▃▆▄▄▅▅▆▄▃▃▃▁▄▁▁▁▃▃▄▃▃▁▁▁▁▃▃▃▄ ▄
  21.1 ms Histogram: frequency by time 22.9 ms <

 Memory estimate: 0 bytes, allocs estimate: 0.

julia> result_lv = zeros(100, 100);

julia> @benchmark min_loops_lv!($m, $mt, $result_lv)
BenchmarkTools.Trial: 6836 samples with 1 evaluation.
 Range (min … max): 664.470 μs … 1.721 ms ┊ GC (min … max): 0.00% … 0.00%
 Time (median): 711.900 μs ┊ GC (median): 0.00%
 Time (mean ± σ): 728.216 μs ± 59.243 μs ┊ GC (mean ± σ): 0.00% ± 0.00%

       ▅ █ ▂                                                 
  ▃▂▁▁▃█▂▂█▆▂▅█▄▃▃▂▂▂▂▂▂▂▂▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁ ▂
  664 μs Histogram: frequency by time 970 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

julia> result_tlv = zeros(100, 100);

julia> @benchmark min_loops_tlv!($m, $mt, $result_tlv)
BenchmarkTools.Trial: 10000 samples with 1 evaluation.
 Range (min … max): 165.561 μs … 404.591 μs ┊ GC (min … max): 0.00% … 0.00%
 Time (median): 194.182 μs ┊ GC (median): 0.00%
 Time (mean ± σ): 201.954 μs ± 29.925 μs ┊ GC (mean ± σ): 0.00% ± 0.00%

    ▆▆ █                                                        
  ▂▁██▂▆█▂▂▃▂▂▂▂▃▄▄▂▂▂▂▂▂▂▂▂▂▄▇▇█▃▃▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▁▁▁▁▁▁▁▁ ▂
  166 μs Histogram: frequency by time 278 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

```