# Matrix chain multiplication with arbitrary number of matrices

**URL:** <https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794>\
**Category:** General Usage\
**Tags:** matrices\
**Created:** [May 13, 2023, 2:29pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794 "2023-05-13T14:29:00Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![Frisus95](https://avatars.discourse-cdn.com/v4/letter/f/77aa72/32.png) [@Frisus95](https://discourse.julialang.org/u/Frisus95)\
**Post date:** [May 13, 2023, 2:29pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/1 "2023-05-13T14:29:01Z")

</div>

Hi everyone,

I want to multiply (or _contract_) an arbitrary number of matrices.

As an example, the _naive_ version with 4 matrices could be written as:

```julia
m1 = rand(10,10)
m2 = rand(10,10)
m3 = rand(10,10)
m4 = rand(10,10)

result = tr(*(m1,m2,m3,m4))

```

But what if the number of matrices must be established at run-time?

Furthermore, if I’m only interested in the trace of the computed matrix, as in the code snipped I wrote above, is there a way to compute just such a trace instead of the entire matrix in the first place?

I’m most probably looking for something very similar to the multi\_dot Python function.

Thank you so much in advance to everyone.

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [May 13, 2023, 2:33pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/2 "2023-05-13T14:33:59Z")

</div>

```julia
M = reduce(*, [m1, m2, m3, m4])
tr(M)

```

---

<div class="post-metadata">

**Author:** ![gdalle](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gdalle/32/27854_2.png) [@gdalle](https://discourse.julialang.org/u/gdalle)\
**Post date:** [May 13, 2023, 4:04pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/3 "2023-05-13T16:04:16Z")

</div>

Do all of these matrices have the same size? Do you care about efficiency?

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [May 13, 2023, 4:08pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/4 "2023-05-13T16:08:51Z")

</div>

> [@Frisus95](#):
>
> is there a way to compute just such a trace instead of the entire matrix in the first place?

You can compute the the trace with one fewer matrix multiplication.

```julia
julia> M = [rand(10,10) for i = 1:4];

julia> tr(prod(M))
630.566913491192

julia> @views dot(M[1]', prod(M[2:end])) # one fewer matrix multiply
630.5669134911918

```

If the matrices are large and you only need to _estimate_ the trace, then potentially there are even faster algorithms that avoid multiplying matrices entirely (by using only matrix–vector products). e.g. Hutchinson’s trace estimation algorithm.

> [@Frisus95](#):
>
> I’m most probably looking for something very similar to the multi\_dot Python function.

If all the matrices have the same size, you can just use `prod` as noted above.

If the matrices have different sizes, in principle there is a [matrix-chain ordering](https://en.wikipedia.org/wiki/Matrix_chain_multiplication) that minimizes the cost. `multi_dot` uses a heuristic algorithm to reduce this; an analogue in Julia is the [MatrixChainMultiply.jl package](https://github.com/AustinPrivett/MatrixChainMultiply.jl), but it needs updating (it dates back to Julia 0.5!).

(In practice there seems to rarely be a need for automatic matrix-chain algorithms — if you are multiplying a bunch of matrices of different sizes, in practice it’s usually only 3–4 matrices and the optimal ordering is easily identified statically by hand. If you are multiplying an arbitrary number of matrices at runtime I’m guessing their sizes are all the same?)

---

<div class="post-metadata">

**Author:** ![mohamed82008](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mohamed82008/32/18171_2.png) [@mohamed82008](https://discourse.julialang.org/u/mohamed82008)\
**Post date:** [May 13, 2023, 4:09pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/5 "2023-05-13T16:09:32Z")

</div>

Write it as a tensor reduction and look at [GitHub - mcabbott/Tullio.jl: ⅀](https://github.com/mcabbott/Tullio.jl).

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [May 13, 2023, 4:10pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/6 "2023-05-13T16:10:39Z")

</div>

> [@mohamed82008](#):
>
> Write it as a tensor reduction and look at [GitHub - mcabbott/Tullio.jl: ⅀](https://github.com/mcabbott/Tullio.jl).

The Tullio README [at one point said](https://github.com/mcabbott/Tullio.jl/blob/5d3105552e87efe6de981fb4d78e23e7df587fbb/README.md?plain=1#L164-L166):

> Chained multiplication is also very slow, because it doesn’t know there’s a better  
> algorithm.

Not sure if this is still true?

---

<div class="post-metadata">

**Author:** ![mohamed82008](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mohamed82008/32/18171_2.png) [@mohamed82008](https://discourse.julialang.org/u/mohamed82008)\
**Post date:** [May 13, 2023, 4:19pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/7 "2023-05-13T16:19:00Z")

</div>

Yes you are right.

---

<div class="post-metadata">

**Author:** ![Frisus95](https://avatars.discourse-cdn.com/v4/letter/f/77aa72/32.png) [@Frisus95](https://discourse.julialang.org/u/Frisus95)\
**Post date:** [May 13, 2023, 5:16pm UTC](https://discourse.julialang.org/t/matrix-chain-multiplication-with-arbitrary-number-of-matrices/98794/8 "2023-05-13T17:16:17Z")

</div>

That’s brilliant. Thank you so much.

P.S. Thanks also to all the others who replied