Merge pull request #315 from pablosanjose/gsliceoutput

GreenFunction indexing refactor

docs/src/advanced/meanfield.md

@@ -41,7 +41,7 @@ julia> h(SA[])
 julia> g = greenfunction(h, GS.Spectrum());
 
 julia> ρ0 = densitymatrix(g[])(0.5, 0) ## density matrix at chemical potential `µ=0.5` and temperature `kBT = 0`  on all sites
-4×4 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+4×4 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
  0.138197+0.0im  0.223607+0.0im  0.223607+0.0im  0.138197+0.0im
  0.223607+0.0im  0.361803+0.0im  0.361803+0.0im  0.223607+0.0im
  0.223607+0.0im  0.361803+0.0im  0.361803+0.0im  0.223607+0.0im

docs/src/tutorial/greenfunctions.md

@@ -213,15 +213,15 @@ julia> g(0.2)
 GreenSolution{Float64,2,0}: Green function at arbitrary positions, but at a fixed energy
 
 julia> g(0.2)[1, 3]
-5×5 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+5×5 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
  -2.56906+0.000123273im  -4.28767+0.00020578im   -4.88512+0.000234514im  -4.28534+0.00020578im    -2.5664+0.000123273im
  -4.28767+0.00020578im   -7.15613+0.00034351im   -8.15346+0.000391475im  -7.15257+0.00034351im    -4.2836+0.000205781im
  -4.88512+0.000234514im  -8.15346+0.000391475im  -9.29002+0.000446138im  -8.14982+0.000391476im  -4.88095+0.000234514im
  -4.28534+0.00020578im   -7.15257+0.00034351im   -8.14982+0.000391476im  -7.14974+0.000343511im  -4.28211+0.000205781im
   -2.5664+0.000123273im   -4.2836+0.000205781im  -4.88095+0.000234514im  -4.28211+0.000205781im  -2.56469+0.000123273im
 
 julia> g(0.2)[siteselector(region = RP.circle(1, (0.5, 0))), 3]
