From 74f7009d0a5793a04b7a98e256638870d448d008 Mon Sep 17 00:00:00 2001 From: Lionel Zoubritzky Date: Fri, 30 Oct 2020 12:44:10 +0100 Subject: [PATCH] Add specialization for 1D; various tests and minor fixes --- src/PeriodicGraphs.jl | 46 +++++++++------ test/runtests.jl | 126 +++++++++++++++++++++++++++++++++++++++--- 2 files changed, 147 insertions(+), 25 deletions(-) diff --git a/src/PeriodicGraphs.jl b/src/PeriodicGraphs.jl index b63e25e..2c30bc7 100644 --- a/src/PeriodicGraphs.jl +++ b/src/PeriodicGraphs.jl @@ -4,6 +4,7 @@ using LinearAlgebra using StaticArrays using LightGraphs export PeriodicVertex, PeriodicEdge, PeriodicGraph, + PeriodicVertex1D, PeriodicEdge1D, PeriodicGraph1D, PeriodicVertex2D, PeriodicEdge2D, PeriodicGraph2D, PeriodicVertex3D, PeriodicEdge3D, PeriodicGraph3D, coordination_sequence, cellgraph, periodiccellgraph, @@ -49,6 +50,7 @@ end isless(x::PeriodicVertex{N}, y::PeriodicVertex{N}) where {N} = cmp(x, y) < 0 ==(x::PeriodicVertex{N}, y::PeriodicVertex{M}) where {N,M} = N == M && iszero(cmp(x, y)) +const PeriodicVertex1D = PeriodicVertex{1} const PeriodicVertex2D = PeriodicVertex{2} const PeriodicVertex3D = PeriodicVertex{3} @@ -71,6 +73,14 @@ function hash_position((x1,x2)::SVector{2,<:Integer}) return b*(x2 + x1^2) + (!b)*(x1 + x2*(x2 + 1) + 1) end +function hash_position((x,)::SVector{1,<:Integer}) + return ZtoN(x) +end + +function hash_position(::SVector{0,<:Integer}) + return 0 +end + """ hash_position(x::PeriodicVertex{N}, n::Integer) where {N} @@ -157,13 +167,11 @@ function show(io::IO, x::PeriodicEdge{N}) where N print(io, '(', x.src, ", ", x.dst.v, ", (", join(x.dst.ofs, ','), ')', ')') end +const PeriodicEdge1D = PeriodicEdge{1} const PeriodicEdge2D = PeriodicEdge{2} const PeriodicEdge3D = PeriodicEdge{3} function LightGraphs.reverse(e::PeriodicEdge{N}) where N - if e.src == e.dst - return e - end return unsafe_edge{N}(e.dst.v, e.src, .-e.dst.ofs) end LightGraphs.src(e::PeriodicEdge) = e.src @@ -199,6 +207,7 @@ struct PeriodicGraph{N} <: AbstractGraph{Int} width::Base.RefValue{Rational{Int}} end +const PeriodicGraph1D = PeriodicGraph{1} const PeriodicGraph2D = PeriodicGraph{2} const PeriodicGraph3D = PeriodicGraph{3} @@ -285,14 +294,18 @@ function iterate(k::KeyString{T}, (next_char, idx)) where T tmp = idx next_char, idx = state end - return (parse(T, SubString(k.x, start:stop)), state) + ret = tryparse(T, SubString(k.x, start:stop)) + if ret isa T + return (ret, state) + end + throw(ArgumentError("Intput string does not represent a graph")) end function iterate(k::KeyString) iterate(k, (' ', k.start[])) end function Base.popfirst!(k::KeyString) next = iterate(k) - next isa Nothing && throw(ArgumentError("collection must be non-empty")) + next isa Nothing && throw(ArgumentError("Intput string does not represent a graph")) char, state = next k.start[] = state isa Nothing ? (ncodeunits(k.x) + 1) : last(state) return char @@ -421,19 +434,19 @@ end """ find_edges(g::PeriodicGraph, s::Int, d::Int) -Return the set of edges of graph `g` such that the source vertex identifier is `s` -and the destination vertex identifier is `d`. +Return the set of PeriodicVertex `v` of graph `g` such that there is an edge +between a source vertex of identifier `s` and `v`, and the identifier of `v` is `d`. """ -function find_edges(g::PeriodicGraph, s::Int, d::Int) +function find_edges(g::PeriodicGraph{N}, s::Int, d::Int) where N ((s < 1) | (s > nv(g))) && return false #=@inbounds=# begin start = g.directedgestart[s] lo, hi = s > d ? (1, start-1) : (start, lastindex(g.nlist[s])) - rng = searchsorted(g.nlist[s], d, lo, hi, Lt((x,y)->isless(x.v, y))) - if s == d.v + rng = searchsorted(g.nlist[s], d, lo, hi, Lt((x,y)->isless(x isa Integer ? x : x.v, y isa Integer ? y : y.v))) + if s == d rng = (2*first(rng) - last(rng) - 1):last(rng) end - return g.nlists[rng] + return g.nlist[s][rng] end end function LightGraphs.has_edge(g::PeriodicGraph, e::PeriodicEdge) @@ -768,7 +781,7 @@ end @noinline __throw_invalid_offsets() = throw(ArgumentError("The size of offsets does not match the number of vertices")) """ - swap_axes!(g::PeriodicGraph, t::AbstractVector{<:Integer}) + swap_axes!(g::PeriodicGraph, t) In-place modifies graph `g` so that the new initial cell corresponds to the previous one with its axes swapped according to the permutation `t`. @@ -776,13 +789,14 @@ one with its axes swapped according to the permutation `t`. Note that the resulting graph is isomorphic to the initial one, only the representation has changed. """ -function swap_axes!(g::PeriodicGraph{N}, t::AbstractVector{<:Integer}) where N +function swap_axes!(g::PeriodicGraph{N}, t) where N length(t) == N || __throw_invalid_axesswap() + tt::Vector{Int} = t isa Vector ? t : collect(t) for i in vertices(g) neighs = g.nlist[i] for j in 1:length(neighs) x = neighs[j] - neigh = PeriodicVertex{N}(x.v, x.ofs[t]) + neigh = PeriodicVertex{N}(x.v, x.ofs[tt]) neighs[j] = neigh end sort!(neighs) @@ -836,7 +850,7 @@ the fields). """ function graph_width!(g::PeriodicGraph{N}) where N distances = floyd_warshall_shortest_paths(cellgraph(g)).dists - extremalpoints = NTuple{N,NTuple{2,Vector{Tuple{Int,Int}}}}((([],[]),([],[]),([],[]))) + extremalpoints = NTuple{N,NTuple{2,Vector{Tuple{Int,Int}}}}([([],[]) for _ in 1:N]) # a, x ∈ extremalpoints[i][j] where i ∈ ⟦1,N⟧ and j ∈ ⟦1,2⟧ means that # vertex x has a neighbor whose offset is a*(-1)^(j-1) along dimension i maxa = 1 @@ -878,7 +892,7 @@ function graph_width!(g::PeriodicGraph{N}) where N g.width[] = width == nv(g)+1 ? Rational(maxa) : width end -function LightGraphs._neighborhood(g::Union{PeriodicGraph{2}, PeriodicGraph{3}}, v::Integer, d::Real, distmx::AbstractMatrix{U}, neighborfn::Function) where U <: Real +function LightGraphs._neighborhood(g::Union{PeriodicGraph{0},PeriodicGraph{1},PeriodicGraph{2},PeriodicGraph{3}}, v::Integer, d::Real, distmx::AbstractMatrix{U}, neighborfn::Function) where U <: Real @assert typeof(neighborfn) === typeof(outneighbors) N = ndims(g) Q = Tuple{PeriodicVertex{N}, U}[] diff --git a/test/runtests.jl b/test/runtests.jl index 3900515..1b6e3c5 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -20,7 +20,7 @@ end @test ndims(PeriodicVertex2D(1)) == 2 @test PeriodicVertex2D(1) != PeriodicVertex3D(1) @test PeriodicVertex{4}(5) == PeriodicVertex{4}(5) - @test PeriodicVertex{1}(1, SVector{1,Int}(3)) == PeriodicVertex{1}(1, [3]) + @test PeriodicVertex1D(1, SVector{1,Int}(3)) == PeriodicVertex1D(1, [3]) @test PeriodicVertex(3, (1,2)) == PeriodicVertex2D(3, (1,2)) == PeriodicVertex(3, SVector{2,Int}(1,2)) @test PeriodicVertex[(1, (-1,0))] == [PeriodicVertex2D(1, (-1,0))] == PeriodicVertex2D[(1, (-1,0))] @test_throws MethodError PeriodicVertex(2, [4]) # Performance trap if the dimension is unknown @@ -45,14 +45,23 @@ end @test_throws LoopException PeriodicEdge(2, PeriodicVertex3D(2)) @test_throws LoopException PeriodicEdge{0}(3, 3, Int[]) @test_throws LoopException PeriodicEdge(1, 1, ()) + try + PeriodicEdge(1, PeriodicVertex2D(1)) + @test false + catch e + @test e isa LoopException + io = IOBuffer() + showerror(io, e) + @test occursin("LoopException", String(take!(io))) + end # There should be no default zero offset to prevent programmers from forgetting the offset - @test_throws MethodError PeriodicEdge{1}(1, 2) + @test_throws MethodError PeriodicEdge1D(1, 2) @test PeriodicEdge2D(1, 2, (0,0)) != PeriodicEdge3D(1, 2, (0,0,0)) @test PeriodicEdge(2, PeriodicVertex2D(1, (1,3))) == PeriodicEdge2D(2, 1, (1,3)) - @test PeriodicEdge(1, 1, SVector{1,Int}(2)) == PeriodicEdge{1}(1, 1, SVector{1,Int}(2)) + @test PeriodicEdge(1, 1, SVector{1,Int}(2)) == PeriodicEdge1D(1, 1, SVector{1,Int}(2)) @test PeriodicEdge(2, 1, ()) == PeriodicEdge{0}((2, 1, SVector{0,Int}())) @test PeriodicEdge((1, 2, (1,0,0))) == PeriodicEdge3D((1, 2, (1,0,0))) - @test PeriodicEdge{1}[(1,2,SVector{1,Int}(3))] == PeriodicEdge[(1,2,(3,))] == [PeriodicEdge(1, 2, (3,))] + @test PeriodicEdge1D[(1,2,SVector{1,Int}(3))] == PeriodicEdge[(1,2,(3,))] == [PeriodicEdge(1, 2, (3,))] @static if VERSION < v"1.6.0-DEV" @test string(PeriodicEdge3D[(1, 2, (1,0,0))]) == "PeriodicEdge{3}[(1, 2, (1,0,0))]" @test string(PeriodicEdge(3,4, (0,0))) == "PeriodicEdge{2}(3, 4, (0,0))" @@ -62,6 +71,7 @@ end end @test reverse(reverse(PeriodicEdge(1, 1, (2,3,1)))) == PeriodicEdge3D(1, 1, (2,3,1)) @test reverse(PeriodicEdge(1, 2, (1,-2))) == PeriodicEdge(2, 1, (-1,2)) + @test reverse(PeriodicEdge(3, 3, (1,0))) == PeriodicEdge(3, 3, (-1,0)) @test src(PeriodicEdge(1, 2, (-1, 0))) == 1 @test dst(PeriodicEdge(2, 2, (1,))) == 2 @test ofs(PeriodicEdge(1, 1, (0,1,3))) == SVector{3,Int}(0,1,3) @@ -92,6 +102,7 @@ end @test ne(gg) == nv(gg) == 0 && isempty(edges(gg)) @test ne(g) == nv(g) == 0 && isempty(edges(g)) @test g == gg + @test_throws ArgumentError swap_axes!(g, [2,1]) gg = swap_axes!(g, Int[1,2,3]) @test ne(gg) == nv(gg) == 0 && isempty(edges(gg)) @test ne(g) == nv(g) == 0 && isempty(edges(g)) @@ -161,6 +172,7 @@ end @test !has_edge(g, 1, 1) @test !has_edge(g, PeriodicEdge3D(1, 1, (0,1,1))) @test add_edge!(g, PeriodicEdge(1, 1, (-1,0,0))) + @test has_edge(g, 1, 1) @test !add_edge!(g, PeriodicEdge(1, 1, (-1,0,0))) @test has_edge(g, 1, 1) @test has_edge(g, PeriodicEdge3D(1, 1, (1,0,0))) @@ -217,12 +229,15 @@ end end @testset "String graph construction" begin + @test_throws MethodError length(PeriodicGraphs.KeyString{Int}("3 1 2 0 0 0")) @test PeriodicGraph("0 1 2 1 3") == PeriodicGraph(PeriodicEdge{0}[(1, 2, ()), (1, 3, ())]) @test PeriodicGraph("3") == PeriodicGraph3D("3") == PeriodicGraph3D() @test ne(PeriodicGraph("3 1 1 1 0 0 1 1 -1 0 0 1 1 3 0 0")) == 2 @test nv(PeriodicGraph("2 11 12 0 0")) == 12 @test_throws ArgumentError PeriodicGraph2D("") @test_throws ArgumentError PeriodicGraph("") + @test_throws ArgumentError PeriodicGraph3D("a") + @test_throws ArgumentError PeriodicGraph2D("2 a") @test_throws DimensionMismatch PeriodicGraph2D("4 1 2 0 0 0 1") @test_throws DimensionMismatch PeriodicGraph3D("2 1 2 0 0 0 1 3 0 1 0") @test_throws ArgumentError PeriodicGraph("1 2 1") @@ -234,6 +249,27 @@ end @test string(PeriodicGraph{5}()) == "5" end +@testset "Basic graph representation and functions" begin + @test zero(PeriodicGraph{6}) == PeriodicGraph{6}() + io = IOBuffer() + show(io, PeriodicGraph("2 1 2 0 0 2 3 1 0")) + g1 = eval(Meta.parse(String(take!(io)))) + @test g1 isa PeriodicGraph2D + @test string(g1) == "2 1 2 0 0 2 3 1 0" + g2 = PeriodicGraph(g1) + g3 = deepcopy(g1) + @test g1 == g2 == g3 + rem_vertex!(g1, 1) + @test g1 != g2 == g3 + @test add_edge!(g2, PeriodicEdge(1, 2, (-1,0))) + @test add_edge!(g2, PeriodicEdge(1, 1, (1,1))) + @test add_edge!(g2, PeriodicEdge(1, 3, (0,0))) + @test find_edges(g2, 1, 1) == PeriodicVertex2D[(1, (-1,-1)), (1, (1,1))] + @test find_edges(g2, 1, 2) == PeriodicVertex2D[(2, (-1,0)), (2, (0,0))] + @test find_edges(g2, 3, 2) == PeriodicVertex2D[(2, (-1,0))] + @test inneighbors(g2, 2) == outneighbors(g2, 2) +end + @testset "Complex graphs modifications" begin edgs = PeriodicEdge3D[(1, 1, (1,0,0)), (1, 2, (0,0,0)), (1, 3, (0,0,1)), (1, 3, (0,1,-1)), (2, 3, (0,0,0)), (2, 4, (0,-1,0))] @@ -249,6 +285,30 @@ end deleteat!(edgs, 4) @test ne(g) == length(edges(g)) == length(edgs) @test collect(edges(g)) == edgs + @test g[[1,2,3,4]] == g + gg = g[[1,4,2,3]] + @test gg != g + @test gg[[1,3,4,2]] == g + gg = g[[3,1]] + @test nv(gg) == 2 + @test collect(edges(gg)) == PeriodicEdge3D[(1, 2, (0,0,-1)), (2, 2, (1,0,0))] + @test_throws ArgumentError g[[3,3]] + @test_throws ArgumentError offset_representatives!(gg, [2]) + gcopy = deepcopy(g) + @test offset_representatives!(g, [0,0,0,0]) == gcopy == g + ggg = offset_representatives!(gg, [(1,-1,1),0]) + @test ggg == gg + @test collect(edges(gg)) == PeriodicEdge3D[(1, 2, (1,-1,0)), (2, 2, (1,0,0))] + gcopy = deepcopy(gg) + offset_representatives!(gg, [(2,1,-3), (2,1,-3)]) + @test gg == gcopy + @test_throws ArgumentError swap_axes!(g, [1,3]) + ggg = swap_axes!(gg, (1,2,3)) + @test gcopy == gg == ggg + ggg = swap_axes!(gg, (2,1,3)) + @test gcopy != gg == ggg + @test collect(edges(gg)) == PeriodicEdge3D[(1, 2, (-1,1,0)), (2, 2, (0,1,0))] + @test edges(gg) < edges(gcopy) # lexicographical ordering on the list of edges end @@ -264,6 +324,29 @@ end @test length(neighbors(g, 2)) == 3 @test has_edge(g, 2, 2) @test neighborhood(g, 1, 1) == pushfirst!(collect(neighbors(g, 1)), PeriodicVertex3D(1, (0, 0, 0))) + @test coordination_sequence(g, 1, 3) == zeros(Int, 3) + @test coordination_sequence(g, 3, 6) == [1, 2, 4, 4, 4, 4] + + g = PeriodicGraph{0}(4) + @test add_edge!(g, PeriodicEdge{0}(1, 2, ())) + @test add_edge!(g, PeriodicEdge{0}(1, 3, ())) + @test add_edge!(g, PeriodicEdge{0}(2, 4, ())) + @test coordination_sequence(g, 1, 10) == [2, 1, 0, 0, 0, 0, 0, 0, 0, 0,] + + g = PeriodicGraph1D(3) + @test add_edge!(g, PeriodicEdge(1, 1, (1,))) + @test add_edge!(g, PeriodicEdge(1, 2, (0,))) + @test add_edge!(g, PeriodicEdge(2, 3, (0,))) + @test add_edge!(g, PeriodicEdge(3, 3, (1,))) + @test coordination_sequence(g, 1, 5) == coordination_sequence(g, 3, 5) == [3, 5, 6, 6, 6] + @test coordination_sequence(g, 2, 5) == [2, 4, 6, 6, 6] + + g = PeriodicGraph{4}(1) + @test add_edge!(g, PeriodicEdge(1, 1, (1,0,0,0))) + @test add_edge!(g, PeriodicEdge(1, 1, (0,1,0,0))) + @test add_edge!(g, PeriodicEdge(1, 1, (0,0,1,0))) + @test add_edge!(g, PeriodicEdge(1, 1, (0,0,0,1))) + @test coordination_sequence(g, 1, 30) == [8i*(i^2+2)/3 for i in 1:30] # sequence A008412 of the OEIS end @testset "Graph reduction" begin @@ -280,7 +363,37 @@ end end @testset "Periodic vertex hash" begin + ## dim = 0 + @test PeriodicGraphs.hash_position(SVector{0,Int}()) == 0 + ## dim = 1 + @test PeriodicGraphs.hash_position(SVector{1,Int}(0)) == 0 + @test PeriodicGraphs.hash_position(SVector{1,Int}(-1)) == 1 + @test PeriodicGraphs.hash_position(SVector{1,Int}(1)) == 2 + @test PeriodicGraphs.hash_position(SVector{1,Int}(-2)) == 3 + @test PeriodicGraphs.hash_position(SVector{1,Int}(2)) == 4 + @test PeriodicGraphs.hash_position(SVector{1,Int}(3)) == 6 + ## dim = 2 n::Int = 4 + seen = falses((2n-1)^2) + for i in -n+1:n-1, j in -n+1:n-1 + x = PeriodicGraphs.hash_position(SVector{2,Int}(i,j))+1 + @test !seen[x] + seen[x] = true + end + @test all(seen) + append!(seen, falses((2n+1)^2-(2n-1)^2)) + for i in (-n,n), j in -n:n + x = PeriodicGraphs.hash_position(SVector{2,Int}(i,j))+1 + @test !seen[x] + seen[x] = true + end + for i in -n+1:n-1, j in (-n,n) + x = PeriodicGraphs.hash_position(SVector{2,Int}(i,j))+1 + @test !seen[x] + seen[x] = true + end + @test all(seen) + ## dim = 3 seen = falses((2n-1)^3) for i in -n+1:n-1, j in -n+1:n-1, k in -n+1:n-1 x = PeriodicGraphs.hash_position(SVector{3,Int}(i,j,k))+1 @@ -305,11 +418,6 @@ end seen[x] = true end - for i in 1:length(seen) - if !seen[i] - @show i - end - end @test all(seen) g = PeriodicGraph3D([PeriodicEdge3D(1, 1, (0, 0, 1)), PeriodicEdge3D(1, 1, (1, 1, 0))])