From 2eb3b82c26ac88a20b38f992f2a0d609b344a538 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 07:42:54 -0500
Subject: [PATCH 01/13] Update permutations.jl

---
 src/permutations.jl | 117 ++++++++++++++++++++------------------------
 1 file changed, 52 insertions(+), 65 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index edab08a..a8847dc 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -15,43 +15,26 @@ struct Permutations{T}
     length::Int
 end
 
-# The following code basically implements `permutations` in terms of `multiset_permutations` as
-#
-#     permutations(a, t::Integer=length(a)) = Iterators.map(
-#         indices -> [a[i] for i in indices],
-#         multiset_permutations(eachindex(a), t))
-#
-# with the difference that we can also define `eltype(::Permutations)`, which is used in some tests.
-
-function Base.iterate(p::Permutations, state=nothing)
-    if state === nothing
-        mp = multiset_permutations(collect(eachindex(p.data)), p.length)
-        it = iterate(mp)
-        if it === nothing return nothing end
-    else
-        mp, mp_state = state
-        it = iterate(mp, mp_state)
-        if it === nothing return nothing end
-    end
-    indices, mp_state = it
-    return [p.data[i] for i in indices], (mp=mp, mp_state=mp_state)
+function Base.iterate(p::Permutations, state::Vector{Int} = collect(eachindex(p.data)))
+    (!isempty(state) && max(state[1], p.length) > length(p.data) || (isempty(state) && p.length > 0)) && return
+    nextpermutation!(p.data, p.length , state)
 end
 
-function Base.length(p::Permutations)
+function Base.length(p::Permutations)::Union{Int, BigInt}
     length(p.data) < p.length && return 0
-    return Int(prod(length(p.data) - p.length + 1:length(p.data)))
+    length(p.data) < 21       && return Int(prod(length(p.data) - p.length + 1:length(p.data)))
+    return prod(big.(length(p.data) - p.length + 1:length(p.data)))
 end
 
-Base.eltype(p::Permutations) = Vector{eltype(p.data)}
+Base.eltype(::Type{Permutations{T}}) where T = Vector{eltype(T)}
 
 Base.IteratorSize(p::Permutations) = Base.HasLength()
 
-
 """
     permutations(a)
 
-Generate all permutations of an indexable object `a` in lexicographic order. Because the number of permutations
-can be very large, this function returns an iterator object.
+Generate all permutations of an indexable object `a` in index-based lexicographic order. 
+Because the number of permutations can be very large, this function returns an iterator object.
 Use `collect(permutations(a))` to get an array of all permutations.
 Only works for `a` with defined length.
 
@@ -89,7 +72,7 @@ If `(t <= 0) || (t > length(a))`, then returns an empty vector of eltype of `a`
 julia> [ (len, permutations(1:3, len)) for len in -1:4 ]
 6-element Vector{Tuple{Int64, Any}}:
  (-1, Vector{Int64}[])
- (0, [Int64[]])
+ (0, [Int64[]])isconcretetype
  (1, [[1], [2], [3]])
  (2, Combinatorics.Permutations{UnitRange{Int64}}(1:3, 2))
  (3, Combinatorics.Permutations{UnitRange{Int64}}(1:3, 3))
@@ -105,18 +88,16 @@ julia> [ (len, collect(permutations(1:3, len))) for len in -1:4 ]
  (4, [])
 ```
 """
-function permutations(a, t::Integer)
-    if t == 0
-        # Correct behavior for a permutation of length 0 is a vector containing a single empty vector
-        return [Vector{eltype(a)}()]
-    elseif t == 1
-        # Easy case, just return each element in its own vector
-        return [[ai] for ai in a]
-    elseif t < 0 || t > length(a)
-        # Correct behavior for a permutation of these lengths is a an empty vector (of the correct type)
-        return Vector{Vector{eltype(a)}}()
+function permutations(a, t::Int)
+    if t < 0
+        t = length(a) + 1
     end
-    return Permutations(a, t)
+    data = eltype(a)[]
+    sizehint!(data, length(a))
+    for i in eachindex(a)
+        @inbounds push!(data, a[i])  # push! flattens `a` even with `CartesianIndex` 
+    end
+    return Permutations(data, t)
 end
 
 """
@@ -153,46 +134,46 @@ julia> derangements("julia") |> collect
  ['a', 'i', 'u', 'l', 'j']
 ```
 """
-derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(t -> t[1] != t[2], zip(a, d)))
-
+derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(a .!= d))
 
-function nextpermutation(m, t, state)
-    perm = [m[state[i]] for i in 1:t]
+function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
+    perm = m[@view state[1:t]]
     n = length(state)
-    if t <= 0
+    if t ≤ 0
         return (perm, [n + 1])
     end
-    s = copy(state)
     if t < n
         j = t + 1
-        while j <= n && s[t] >= s[j]
+        @inbounds while j ≤ n && state[t] ≥ state[j]
             j += 1
         end
     end
-    if t < n && j <= n
-        s[t], s[j] = s[j], s[t]
+    @inbounds if t < n && j <= n
+        state[t], state[j] = state[j], state[t]
     else
         if t < n
-            reverse!(s, t + 1)
+            reverse!(state, t + 1)
         end
         i = t - 1
-        while i >= 1 && s[i] >= s[i+1]
+        while i ≥ 1 && state[i] ≥ state[i+1]
             i -= 1
         end
         if i > 0
             j = n
