distinguish between null and unspecified selectors (#254)

more fastpaths and tests

src/apply.jl

@@ -18,7 +18,11 @@ function apply(s::SiteSelector, lat::Lattice{T,E,L}, cells...) where {T,E,L}
     cells = recursive_push!(SVector{L,Int}[], sanitize_cells(s.cells, Val(L)), cells...)
     unique!(sort!(sublats))
     unique!(sort!(cells))
-    return AppliedSiteSelector{T,E,L}(lat, region, sublats, cells)
+    # isnull: to distinguish in a type-stable way between s.cells === missing and no-selected-cells
+    # and the same for sublats
+    isnull = (s.cells !== missing && isempty(cells)) ||
+        (s.sublats !== missing && isempty(sublats))
+    return AppliedSiteSelector{T,E,L}(lat, region, sublats, cells, isnull)
 end
 
 function apply(s::HopSelector, lat::Lattice{T,E,L}, cells...) where {T,E,L}
@@ -36,7 +40,9 @@ function apply(s::HopSelector, lat::Lattice{T,E,L}, cells...) where {T,E,L}
         sublats .= reverse.(sublats)
         dcells .*= -1
     end
-    return AppliedHopSelector{T,E,L}(lat, region, sublats, dcells, (rmin, rmax))
+    isnull = (s.dcells !== missing && isempty(dcells)) ||
+        (s.sublats !== missing && isempty(sublats))
+    return AppliedHopSelector{T,E,L}(lat, region, sublats, dcells, (rmin, rmax), isnull)
 end
 
 sublatindex_or_zero(lat, ::Missing) = missing
@@ -66,7 +72,8 @@ applyrange(r::Real, lat) = r
 padrange(r::Real, m) = isfinite(r) ? float(r) + m * sqrt(eps(float(r))) : float(r)
 
 applied_region(r, ::Missing) = true
-applied_region((r, dr)::Tuple{SVector,SVector}, region::Function) = ifelse(region(r, dr), true, false)
+applied_region((r, dr)::Tuple{SVector,SVector}, region::Function) =
+    ifelse(region(r, dr), true, false)
 applied_region(r::SVector, region::Function) = ifelse(region(r), true, false)
 
 recursive_apply(f, t::Tuple) = recursive_apply.(f, t)
@@ -181,6 +188,7 @@ end
 function pointers(h::Hamiltonian{T,E,L,B}, s::AppliedSiteSelector{T,E,L}, shifts) where {T,E,L,B}
     isempty(cells(s)) || argerror("Cannot constrain cells in an onsite modifier, cell periodicity is assumed.")
     ptrs = Tuple{Int,SVector{E,T},CellSitePos{T,E,L,B},Int}[]
+    isnull(s) && return ptrs
     lat = lattice(h)
     har0 = first(harmonics(h))
     dn0 = zerocell(lat)
@@ -204,6 +212,7 @@ end
 function pointers(h::Hamiltonian{T,E,L,B}, s::AppliedHopSelector{T,E,L}, shifts) where {T,E,L,B}
     hars = harmonics(h)
     ptrs = [Tuple{Int,SVector{E,T},SVector{E,T},CellSitePos{T,E,L,B},CellSitePos{T,E,L,B},Tuple{Int,Int}}[] for _ in hars]
+    isnull(s) && return ptrs
     lat = lattice(h)
     dn0 = zerocell(lat)
     norbs = norbitals(h)

src/selectors.jl

@@ -24,12 +24,12 @@ neighbors(n::Int, lat::Lattice) = nrange(n, lat)
 #region
 
 function Base.in((s, r)::Tuple{Int,SVector{E,T}}, sel::AppliedSiteSelector{T,E}) where {T,E}
-    return inregion(r, sel) &&
+    return !isnull(sel) && inregion(r, sel) &&
            insublats(s, sel)
 end
 
 function Base.in((s, r, cell)::Tuple{Int,SVector{E,T},SVector{L,Int}}, sel::AppliedSiteSelector{T,E,L}) where {T,E,L}
-    return incells(cell, sel) &&
+    return !isnull(sel) && incells(cell, sel) &&
            inregion(r, sel) &&
            insublats(s, sel)
 end
@@ -51,7 +51,7 @@ end
 # end
 
 function Base.in(((sj, si), (r, dr), dcell)::Tuple{Pair,Tuple,SVector}, sel::AppliedHopSelector)
-    return !isonsite(dr) &&
+    return !isnull(sel) && !isonsite(dr) &&
             indcells(dcell, sel) &&
             insublats(sj => si, sel) &&
             iswithinrange(dr, sel) &&
@@ -70,6 +70,7 @@ isonsite(dr) = iszero(dr)
 # foreach_cell(f,...) should be called with a boolean function f that returns whether the
 # cell should be mark as accepted when BoxIterated
 function foreach_cell(f, sel::AppliedSiteSelector)
+    isnull(sel) && return nothing
     lat = lattice(sel)
     cells_list = cells(sel)
     if isempty(cells_list)      # no cells specified
@@ -92,6 +93,7 @@ end
 
 
 function foreach_cell(f, sel::AppliedHopSelector)
+    isnull(sel) && return nothing
     lat = lattice(sel)
     dcells_list = dcells(sel)
     if isempty(dcells_list) # no dcells specified
@@ -109,6 +111,7 @@ function foreach_cell(f, sel::AppliedHopSelector)
 end
 
 function foreach_site(f, sel::AppliedSiteSelector, cell::SVector)
+    isnull(sel) && return nothing
     lat = lattice(sel)
     for s in sublats(lat)
         insublats(s, sel) || continue
@@ -122,6 +125,7 @@ function foreach_site(f, sel::AppliedSiteSelector, cell::SVector)
 end
 
 function foreach_hop(f, sel::AppliedHopSelector, kdtrees::Vector{<:KDTree}, ni::SVector = zerocell(lattice(sel)))
+    isnull(sel) && return nothing
     lat = lattice(sel)
     _, rmax = sel.range
     # source cell at origin

src/slices.jl

@@ -14,18 +14,20 @@ Base.getindex(lat::Lattice, ::SiteSelectorAll) = lat[siteselector(; cells = zero
 function Base.getindex(lat::Lattice, as::AppliedSiteSelector)
     L = latdim(lat)
     csites = CellSites{L,Vector{Int}}[]
-    sinds = Int[]
-    foreach_cell(as) do cell
-        isempty(sinds) || (sinds = Int[])
-        cs = CellSites(cell, sinds)
-        foreach_site(as, cell) do s, i, r
-            push!(siteindices(cs), i)
-        end
-        if isempty(cs)
-            return false
-        else
-            push!(csites, cs)
-            return true
+    if !isnull(as)
+        sinds = Int[]
+        foreach_cell(as) do cell
+            isempty(sinds) || (sinds = Int[])
+            cs = CellSites(cell, sinds)
+            foreach_site(as, cell) do s, i, r
+                push!(siteindices(cs), i)
+            end
+            if isempty(cs)
+                return false
+            else
+                push!(csites, cs)
+                return true
+            end
         end
     end
     cellsdict = CellSitesDict{L}(cell.(csites), csites)
@@ -60,6 +62,7 @@ function Base.getindex(latslice::SiteSlice, as::AppliedSiteSelector)
     lat = parent(latslice)
     L = latdim(lat)
     cellsdict´ = CellSitesDict{L}()
+    isnull(as) && return SiteSlice(lat, cellsdict´)
     sinds = Int[]
     for subcell in cellsdict(latslice)
         dn = cell(subcell)
@@ -82,6 +85,7 @@ function Base.getindex(latslice::OrbitalSliceGrouped, as::AppliedSiteSelector)
     lat = parent(latslice)
     L = latdim(lat)
     cellsdict´ = CellOrbitalsGroupedDict{L}()
+    isnull(as) && return OrbitalSliceGrouped(lat, cellsdict´)
     oinds = Int[]
     ogroups = Dictionary{Int,UnitRange{Int}}()
     for subcell in cellsdict(latslice)

src/supercell.jl

@@ -88,6 +88,7 @@ function supercell_sitelist!!(sitelist, masklist, smatperp, seedperp, lat, appli
     masklist´ = copy(masklist)
     empty!(masklist)
     empty!(sitelist)
+    isnull(applied_selector) && return nothing
     cs = cells(applied_selector)
     csbool = zeros(Bool, length(cs))
     iter = BoxIterator(seedperp)

src/types.jl

@@ -194,6 +194,7 @@ struct AppliedSiteSelector{T,E,L}
     region::FunctionWrapper{Bool,Tuple{SVector{E,T}}}
     sublats::Vector{Int}
     cells::Vector{SVector{L,Int}}
+    isnull::Bool    # if isnull, the selector selects nothing, regardless of other fields
 end
 
 struct HopSelector{F,S,D,R}
@@ -210,6 +211,7 @@ struct AppliedHopSelector{T,E,L}
     sublats::Vector{Pair{Int,Int}}
     dcells::Vector{SVector{L,Int}}
     range::Tuple{T,T}
+    isnull::Bool    # if isnull, the selector selects nothing, regardless of other fields
 end
 
 struct Neighbors
@@ -252,6 +254,9 @@ iswithinrange(dr, (rmin, rmax)::Tuple{Real,Real}) =  ifelse(sign(rmin)*rmin^2 <=
 isbelowrange(dr, s::AppliedHopSelector) = isbelowrange(dr, s.range)
 isbelowrange(dr, (rmin, rmax)::Tuple{Real,Real}) =  ifelse(dr'dr < rmin^2, true, false)
 
+isnull(s::AppliedSiteSelector) = s.isnull
+isnull(s::AppliedHopSelector) = s.isnull
+
 Base.adjoint(s::HopSelector) = HopSelector(s.region, s.sublats, s.dcells, s.range, !s.adjoint)
 
 Base.NamedTuple(s::SiteSelector) =
@@ -1295,12 +1300,14 @@ Base.size(h::Harmonic, i...) = size(matrix(h), i...)
 
 Base.isless(h::Harmonic, h´::Harmonic) = sum(abs2, dcell(h)) < sum(abs2, dcell(h´))
 
-Base.zero(h::Harmonic{<:Any,<:Any,B}) where B = Harmonic(zero(dcell(h)), zero(matrix(h)))
+Base.zero(h::Harmonic{<:Any,<:Any,B}) where {B} = Harmonic(zero(dcell(h)), zero(matrix(h)))
 
 Base.copy(h::Harmonic) = Harmonic(dcell(h), copy(matrix(h)))
 
 Base.:(==)(h::Harmonic, h´::Harmonic) = h.dn == h´.dn && unflat(h.h) == unflat(h´.h)
 
+Base.iszero(h::Harmonic) = iszero(flat(h))
+
 #endregion
 #endregion
 
@@ -1427,6 +1434,7 @@ end
 
 Base.size(h::Hamiltonian, i...) = size(bloch(h), i...)
 Base.axes(h::Hamiltonian, i...) = axes(bloch(h), i...)
+Base.iszero(h::Hamiltonian) = all(iszero, harmonics(h))
 
 Base.copy(h::Hamiltonian) = Hamiltonian(
     copy(lattice(h)), copy(blockstructure(h)), copy.(harmonics(h)), copy(bloch(h)))

test/test_hamiltonian.jl

@@ -341,7 +341,11 @@ end
     # non-spatial models
     h = LP.linear() |> @hopping((i,j) --> ind(i) + ind(j)) + @onsite((i; k = 1) --> pos(i)[k])
     @test ishermitian(h())
-
+    # null selectors
+    h0 = LP.square() |> onsite(0) + hopping(0) |> supercell(3) |> @onsite!((t, r) -> 1; sublats = Symbol[])
+    @test iszero(h0())
+    h0 = LP.square() |> hopping(0) |> supercell(3) |> @hopping!((t, r, dr) -> 1; dcells = SVector{2,Int}[])
+    @test iszero(h0())
 end
 
 

test/test_lattice.jl

@@ -183,4 +183,12 @@ end
     ls = lat[cellsites(SA[1,0], 2)]
     @test ls isa Quantica.SiteSlice
     @test nsites(ls) == 1
+    # test the difference between a null selector and an unspecified one
+    ls = lat[cells = SVector{2,Int}[]]
+    @test isempty(ls)
+    ls = lat[cells = missing]
+    @test !isempty(ls)
+    ls = lat[region = RP.circle(4)]
+    ls´ = ls[cells = n -> !isreal(n[1])]
+    @test isempty(ls´)
 end