diff --git a/Project.toml b/Project.toml
index 6e4c40b0..157b96e6 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "SparseMatrixColorings"
 uuid = "0a514795-09f3-496d-8182-132a7b665d35"
 authors = ["Guillaume Dalle", "Alexis Montoison"]
-version = "0.4.14"
+version = "0.4.15"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
@@ -12,13 +12,16 @@ Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 
 [weakdeps]
+CliqueTrees = "60701a23-6482-424a-84db-faee86b9b1f8"
 Colors = "5ae59095-9a9b-59fe-a467-6f913c188581"
 
 [extensions]
+SparseMatrixColoringsCliqueTreesExt = "CliqueTrees"
 SparseMatrixColoringsColorsExt = "Colors"
 
 [compat]
 ADTypes = "1.2.1"
+CliqueTrees = "0.5.2"
 Colors = "0.12.11, 0.13"
 DataStructures = "0.18"
 DocStringExtensions = "0.8,0.9"
diff --git a/docs/Project.toml b/docs/Project.toml
index fc5d3c54..61f164de 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -6,3 +6,6 @@ DocumenterInterLinks = "d12716ef-a0f6-4df4-a9f1-a5a34e75c656"
 Images = "916415d5-f1e6-5110-898d-aaa5f9f070e0"
 SparseMatrixColorings = "0a514795-09f3-496d-8182-132a7b665d35"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
+
+[sources]
+SparseMatrixColorings = {path=".."}
diff --git a/docs/src/api.md b/docs/src/api.md
index d2ca870e..c121d1e3 100644
--- a/docs/src/api.md
+++ b/docs/src/api.md
@@ -53,4 +53,5 @@ SmallestLast
 IncidenceDegree
 DynamicLargestFirst
 DynamicDegreeBasedOrder