-            while j > i && s[i] >= s[j]
+            while j > i && state[i] ≥ state[j]
                 j -= 1
             end
-            s[i], s[j] = s[j], s[i]
-            reverse!(s, i + 1)
+            state[i], state[j] = state[j], state[i]
+            reverse!(state, i + 1)
         else
-            s[1] = n + 1
+            state[1] = n + 1
         end
     end
-    return (perm, s)
+    return (perm, state)
 end
 
+nextpermutation(m::Vector, t::Int, state::Vector{Int}) = nextpermutation!(m, t, copy(state))
+
 struct MultiSetPermutations{T}
     m::T
     f::Vector{Int}
@@ -202,15 +183,15 @@ end
 
 Base.eltype(::Type{MultiSetPermutations{T}}) where {T} = Vector{eltype(T)}
 
-function Base.length(c::MultiSetPermutations)
+function Base.length(c::MultiSetPermutations)::Int
     t = c.t
     if t > length(c.ref)
         return 0
     end
     if t > 20
-        g = [factorial(big(i)) for i in 0:t]
+        g = factorial.(big.(0:t))
     else
-        g = [factorial(i) for i in 0:t]
+        g = factorial.(0:t)
     end
     p = [g[t+1]; zeros(Float64, t)]
     for i in 1:length(c.f)
@@ -256,7 +237,7 @@ julia> collect(permutations([1,1,1], 2))
  [1, 1]
  [1, 1]
 
-julia> collect(multiset_permutations([1,1,1], 2))
+julia> co1llect(multiset_permutations([1,1,1], 2))
 1-element Vector{Vector{Int64}}:
  [1, 1]
 
@@ -268,23 +249,29 @@ julia> collect(multiset_permutations([1,1,2], 3))
 ```
 """
 function multiset_permutations(a, t::Integer)
-    m = unique(a)
-    f = [sum(c == x for c in a)::Int for x in m]
+    counts = Dict{eltype(a), Int}()
+    m = eltype(a)[]
+    @inbounds for i in eachindex(a)
+        n = get(counts, a[i], 0) + 1
+        counts[a[i]] = n
+        isone(n) && push!(m, a[i])
+    end
+    f = [counts[key] for key in m]
     multiset_permutations(m, f, t)
 end
 
-function multiset_permutations(m, f::Vector{<:Integer}, t::Integer)
+function multiset_permutations(m::Vector, f::Vector{<:Integer}, t::Integer)
     length(m) == length(f) || error("Lengths of m and f are not the same.")
-    ref = length(f) > 0 ? vcat([[i for j in 1:f[i]] for i in 1:length(f)]...) : Int[]
+    ref = length(f) > 0 ? vcat(fill.(1:length(f), f)...) : Int[]
     if t < 0
         t = length(ref) + 1
     end
     MultiSetPermutations(m, f, t, ref)
 end
 
-function Base.iterate(p::MultiSetPermutations, s=p.ref)
-    (!isempty(s) && max(s[1], p.t) > length(p.ref) || (isempty(s) && p.t > 0)) && return
-    nextpermutation(p.m, p.t, s)
+function Base.iterate(p::MultiSetPermutations, state::Vector{Int} = copy(p.ref))
+    (!isempty(state) && max(state[1], p.t) > length(p.ref) || (isempty(state) && p.t > 0)) && return
+    nextpermutation!(p.m, p.t, state)
 end
 
 

From 909114441fd864bfba1ad3372d3d8d32f92fe895 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 09:05:41 -0500
Subject: [PATCH 02/13] Update permutations.jl

Update docstrings
---
 src/permutations.jl | 18 +++++++++---------
 1 file changed, 9 insertions(+), 9 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index a8847dc..34a13fd 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -41,7 +41,7 @@ Only works for `a` with defined length.
 # Examples
 ```jldoctest
 julia> permutations(1:2)
-Combinatorics.Permutations{UnitRange{Int64}}(1:2, 2)
+Combinatorics.Permutations{Vector{Int64}}([1, 2], 2)
 
 julia> collect(permutations(1:2))
 2-element Vector{Vector{Int64}}:
