I am trying to apply an exponential fit to my data to determine the point at which the value drops by 1/e. When plotted, the fit seems to favor smaller values and does not portray the true relationship.
import numpy as np
import matplotlib
matplotlib.use("TkAgg") # need to set the TkAgg backend explicitly otherwise it introduced a low-level error
from matplotlib import pyplot as plt
import scipy as sc
def autoCorrelation(sample, longTime, temp, plotTau = False ):
# compute empirical autocovariance with lag tau averaged over time longTime
sample.takeTimeStep(timesteps=1500) # 1500 timesteps to let sample reach equilibrium
M = np.zeros(longTime)
for tau in range(longTime):
M[tau] = sample.calcMagnetisation()
sample.takeTimeStep()
M_ave = np.average(M) #time - average
M = (M - M_ave)
autocorrelation = np.correlate(M, M, mode='full')
autocorrelation /= autocorrelation.max() # normalise such that max autocorrelation is 1
autocorrelationArray = autocorrelation[int(len(autocorrelation)/2):]
x = np.arange(0, len(autocorrelationArray), 1)
# apply exponential fit
def exponenial(x, a, b):
return a * np.exp(-b * x)
popt, pcov = curve_fit(exponenial, x, np.absolute(autocorrelationArray)) # array, 2d array
yy = exponenial(x, *popt)
plt.plot(x, np.absolute(autocorrelationArray), 'o', x, yy)
plt.title('Exponential Fit of Magnetisation Autocorrelation against Time for Temperature = ' + str(T) + ' J/k')
plt.xlabel('Time / Number of Iterations ')
plt.ylabel('Magnetisation Autocorrelation')
plt.show()
# prints tau_e value b from exponential a * np.exp(-b * x)
print('tau_e is ' + str(1/popt[1])) # units converted to time steps by taking reciprocal
if __name__ == '__main__':
#plot autocorrelation against time
longTime = 100
temp = [1, 2, 2.3, 2.6, 3, 4]
for T in temp:
magnet = Ising(30, T) # (N, temp)
autoCorrelation(magnet, longTime, temp)
Note: Ising is a class in another .py file containing the functions takeTimeStep and calcMagnetisation.
Expect greater values of tau_e
