Factor ZigZag X2X4.py

From Werner KRAUTH

(Difference between revisions)
Jump to: navigation, search
Revision as of 07:03, 11 June 2024
Werner (Talk | contribs)

← Previous diff
Revision as of 07:04, 11 June 2024
Werner (Talk | contribs)

Next diff →
Line 28: Line 28:
time_ev = new_time_ev time_ev = new_time_ev
sigma *= -1 sigma *= -1
-  
-==Further information== 
-==References== 
plt.title('Factor-Zig-Zag algorithm, anharmonic oscillator') plt.title('Factor-Zig-Zag algorithm, anharmonic oscillator')
plt.xlabel('$x$') plt.xlabel('$x$')
Line 44: Line 41:
plt.legend(loc='upper right') plt.legend(loc='upper right')
plt.show() plt.show()
 +==Further information==
 +==References==

Revision as of 07:04, 11 June 2024

Contents

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 random
import matplotlib.pyplot as plt
def u(x):
    return x ** 2 / 2.0 + x ** 4 / 4.0

x = 0.0
time_ev = 0.0
sigma = 1 if x <= 0.0 else -1

data = []
n_samples = 10 ** 6
while len(data) < n_samples:
    delta_u2 = -math.log(random.uniform(0.0, 1.0))
    delta_u4 = -math.log(random.uniform(0.0, 1.0))
    new_x2 = math.sqrt(2.0 * delta_u2)
    new_x4 = (4.0 * delta_u4) ** 0.25
    new_x = sigma * min(abs(new_x2), abs(new_x4))
    new_time_ev = time_ev + abs(new_x - x)
    for t in range(math.ceil(time_ev), math.floor(new_time_ev) + 1):
        data.append(x + sigma * (t - time_ev))
    x = new_x
    time_ev = new_time_ev
    sigma *= -1
plt.title('Factor-Zig-Zag algorithm, anharmonic oscillator')
plt.xlabel('$x$')
plt.ylabel('$\pi(x)$')
plt.hist(data, bins=100, density=True, label='data')
XValues = []
YValues = []
for i in range(-1000,1000):
    x = i / 400.0
    XValues.append(x)
    YValues.append(math.exp(- u(x)) / 1.93525)
plt.plot(XValues, YValues, label='theory')
plt.legend(loc='upper right')
plt.show()

Further information

References

Personal tools