-2×5 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+2×5 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
  0.0749214+3.15744e-8im   0.124325+5.27948e-8im   0.141366+6.01987e-8im   0.124325+5.27948e-8im  0.0749214+3.15744e-8im
  -0.374862+2.15287e-5im  -0.625946+3.5938e-5im   -0.712983+4.09561e-5im  -0.624747+3.59379e-5im   -0.37348+2.15285e-5im
 ```

docs/src/tutorial/observables.md

@@ -42,7 +42,7 @@ DensityMatrix: density matrix on specified sites using solver of type DensityMat
 
 julia> @time ρ(4)
   6.150548 seconds (57.84 k allocations: 5.670 GiB, 1.12% gc time)
-5×5 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+5×5 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
           0.5+0.0im          -7.34893e-10-3.94035e-15im  0.204478+1.9366e-14im   -7.34889e-10-1.44892e-15im  -5.70089e-10+5.48867e-15im
  -7.34893e-10+3.94035e-15im           0.5+0.0im          0.200693-2.6646e-14im   -5.70089e-10-1.95251e-15im  -7.34891e-10-2.13804e-15im
      0.204478-1.9366e-14im       0.200693+2.6646e-14im        0.5+0.0im              0.200693+3.55692e-14im      0.204779-4.27255e-14im
@@ -62,7 +62,7 @@ DensityMatrix: density matrix on specified sites with solver of type DensityMatr
 
 julia> @time ρ(4)
   0.001659 seconds (9 allocations: 430.906 KiB)
-5×5 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+5×5 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
           0.5+0.0im  -2.21437e-15+0.0im  0.204478+0.0im   2.67668e-15+0.0im   3.49438e-16+0.0im
  -2.21437e-15+0.0im           0.5+0.0im  0.200693+0.0im  -1.40057e-15+0.0im  -2.92995e-15+0.0im
      0.204478+0.0im      0.200693+0.0im       0.5+0.0im      0.200693+0.0im      0.204779+0.0im
@@ -81,7 +81,7 @@ DensityMatrix: density matrix on specified sites with solver of type DensityMatr
 
 julia> @time ρ(4)
   0.006580 seconds (3 allocations: 1.156 KiB)
-5×5 OrbitalSliceMatrix{Matrix{ComplexF64}}:
+5×5 OrbitalSliceMatrix{ComplexF64,Matrix{ComplexF64}}:
          0.5+0.0im  2.15097e-17+0.0im   0.20456+0.0im  2.15097e-17+0.0im   3.9251e-17+0.0im
  2.15097e-17+0.0im          0.5+0.0im  0.200631+0.0im  1.05873e-16+0.0im  1.70531e-18+0.0im
      0.20456+0.0im     0.200631+0.0im       0.5+0.0im     0.200631+0.0im      0.20482+0.0im

ext/QuanticaMakieExt/plotlattice.jl

@@ -150,7 +150,7 @@ function _hoppingprimitives(ls::LatticeSlice{<:Any,E}, ls´, h, opts, siteradii)
     return hp, hp´
 end
 
-maybe_evaluate_observable(o::Quantica.IndexableObservable, ls) = o[ls]
+maybe_evaluate_observable(o::Quantica.IndexableObservable, ls) = parent(o[ls])
 maybe_evaluate_observable(x, ls) = x
 
 maybe_getindex(v::AbstractVector{<:Number}, i) = v[i]
@@ -189,9 +189,9 @@ push_siteopacity!(sp, siteopacity::Real, shellopacity, i, r, is_shell) = push!(s
 push_siteopacity!(sp, siteopacity::Function, shellopacity, i, r, is_shell) = push!(sp.opacities, is_shell ? siteopacity(i, r) : shellopacity)
 push_siteopacity!(sp, siteopacity, shellopacity, i, r, is_shell) = argerror("Unrecognized siteradius")
 
-push_siteradius!(sp, siteradius::Real, i, r) = push!(sp.radii, siteradius)
+push_siteradius!(sp, siteradius::Number, i, r) = push!(sp.radii, siteradius)
 push_siteradius!(sp, siteradius::Function, i, r) = push!(sp.radii, siteradius(i, r))
-push_siteradius!(sp, siteradius, i, r) = argerror("Unrecognized siteradius")
+push_siteradius!(sp, siteradius, i, r) = argerror("Unrecognized siteradius of type $(typeof(siteradius))")
 
 push_sitetooltip!(sp, i, r, mat) = push!(sp.tooltips, matrixstring(i, mat))
 push_sitetooltip!(sp, i, r) = push!(sp.tooltips, positionstring(i, r))
@@ -230,9 +230,9 @@ push_hopopacity!(hp, hopopacity::Real, shellopacity, ij, rdr, is_shell) = push!(
 push_hopopacity!(hp, hopopacity::Function, shellopacity, ij, rdr, is_shell) = push!(hp.opacities, hopopacity(ij, rdr))
 push_hopopacity!(hp, hopopacity, shellopacity, ij, rdr, is_shell) = argerror("Unrecognized hopradius")
 
-push_hopradius!(hp, hopradius::Real, ij, rdr) = push!(hp.radii, hopradius)
+push_hopradius!(hp, hopradius::Number, ij, rdr) = push!(hp.radii, hopradius)
 push_hopradius!(hp, hopradius::Function, ij, rdr) = push!(hp.radii, hopradius(ij, rdr))
-push_hopradius!(hp, hopradius, ij, rdr) = argerror("Unrecognized hopradius")
+push_hopradius!(hp, hopradius, ij, rdr) = argerror("Unrecognized hopradius of type $(typeof(hopradius))")
 
 push_hoptooltip!(hp, (i, j), mat) = push!(hp.tooltips, matrixstring(i, j, mat))
 

src/Quantica.jl

@@ -25,14 +25,15 @@ using DelimitedFiles
 export sublat, bravais_matrix, lattice, sites, supercell, hamiltonian,
        hopping, onsite, @onsite, @hopping, @onsite!, @hopping!, pos, ind, cell,
        position,
-       plusadjoint, neighbors, siteselector, hopselector, diagonal,
+       plusadjoint, neighbors, siteselector, hopselector, diagonal, sitepairs,
        unflat, torus, transform, translate, combine,
        spectrum, energies, states, bands, subdiv,
        greenfunction, selfenergy, attach,
        plotlattice, plotlattice!, plotbands, plotbands!, qplot, qplot!, qplotdefaults,
        conductance, josephson, ldos, current, transmission, densitymatrix,
-       OrbitalSliceArray, OrbitalSliceVector, OrbitalSliceMatrix, orbaxes,
-       serializer, serialize, serialize!, deserialize!, deserialize
+       OrbitalSliceArray, OrbitalSliceVector, OrbitalSliceMatrix, orbaxes, siteindexdict,
+       serializer, serialize, serialize!, deserialize!, deserialize,
+       meanfield
 
 export LatticePresets, LP, RegionPresets, RP, HamiltonianPresets, HP, ExternalPresets, EP
 export EigenSolvers, ES, GreenSolvers, GS
@@ -63,6 +64,7 @@ include("mesh.jl")
 include("bands.jl")
 include("greenfunction.jl")
 include("observables.jl")
+include("meanfield.jl")
 # Plumbing
 include("apply.jl")
 include("show.jl")

src/apply.jl

@@ -218,7 +218,7 @@ function push_pointers!(ptrs, h, har0, s::AppliedSiteSelector, shifts = missing,
     lat = lattice(h)
     umat = unflat(har0)
     rows = rowvals(umat)
-    norbs = norbitals(h)
+    norbs = norbitals(h, :)
     for scol in sublats(lat), col in siterange(lat, scol)
         isinblock(col, block) || continue
         for p in nzrange(umat, col)
@@ -243,7 +243,7 @@ function push_pointers!(ptrs, h, har, s::AppliedHopSelector, shifts = missing, b
     lat = lattice(h)
     dn0 = zerocell(lat)
     dn = dcell(har)
-    norbs = norbitals(h)
+    norbs = norbitals(h, :)
     umat = unflat(har)
     rows = rowvals(umat)
     for scol in sublats(lat), col in siterange(lat, scol), p in nzrange(umat, col)