SU2xSU2 package
SU2xSU2.SU2_mat_routines module
- SU2xSU2.SU2_mat_routines.alpha_to_a(alpha)[source]
A lattice of SU(2) matrices may be defined through the exponential map with parameters alpha. The exponential map can be explicitly evaluated, resulting in a linear combination of the Pauli matrices and unity: U = exp(i alpha_i sigma_i) = a_0*1 + i*a_i sigma_i This function finds the coefficients
a
based onalpha
.- Parameters:
alpha ((L,L,3) array) – parameters when representing a SU(2) group element via the exponential map at every lattice site
- Returns:
a – parameters of matrices at each lattice site when explicitly evaluating the exponential map
- Return type:
(L,L,4) array
- SU2xSU2.SU2_mat_routines.det(a)[source]
The determinant of an SU(2) matrix is given by the squared length of the parameter vector.
- Parameters:
a ((L,L,4) array) – parameters of the SU(2) valued lattice
- Returns:
determinant – determinants of the SU(2) valued lattice
- Return type:
(L,L) array
- SU2xSU2.SU2_mat_routines.dot(a, b)[source]
Computes the elementwise matrix product between two lattices of SU(2) matrices with parameter vectors
a
andb
.- Parameters:
a ((L,L,4) array) – parameters of first SU(2) valued lattice
b ((L,L,4) array) – parameters of second SU(2) valued lattice
- Returns:
c – parameters of SU(2) valued lattice resulting from the elementwise matrix products
- Return type:
(L,L,4) array
- SU2xSU2.SU2_mat_routines.hc(a)[source]
Returns the parameter vector of the hermitian conjugate at each lattice site.
- Parameters:
a ((L,L,4) array) – parameters of the SU(2) valued lattice
- Returns:
a_hc – parameters of hermitian conjugate SU(2) valued lattice
- Return type:
(L,L,4) array
- SU2xSU2.SU2_mat_routines.make_mats(a)[source]
Constructs explicit matrices corresponding to parameter vector
a
.- Parameters:
a ((L,L,4) array) – parameters of the SU(2) valued lattice
- Returns:
mats – np.matrix instance at every site
- Return type:
(L,L) object array
- SU2xSU2.SU2_mat_routines.norm2(a)[source]
Returns squared norm of the parameter vector
a
.- Parameters:
a ((L,L,4) array) – parameters of the SU(2) valued lattice
- Returns:
norm_sq – containing the norm at each site
- Return type:
(L,L) array
- SU2xSU2.SU2_mat_routines.renorm(a)[source]
Renormalises matrix to have det = 1
- Parameters:
a ((L,L,4) array) – parameters of the SU(2) valued lattice
- Returns:
renormed – renormalised parameters of the SU(2) valued lattice
- Return type:
(L,L,4) array
- SU2xSU2.SU2_mat_routines.sum(a, b)[source]
Computes the elementwise sum of two SU(2) valued lattices A and B with parameters
a
andb
. Let C = A + B, i.e. c = a + b. Note that the sum of two SU(2) matrices is proportional to an SU(2) matrix with proportionality constantk
, meaning D = C/k = 1/k (A + B) is in SU(2). To only having to perform manipulations on SU(2) matrices, the parametersd
of the SU(2) valued lattice D and the constantk
are returned such that their product gives the parameter vectors of C, the sum of lattice A and B.- Parameters:
a ((L,L,4) array) – parameters of first SU(2) valued lattice
b ((L,L,4) array) – parameters of second SU(2) valued lattice
- Returns:
d ((L,L,4) array) – parameters of SU(2) valued lattice proportional to a+b
k ((L,L,1) array) – proportionality constant between d and a+b
SU2xSU2.SU2xSU2 module
- class SU2xSU2.SU2xSU2.SU2xSU2(L, a, ell, eps, beta, mass=0.1)[source]
Bases:
object
Each instance describes a realization of the SU(2)xSU(2) model on a square lattice whose dynamics can be simulated using the Fourier Accelerated Hybrid Monte Carlo algorithm.
- Ham(phi, pi)[source]
Computes the Hamiltonian for a lattice configuration
phi
and momentum configurationpi
. The kinetic term is chosen such that the momenta follow a standard Gaussian.- Parameters:
phi ((L,L,4) array) – parameters of the SU(2) valued field
pi ((L,L,3) array) – parameter values of the conjugate momenta at each lattice site
- Returns:
H – the Hamiltonian as the sum of the action and the kinetic term
- Return type:
float
- Ham_FA(phi, pi)[source]
Analogous to
SU2xSU2.SU2xSU2.Ham()
but computes the modified Hamiltonian used to accelerate the dynamics.- Parameters:
phi ((L,L,4) array) – parameters of the SU(2) valued field
pi ((L,L,3) array) – parameter values of the conjugate momenta at each lattice site
- Returns:
H – the modified Hamiltonian as the sum of the action and the modified kinetic term
- Return type:
float
See also
- __init__(L, a, ell, eps, beta, mass=0.1)[source]
- Parameters:
L (int) – Number of lattice sites along one dimension. Must be even for implementation of Fourier acceleration to work properly
a (float) – Lattice spacing
ell (int) – Number of steps to integrate Hamilton’s equation, each of size eps
eps (float) – Step size for integrating Hamilton’s equations
beta (float) – Model parameter defined in terms of the nearest neighbor coupling parameter g via beta = 1/(2g^2)
mass (float (optional)) – Mass parameter used in Fourier Acceleration. The physics is independent of the parameter value which only effects the simulation performance. The most efficient choice depends on beta.
- action(phi)[source]
Computes the action for lattice configuration
phi
.- Parameters:
phi ((L,L,4) array) – parameters of the SU(2) valued field at each lattice site
- Returns:
S – the action
- Return type:
float
- exp_update(pi_dt)[source]
Computes the update matrix of the field for the transition t -> t + dt. The update matrix is given by the exponential of the conjugate momenta pi times i dt. As the momentum field is a linear combination of the SU(2) generators, the update matrix is itself an SU(2) element, described by the parameter vector of pi scaled by dt.
- Parameters:
pi_dt ((L,L,4) array) – parameters of the momenta times the integration step size
- Returns:
update – parameters of the update matrices
- Return type:
(L,L,4) array
- internal_energy_density(phi)[source]
Computes the internal energy (action) per site of a lattice configuration.
- Parameters:
phi ((L,L,4) array) – lattice configuration to perform computation on
- Returns:
e – action per site
- Return type:
float
- kernel_inv_F()[source]
Finds inverse of the action kernel computed in the Fourier space, here referred to as
A
.- Returns:
A – inverse action kernel in Fourier space
- Return type:
(L,L) array
- leapfrog(phi_old, pi_old)[source]
Performs the leapfrog integration of Hamilton’s equations for
self.ell
steps, each of sizeself.eps
. The passed arguments define the initial conditions.- Parameters:
phi_old ((L,L,4) array) – last accepted sample in the chain of the SU(2) valued field
pi_old ((L,L,3) array) – parameters of conjugate momenta associated with phi_old
- Returns:
phi_cur ((L,L,4) array) – SU(2) matrix parameters after simulating dynamics
pi_cur ((L,L,3) array) – momenta parameters after simulating dynamics
- leapfrog_FA(phi_old, pi_old)[source]
Leapfrog integrator analogous to
SU2xSU2.SU2xSU2.leapfrog()
but using the modified EoMs.See also
- make_NN_mask()[source]
Makes mask to apply to SU(2) valued field phi or momentum field pi which then gives the matrix parameter values of the nearest neighbors (NN) for each lattice site. Hence
phi[self.NN_mask]
is of shape (L,L,#neighbors,#parameters) i.e (L,L,4,4).- Returns:
NN_mask – tuple of two (L,L,4,1) arrays, each giving the row and column coordinate for all nearest neighbors
- Return type:
tuple
- pi_dot(phi)[source]
Computes the derivative of
pi
with respect to the Hamiltonian time. From Hamilton’s equations, this is given as i times the derivative of the action wrt.phi
.pi
and its time derivative are linear combinations of the Pauli matrices and hence described by 3 real parameters alpha.- Parameters:
phi ((L,L,4) array) – parameters of the SU(2) valued field
- Returns:
pi_t – parameters of the time derivative of the conjugate momenta
- Return type:
(L,L,3) array
- pi_samples()[source]
Returns real space sample of momenta according to the distribution based on the modified kinetic term in the modified Hamiltonian. The sampling is easiest done in Fourier space and in terms of a real and hermitian object
PI
from which the momentum samples can be reconstructed (both in Fourier space) The size of the lattice along one dimension L is assumed to be even.- Returns:
pi – parameters for the sample of the conjugate momenta in real space
- Return type:
(L,L,3) array
- prod_A_pi(pi_F)[source]
Computes the elementwise product of the inverse action kernel in Fourier space (denoted as ‘A’) and the momenta in Fourier space.
- Parameters:
pi_F ((L,L,3) array) – parameter vector of momenta in Fourier space
- Returns:
prod – parameter vector of momenta in Fourier space, each site being weighted by the inverse Fourier space kernel
- Return type:
(L,L,3)
- run_HMC(M, thin_freq, burnin_frac, accel=True, measurements=[], chain_paths=[], saving_bool=True, partial_save=5000, starting_config_path='', RNG_state_path='', chain_state_dir='data/chain_state/', renorm_freq=10000)[source]
Performs the Hybrid Monte Carlo simulation, by default with Fourier Acceleration.
A total of
M
trajectories will be simulated and measurements are taken everythin_freq
steps after the thermalisation period defined as the firstM*burnin_frac
samples. Due to accumulating rounding errors, the molecular dynamics cause field to leave the group manifold slowly. This is counteracted by enforcing the unitarity constraint of the matrix parameters everyrenorm_freq
-th step. The Monte Carlo chain is fully defined (and thus reproducible) by the model and simulation parameters as well as the initial configuration of the chain and the state of the random number generator (RNG). By using the last configuration of a previous chain and the associated RNG state, one can continue the chain seamlessly in a new simulation.- Parameters:
M (int) – number of HMC trajectories and thus total number of generated samples
thin_freq (int) – frequency at which measurements will be taken
burin_frac (float) – fraction of total HMC samples which are rejected as burn in
accel (bool (optional)) – By default True, indicating to use Fourier Acceleration
measurements (list of callables (optional)) – can select from: internal_energy_density, susceptibility, ww_correlation_func to measure the internal energy density, susceptibility, and wall-to-wall correlation respectively
chain_paths (list of str (optional)) – Only required if
saving_bool=True
, otherwise can be left empty. Listing the file paths relative to current working directory to store the measurements. The data will always be saved as a .npy file, allowing to omit the file extension.saving_bool (bool (optional)) – save measurement data
partial_save (int (optional)) – after how many steps preliminary measurements and chain state is saved to disk. Requires
saving_bool=True
.starting_config_path (str (optional)) – path to configuration to initialize the chain (.npy file). If not passed a disordered (i.e hot) start will be used.
RNG_state_path (str (optional)) – relative path to a .obj file containing the internal state of the random number generator from a previous run. When using the final configuration and the associated RNG state of a previous run as the starting configuration for a new one, the chain is seamlessly continued.
chain_state_dir (str (optional)) – path to directory in which the RNG state and the last configuration will be saved
renorm_freq (int (optional)) – after how many trajectories the SU(2) valued fields are projected back to the group manifold. Set to
None
to never renormalize
- susceptibility(phi)[source]
Computes the susceptibility (the average point to point correlation) for the configuration
phi
which can be obtained by summing up the correlation function.- Parameters:
phi ((L,L,4) array) – configuration to perform computation on
- Returns:
chi – susceptibility of phi
- Return type:
float
- ww_correlation_func(phi)[source]
Computes the non-normalized wall to wall correlation function for the lattice configuration
phi
. Observing that the correlation of two variables can be written as a convolution, one can apply the cross correlation theorem to efficiently compute the latter using FFTs.- Parameters:
phi ((L,L,4) array) – configuration to perform computation on
- Returns:
ww_cor – wall to wall correlation evaluated for wall separation in interval [0, L), not normalized to value at zero separation.
- Return type:
(L,) array
SU2xSU2.analysis module
- SU2xSU2.analysis.acceleration_mass_search(num_traj, burnin_frac, beta, L, corlength, masses, chain_dir='data/acceleration_mass/', data_path='data/accelertion_mass.txt', plot_path='plots/acceleration_mass.pdf')[source]
Performs a grid search of the mass acceleration parameter for the passed value pair of
beta
,L
. To assess the quality of the canonical choice of M=1/correlation length against other values of the parameter, a plot is produced comparing the simulation cost agaist values of the acceleration mass. The plot is stored atplot_path
and the associated data is stored atdata_path
. Computing the cost function requires measurming the susceptibility integrated autocorrelation time. The susceptibility chain is stored in the directorychain_dir
with the file being labeled by the value ofbeta
andL
.- Parameters:
num_traj (int) – number of trajectories in each simulation using a different acceleration mass
burnin_frac (float) – fraction of trajectories discarded as burn in
beta (float) – model parameter beta for which the simulations are performed
L (int) – lattice size along one direction. Must be even.
corlength (float) – correlation length infered for simulations at beta
masses (array) – values of the acceleration parameter to include in the grid search
chain_dir (str (optional)) – path to directory containing susceptibility measurement chain.
data_path (str (optional)) – path at which the collected data will be stored: masses, cost function and its error, simulation time, acceptance rate, susceptbility IAT and its error. Must be a text file.
plot_path (str (optional)) – path at which the plot of the cost fucntion vs the acceleration mass is stored
- SU2xSU2.analysis.critical_slowingdown(num_traj, burnin_frac, corlengthdata_path='data/corlength.txt', chi_chain_dir='data/slowdown/rawchains/', chi_data_dir='data/slowdown/', slowdownplot_path='plots/slowdown.pdf', costfuncplot_path='plots/costfunc.pdf', xticks=[2, 5, 10, 20, 50, 100], yticks=[1, 5, 10, 25])[source]
Computes the integrated autocorrelation time (IAT) of the susceptibility using beta, L form previous simulations where the correlation length was measured and stored. For example, the data produced in
SU2xSU2.analysis.mass_lambda()
and stored at its argumentcorlengthdata_path
can be used. Its is assumed that the data has a row wise structure of Ls, betas, correlation length. The computation is done for accelerated and unaccelerated dynamics, with the measurement chains being stored in directoriesaccel/
,unaccel/
atchi_chain_dir
. Files are labeled by the used value of beta and L. The ensemble average results are stored as text files atchi_data_dir
while the associated plot of IAT vs correlation length is stored atslowdownplot_path
. Power laws are fitted for either acceleration choice, allowing to quantify the degree of critical slowing down through the fitted value of the critical exponent ‘z’. A further plot, showing the simulation cost of either acceleration choice, is stored atcostfuncplot_path
.- Parameters:
num_traj (int) – number of trajectories in each simulation
burnin_frac (float) – fraction of trajectories discarded as burn in
corlengthdata_path (str (optional)) – path to the Ls, betas, correlation length data (must be .txt). The default causes to use the result from
SU2xSU2.analysis.mass_lambda()
when its argument with the same name is also left as default.chi_chain_dir (str (optional)) – path of directory where two directories (named accel and unaccel) will be created to store the measurement chains of the susceptibility.
chi_data_dir (str (optional)) – path of directory where two directories (named accel and unaccel) will be created to store L, beta, the IAT, its error, the ensemble averaged susceptibility, its error, the simulation time, and the acceptance rate (row wise).
slowdownplot_path (str (optional)) – path to plot showing the power law scaling of the susceptibility IAT with the correlation length
costfuncplot_path (str (optional)) – path to save the simulation cost function plot at
xticks (list (optional)) – list of ints or floats specifying the tick labels for the correlation length in the two plots
yticks (list (optional)) – list of ints or floats specifying the tick labels for the ratio unaccel cost function / accel cost function in the cost function plot
- SU2xSU2.analysis.get_avg_error(data, get_IAT=False)[source]
Performs the ensemble average based on a sequence of measurements of a single observable and estimates its error as the standard error on the mean (SEM), corrected by the square root of the observable’s integrated autocorrelation time. Optionally also returns integrated autocorrelation time and its error. Each element of data is refereed to as a data point and can either be a float or a 1D array of length N. The returned average and error will be of the same shape as the data point.
- Parameters:
data ((M,) or (M,N) array) – set of measurements made during the chain
get_IAT (bool (optional)) – False by default. Set to True to return IAT and its error
- Returns:
avg (float or (N, )) – ensemble average of data
error (float or (N, )) – autocorrelation corrected SEM
IAT (float (if
get_IAT==True
)) – integrated autocorrelation time (IAT)IAT_err (float (if
get_IAT==True
)) – error of IAT
- SU2xSU2.analysis.get_corlength(ww_cor, ww_cor_err, data_save_path)[source]
Infers the correlation length based on the ensemble averaged wall to wall correlation function and its error. The correlation function data is processed by normalizing and averaging it about its symmetry axis at L/2 which will be saved at
data_save_path
.- Parameters:
ww_cor ((L,) array) – ensemble average of correlation function on a lattice of length L
ww_cor_err ((L,) array) – associated error for all possible wall separations
data_save_path (str) – relative path to cwd at which the wall separations, the processed correlation function and its error will be stored row wise as an .npy file (extension not needed in path)
- Returns:
cor_length (float) – fitted correlation length in units of the lattice spacing
cor_length_err (float) – error in the fitted correlation length
reduced_chi2 (float) – chi-square per degree of freedom as a goodness of fit proxy
- SU2xSU2.analysis.internal_energy_coupling_exp(betas, Ls, num_traj, burnin_frac, accel=True, chaindata_pathbase='data/energy_density/', simdata_path='data/energy_density.txt', plot_path='plots/coupling_expansion.pdf')[source]
Computes and stores the internal energy per site for the passed value pairs of
betas
,Ls
. As a density is computed, the value of the lattice size is not crucial and finite size effects are often negligible. For each value pair,num_traj
trajectories are simulated, withburnin_frac
specifying the fraction of these rejected as burn in. By default the Fourier accelerated Hybrid Monte Carlo algorithm is used. For each simulation, the chain of internal energy measurements is stored atchaindata_pathbase
and files are labeled by the used value of beta and L. The ensemble average values of the internal energy densities are stored atsimdata_path
as a text file. A plot, stored atplot_path
, is produced comparing the simulation result with the weak and strong coupling expansions (w.c. and s.c.).- Parameters:
betas ((n,) array) – values of the model parameter beta for which simulations are performed in ascending order
Ls (int or (n,) array) – size of the lattice along one dimension for each simulation at different beta. When an integer is passed, the size will be assumed for all values of beta.
num_traj (int) – number of trajectories in each simulation
burnin_frac (float) – fraction of trajectories discarded as burn in
accel (bool (optional)) – using Fourier Acceleration by default
chaindata_pathbase (str) – path of the directory where the internal energy measurement chains for all value pairs of beta, L will be saved
simdata_path (str (optional)) – path of .txt file (relative to the current working directory) to store the ensemble averaged simulation results for the internal energy density and its error (with file extension)
plot_path (str) – path of final plot file (with file extension)
- SU2xSU2.analysis.mass_lambda(betas, Ls, num_traj, burnin_frac, accel=True, corlengthdata_path='data/corlength.txt', corfuncs_chain_dir='data/corfuncs/rawchains/', corfuncs_dir='data/corfuncs/', corfuncs_plot_dir='plots/corfuncs/', plot_path='plots/asymptotic_scaling.pdf')[source]
Computes the mass over lambda parameter ratio for for the passed value pairs of
betas
,Ls
(lattice size along one dimension). The lattice size must be chosen sufficiently large (meaning a multiple of the correlation length) to avoid finite size effects. In general, the required lattice size increases quickly with beta.Computing the mass over lambda ratio requires the correlation length for each beta,L value pair and will be stored in the text file at
corlengthdata_path
. For each of the considered beta,L value pairs, the raw measurement chain for the correlation function is saved in the directorycorfuncs_chain_dir
, while the the processed correlation function data (normalized as well as ensemble averaged and averaged across the symmetry axis at L/2) is stored atcorfuncs_dir
. Plots of the correlation functions, with the fitted analytical expectation are stored atcorfuncs_plot_dir
. In all cases the file names are given by considered values of beta,L.A plot of mass over lambda ratio is produced and stored at
plot_path
, allowing to assess the convergence of the simulation data to the continuum mass gap prediction as beta gets large.- Parameters:
betas ((n,) array) – values of the model parameter beta for which simulations are performed in ascending order
Ls ((n,) array) – size of the lattice along one dimension for each simulation at different beta. Must be even.
num_traj (int) – number of trajectories in each simulation
burnin_frac (float) – fraction of trajectories discarded as burn in
accel (bool (optional)) – using Fourier Acceleration by default
couplingdata_path (str (optional)) – path of .txt file (relative to the current working directory and with file extension) to store the correlation length, its error and the associated chi-squared value for all betas,Ls.
chaindata_pathbase (str (optional)) – path of the directory where the correlation function measurement chains for all value pairs of beta,L will be saved
corfuncs_dir (str (optional)) – path of the directory where the processed correlation function data for all value pairs of beta,L will be saved
corfuncs_plot_dir (str (optional)) – path of the directory where the plots of the processed correlation function for all value pairs of beta,L will be saved
plot_path (str (optional)) – path (with file extension) of the plot showing the numerical mass over lambda ratio against beta.
SU2xSU2.calibrate_paras module
- SU2xSU2.calibrate_paras.calibrate(model_paras, accel=False, sim_paras=None)[source]
For a model, specified by the dictionary
model_paras
, this function calibrates the number of leapfrog integration steps and their size under the constraint that the trajectory length is 1 and that the acceptance rate is with in desireable range between 60 and 75%.The calibration is done by performing short simulation runs (500 trajectories with 50% burn in unless overridden by passing
sim_paras
), extracting the acceptance rate and, if the acceptance rate is outside this range, adjusting the number of steps according to the difference to the ideal acceptance rate of 65% for the next run. The number of steps is inferred from constraining the trajectory length to unity. It is not guaranteed that the calibration is successful for all possible model specification given the fixed trajectory length. Hence the calibration is limited to 10 iterations. An indicator that longer trajectories are required is when the calibration algorithm tries to reduce the number of steps below one.Recommended to use the returned calibrated
model_paras
to define the model for the production run.- Parameters:
model_paras (dict) –
{L, a, ell, eps, beta}
denoting lattice size, lattice spacing, number of integration steps, integration step size and the SU(2)xSU(2) model parameter beta respectively The values of ell, eps are used as guesses to start the calibration and their product must be 1.accel (bool, optional) – use Fourier Acceleration or not
sim_paras (dict, optional) –
{M, thin_freq, burnin_frac, accel=True, measurements=[], chain_paths=[], saving_bool=True, partial_save=5000, starting_config=None, RGN_state=None, renorm_freq=10000}
Specifying the simulation parameters for the calibration run by callingSU2xSU2.SU2xSU2.run_HMC()
. Consider the associated docstring for definitions of these parameters.
- Returns:
model_paras – calibrated model parameters and of the same form as model_paras
- Return type:
dict
SU2xSU2.correlations module
- SU2xSU2.correlations.auto_window(IATs, c)[source]
Windowing procedure of Caracciolo and Sokal, Dynamic critical exponent of some Monte Carlo algorithms for the self-avoiding walk, J. Phys. A: Math. Gen. 19 (1986) 797 to truncate the sum for the IAT.
- Parameters:
IATs (array) – integrated autocorrelation time with increasingly late termination of the sum
c (float) – defines the window width. For correlations that decay exponentially, a value of 4 or 5 is conventional
- Returns:
idx – index for array IATs, representing the IAT estimate from the windowing procedure
- Return type:
int
- SU2xSU2.correlations.autocorr_func_1d(x)[source]
Computes the autocorrelation of a 1D array
x
using FFT and the cross correlation theorem. As FFTs yield circular convolutions and work most efficiently when the number of elements is a power of 2, pad the data with zeros to the next power of 2.- Parameters:
x ((N,) array) – 1 dimensional data series of N elements to find the autocorrelation function for
- Returns:
acf – the autocorrelation function
- Return type:
(N,) array
- SU2xSU2.correlations.autocorrelator(data, c=4.0)[source]
Alias for
SU2xSU2.correlations.autocorrelator_repeats()
when the data has been observed only once, meaning data is of shape (M,).See also
correlations.autocorrelator_repeats
- SU2xSU2.correlations.autocorrelator_repeats(data, c=4.0)[source]
Based on the implementation in emcee. Computes the autocorrelation function and integrated autocorrelation time (IAT) for passed
data
which is a 2D array such that each row represents a sample of observations. The correlations are computed between rows of the data by finding the 1D autocorrelation function for each column of the data. The overall ACF between rows is estimated as the average of correlations across columns. This also allows to get an estimate of the error of the ACF. An alternative is to average along rows first and to estimate the ACF as the autocorrelation of the final column (Goodman, Weare 2010).- Parameters:
data ((M,L) array) – each row contains L samples of the same observable while different rows correspond to different positions in the chain of M observations. The correlation is computed across axis 1 i.e. gives the AFC for samples from different positions in the chain.
c (float, optional) – parameter used in determining when to truncate the sum for the IAT
- Returns:
ts ((M,) array) – array of considered separations between two samples
ACF ((M,) array) – autocorrelation function
ACF_err ((M,) array) – error of the autocorrelation function
IAT (float) – integrated autocorrelation time, showing how many rows in the data lie between uncorrelated samples
IAT_err (float) – error of the autocorrelation time
delta_t (float) – time needed to compute the ACF and the IAT
- SU2xSU2.correlations.corr_func_1D(x, y)[source]
Computes the correlation between the equal length 1D arrays
x
,y
using the cross correlation theorem. FFTs yield circular convolutions, such thatx
andy
are assumed to have periodic boundary conditions. When the data is not circular, need to pad it with zeros as done inSU2xSU2.correlations.autocorr_func_1d()
.- Parameters:
x ((N,) array) – first data series to correlate
y ((N,) array) – second data series to correlate
- Returns:
cf – correlation function between the data in x and y and of the same length as either array
- Return type:
(N,) array
- SU2xSU2.correlations.correlator(xs, ys)[source]
Alias for
SU2xSU2.correlations.correlator_repeats()
when thexs
andys
have been observed only once, i.e. are both of shape (L,) This implies that the error of the correlation function cannot be estimated.See also
correlations.correlator_repeats
- SU2xSU2.correlations.correlator_repeats(xs, ys)[source]
Computes the correlation function (CF) using FFTs between two equally sized 1D arrays for which several measurements exists. Each column presents a new observation and correlations are computed along axis 0. The final CF is the average along axis 1.
- Parameters:
xs ((L,M) array) – M measurements of first data vector of length L
ys ((L,M) array) – M measurements of second data vector of length L
- Returns:
CF ((L,) array) – average correlation function based on the M measurements
CF_err ((L,) array) – uncertainty in the CF, estimated as the standard error on the mean corrected by the square root of the integrated autocorrelation time
SU2xSU2.plotting module
- SU2xSU2.plotting.correlation_func_plot(data_path, plot_path, fit_upper=None, show_plot=True, ybottom=None, xright=None)[source]
Produces a plot (with a logarithmic y axis) of the correlation function data (stored at
data_path
) and the fitted, analytical expectation. The plot is saved atplot_path
with the fitting range and other plot parameters can be adjusted through the arguments.Allows manual adjustment of fitting the correlation length to the processed correlation function data (normalized and mirror averaged). A plot with the new fitting is produced and the inferred correlation length, its error and the associated chi2 are printed. These can then be manually added to (for example) a data/corlen_data.txt file.
- Parameters:
data_path (str) – file path to the correlation function data that will be plotted. The file is assumed to be .npy and contain the rows separation [0,L/2], correlation function, correlation function error respectively
plot_path (str) – file path (including file extension) to store the plot at
fit_upper (int (optional)) – largest separation (in units of the lattice spacing) to be included in the fit. If left as
None
, only all non-zero values of the correlation function will be included in the fitshow_plot (bool (optional)) – shows produced plot before saving it
ybottom (float (optional)) – lower limit for y axis
xright (float (optional)) – upper limit for x axis
- SU2xSU2.plotting.effective_mass_plot(data_path, xright=None, ytop=None, ybottom=None, show_plot=True, plot_path='plots/effective_mass.pdf')[source]
Produces an effective mass plot using the processed correlation function data stored at
data_path
. The data is expected to contain the separation, the correlation function value and its error row wise. The effective mass will be computed based on the assumption that the correlation function follows the shape of a cosh as analytically expected due to periodic boundary conditions. As the effective mass becomes noisy for large separations, the plot range can be adjusted using the remaining keywords.- Parameters:
data_path (str) – path to .npy file containing the averaged (ensemble and across the symmetry axis at L/2) and normalized correlation function
xright (float) – upper limit of the separation shown in the plot
ytop (float) – upper limit for the effective mass
ybottom (float) – lower limit for the effective mass
- SU2xSU2.plotting.plot_chain(data_chain_path, ylabel, start_idx=0)[source]
Plots raw measurement chain against computer time to gauge when convergence has been achieved. Observables of particular interest are the susceptibility (as sufficient burn is needed for the critical slowing down analysis) and the correlation function at a fixed separation (ideally chosen close to the correlation length as this is a slowly converging quantity and thus gives a lower bound for burn in). If the stored chain has already had some burn in removed, the starting value of the computer time can be adjusted by
start_idx
.- Parameters:
data_chain_path (str) – path to measurement chain file (must be .npy)
ylabel (str) – label for the plotted measurement. When Latex commands are contained, pass as raw string i.e. r’$chi$’
start_idx (int (optional)) – amount of burn in already removed from the measurement chain