@@ -70,13 +70,13 @@ If `(t <= 0) || (t > length(a))`, then returns an empty vector of eltype of `a`
 # Examples
 ```jldoctest
 julia> [ (len, permutations(1:3, len)) for len in -1:4 ]
-6-element Vector{Tuple{Int64, Any}}:
- (-1, Vector{Int64}[])
- (0, [Int64[]])isconcretetype
- (1, [[1], [2], [3]])
- (2, Combinatorics.Permutations{UnitRange{Int64}}(1:3, 2))
- (3, Combinatorics.Permutations{UnitRange{Int64}}(1:3, 3))
- (4, Vector{Int64}[])
+6-element Vector{Tuple{Int64, Combinatorics.Permutations{Vector{Int64}}}}:
+ (-1, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 4))
+ (0, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 0))
+ (1, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 1))
+ (2, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 2))
+ (3, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 3))
+ (4, Combinatorics.Permutations{Vector{Int64}}([1, 2, 3], 4))
 
 julia> [ (len, collect(permutations(1:3, len))) for len in -1:4 ]
 6-element Vector{Tuple{Int64, Vector{Vector{Int64}}}}:
@@ -237,7 +237,7 @@ julia> collect(permutations([1,1,1], 2))
  [1, 1]
  [1, 1]
 
-julia> co1llect(multiset_permutations([1,1,1], 2))
+julia> collect(multiset_permutations([1,1,1], 2))
 1-element Vector{Vector{Int64}}:
  [1, 1]
 

From 2a634d713b938f507a75d3f421bfaa3127992357 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 09:17:54 -0500
Subject: [PATCH 03/13] Update permutations.jl

Update for string comparison in `derangements`
---
 src/permutations.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 34a13fd..308549f 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -134,7 +134,7 @@ julia> derangements("julia") |> collect
  ['a', 'i', 'u', 'l', 'j']
 ```
 """
-derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(a .!= d))
+derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(collect(a) .≠ collect(d)))
 
 function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
     perm = m[@view state[1:t]]

From 6139b8cfcf66385e30980a8e0410558784f7aa2f Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 09:32:24 -0500
Subject: [PATCH 04/13] Update permutations.jl

---
 src/permutations.jl | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 308549f..53df312 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -23,7 +23,7 @@ end
 function Base.length(p::Permutations)::Union{Int, BigInt}
     length(p.data) < p.length && return 0
     length(p.data) < 21       && return Int(prod(length(p.data) - p.length + 1:length(p.data)))
-    return prod(big.(length(p.data) - p.length + 1:length(p.data)))
+    return prod(BigInt(length(p.data)) - p.length + 1:length(p.data))
 end
 
 Base.eltype(::Type{Permutations{T}}) where T = Vector{eltype(T)}
@@ -134,7 +134,7 @@ julia> derangements("julia") |> collect
  ['a', 'i', 'u', 'l', 'j']
 ```
 """
-derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(collect(a) .≠ collect(d)))
+derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(collect(a) .≠ d))
 
 function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
     perm = m[@view state[1:t]]
@@ -172,7 +172,7 @@ function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
     return (perm, state)
 end
 
-nextpermutation(m::Vector, t::Int, state::Vector{Int}) = nextpermutation!(m, t, copy(state))
+nextpermutation(m::Vector, t::Int, state::Vector{Int}) = nextpermutation!(m, t, collect(state))
 
 struct MultiSetPermutations{T}
     m::T

From 819f2eb861d15ca049bb9859b378d2a4dce0318e Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 09:39:21 -0500
Subject: [PATCH 05/13] Update permutations.jl

Delete non-mutating `nextpermutation()`
---
 src/permutations.jl | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 53df312..ab8f422 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -23,7 +23,7 @@ end
 function Base.length(p::Permutations)::Union{Int, BigInt}
     length(p.data) < p.length && return 0
     length(p.data) < 21       && return Int(prod(length(p.data) - p.length + 1:length(p.data)))
-    return prod(BigInt(length(p.data)) - p.length + 1:length(p.data))
+    return prod(big(length(p.data)) - p.length + 1:big(length(p.data)))
 end
 
 Base.eltype(::Type{Permutations{T}}) where T = Vector{eltype(T)}
@@ -172,7 +172,6 @@ function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
     return (perm, state)
 end
 
-nextpermutation(m::Vector, t::Int, state::Vector{Int}) = nextpermutation!(m, t, collect(state))
 
 struct MultiSetPermutations{T}
     m::T

From ba1fdb141eb0497e5ae1d61e4b54fde034167566 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Thu, 18 Dec 2025 11:14:19 -0500
Subject: [PATCH 06/13] Update permutations.jl

---
 src/permutations.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index ab8f422..2c0c40c 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -23,7 +23,7 @@ end
 function Base.length(p::Permutations)::Union{Int, BigInt}
     length(p.data) < p.length && return 0
     length(p.data) < 21       && return Int(prod(length(p.data) - p.length + 1:length(p.data)))
-    return prod(big(length(p.data)) - p.length + 1:big(length(p.data)))
+    return prod(big(length(p.data) - p.length + 1):big(length(p.data)))
 end
 
 Base.eltype(::Type{Permutations{T}}) where T = Vector{eltype(T)}

From d4898d8e32d22dd47c3677093f7b555cdc4e6b0a Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Fri, 19 Dec 2025 11:09:46 -0500
Subject: [PATCH 07/13] Update permutations.jl

Provide a derangement-specific implementation. Improving performance and providing further functionality.
---
 src/permutations.jl | 133 ++++++++++++++++++++++++++++++++------------
 1 file changed, 98 insertions(+), 35 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 2c0c40c..a2a0f8c 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -100,41 +100,6 @@ function permutations(a, t::Int)
     return Permutations(data, t)
 end
 
-"""
-    derangements(a)
-
-Generate all derangements of an indexable object `a` in lexicographic order.
-Because the number of derangements can be very large, this function returns an iterator object.
-Use `collect(derangements(a))` to get an array of all derangements.
-Only works for `a` with defined length.
-
-# Examples
-```jldoctest
-julia> derangements("julia") |> collect
-44-element Vector{Vector{Char}}:
- ['u', 'j', 'i', 'a', 'l']
- ['u', 'j', 'a', 'l', 'i']
- ['u', 'l', 'j', 'a', 'i']
- ['u', 'l', 'i', 'a', 'j']
- ['u', 'l', 'a', 'j', 'i']
- ['u', 'i', 'j', 'a', 'l']
- ['u', 'i', 'a', 'j', 'l']
- ['u', 'i', 'a', 'l', 'j']
- ['u', 'a', 'j', 'l', 'i']
- ['u', 'a', 'i', 'j', 'l']
- ⋮
- ['a', 'j', 'i', 'l', 'u']
- ['a', 'l', 'j', 'u', 'i']
- ['a', 'l', 'u', 'j', 'i']
- ['a', 'l', 'i', 'j', 'u']
- ['a', 'l', 'i', 'u', 'j']
- ['a', 'i', 'j', 'u', 'l']
- ['a', 'i', 'j', 'l', 'u']
- ['a', 'i', 'u', 'j', 'l']
- ['a', 'i', 'u', 'l', 'j']
-```
-"""
-derangements(a) = (d for d in multiset_permutations(a, length(a)) if all(collect(a) .≠ d))
 
 function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
     perm = m[@view state[1:t]]