+PerfectEliminationOrder
 ```
diff --git a/ext/SparseMatrixColoringsCliqueTreesExt.jl b/ext/SparseMatrixColoringsCliqueTreesExt.jl
new file mode 100644
index 00000000..5dde067f
--- /dev/null
+++ b/ext/SparseMatrixColoringsCliqueTreesExt.jl
@@ -0,0 +1,23 @@
+module SparseMatrixColoringsCliqueTreesExt
+
+import CliqueTrees: permutation, EliminationAlgorithm, MCS
+import SparseArrays: SparseMatrixCSC, rowvals, nnz
+import SparseMatrixColorings:
+    AdjacencyGraph, BipartiteGraph, PerfectEliminationOrder, pattern, vertices
+
+PerfectEliminationOrder() = PerfectEliminationOrder(MCS())
+
+function vertices(g::AdjacencyGraph{T}, order::PerfectEliminationOrder) where {T}
+    S = pattern(g)
+
+    # construct matrix with sparsity pattern S
+    M = SparseMatrixCSC{Bool,T}(size(S)..., S.colptr, rowvals(S), ones(Bool, nnz(S)))
+
+    # construct a perfect elimination order
+    # self-loops are ignored
+    order, _ = permutation(M; alg=order.elimination_algorithm)
+
+    return reverse!(order)
+end
+
+end # module
diff --git a/src/SparseMatrixColorings.jl b/src/SparseMatrixColorings.jl
index 7d879c9a..35b9857c 100644
--- a/src/SparseMatrixColorings.jl
+++ b/src/SparseMatrixColorings.jl
@@ -57,6 +57,7 @@ include("show_colors.jl")
 
 export NaturalOrder, RandomOrder, LargestFirst
 export DynamicDegreeBasedOrder, SmallestLast, IncidenceDegree, DynamicLargestFirst
+export PerfectEliminationOrder
 export ColoringProblem, GreedyColoringAlgorithm, AbstractColoringResult
 export ConstantColoringAlgorithm
 export coloring, fast_coloring
diff --git a/src/order.jl b/src/order.jl
index 0635486b..319279fb 100644
--- a/src/order.jl
+++ b/src/order.jl
@@ -14,6 +14,7 @@ Depending on how the vertices are ordered, the number of colors necessary may va
 - [`IncidenceDegree`](@ref) (experimental)
 - [`SmallestLast`](@ref) (experimental)
 - [`DynamicLargestFirst`](@ref) (experimental)
+- [`PerfectEliminationOrder`](@ref)
 """
 abstract type AbstractOrder end
 
@@ -301,3 +302,24 @@ Instance of [`AbstractOrder`](@ref) which sorts vertices from lowest to highest
 - [`DynamicDegreeBasedOrder`](@ref)
 """
 const DynamicLargestFirst = DynamicDegreeBasedOrder{:forward,:low2high}
+
+"""
+    PerfectEliminationOrder(elimination_algorithm=CliqueTrees.MCS())
+
+Instance of [`AbstractOrder`](@ref) which computes a perfect elimination ordering when the underlying graph is [chordal](https://en.wikipedia.org/wiki/Chordal_graph). For non-chordal graphs, it computes a suboptimal ordering.
+
+The `elimination_algorithm` must be an instance of `CliqueTrees.EliminationAlgorithm`.
+
+!!! warning
+    This order can only be applied for symmetric or bidirectional coloring problems. Furthermore, its theoretical guarantees only hold for decompression by substitution.
+
+!!! danger
+    This order is implemented as a package extension and requires loading [CliqueTrees.jl](https://github.com/AlgebraicJulia/CliqueTrees.jl).
+
+# References
+
+- [Simple Linear-Time Algorithms to Test Chordality of Graphs, Test Acyclicity of Hypergraphs, and Selectively Reduce Acyclic Hypergraphs](https://epubs.siam.org/doi/10.1137/0213035), Tarjan and Yannakakis (1984)
+"""
+struct PerfectEliminationOrder{E} <: AbstractOrder
+    elimination_algorithm::E
+end
diff --git a/test/Project.toml b/test/Project.toml
index 05d0dedf..00ee7348 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -7,6 +7,7 @@ BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e"
 BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0"
 CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
 Chairmarks = "0ca39b1e-fe0b-4e98-acfc-b1656634c4de"
+CliqueTrees = "60701a23-6482-424a-84db-faee86b9b1f8"
 Colors = "5ae59095-9a9b-59fe-a467-6f913c188581"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
diff --git a/test/order.jl b/test/order.jl
index 8071603f..347e4b28 100644
--- a/test/order.jl
+++ b/test/order.jl
@@ -1,3 +1,5 @@
+using CliqueTrees: CliqueTrees
+using BandedMatrices
 using LinearAlgebra
 using SparseArrays
 using SparseMatrixColorings
@@ -7,6 +9,7 @@ using SparseMatrixColorings:
     LargestFirst,
     NaturalOrder,
     RandomOrder,
+    PerfectEliminationOrder,
     degree_dist2,
     nb_vertices,
     valid_dynamic_order,
@@ -110,3 +113,26 @@ end;
         end
     end
 end;
+
+@testset "PerfectEliminationOrder" begin
+    problem = ColoringProblem(; structure=:symmetric, partition=:column)
+    substitution_algo = GreedyColoringAlgorithm(
+        PerfectEliminationOrder(); decompression=:substitution
+    )
+
+    # band graphs
+    for (n, m) in ((800, 80), (400, 40), (200, 20), (100, 10))
+        perm = randperm(rng, n)
+        matrix = permute!(sparse(Symmetric(brand(n, n, m, 0), :L)), perm, perm)
+        π = vertices(AdjacencyGraph(matrix), PerfectEliminationOrder())
+        @test isperm(π)
+        @test ncolors(coloring(matrix, problem, substitution_algo)) == m + 1
+    end
+
+    # random graphs
+    for (n, p) in Iterators.product(20:20:100, 0.0:0.1:0.2)
+        matrix = sparse(Symmetric(sprand(rng, Bool, n, n, p)))
+        π = vertices(AdjacencyGraph(matrix), PerfectEliminationOrder())
+        @test isperm(π)
+    end
+end