added hotrg

hotrg.py

+import torch 
+
+def gauge(M, Dcut, side):
+    if side=='l':
+        D2 = M.shape[0]
+        M = M.view(D2, -1)
+    elif side=='r':
+        D2 = M.shape[3]
+        M = M.permute(3,0,1,2).view(D2, -1)
+    M2 = M@M.t()
+
+    w, v = torch.symeig(M2, eigenvectors=True)
+    Dnew = min(D2, Dcut)
+    truncation = w[-Dnew:].sum() / w.sum()
+    v = v[:,-Dnew:].view(D2, Dnew)
+    return v, truncation
+
+def hotrg(T, D, Dcut, no_iter):
+    lnZ = 0.0
+    truncation_error = 0.0
+    for n in range(no_iter):
+        f = torch.norm(T)
+        lnZ += 2**(no_iter-n)*torch.log(f)
+        T = T / f
+
+        #contract up-down
+        Dh, Dv = T.shape[0], T.shape[1] # horizontal and vertical bond dimension of the tensor
+        M = torch.einsum('ludr,adcb->luracb', (T, T)).permute(0,3, 1, 4, 2,5).contiguous()
+        M = M.view(Dh**2, Dv, Dv, Dh**2)
+        vl, el = gauge(M, Dcut,'l')
+        vr, er = gauge(M, Dcut,'r')
+        U = vl if (el < er) else vr
+        truncation_error += min(el, er)
+        T = torch.einsum('li,ludr,rj->iudj', (U, M, U))
+       
+        #rotate tensor 
+        T = T.permute(1, 3, 0, 2).contiguous()
+
+    trace = 0.0 
+    for x in range(T.shape[0]):
+        for y in range(T.shape[1]):
+            trace += T[x, y, x, y]
+    lnZ += torch.log(trace)
+    return lnZ, truncation_error
+
+if __name__=='__main__':
+    K = torch.tensor([0.44])
+    Dcut = 10
+    n = 30
+
+    #Boltzmann factor on a bond M=LR^T
+    M = torch.stack([torch.cat([torch.exp(K),  torch.exp(-K)]), 
+                     torch.cat([torch.exp(-K), torch.exp(K)])
+                    ])
+    U, S, V = torch.svd(M)
+    L = U*torch.sqrt(S) 
+    R = V*torch.sqrt(S)
+
+    #            L 
+    #            |
+    # T =  R^{T}-o-L
+    #            |
+    #            R^{T}
+
+    T = torch.einsum('ai,aj,ak,al->ijkl', (L, L, R, R)) 
+    lnZ, _ = hotrg(T, 2, Dcut, n)
+    print (lnZ.item()/2**n)

trg.py

 import torch
 
-def levin_nave_trg(T, D, Dcut, no_iter, device='cpu'):
+def levin_nave_trg(T, D, Dcut, no_iter):
     
     lnZ = 0.0
     truncation_error = 0.0
@@ -37,4 +37,25 @@ def levin_nave_trg(T, D, Dcut, no_iter, device='cpu'):
 
     return lnZ, truncation_error
 
-
+if __name__=='__main__':
+    K = torch.tensor([0.44])
+    Dcut = 20
+    n = 20
+
+    #Boltzmann factor on a bond M=LR^T
+    M = torch.stack([torch.cat([torch.exp(K),  torch.exp(-K)]), 
+                     torch.cat([torch.exp(-K), torch.exp(K)])
+                    ])
+    U, S, V = torch.svd(M)
+    L = U*torch.sqrt(S) 
+    R = V*torch.sqrt(S)
+
+    #            L 
+    #            |
+    # T =  R^{T}-o-L
+    #            |
+    #            R^{T}
+
+    T = torch.einsum('ai,aj,ak,al->ijkl', (L, L, R, R)) 
+    lnZ, _ = levin_nave_trg(T, 2, Dcut, n)
+    print (lnZ.item()/2**n)