@@ -274,6 +239,104 @@ function Base.iterate(p::MultiSetPermutations, state::Vector{Int} = copy(p.ref))
 end
 
 
+#Derangements
+
+struct Derangements{T}
+    data::T
+    order::T
+    counts::Vector{Int}
+    t::Int
+end
+
+function derangements(a, t::Int=length(a))
+    data, order, counts = eltype(a)[], eltype(a)[], Dict{eltype(a), Int}()
+    sizehint!(data, length(a))
+    for i in eachindex(a)
+        n = get(counts, a[i], 0) + 1
+        counts[a[i]] = n
+        push!(data, a[i])
+        isone(n) && push!(order, a[i])
+    end
+    Derangements(data, order, [counts[key] for key in order], t)
+end
+
+Base.eltype(::Type{Derangements{T}}) where {T} = Vector{eltype(T)}
+
+Base.IteratorSize(::Derangements) = Base.SizeUnknown()
+
+function Base.iterate(d::Derangements)
+    2maximum(d.counts) > length(d.data) && return
+    nextderangement(d, ones(Int, length(d.data)), copy(d.counts), 1, ones(Int, length(d.data)))
+end
+
+function Base.iterate(d::Derangements, state)
+    derangement, state = nextderangement(d, state...)
+    all(isone, last(state)) ? nothing : (derangement, state)
+end
+
+"""
+    derangements(a)
+
+Generate all derangements of an indexable object `a` in index-based lexicographic order.
+Because the number of derangements can be very large, this function returns an iterator object.
+Use `collect(derangements(a))` to get an array of all derangements.
+Only works for `a` with defined length.
+
+# Examples
+```jldoctest
+julia> derangements("julia") |> collect
+44-element Vector{Vector{Char}}:
+ ['u', 'j', 'i', 'a', 'l']
+ ['u', 'j', 'a', 'l', 'i']
+ ['u', 'l', 'j', 'a', 'i']
+ ['u', 'l', 'i', 'a', 'j']
+ ['u', 'l', 'a', 'j', 'i']
+ ['u', 'i', 'j', 'a', 'l']
+ ['u', 'i', 'a', 'j', 'l']
+ ['u', 'i', 'a', 'l', 'j']
+ ['u', 'a', 'j', 'l', 'i']
+ ['u', 'a', 'i', 'j', 'l']
+ ⋮
+ ['a', 'j', 'i', 'l', 'u']
+ ['a', 'l', 'j', 'u', 'i']
+ ['a', 'l', 'u', 'j', 'i']
+ ['a', 'l', 'i', 'j', 'u']
+ ['a', 'l', 'i', 'u', 'j']
+ ['a', 'i', 'j', 'u', 'l']
+ ['a', 'i', 'j', 'l', 'u']
+ ['a', 'i', 'u', 'j', 'l']
+ ['a', 'i', 'u', 'l', 'j']
+```
+"""
+derangements(a) = derangements(a, length(a))
+
+function nextderangement(d::Derangements, perm::Vector{Int}, counts::Vector{Int}, idx::Int, iterstate::Vector{Int})
+    ordlen = length(d.order)
+    while true
+        depth = idx
+        @inbounds for i in iterstate[idx]:ordlen
+            iterstate[idx] = i + 1
+            if counts[i] ≥ 1 && d.order[i] ≠ d.data[idx]
+                perm[idx] = i
+                counts[i] -= 1
+                idx += 1
+                break
+            end
+        end
+        @inbounds if idx > d.t
+            idx -= 1
+            counts[perm[idx]] += 1
+            return d.order[@view perm[1:d.t]], (perm, counts, idx, iterstate)
+        elseif iterstate[idx] == ordlen + 1 && depth == idx
+            iterstate[idx] = 1
+            idx -= 1
+            iszero(idx) && return d.order[@view perm[1:d.t]], (perm, counts, idx, iterstate)
+            counts[perm[idx]] += 1
+        end
+    end
+end
+
+
 """
     nthperm!(a, k)
 

From c84b4cb2dd7cd2590dc95a0a4cc2a7af8a5baa0b Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Fri, 19 Dec 2025 12:03:19 -0500
Subject: [PATCH 08/13] Update derangements

---
 src/permutations.jl | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index a2a0f8c..59313d9 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -265,11 +265,13 @@ Base.eltype(::Type{Derangements{T}}) where {T} = Vector{eltype(T)}
 Base.IteratorSize(::Derangements) = Base.SizeUnknown()
 
 function Base.iterate(d::Derangements)
-    2maximum(d.counts) > length(d.data) && return
+    (isempty(d.data) || iszero(d.t)) && return eltype(d)[], nothing
+    (d.t > length(d.data) || d.t < 0 || 2maximum(d.counts) > length(d.data)) && return
     nextderangement(d, ones(Int, length(d.data)), copy(d.counts), 1, ones(Int, length(d.data)))
 end
 
 function Base.iterate(d::Derangements, state)
+    isnothing(state) && return nothing
     derangement, state = nextderangement(d, state...)
     all(isone, last(state)) ? nothing : (derangement, state)
 end
@@ -277,7 +279,7 @@ end
 """
     derangements(a)
 
