# Constructing a permutation matrix

**URL:** <https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234>\
**Category:** General Usage\
**Tags:** linearalgebra\
**Created:** [May 21, 2025, 11:04pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234 "2025-05-21T23:04:23Z")\
**Posts on this page:** 16\
**Page:** 1

<div class="post-metadata">

**Author:** ![cnaak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cnaak/32/10287_2.png) [@cnaak](https://discourse.julialang.org/u/cnaak)\
**Post date:** [May 21, 2025, 11:04pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/1 "2025-05-21T23:04:23Z")

</div>

Approaching the problem from a slightly different angle, one can first transform a permutation vector `p` into a permutation matrix, for instance:

```julia
julia> function MP(p)
           ℓ = length(p)
           vcat([begin t=reshape(repeat([0], ℓ), 1, ℓ); t[i]=1; t; end for i in p ]...)
       end
MP (generic function with 1 method)

```

I’m pretty sure the implementation can be improved / this is what first worked for me.

`MP` assembles the corresponding (left multiplicative) permutation matrix. Therefore:

```julia
julia> MP([2, 1, 4, 3])
4×4 Matrix{Int64}:
 0 1 0 0
 1 0 0 0
 0 0 0 1
 0 0 1 0

julia> A = [1 1 1 1; 2 2 2 2; 3 3 3 3; 4 4 4 4]
4×4 Matrix{Int64}:
 1 1 1 1
 2 2 2 2
 3 3 3 3
 4 4 4 4

julia> MP([2, 1, 4, 3]) * A
4×4 Matrix{Int64}:
 2 2 2 2
 1 1 1 1
 4 4 4 4
 3 3 3 3

```

---

<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 21, 2025, 11:19pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/2 "2025-05-21T23:19:03Z")

</div>

