Top to random TVD.py

From Werner KRAUTH

Revision as of 13:28, 6 June 2024 Werner (Talk | contribs) ← Previous diff		Current revision Werner (Talk | contribs)
Line 1:		Line 1:
		+	==Context==
		+	This page is part of my [[BegRohu_Lectures_2024|2024 Beg Rohu Lectures]] on "The second Markov chain revolution" at the [https://www.ipht.fr/Meetings/BegRohu2024/index.html Summer School] "Concepts and Methods of Statistical Physics" (3 - 15 June 2024).
		+
		+	==Python program==
	import math		import math
	import copy		import copy

---

Current revision

Context

This page is part of my 2024 Beg Rohu Lectures on "The second Markov chain revolution" at the Summer School "Concepts and Methods of Statistical Physics" (3 - 15 June 2024).

Python program

import math
import copy
import matplotlib.pyplot as plt

def factorial(n):
    return 1 if n == 0 else (0 if n == 0 else factorial(n - 1) * n)

for N in [3, 4, 5, 6, 7, 8, 9, 10]:
    print(N, 'N')
    XValues = []
    YValues = []
#
# lifted TASEP: sample a random initial active particle
#
    Pit = {}
    Conf = [k for k in range(N)]
    Pit[tuple(Conf)] = 1.0
    NConfs = factorial(N)
    weight = 1.0 / NConfs
    DiameterFlag = True
    for step in range(40000):
#    print(Pit)
        TVD = (NConfs - len(Pit)) * weight / 2.0
        for k in Pit:
            TVD += abs(Pit[k] - weight) / 2.0
#        print('============', N, step, len(Pit) / NConfs, TVD, '==============')
        if len(Pit) == NConfs and DiameterFlag == True:
            DiameterFlag = False
            print(N, step, 'N, Diameter / N ** 2')
        XValues.append(step / N / math.log(N))
        YValues.append(TVD)
        if TVD < 0.02: break
        NewPit = {}
        for OldConf, OldPi in Pit.items():
            OldConf = list(OldConf)
            a = OldConf.pop(0)
            for k in range(N):
                Conf = OldConf[0:k] + [a] + OldConf[k:N-1]
                if tuple(Conf) in NewPit.keys():
                    NewPit[tuple(Conf)] += OldPi / N
                else:
                    NewPit[tuple(Conf)] = OldPi / N
        Pit = copy.deepcopy(NewPit)
    plt.plot(XValues, YValues, label = '$N= $' + str(N))
plt.title('TVD of Top-to-random shuffle')
plt.xlabel('$t /( N \log N)$ (rescaled time)')
plt.ylabel('TVD')
plt.legend(loc='upper center')
plt.show()