-Generate all derangements of an indexable object `a` in index-based lexicographic order.
+Generate all derangements of an indexable object `a` in lexicographic order.
 Because the number of derangements can be very large, this function returns an iterator object.
 Use `collect(derangements(a))` to get an array of all derangements.
 Only works for `a` with defined length.

From c1d6350cd5ae0e5c8d0979a28d82be47f813aaaf Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Fri, 19 Dec 2025 12:36:38 -0500
Subject: [PATCH 09/13] Replace isnothing for Julia 1.0 comaptibility

---
 src/permutations.jl | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 59313d9..0a45c9f 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -113,7 +113,7 @@ function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
             j += 1
         end
     end
-    @inbounds if t < n && j <= n
+    @inbounds if t < n && j ≤ n
         state[t], state[j] = state[j], state[t]
     else
         if t < n
@@ -271,7 +271,7 @@ function Base.iterate(d::Derangements)
 end
 
 function Base.iterate(d::Derangements, state)
-    isnothing(state) && return nothing
+    state === nothing && return nothing
     derangement, state = nextderangement(d, state...)
     all(isone, last(state)) ? nothing : (derangement, state)
 end

From 287b265352bcb25045bc9964948b3a42d1e5bb14 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Sat, 20 Dec 2025 14:16:32 -0500
Subject: [PATCH 10/13] add DerangementsIterState type

---
 src/permutations.jl | 67 +++++++++++++++++++++++++++------------------
 1 file changed, 41 insertions(+), 26 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 0a45c9f..610211f 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -147,7 +147,7 @@ end
 
 Base.eltype(::Type{MultiSetPermutations{T}}) where {T} = Vector{eltype(T)}
 
-function Base.length(c::MultiSetPermutations)::Int
+function Base.length(c::MultiSetPermutations)
     t = c.t
     if t > length(c.ref)
         return 0
@@ -174,7 +174,7 @@ function Base.length(c::MultiSetPermutations)::Int
             end
         end
     end
-    return round(Int, p[t+1])
+    return round(p[t+1] > typemax(Int) ? BigInt : Int, p[t+1])
 end
 
 
@@ -248,7 +248,18 @@ struct Derangements{T}
     t::Int
 end
 
-function derangements(a, t::Int=length(a))
+mutable struct DerangementsIterState
+    idx::Int
+    iterstate::Vector{Int}
+    perm::Vector{Int}
+    counts::Vector{Int}
+end
+
+function DerangementsIterState(d::Derangements)
+    DerangementsIterState(1, ones(Int, length(d.data)), ones(Int, length(d.data)), copy(d.counts))
+end
+
+function derangements(a, t::Int)
     data, order, counts = eltype(a)[], eltype(a)[], Dict{eltype(a), Int}()
     sizehint!(data, length(a))
     for i in eachindex(a)
@@ -265,15 +276,19 @@ Base.eltype(::Type{Derangements{T}}) where {T} = Vector{eltype(T)}
 Base.IteratorSize(::Derangements) = Base.SizeUnknown()
 
 function Base.iterate(d::Derangements)
-    (isempty(d.data) || iszero(d.t)) && return eltype(d)[], nothing
+    state = DerangementsIterState(d)
+    if isempty(d.data) || iszero(d.t)
+        state.idx = 0
+        return eltype(d)[], state
+    end
     (d.t > length(d.data) || d.t < 0 || 2maximum(d.counts) > length(d.data)) && return
-    nextderangement(d, ones(Int, length(d.data)), copy(d.counts), 1, ones(Int, length(d.data)))
+    nextderangement(d, state)
 end
 
-function Base.iterate(d::Derangements, state)
-    state === nothing && return nothing
-    derangement, state = nextderangement(d, state...)
-    all(isone, last(state)) ? nothing : (derangement, state)
+function Base.iterate(d::Derangements, state::DerangementsIterState)
+    iszero(state.idx) && return nothing
+    derangement, state = nextderangement(d, state)
+    iszero(state.idx) ? nothing : (derangement, state)
 end
 
 """
@@ -312,28 +327,28 @@ julia> derangements("julia") |> collect
 """
 derangements(a) = derangements(a, length(a))
 
-function nextderangement(d::Derangements, perm::Vector{Int}, counts::Vector{Int}, idx::Int, iterstate::Vector{Int})
+function nextderangement(d::Derangements, state::DerangementsIterState)
     ordlen = length(d.order)
     while true
-        depth = idx
-        @inbounds for i in iterstate[idx]:ordlen
-            iterstate[idx] = i + 1
-            if counts[i] ≥ 1 && d.order[i] ≠ d.data[idx]
-                perm[idx] = i
-                counts[i] -= 1
-                idx += 1
+        depth = state.idx
+        @inbounds for i in state.iterstate[state.idx]:ordlen
+            state.iterstate[state.idx] = i + 1
+            if state.counts[i] ≥ 1 && d.order[i] ≠ d.data[state.idx]
+                state.perm[state.idx] = i
+                state.counts[i] -= 1
+                state.idx += 1
                 break
             end
         end
