From b02ac863609eb5892723ec1ff579a2c4a3dfd898 Mon Sep 17 00:00:00 2001
From: Pablo San-Jose <lekand@gmail.com>
Date: Thu, 21 Nov 2024 22:19:26 +0100
Subject: [PATCH] don't drop structural zeros in Fock potential

---
 src/meanfield.jl           | 7 ++++---
 test/test_greenfunction.jl | 1 -
 2 files changed, 4 insertions(+), 4 deletions(-)

diff --git a/src/meanfield.jl b/src/meanfield.jl
index 3cab3f15..bdff0887 100644
--- a/src/meanfield.jl
+++ b/src/meanfield.jl
@@ -48,7 +48,7 @@ function meanfield(g::GreenFunction{T,E}, args...;
     hFock = lat |> hopping((r, dr) -> iszero(dr) ? Uf : T(Vf(dr)); selector..., includeonsite = true)
     hHartree = (Uf == U && Vh === Vf) ? hFock :
         lat |> hopping((r, dr) -> iszero(dr) ? U : T(Vh(dr)); selector..., includeonsite = true)
-
+    # this drops zeros
     potHartree = T.(sum(unflat, harmonics(hHartree)))
 
     oaxes = orbaxes(call!_output(gs))
@@ -61,7 +61,8 @@ function meanfield(g::GreenFunction{T,E}, args...;
 
     # this is important for the fast orbrange-based implementation of MeanField evaluation
     check_cell_order(hFock_slice, rho)
-    potFock = T.(parent(hFock_slice))
+    # this does not drop zeros. Important to keep diagonal zeros
+    potFock = convert(SparseMatrixCSC{T,Int}, parent(hFock_slice))
 
     encoder, decoder = nambu ? NambuEncoderDecoder(is_nambu_rotated´) : (identity, identity)
     S = typeof(encoder(zero(Q)))
@@ -208,7 +209,7 @@ function maybe_nambufy_traces!(traces, m::MeanField{B}) where {B}
 end
 
 diag_real_tr_rho_Q(ρ, Q) =
-    [real(unsafe_trace_prod(Q, view(ρ, rng, rng))) for rng in siteindexdict(orbaxes(ρ, 2))]
+    [real(unsafe_trace_prod(Q, view(parent(ρ), rng, rng))) for rng in siteindexdict(orbaxes(ρ, 2))]
 
 hole_id(::Type{<:SMatrix{2,2}}) = SA[0 0; 0 1]
 hole_id(::Type{<:SMatrix{4,4}}) = SA[0 0 0 0; 0 0 0 0; 0 0 1 0; 0 0 0 1]
diff --git a/test/test_greenfunction.jl b/test/test_greenfunction.jl
index 2e61df47..a6c4eb6f 100644
--- a/test/test_greenfunction.jl
+++ b/test/test_greenfunction.jl
@@ -601,7 +601,6 @@ end
     @test_throws ArgumentError meanfield(g)
     g = LP.honeycomb() |> hamiltonian(hopping(I), orbitals = (2,1)) |> greenfunction
     @test_throws ArgumentError meanfield(g)
-
 end
 
 @testset begin "OrbitalSliceArray serialization"