TVDTemperingLift.py

From Werner KRAUTH

(Difference between revisions)
Jump to: navigation, search
Revision as of 08:31, 10 September 2022
Werner (Talk | contribs)

← Previous diff
Current revision
Werner (Talk | contribs)

Line 1: Line 1:
 +This is the final program TVDTemperingLift.py to be discussed in my fourth lecture at the summerschool .... at the International Center for Theoretical Physics in Trieste, in September 2022.
 +
 +In previous lecture, we had started the discussion of the diffusion on the path graph P_n, in other words the one-dimensional linear lattice with n vertices and without periodic boundary conditions.
 +
 +
import random import random
import pylab import pylab

Current revision

This is the final program TVDTemperingLift.py to be discussed in my fourth lecture at the summerschool .... at the International Center for Theoretical Physics in Trieste, in September 2022.

In previous lecture, we had started the discussion of the diffusion on the path graph P_n, in other words the one-dimensional linear lattice with n vertices and without periodic boundary conditions.


import random
import pylab
import numpy as np
for n in [10, 20, 40, 80, 160, 320]:
    ReplicaChange = 0.1
    const = 4.0 / n ** 2
    PiStat = {}
    Table = []
   
    for x in range(1, n + 1):
        Table.append((x, -1, 0))
        Table.append((x,  1, 0))
        Table.append((x, -1, 1))
        Table.append((x,  1, 1))
#
#   factor of 1/4 because the total must be normalized, with two liftings and
#   two replicas  
#
        PiStat[(x, -1,  0)] = 1.0 / float(n) / 4.0
        PiStat[(x,  1,  0)] = 1.0 / float(n) / 4.0
        PiStat[(x, -1, 1)] = const * abs( (n + 1) / 2 - x) / 4.0
        PiStat[(x,  1, 1)] = const * abs( (n + 1) / 2 - x) / 4.0
    PiStat[(0, -1, 0)] = 0.0 
    PiStat[(0,  1, 0)] = 0.0 
    PiStat[(0, -1, 1)] = 0.0
    PiStat[(0,  1, 1)] = 0.0
    PiStat[(n + 1, -1,  0)] = 0.0
    PiStat[(n + 1,  1,  0)] = 0.0
    PiStat[(n + 1, -1, 1)] = 0.0
    PiStat[(n + 1,  1, 1)] = 0.0
    Position = (1,  1, 0)
    PTrans   = np.zeros((4 * n, 4 * n))
    Pi = np.zeros([4 * n])
    for x in range(1, n + 1):
        for Dir in [-1, 1]:
            for Rep in [0, 1]: 
                i = Table.index((x, Dir, Rep))
                Pi[i] = PiStat[(x, Dir, Rep)]
                k = Table.index((x, -Dir, Rep))
                if PiStat[(x + Dir, Dir, Rep)] > 0.0:
                    j = Table.index((x + Dir, Dir, Rep))
                    PTrans[i, j] = min(1.0, PiStat[(x + Dir, Dir,  Rep)] / PiStat[(x, Dir, Rep)])
                    PTrans[i, k] = 1.0 - PTrans[i, j]
                else:
                    PTrans[i, k] = 1.0  

    PReplica = np.zeros((4 * n, 4 * n)) 
    for x in range(1, n + 1):
        for Dir in [-1, 1]:
            i = Table.index((x, Dir, 0))
            j = Table.index((x, Dir, 1))
            PReplica[i, j] = ReplicaChange * min(1.0, PiStat[(x, Dir, 1)] / PiStat[(x, Dir, 0)])
            PReplica[i, i] = 1.0 - PReplica[i, j]
            PReplica[j, i] = ReplicaChange * min(1.0, PiStat[(x, Dir,  0)] / PiStat[(x, Dir, 1)])
            PReplica[j, j] = 1.0 - PReplica[j, i]
    PResampling = np.zeros((4 * n, 4 * n))
    for x in range(1, n + 1):
        for Dir in [-1, 1]: 
            i = Table.index((x, Dir, 0))
            j = Table.index((x, -Dir, 0))
            k = Table.index((x, Dir, 1))
            l = Table.index((x, -Dir, 1))
            PResampling[i, j] = 1.0 / (4.0 * n)
            PResampling[i, i] =  1 - 1. / (4.0 * n)
            PResampling[k, l] = 1.0 / (4.0 * n)
            PResampling[k, k] =  1 - 1. / (4.0 * n)

    P = PTrans @ PReplica @ PResampling
    Pit = np.zeros([4 * n])
    Pit[0] = 1.0
    xvalues = []
    yvalues = []
    iter = 0
    while True:
        iter += 1
        Pit = Pit @ P
        TVD = sum(np.absolute(Pi - Pit) / 2.0)
        xvalues.append(iter / float(n))
        yvalues.append(TVD)
        if TVD < 0.1: break
    pylab.plot(xvalues,yvalues, label='$n =$ '+str(n))
pylab.legend(loc='upper right')
pylab.xlabel("$t/ n$ (rescaled time) ")
pylab.ylabel("TVD")
pylab.title("TVD lift tempering on the path graph of $n$ sites")
pylab.show()
Personal tools