-        @inbounds if idx > d.t
-            idx -= 1
-            counts[perm[idx]] += 1
-            return d.order[@view perm[1:d.t]], (perm, counts, idx, iterstate)
-        elseif iterstate[idx] == ordlen + 1 && depth == idx
-            iterstate[idx] = 1
-            idx -= 1
-            iszero(idx) && return d.order[@view perm[1:d.t]], (perm, counts, idx, iterstate)
-            counts[perm[idx]] += 1
+        @inbounds if state.idx > d.t
+            state.idx -= 1
+            state.counts[state.perm[state.idx]] += 1
+            return d.order[@view state.perm[1:d.t]], state
+        elseif state.iterstate[state.idx] == ordlen + 1 && depth == state.idx
+            state.iterstate[state.idx] = 1
+            state.idx -= 1
+            iszero(state.idx) && return d.order[@view state.perm[1:d.t]], state
+            state.counts[state.perm[state.idx]] += 1
         end
     end
 end

From 51426e89c50f431bd6cfd09403fd8625544a99f9 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Sun, 21 Dec 2025 12:10:33 -0500
Subject: [PATCH 11/13] Further optimization of derangements

Iteration has been streamlined a bit more and some comments have been added to help with future maintenance.
---
 src/permutations.jl | 101 +++++++++++++++++++++++++++++++-------------
 1 file changed, 72 insertions(+), 29 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 610211f..bb3b1e3 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -95,7 +95,7 @@ function permutations(a, t::Int)
     data = eltype(a)[]
     sizehint!(data, length(a))
     for i in eachindex(a)
-        @inbounds push!(data, a[i])  # push! flattens `a` even with `CartesianIndex` 
+        @inbounds push!(data, a[i])  # `push!` to a `Vector` flattens `a`
     end
     return Permutations(data, t)
 end
@@ -221,10 +221,10 @@ function multiset_permutations(a, t::Integer)
         isone(n) && push!(m, a[i])
     end
     f = [counts[key] for key in m]
-    multiset_permutations(m, f, t)
+    MultiSetPermutations(m, f, t)
 end
 
-function multiset_permutations(m::Vector, f::Vector{<:Integer}, t::Integer)
+function MultiSetPermutations(m::Vector, f::Vector{<:Integer}, t::Integer)
     length(m) == length(f) || error("Lengths of m and f are not the same.")
     ref = length(f) > 0 ? vcat(fill.(1:length(f), f)...) : Int[]
     if t < 0
@@ -259,42 +259,74 @@ function DerangementsIterState(d::Derangements)
     DerangementsIterState(1, ones(Int, length(d.data)), ones(Int, length(d.data)), copy(d.counts))
 end
 
-function derangements(a, t::Int)
-    data, order, counts = eltype(a)[], eltype(a)[], Dict{eltype(a), Int}()
-    sizehint!(data, length(a))
-    for i in eachindex(a)
-        n = get(counts, a[i], 0) + 1
-        counts[a[i]] = n
-        push!(data, a[i])
-        isone(n) && push!(order, a[i])
-    end
-    Derangements(data, order, [counts[key] for key in order], t)
-end
-
 Base.eltype(::Type{Derangements{T}}) where {T} = Vector{eltype(T)}
 
 Base.IteratorSize(::Derangements) = Base.SizeUnknown()
 
 function Base.iterate(d::Derangements)
-    state = DerangementsIterState(d)
-    if isempty(d.data) || iszero(d.t)
-        state.idx = 0
-        return eltype(d)[], state
-    end
-    (d.t > length(d.data) || d.t < 0 || 2maximum(d.counts) > length(d.data)) && return
+    state = (isempty(d.data) || iszero(d.t)) ? DerangementsIterState(0, Int[], Int[], Int[]) : DerangementsIterState(d)
+    (d.t > length(d.data) || d.t < 0 || 2maximum([d.counts; 0]) > length(d.data)) && return
     nextderangement(d, state)
 end
 
 function Base.iterate(d::Derangements, state::DerangementsIterState)
-    iszero(state.idx) && return nothing
     derangement, state = nextderangement(d, state)
     iszero(state.idx) ? nothing : (derangement, state)
 end
 
