From caa3fe0fa090b52bb003492cbcf71e8a6062b7cc Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Brieuc=20Le=20D=C3=A9?= <ledebrieuc@gmail.com>
Date: Tue, 30 Jul 2024 16:59:50 +0200
Subject: [PATCH] Add measure2siteoperator tree w/o range

---
 src/treeMeasure.jl | 67 ++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 67 insertions(+)

diff --git a/src/treeMeasure.jl b/src/treeMeasure.jl
index a6a5ea8..ec99880 100644
--- a/src/treeMeasure.jl
+++ b/src/treeMeasure.jl
@@ -118,6 +118,73 @@ function measure2siteoperator(net::TreeNetwork, O1, O2, sites::Tuple{Int,Int})
     end    
 end
 
+"""
+    measure2siteoperator(net::TreeNetwork, O1, O2, site1::Int64, site2::Int64})
+
+For a Tree, compute the local expectation value of two one-site operators O1 and O2 only for site1 and site2.
+
+"""
+function measure2siteoperator(net::TreeNetwork, O1, O2, site1::Int64,site2::Int64)
+
+    sites=(site1,site2)
+    herm_cis = ishermitian(O1*O2)
+    herm_trans = ishermitian(O1) && ishermitian(O2)
+    conjpair = O1 == O2'
+
+    if site1==site2
+        return measure1siteoperator(A, O1*O2, site1)
+    end
+
+    ρ1 = ones(ComplexF64, 1, 1)
+    ρ2 = ones(ComplexF64, 1, 1)
+    ρ = ones(ComplexF64, 1, 1)
+    firstsite = sites[1]
+    lastsite = sites[2]
+    N = length(path(net, firstsite, lastsite))
+    T = (herm_cis && herm_trans) ? Float64 : ComplexF64
+    expval = Array{T,2}(undef, 2,2)
+
+    if !in(findheadnode(net), path(net, firstsite, lastsite)[2:end-1])
+        if length(path(net, firstsite)) < length(path(net, lastsite))
+            nearestsite = firstsite
+            farthestsite = lastsite
+        else
+            nearestsite = lastsite
+            farthestsite = firstsite
+        end
+        for (A, id, dir) in Path(net, nearestsite)[1:end-1]
+            ρ = rhoAAstar(ρ, A, 1, dir)
+        end
+        for (i, (A, id, dir)) in enumerate(Path(net, nearestsite, farthestsite)[1:end-1])
+	    if id == nearestsite
+               v = rhoAOAstar(ρ, A, O1*O2, 1, nothing)
+               expval[1,1] = v
+	    end
+            ρ1 = rhoAOAstar(ρ, A, O1, 1, dir)
+            ρ2 = rhoAOAstar(ρ, A, O2, 1, dir)
+            for (j, (B, jd, djr)) in enumerate(Path(net, id, farthestsite)[2:end])
+		if id==nearestsite
+                   v = rhoAOAstar(ρ1, B, O2, 1, nothing)
+                   expval[1,2] = v
+                   v = rhoAOAstar(ρ2, B, O1, 1, nothing)
+                   expval[2,1] = v
+		end
+                djr != nothing && (ρ1 = rhoAAstar(ρ1, B, 1, djr))
+                djr != nothing && (ρ2 = rhoAAstar(ρ2, B, 1, djr))
+            end
+            ρ = rhoAAstar(ρ, A, 1, dir)
+        end
+        v = rhoAOAstar(ρ, net[farthestsite], O1*O2, 1, nothing)
+        expval[2,2] = v
+        if nearestsite == lastsite
+            expval = reverse(reverse(expval, dims=1),dims=2)
+        end
+        return expval
+    else
+        throw("missing method")
+    end
+end
+
 function measure2siteoperator_herm(net::TreeNetwork, O1, O2, sites::Tuple{Int,Int})
 
     herm_cis = ishermitian(O1*O2)