Replace Set with OrderedSets in grid collections (#834)

This patch changes the use of `Set` to `OrderedSet` in the grid. This
means that e.g. loops over these sets follow the original specified
order. This is important for performance since it can reduce
cache-misses, for example. Fixes #631. Closes #654.

Co-authored-by: Fredrik Ekre <ekrefredrik@gmail.com>
Co-authored-by: Dennis Ogiermann <termi-official@users.noreply.github.com>

Project.toml

@@ -6,6 +6,7 @@ version = "0.3.14"
 EnumX = "4e289a0a-7415-4d19-859d-a7e5c4648b56"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
+OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
 Preferences = "21216c6a-2e73-6563-6e65-726566657250"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
@@ -26,6 +27,7 @@ BlockArrays = "0.16"
 EnumX = "1"
 Metis = "1.3"
 NearestNeighbors = "0.4"
+OrderedCollections = "1"
 Preferences = "1"
 Reexport = "1"
 Tensors = "1.14"

docs/Manifest.toml

@@ -429,7 +429,7 @@ uuid = "29a986be-02c6-4525-aec4-84b980013641"
 version = "2.0.0"
 
 [[deps.Ferrite]]
-deps = ["EnumX", "LinearAlgebra", "NearestNeighbors", "Preferences", "Reexport", "SparseArrays", "StaticArrays", "Tensors", "WriteVTK"]
+deps = ["EnumX", "LinearAlgebra", "NearestNeighbors", "OrderedCollections", "Preferences", "Reexport", "SparseArrays", "StaticArrays", "Tensors", "WriteVTK"]
 path = ".."
 uuid = "c061ca5d-56c9-439f-9c0e-210fe06d3992"
 version = "0.3.14"
@@ -450,9 +450,9 @@ version = "1.0.1"
 
 [[deps.FerriteMeshParser]]
 deps = ["Ferrite"]
-git-tree-sha1 = "8b948577bc4066e9c8693438fd511309c7383761"
+git-tree-sha1 = "47a1215bc76454082e265cb5193cffc9700bccb7"
 uuid = "0f8c756f-80dd-4a75-85c6-b0a5ab9d4620"
-version = "0.1.7"
+version = "0.1.8"
 
 [[deps.FileIO]]
 deps = ["Pkg", "Requires", "UUIDs"]

docs/src/literate-tutorials/porous_media.jl

@@ -99,7 +99,7 @@
 using Ferrite, FerriteMeshParser, Tensors
 
 # Temporary overload
-function FerriteMeshParser.create_faceset(grid::Grid, nodeset::Set{Int}, cellset=1:getncells(grid))
+function FerriteMeshParser.create_faceset(grid::Grid, nodeset::Ferrite.OrderedSet{Int}, cellset=1:getncells(grid))
     faceset = sizehint!(Set{FaceIndex}(), length(nodeset))
     for (cellid, cell) in enumerate(getcells(grid))
         cellid ∈ cellset || continue

src/Dofs/ConstraintHandler.jl

@@ -33,14 +33,14 @@ which applies the condition via [`apply!`](@ref) and/or [`apply_zero!`](@ref).
 """
 struct Dirichlet # <: Constraint
     f::Function # f(x) or f(x,t) -> value(s)
-    facets::Union{Set{Int},Set{FacetIndex},Set{FaceIndex},Set{EdgeIndex},Set{VertexIndex}}
+    facets::OrderedSet{T} where T <: Union{Int, FacetIndex, FaceIndex, EdgeIndex, VertexIndex}
     field_name::Symbol
     components::Vector{Int} # components of the field
     local_facet_dofs::Vector{Int}
     local_facet_dofs_offset::Vector{Int}
 end
-function Dirichlet(field_name::Symbol, facets::Set, f::Function, components=nothing)
-    return Dirichlet(f, facets, field_name, __to_components(components), Int[], Int[])
+function Dirichlet(field_name::Symbol, facets::AbstractVecOrSet, f::Function, components=nothing)
+    return Dirichlet(f, convert_to_orderedset(facets), field_name, __to_components(components), Int[], Int[])
 end
 
 # components=nothing is default and means that all components should be constrained
@@ -293,7 +293,7 @@ function add_prescribed_dof!(ch::ConstraintHandler, constrained_dof::Int, inhomo
 end
 
 # Dirichlet on (facet|face|edge|vertex)set
-function _add!(ch::ConstraintHandler, dbc::Dirichlet, bcfacets::Set{Index}, interpolation::Interpolation, field_dim::Int, offset::Int, bcvalue::BCValues, _) where {Index<:BoundaryIndex}
+function _add!(ch::ConstraintHandler, dbc::Dirichlet, bcfacets::AbstractVecOrSet{Index}, interpolation::Interpolation, field_dim::Int, offset::Int, bcvalue::BCValues, _) where {Index<:BoundaryIndex}
     local_facet_dofs, local_facet_dofs_offset =
         _local_facet_dofs_for_bc(interpolation, field_dim, dbc.components, offset, dirichlet_boundarydof_indices(eltype(bcfacets)))
     copy!(dbc.local_facet_dofs, local_facet_dofs)
@@ -335,7 +335,7 @@ function _local_facet_dofs_for_bc(interpolation, field_dim, components, offset,
     return local_facet_dofs, local_facet_dofs_offset
 end
 
-function _add!(ch::ConstraintHandler, dbc::Dirichlet, bcnodes::Set{Int}, interpolation::Interpolation, field_dim::Int, offset::Int, bcvalue::BCValues, cellset::Set{Int}=Set{Int}(1:getncells(get_grid(ch.dh))))
+function _add!(ch::ConstraintHandler, dbc::Dirichlet, bcnodes::AbstractVecOrSet{Int}, interpolation::Interpolation, field_dim::Int, offset::Int, bcvalue::BCValues, cellset::AbstractVecOrSet{Int}=OrderedSet{Int}(1:getncells(get_grid(ch.dh))))
     grid = get_grid(ch.dh)
     if interpolation !== default_interpolation(getcelltype(grid, first(cellset)))
         @warn("adding constraint to nodeset is not recommended for sub/super-parametric approximations.")
@@ -427,7 +427,7 @@ function update!(ch::ConstraintHandler, time::Real=0.0)
 end
 
 # for facets, vertices, faces and edges
-function _update!(inhomogeneities::Vector{T}, f::Function, boundary_entities::Set{<:BoundaryIndex}, field::Symbol, local_facet_dofs::Vector{Int}, local_facet_dofs_offset::Vector{Int},
+function _update!(inhomogeneities::Vector{T}, f::Function, boundary_entities::AbstractVecOrSet{<:BoundaryIndex}, field::Symbol, local_facet_dofs::Vector{Int}, local_facet_dofs_offset::Vector{Int},
                   components::Vector{Int}, dh::AbstractDofHandler, boundaryvalues::BCValues,
                   dofmapping::Dict{Int,Int}, dofcoefficients::Vector{Union{Nothing,DofCoefficients{T}}}, time::Real) where {T}
 
@@ -464,7 +464,7 @@ function _update!(inhomogeneities::Vector{T}, f::Function, boundary_entities::Se
 end
 
 # for nodes
-function _update!(inhomogeneities::Vector{T}, f::Function, ::Set{Int}, field::Symbol, nodeidxs::Vector{Int}, globaldofs::Vector{Int},
+function _update!(inhomogeneities::Vector{T}, f::Function, ::AbstractVecOrSet{Int}, field::Symbol, nodeidxs::Vector{Int}, globaldofs::Vector{Int},
                   components::Vector{Int}, dh::AbstractDofHandler, facetvalues::BCValues,
                   dofmapping::Dict{Int,Int}, dofcoefficients::Vector{Union{Nothing,DofCoefficients{T}}}, time::Real) where T
     counter = 1
@@ -859,7 +859,7 @@ end
 
 function filter_dbc_set(grid::AbstractGrid, fhset::AbstractSet{Int}, dbcset::AbstractSet{Int})
     ret = empty(dbcset)
-    nodes_in_fhset = Set{Int}()
+    nodes_in_fhset = OrderedSet{Int}()
     for cc in CellIterator(grid, fhset, UpdateFlags(; nodes=true, coords=false))
         union!(nodes_in_fhset, cc.nodes)
     end
@@ -1060,7 +1060,7 @@ function _add!(ch::ConstraintHandler, pdbc::PeriodicDirichlet, interpolation::In
         points = construct_cornerish(min_x, max_x)
         tree = KDTree(Tx[get_node_coordinate(grid, i) for i in all_node_idxs_v])
         idxs, _ = NearestNeighbors.nn(tree, points)
-        corner_set = Set{Int}(all_node_idxs_v[i] for i in idxs)
+        corner_set = OrderedSet{Int}(all_node_idxs_v[i] for i in idxs)
 
         dbc = Dirichlet(pdbc.field_name, corner_set,
             pdbc.func === nothing ? (x, _) -> pdbc.components * eltype(x)(0) : pdbc.func,
@@ -1087,7 +1087,7 @@ function _add!(ch::ConstraintHandler, pdbc::PeriodicDirichlet, interpolation::In
     if pdbc.func !== nothing
         # Create another temp constraint handler if we need to compute inhomogeneities
         chtmp2 = ConstraintHandler(ch.dh)
-        all_facets = Set{FacetIndex}()
+        all_facets = OrderedSet{FacetIndex}()
         union!(all_facets, (x.mirror for x in facet_map))
         union!(all_facets, (x.image for x in facet_map))
         dbc_all = Dirichlet(pdbc.field_name, all_facets, pdbc.func, pdbc.components)
@@ -1256,7 +1256,7 @@ dictionary which maps each mirror facet to a image facet. The result can then be
 [`PeriodicDirichlet`](@ref).
 
 `mset` and `iset` can be given as a `String` (an existing facet set in the grid) or as a
-`Set{FacetIndex}` directly.
+`AbstractSet{FacetIndex}` directly.
 
 By default this function looks for a matching facet in the directions of the coordinate
 system. For other types of periodicities the `transform` function can be used. The
@@ -1268,12 +1268,12 @@ between a image-facet and mirror-facet, for them to be considered matched.
 
 See also: [`collect_periodic_facets!`](@ref), [`PeriodicDirichlet`](@ref).
 """
-function collect_periodic_facets(grid::Grid, mset::Union{Set{FacetIndex},String}, iset::Union{Set{FacetIndex},String}, transform::Union{Function,Nothing}=nothing; tol::Float64=1e-12)
+function collect_periodic_facets(grid::Grid, mset::Union{AbstractSet{FacetIndex},String}, iset::Union{AbstractSet{FacetIndex},String}, transform::Union{Function,Nothing}=nothing; tol::Float64=1e-12)
     return collect_periodic_facets!(PeriodicFacetPair[], grid, mset, iset, transform; tol)
 end
 
 """
-    collect_periodic_facets(grid::Grid, all_facets::Union{Set{FacetIndex},String,Nothing}=nothing; tol=1e-12)
+    collect_periodic_facets(grid::Grid, all_facets::Union{AbstractSet{FacetIndex},String,Nothing}=nothing; tol=1e-12)
 
 Split all facets in `all_facets` into image and mirror sets. For each matching pair, the facet
 located further along the vector `(1, 1, 1)` becomes the image facet.
@@ -1283,7 +1283,7 @@ have a neighbor) is used.
 
 See also: [`collect_periodic_facets!`](@ref), [`PeriodicDirichlet`](@ref).
 """
-function collect_periodic_facets(grid::Grid, all_facets::Union{Set{FacetIndex},String,Nothing}=nothing; tol::Float64=1e-12)
+function collect_periodic_facets(grid::Grid, all_facets::Union{AbstractSet{FacetIndex},String,Nothing}=nothing; tol::Float64=1e-12)
     return collect_periodic_facets!(PeriodicFacetPair[], grid, all_facets; tol)
 end
 
@@ -1293,7 +1293,7 @@ end
 
 Same as [`collect_periodic_facets`](@ref) but adds all matches to the existing `facet_map`.
 """
-function collect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset::Union{Set{FacetIndex},String}, iset::Union{Set{FacetIndex},String}, transform::Union{Function,Nothing}=nothing; tol::Float64=1e-12)
+function collect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset::Union{AbstractSet{FacetIndex},String}, iset::Union{AbstractSet{FacetIndex},String}, transform::Union{Function,Nothing}=nothing; tol::Float64=1e-12)
     mset = __to_facetset(grid, mset)
     iset = __to_facetset(grid, iset)
     if transform === nothing
@@ -1306,15 +1306,15 @@ function collect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Gr
     return facet_map
 end
 
-function collect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, facetset::Union{Set{FacetIndex},String,Nothing}; tol::Float64=1e-12)
+function collect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, facetset::Union{AbstractSet{FacetIndex},String,Nothing}; tol::Float64=1e-12)
     facetset = facetset === nothing ? __collect_boundary_facets(grid) : copy(__to_facetset(grid, facetset))
     if mod(length(facetset), 2) != 0
         error("uneven number of facets")
     end
     return __collect_periodic_facets_bruteforce!(facet_map, grid, facetset, facetset, #=known_order=#false, tol)
 end
 
-__to_facetset(_, set::Set{FacetIndex}) = set
+__to_facetset(_, set::AbstractSet{FacetIndex}) = set
 __to_facetset(grid, set::String) = getfacetset(grid, set)
 function __collect_boundary_facets(grid::Grid)
     candidates = Dict{Tuple, FacetIndex}()
@@ -1328,7 +1328,7 @@ function __collect_boundary_facets(grid::Grid)
             end
         end
     end
-    return Set{FacetIndex}(values(candidates))
+    return OrderedSet{FacetIndex}(values(candidates))
 end
 
 function __collect_periodic_facets_tree!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset::Vector{FacetIndex}, iset::Vector{FacetIndex}, transformation::F, tol::Float64) where F <: Function
@@ -1373,7 +1373,7 @@ function __collect_periodic_facets_tree!(facet_map::Vector{PeriodicFacetPair}, g
 end
 
 # This method empties mset and iset
-function __collect_periodic_facets_bruteforce!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset::Set{FacetIndex}, iset::Set{FacetIndex}, known_order::Bool, tol::Float64)
+function __collect_periodic_facets_bruteforce!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset::AbstractSet{FacetIndex}, iset::AbstractSet{FacetIndex}, known_order::Bool, tol::Float64)
     if length(mset) != length(iset)
         error("different facets in mirror and image")
     end

src/Dofs/DofHandler.jl

@@ -15,7 +15,7 @@ get_grid(dh::AbstractDofHandler)
 struct SubDofHandler{DH} <: AbstractDofHandler
     # From constructor
     dh::DH
-    cellset::Set{Int}
+    cellset::OrderedSet{Int}
     # Populated in add!
     field_names::Vector{Symbol}
     field_interpolations::Vector{Interpolation}
@@ -26,7 +26,7 @@ struct SubDofHandler{DH} <: AbstractDofHandler
 end
 
 """
-    SubDofHandler(dh::AbstractDofHandler, cellset::Set{Int})
+    SubDofHandler(dh::AbstractDofHandler, cellset::AbstractVecOrSet{Int})
 
 Create an `sdh::SubDofHandler` from the parent `dh`, pertaining to the 
 cells in `cellset`. This allows you to add fields to parts of the domain, or using 
@@ -53,7 +53,7 @@ add!(sdh_quad, :u, ip_quad)
 close!(dh) # Finalize by closing the parent 
 ```
 """
-function SubDofHandler(dh::DH, cellset) where {DH <: AbstractDofHandler}
+function SubDofHandler(dh::DH, cellset::AbstractVecOrSet{Int}) where {DH <: AbstractDofHandler}
     # TODO: Should be an inner constructor.
     isclosed(dh) && error("DofHandler already closed")
     # Compute the celltype and make sure all elements have the same one
@@ -68,7 +68,7 @@ function SubDofHandler(dh::DH, cellset) where {DH <: AbstractDofHandler}
         end
     end
     # Construct and insert into the parent dh
-    sdh = SubDofHandler{typeof(dh)}(dh, cellset, Symbol[], Interpolation[], Int[], ScalarWrapper(-1))
+    sdh = SubDofHandler{typeof(dh)}(dh, convert_to_orderedset(cellset), Symbol[], Interpolation[], Int[], ScalarWrapper(-1))
     push!(dh.subdofhandlers, sdh)
     return sdh
 end
@@ -298,7 +298,7 @@ function add!(dh::DofHandler, name::Symbol, ip::Interpolation)
     @assert isconcretetype(celltype)
     if isempty(dh.subdofhandlers)
         # Create a new SubDofHandler for all cells
-        sdh = SubDofHandler(dh, Set(1:getncells(get_grid(dh))))
+        sdh = SubDofHandler(dh, OrderedSet(1:getncells(get_grid(dh))))
     elseif length(dh.subdofhandlers) == 1
         # Add to existing SubDofHandler (if it covers all cells)
         sdh = dh.subdofhandlers[1]
@@ -437,8 +437,7 @@ function _close_subdofhandler!(dh::DofHandler{sdim}, sdh::SubDofHandler, sdh_ind
     global_fidxs = Int[findfirst(gname -> gname === lname, dh.field_names) for lname in sdh.field_names]
 
     # loop over all the cells, and distribute dofs for all the fields
-    # TODO: Remove BitSet construction when SubDofHandler ensures sorted collections
-    for ci in BitSet(sdh.cellset)
+    for ci in sdh.cellset
         @debug println("Creating dofs for cell #$ci")
 
         # TODO: _check_cellset_intersections can be removed in favor of this assertion

src/Ferrite.jl

@@ -12,6 +12,8 @@ using LinearAlgebra:
     pinv, tr
 using NearestNeighbors:
     NearestNeighbors, KDTree, knn
+using OrderedCollections:
+    OrderedSet
 using SparseArrays:
     SparseArrays, SparseMatrixCSC, nonzeros, nzrange, rowvals, sparse, spzeros
 using StaticArrays:
@@ -92,6 +94,9 @@ struct FacetIndex <: BoundaryIndex
     idx::Tuple{Int,Int} # cell and side
 end
 
+const AbstractVecOrSet{T} = Union{AbstractSet{T}, AbstractVector{T}}
+const IntegerCollection = AbstractVecOrSet{<:Integer}
+
 include("utils.jl")
 
 # Matrix/Vector utilities

src/Grid/grid.jl

@@ -292,30 +292,34 @@ There are multiple helper structures to apply boundary conditions or define subd
 - `nodesets::Dict{String,Set{Int}}`: maps a `String` key to a `Set` of global node ids
 - `facetsets::Dict{String,Set{FacetIndex}}`: maps a `String` to a `Set` of `Set{FacetIndex} (global_cell_id, local_facet_id)`
 - `vertexsets::Dict{String,Set{VertexIndex}}`: maps a `String` key to a `Set` of local vertex ids
+- `cellsets::Dict{String,OrderedSet{Int}}`: maps a `String` key to an `OrderedSet` of cell ids
+- `nodesets::Dict{String,OrderedSet{Int}}`: maps a `String` key to an `OrderedSet` of global node ids
+- `facetsets::Dict{String,OrderedSet{FacetIndex}}`: maps a `String` to an `OrderedSet` of `FacetIndex`
+- `vertexsets::Dict{String,OrderedSet{VertexIndex}}`: maps a `String` key to an `OrderedSet` of `VertexIndex`
 """
 mutable struct Grid{dim,C<:AbstractCell,T<:Real} <: AbstractGrid{dim}
     cells::Vector{C}
     nodes::Vector{Node{dim,T}}
     # Sets
-    cellsets::Dict{String,Set{Int}}
-    nodesets::Dict{String,Set{Int}}
-    facetsets::Dict{String,Set{FacetIndex}}
-    vertexsets::Dict{String,Set{VertexIndex}}
+    cellsets::Dict{String,OrderedSet{Int}}
+    nodesets::Dict{String,OrderedSet{Int}}
+    facetsets::Dict{String,OrderedSet{FacetIndex}}
+    vertexsets::Dict{String,OrderedSet{VertexIndex}}
 end
 
 function Grid(cells::Vector{C},
               nodes::Vector{Node{dim,T}};
-              cellsets::Dict{String,Set{Int}}=Dict{String,Set{Int}}(),
-              nodesets::Dict{String,Set{Int}}=Dict{String,Set{Int}}(),
-              facetsets::Dict{String,Set{FacetIndex}}=Dict{String,Set{FacetIndex}}(),
-              facesets = nothing,
-              vertexsets::Dict{String,Set{VertexIndex}}=Dict{String,Set{VertexIndex}}(),
+              cellsets::Dict{String, <:AbstractVecOrSet{Int}}=Dict{String,OrderedSet{Int}}(),
+              nodesets::Dict{String, <:AbstractVecOrSet{Int}}=Dict{String,OrderedSet{Int}}(),
+              facetsets::Dict{String, <:AbstractVecOrSet{FacetIndex}}=Dict{String,OrderedSet{FacetIndex}}(),
+              facesets=nothing, # deprecated
+              vertexsets::Dict{String, <:AbstractVecOrSet{VertexIndex}}=Dict{String,OrderedSet{VertexIndex}}(),
               boundary_matrix = nothing) where {dim,C,T}
     if facesets !== nothing 
         if isempty(facetsets)
             @warn "facesets in Grid is deprecated, use facetsets instead" maxlog=1
             for (key, set) in facesets
-                facetsets[key] = Set(FacetIndex(cellnr, facenr) for (cellnr, facenr) in set)
+                facetsets[key] = OrderedSet(FacetIndex(cellnr, facenr) for (cellnr, facenr) in set)
             end
         else
             error("facesets are deprecated, use only facetsets")
@@ -324,7 +328,14 @@ function Grid(cells::Vector{C},
     if boundary_matrix !== nothing
         error("`boundary_matrix` is not part of the Grid anymore and thus not a supported keyword argument.")
     end
-    return Grid(cells, nodes, cellsets, nodesets, facetsets, vertexsets)
+    return Grid(
+        cells,
+        nodes,
+        convert_to_orderedsets(cellsets),
+        convert_to_orderedsets(nodesets),
+        convert_to_orderedsets(facetsets),
+        convert_to_orderedsets(vertexsets),
+    )
 end
 
 ##########################
@@ -416,7 +427,7 @@ end
 """
     getcellset(grid::AbstractGrid, setname::String)
 
-Returns all cells as cellid in a `Set` of a given `setname`.
+Returns all cells as cellid in the set with name `setname`.
 """
 @inline getcellset(grid::AbstractGrid, setname::String) = grid.cellsets[setname]
 """
@@ -429,7 +440,7 @@ Returns all cellsets of the `grid`.
 """
     getnodeset(grid::AbstractGrid, setname::String)
 
-Returns all nodes as nodeid in a `Set` of a given `setname`.
+Returns all nodes as nodeid in the set with name `setname`.
 """
 @inline getnodeset(grid::AbstractGrid, setname::String) = grid.nodesets[setname]
 """
@@ -442,7 +453,7 @@ Returns all nodesets of the `grid`.
 """
     getfacetset(grid::AbstractGrid, setname::String)
 
-Returns all facets as `FacetIndex` in a `Set` of a given `setname`.
+Returns all faces as `FacetIndex` in the set with name `setname`.
 """
 @inline getfacetset(grid::AbstractGrid, setname::String) = grid.facetsets[setname]
 """
@@ -456,7 +467,7 @@ Returns all facet sets of the `grid`.
 """
     getvertexset(grid::AbstractGrid, setname::String)
 
-Returns all vertices as `VertexIndex` in a `Set` of a given `setname`.
+Returns all vertices as `VertexIndex` in the set with name `setname`.
 """
 @inline getvertexset(grid::AbstractGrid, setname::String) = grid.vertexsets[setname]
 """
@@ -535,7 +546,7 @@ function Base.show(io::IO, ::MIME"text/plain", grid::Grid)
     if isconcretetype(eltype(grid.cells))
         typestrs = [repr(eltype(grid.cells))]
     else
-        typestrs = sort!(repr.(Set(typeof(x) for x in grid.cells)))
+        typestrs = sort!(repr.(OrderedSet(typeof(x) for x in grid.cells)))
     end
     join(io, typestrs, '/')
     print(io, " cells and $(getnnodes(grid)) nodes")