> [@cnaak](#):
>
> Approaching the problem from a slightly different angle, one can first transform a permutation vector `p` into a permutation matrix, for instance:

Using an explicit multiplication PA of dense matrices to permute the rows of A transforms the cost from O(n^2) to O(n^3) for an n \times n matrix.

But if you do want a permutation matrix, it’s much easier (and faster) to simply do

```julia
using LinearAlgebra
MP(p) = I(length(p))[p,:]

```

---

<div class="post-metadata">

**Author:** ![cnaak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cnaak/32/10287_2.png) [@cnaak](https://discourse.julialang.org/u/cnaak)\
**Post date:** [May 22, 2025, 12:41am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/3 "2025-05-22T00:41:24Z")

</div>

You’re spot-on on the algorithm complexity increase in general (for matrices without assumptions). I wonder whether there is or could there be a sort of “one-hot-sparse” matrix type (i.e., line/column permutations of identity matrices) with specialized operations that would bring down the complexity of matrix multiplications, since there would be logical ones and zeroes and trivial operations…

---

<div class="post-metadata">

**Author:** ![nsajko](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nsajko/32/221187_2.png) [@nsajko](https://discourse.julialang.org/u/nsajko)\
**Post date:** [May 22, 2025, 6:53am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/4 "2025-05-22T06:53:32Z")

</div>

> [@cnaak](#):
>
> I wonder whether there is or could there be a sort of “one-hot-sparse” matrix type

Sure. StaticPermutations.jl is relevant, however it’s in the type domain (compile-time computation), which might or might not be what you want. In `Base` there is already `PermutedDimsArray`, which I guess can be thought of as a lazy multiplication with a permutation matrix. There’s also Permutations.jl.

---

<div class="post-metadata">

**Author:** ![abraemer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/abraemer/32/51403_2.png) [@abraemer](https://discourse.julialang.org/u/abraemer)\
**Post date:** [May 22, 2025, 8:51am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/5 "2025-05-22T08:51:02Z")

</div>

In a sense the “sparse” representation of your permutation matrix is just the permutation vector. So you could just wrap that into a type that overloads multiplication to instead perform the permutation of rows/columns.

---

<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 22, 2025, 11:58am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/6 "2025-05-22T11:58:46Z")

</div>

> [@cnaak](#):
>
> I wonder whether there is or could there be a sort of “one-hot-sparse” matrix type

You could, of course, use SparseArrays.jl for this.

But representing the linear operator implicitly by the permutation vector is even more efficient. You could wrap it in a matrix-like type with LinearMaps.jl or with a custom `AbstractArray` type if you wanted to, though you’d need some compelling motivation to go to the trouble.

---

<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 22, 2025, 12:03pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/7 "2025-05-22T12:03:52Z")

</div>

The trouble is relatively minimal, at least to get the basics working:

```julia
struct PermutationMatrix <: AbstractMatrix{Bool}
    p::Vector{Int}
end

Base.size(P::PermutationMatrix) = (length(P.p), length(P.p))
Base.getindex(P::PermutationMatrix, i::Integer, j::Integer) = i == P.p[j]

Base.:*(P::PermutationMatrix, A::AbstractMatrix) = A[P.p, :]
Base.:*(A::AbstractMatrix, P::PermutationMatrix) = A[:, P.p]

```

```julia-repl
julia> P = PermutationMatrix([2, 1, 4, 3])
4×4 PermutationMatrix:
 0 1 0 0
 1 0 0 0
 0 0 0 1
 0 0 1 0

julia> A = Matrix(reshape(1:16, 4, 4))
4×4 Matrix{Int64}:
 1 5 9 13
 2 6 10 14
 3 7 11 15
 4 8 12 16

julia> P * A
4×4 Matrix{Int64}:
 2 6 10 14
 1 5 9 13
 4 8 12 16
 3 7 11 15

julia> A * P
4×4 Matrix{Int64}:
 5 1 13 9
 6 2 14 10
 7 3 15 11

```

---

<div class="post-metadata">

**Author:** ![mikmoore](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mikmoore/32/31109_2.png) [@mikmoore](https://discourse.julialang.org/u/mikmoore)\
**Post date:** [May 22, 2025, 2:11pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/8 "2025-05-22T14:11:49Z")

</div>

> [@gdalle](#):
>
> ```julia
> Base.:*(P::PermutationMatrix, A::AbstractMatrix) = A[P.p, :]
> Base.:*(A::AbstractMatrix, P::PermutationMatrix) = A[:, P.p]
> 
> ```

Although the permutation should be inverse-applied (i.e., with its transpose) on one side or the other if it’s actually to be treated as a permutation matrix, which would make this a bit more complicated.

---

<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 22, 2025, 6:37pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/9 "2025-05-22T18:37:16Z")

</div>

> [@mikmoore](#):
>
> Although the permutation should be inverse-applied (i.e., with its transpose) on one side or the other if it’s actually to be treated as a permutation matrix, which would make this a bit more complicated.

Yes, the `getindex` should be

```julia
Base.getindex(P::PermutationMatrix, i::Integer, j::Integer) = P.p[i] == j

```

and the right-multiplication should be:

```julia
function Base.:*(A::AbstractMatrix, P::PermutationMatrix)
    AP = similar(A)
    AP[:, P.p] = A # computes AP = (PᵀAᵀ)ᵀ = (P⁻¹Aᵀ)ᵀ
    return AP
end

```

(The test permutation `P = PermutationMatrix([2, 1, 4, 3])` happens to be symmetric, which hides these bugs. Try e.g. `PermutationMatrix([2, 4, 3, 1])` instead.)

If you wanted to be fancy, you could define all sorts of specialized routines, most obviously for `inv(P)` and `P \ x` (`ldiv!`), but also for things like eigenvectors and eigenvalues (which can be computed from the cycles of the permutation, and moreover you can then multiply by the eigenvectors using FFTs, giving you fast algorithms for things like matrix exponentials and logarithms of permutations etcetera), and more.

It’s a fair amount of work to get something truly full-featured, but it might be fun for someone. Still I’m not sure of the practical motivation.

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [May 22, 2025, 6:43pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/10 "2025-05-22T18:43:24Z")

</div>

This package:

> **[GitHub - jlapeyre/PermutationsA.jl: Permutation types and methods](https://github.com/jlapeyre/PermutationsA.jl)**
>
> Permutation types and methods

is neglected, but has the same ideas (matrix form, cycle counts, sparse form). Perhaps it can be alived.

---

<div class="post-metadata">

**Author:** ![cnaak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cnaak/32/10287_2.png) [@cnaak](https://discourse.julialang.org/u/cnaak)\
**Post date:** [May 22, 2025, 10:47pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/11 "2025-05-22T22:47:52Z")

</div>

From [this article](https://en.wikipedia.org/wiki/Permutation_matrix), it seems to me that the proposed implementation of @gdalle is already correct, isn’t it?

---

<div class="post-metadata">

**Author:** ![cnaak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cnaak/32/10287_2.png) [@cnaak](https://discourse.julialang.org/u/cnaak)\
**Post date:** [May 22, 2025, 10:51pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/12 "2025-05-22T22:51:53Z")

</div>

Great implementation! Concise, sweet, and readable!!!

---

<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 22, 2025, 11:21pm UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/13 "2025-05-22T23:21:16Z")

</div>

> [@cnaak](#):
>
> From [this article](https://en.wikipedia.org/wiki/Permutation_matrix), it seems to me that the proposed implementation of @gdalle is already correct, isn’t it?

Nope. It implements the transposed `getindex` operation, so it is equivalent to P^T. Even if you fix that so that it is the right matrix, it implements `A * P` as `A[:, P.p]`, which is then equivalent to AP^T. (Even if you reinterpret what you mean by the “permutation vector” `p` to make his `P` correct, his `A * P` and `P * A` were inconsistent.) You can’t tell for your example P, since that matrix happened to be symmetric.

That’s why I gave the corrected implementation [in my post](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/9).

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [May 23, 2025, 12:11am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/14 "2025-05-23T00:11:24Z")

</div>

Also, this seems like a good problem to be lazy about, by defining a `PermutedMatrix` in addition to `PermutationMatrix`.  
Suppose,

```julia
P = PermutationMatrix([1,3,4,2])
A = rand(4,4)
B = rand(4,4)
M = P*A // would be a PermutedMatrix
R = B*M // == B*(P*A) == (B*P)*A
        // and performing permutation on columns
        // is easier and faster.
        // `getindex` on PermutedMatrix would go
        // through the permutation

```

---

<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 23, 2025, 12:13am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/15 "2025-05-23T00:13:04Z")

</div>

> [@Dan](#):
>
> Also, this seems like a good problem to be lazy about,

Hard to say without more context. Storing the permuted matrix lazily will make many operations on it much slower.

---

<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 23, 2025, 6:47am UTC](https://discourse.julialang.org/t/constructing-a-permutation-matrix/129234/16 "2025-05-23T06:47:39Z")

</div>

My bad, I wrote it in a hurry based on C\_\pi from the Wikipedia article but I did screw up one of the multiplications, thanks for fixing the issues.  
Ultimately it also depends how you define the permutation matrix (row or column wise).