+"""
+    derangements(a, t)
+
+Generate all derangements of an indexable object `a` of length `t` in index-based lexicographic order.
+Because the number of derangements can be very large, this function returns an iterator object.
+Use `collect(derangements(a))` to get an array of all derangements.
+Only works for `a` with defined length.
+
+# Examples
+```jldoctest
+julia> derangements(1:4, 4) |> collect
+9-element Vector{Vector{Int64}}:
+ [2, 1, 4, 3]
+ [2, 3, 4, 1]
+ [2, 4, 1, 3]
+ [3, 1, 4, 2]
+ [3, 4, 1, 2]
+ [3, 4, 2, 1]
+ [4, 1, 2, 3]
+ [4, 3, 1, 2]
+ [4, 3, 2, 1]
+
+julia> derangements(1:4, 3) |> collect
+11-element Vector{Vector{Int64}}:
+ [2, 1, 4]
+ [2, 3, 1]
+ [2, 3, 4]
+ [2, 4, 1]
+ [3, 1, 2]
+ [3, 1, 4]
+ [3, 4, 1]
+ [3, 4, 2]
+ [4, 1, 2]
+ [4, 3, 1]
+ [4, 3, 2]
+```
+"""
+function derangements(a, t::Int)
+    data, order, counts = eltype(a)[], eltype(a)[], Dict{eltype(a), Int}()
+    sizehint!(data, length(a))
+    for i in eachindex(a)
+        n = get(counts, a[i], 0) + 1
+        counts[a[i]] = n
+        push!(data, a[i])
+        isone(n) && push!(order, a[i])
+    end
+    Derangements(data, order, [counts[key] for key in order], t)
+end
+
 """
     derangements(a)
 
-Generate all derangements of an indexable object `a` in lexicographic order.
+Generate all derangements of an indexable object `a` in index-based lexicographic order.
 Because the number of derangements can be very large, this function returns an iterator object.
 Use `collect(derangements(a))` to get an array of all derangements.
 Only works for `a` with defined length.
@@ -327,33 +359,44 @@ julia> derangements("julia") |> collect
 """
 derangements(a) = derangements(a, length(a))
 
+############################################################################################
+# `nextderangement` is a iterative translation of a pruned-dfs with backtracking algoritm. #
+# The iteration state is kept in `state` where: `idx` is the depth, `iterstate` is the     #
+# position of the for loop at each depth, `perm` is the sort permutation used to slice     #
+# `d.order`, and `counts[i]` is the remaining number of elements from `d.order[i]` that    #
+# still need to be accounted for in the derangment.                                        #
+# The source of the original algorithm is at https://arxiv.org/pdf/1009.4214               #
+############################################################################################
+
 function nextderangement(d::Derangements, state::DerangementsIterState)
     ordlen = length(d.order)
-    while true
+    while 0 < state.idx
         depth = state.idx
         @inbounds for i in state.iterstate[state.idx]:ordlen
             state.iterstate[state.idx] = i + 1
-            if state.counts[i] ≥ 1 && d.order[i] ≠ d.data[state.idx]
+            if state.counts[i] ≥ 1 && d.order[i] ≠ d.data[state.idx] # If true, candidate element for position idx has been found
                 state.perm[state.idx] = i
                 state.counts[i] -= 1
                 state.idx += 1
                 break
             end
         end
-        @inbounds if state.idx > d.t
+        @inbounds if state.idx > d.t # If true, a derangement of length `t` has been found
             state.idx -= 1
             state.counts[state.perm[state.idx]] += 1
             return d.order[@view state.perm[1:d.t]], state
-        elseif state.iterstate[state.idx] == ordlen + 1 && depth == state.idx
+        elseif state.iterstate[state.idx] == ordlen + 1 && depth == state.idx # If true, the for loop at this depth has terminated
             state.iterstate[state.idx] = 1
             state.idx -= 1
-            iszero(state.idx) && return d.order[@view state.perm[1:d.t]], state
-            state.counts[state.perm[state.idx]] += 1
+            !iszero(state.idx) && (state.counts[state.perm[state.idx]] += 1) 
         end
     end
+    eltype(d.data)[], state # Termination of algorithm
 end
 
 
+# nthperm
+
 """
     nthperm!(a, k)
 

From 153f827add4df70d8fd59e75316fab62abcf50de Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Tue, 23 Dec 2025 10:43:28 -0500
Subject: [PATCH 12/13] Update permutations.jl

---
 src/permutations.jl | 33 ++++++++++++++++++---------------
 1 file changed, 18 insertions(+), 15 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index bb3b1e3..8fe8642 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -250,7 +250,7 @@ end
 
 mutable struct DerangementsIterState
     idx::Int
-    iterstate::Vector{Int}
+    inner::Vector{Int}
     perm::Vector{Int}
     counts::Vector{Int}
 end
@@ -361,36 +361,39 @@ derangements(a) = derangements(a, length(a))
 
 ############################################################################################
 # `nextderangement` is a iterative translation of a pruned-dfs with backtracking algoritm. #
-# The iteration state is kept in `state` where: `idx` is the depth, `iterstate` is the     #
-# position of the for loop at each depth, `perm` is the sort permutation used to slice     #
-# `d.order`, and `counts[i]` is the remaining number of elements from `d.order[i]` that    #
-# still need to be accounted for in the derangment.                                        #
-# The source of the original algorithm is at https://arxiv.org/pdf/1009.4214               #
+# The iteration state is kept in `state` where: `idx` is the depth, `inner` is the         #
+# position of the inner loop at each depth, `perm` has the derangement indices used to     #
+# slice `d.order`, and `counts[i]` is the remaining number of elements from `d.order[i]`   #
+# that still need to be included in the derangment.                                        #
+# The source of the original translated algorithm is at https://arxiv.org/pdf/1009.4214    #
 ############################################################################################
 
 function nextderangement(d::Derangements, state::DerangementsIterState)
     ordlen = length(d.order)
+    
     while 0 < state.idx
-        depth = state.idx
-        @inbounds for i in state.iterstate[state.idx]:ordlen
-            state.iterstate[state.idx] = i + 1
-            if state.counts[i] ≥ 1 && d.order[i] ≠ d.data[state.idx] # If true, candidate element for position idx has been found
-                state.perm[state.idx] = i
-                state.counts[i] -= 1
+        i = state.idx
+        @inbounds while i == state.idx && state.inner[i] ≤ ordlen
+            if state.counts[state.inner[i]] ≥ 1 && d.order[state.inner[i]] ≠ d.data[i] # If true, candidate element for index idx has been found
+                state.perm[i] = state.inner[i]
+                state.counts[state.inner[i]] -= 1
                 state.idx += 1
-                break
             end
+            
+            state.inner[i] += 1
         end
+        
         @inbounds if state.idx > d.t # If true, a derangement of length `t` has been found
             state.idx -= 1
             state.counts[state.perm[state.idx]] += 1
             return d.order[@view state.perm[1:d.t]], state
-        elseif state.iterstate[state.idx] == ordlen + 1 && depth == state.idx # If true, the for loop at this depth has terminated
-            state.iterstate[state.idx] = 1
+        elseif state.inner[state.idx] == ordlen + 1 && i == state.idx # If true, the inner loop at this depth has terminated
+            state.inner[state.idx] = 1
             state.idx -= 1
             !iszero(state.idx) && (state.counts[state.perm[state.idx]] += 1) 
         end
     end
+    
     eltype(d.data)[], state # Termination of algorithm
 end
 

From b9479c40d7b300cdc997c680a136d6fac71b4462 Mon Sep 17 00:00:00 2001
From: depial <91621102+depial@users.noreply.github.com>
Date: Tue, 27 Jan 2026 16:35:30 -0500
Subject: [PATCH 13/13] Update api for Derangements

Limits constructors for the types `Derangements` and `DerangementsIterState` to one inner constructor a piece.
---
 src/permutations.jl | 48 +++++++++++++++++++++++++--------------------
 1 file changed, 27 insertions(+), 21 deletions(-)

diff --git a/src/permutations.jl b/src/permutations.jl
index 8fe8642..3bb64bf 100644
--- a/src/permutations.jl
+++ b/src/permutations.jl
@@ -88,7 +88,7 @@ julia> [ (len, collect(permutations(1:3, len))) for len in -1:4 ]
  (4, [])
 ```
 """
