import random import math import matplotlib.pyplot as plt from tqdm import tqdm niters = [int(1e3),int(1e5),int(1e7)] # niter = 1000 for niter in tqdm(niters): step_size = 0.5 x = 0 naccept = 0 xs = [] for i in range(niter): backup_x = x action_init = 0.5*x**2 dx = random.random() dx = (dx-0.5)*step_size*2 x = x+dx action_fin = 0.5*x**2 metropolis = random.random() if (math.exp(action_init-action_fin) > metropolis): # accept naccept = naccept+1 else: x = backup_x # print(f"{x},{naccept}/{i}") xs.append(x) plt.figure() plt.hist(xs,100) plt.show()