-function permutations(a, t::Int)
+function permutations(a, t::Integer)
     if t < 0
         t = length(a) + 1
     end
@@ -138,6 +138,8 @@ function nextpermutation!(m::Vector, t::Int, state::Vector{Int})
 end
 
 
+# Multiset Permutations
+
 struct MultiSetPermutations{T}
     m::T
     f::Vector{Int}
@@ -239,13 +241,26 @@ function Base.iterate(p::MultiSetPermutations, state::Vector{Int} = copy(p.ref))
 end
 
 
-#Derangements
+# Derangements
 
 struct Derangements{T}
     data::T
     order::T
     counts::Vector{Int}
     t::Int
+    function Derangements{T}(a, t::Integer) where T<:Vector
+        order, counts = eltype(T)[], Dict{eltype(T), Int}()
+        data = sizehint!(eltype(T)[], length(a))
+        
+        for i in eachindex(a)
+            n = get(counts, a[i], 0) + 1
+            counts[a[i]] = n
+            push!(data, a[i])
+            isone(n) && push!(order, a[i])
+        end
+        
+        new{T}(data, order, [counts[key] for key in order], t)
+    end
 end
 
 mutable struct DerangementsIterState
@@ -253,20 +268,21 @@ mutable struct DerangementsIterState
     inner::Vector{Int}
     perm::Vector{Int}
     counts::Vector{Int}
+    function DerangementsIterState(d::Derangements)
+        if isempty(d.data) || iszero(d.t)
+            return new(0, Int[], Int[], Int[])
+        end
+        new(1, ones(Int, length(d.data)), ones(Int, length(d.data)), copy(d.counts))
+    end
 end
 
-function DerangementsIterState(d::Derangements)
-    DerangementsIterState(1, ones(Int, length(d.data)), ones(Int, length(d.data)), copy(d.counts))
-end
-
-Base.eltype(::Type{Derangements{T}}) where {T} = Vector{eltype(T)}
+Base.eltype(::Type{Derangements{T}}) where T<:Vector = T
 
-Base.IteratorSize(::Derangements) = Base.SizeUnknown()
+Base.IteratorSize(::Type{Derangements{T}}) where T<:Vector = Base.SizeUnknown()
 
 function Base.iterate(d::Derangements)
-    state = (isempty(d.data) || iszero(d.t)) ? DerangementsIterState(0, Int[], Int[], Int[]) : DerangementsIterState(d)
     (d.t > length(d.data) || d.t < 0 || 2maximum([d.counts; 0]) > length(d.data)) && return
-    nextderangement(d, state)
+    nextderangement(d, DerangementsIterState(d))
 end
 
 function Base.iterate(d::Derangements, state::DerangementsIterState)
@@ -311,17 +327,7 @@ julia> derangements(1:4, 3) |> collect
  [4, 3, 2]
 ```
 """
-function derangements(a, t::Int)
-    data, order, counts = eltype(a)[], eltype(a)[], Dict{eltype(a), Int}()
-    sizehint!(data, length(a))
-    for i in eachindex(a)
-        n = get(counts, a[i], 0) + 1
-        counts[a[i]] = n
-        push!(data, a[i])
-        isone(n) && push!(order, a[i])
-    end
-    Derangements(data, order, [counts[key] for key in order], t)
-end
+derangements(a, t::Integer) = Derangements{Vector{eltype(a)}}(a, t)
 
 """
     derangements(a)