docstrings updated

This commit is contained in:
Fabian Joswig 2021-11-07 21:44:22 +00:00
parent a23a97aed1
commit effccb1cc8
8 changed files with 173 additions and 159 deletions

View file

@ -226,8 +226,8 @@ class Corr:
def roll(self, dt):
"""Periodically shift the correlator by dt timeslices
Attributes:
-----------
Parameters
----------
dt : int
number of timeslices
"""
@ -264,9 +264,6 @@ class Corr:
weight : Obs
Reweighting factor. An Observable that has to be defined on a superset of the
configurations in obs[i].idl for all i.
Keyword arguments
-----------------
all_configs : bool
if True, the reweighted observables are normalized by the average of
the reweighting factor on all configurations in weight.idl and not
@ -283,8 +280,8 @@ class Corr:
def T_symmetry(self, partner, parity=+1):
"""Return the time symmetry average of the correlator and its partner
Attributes:
-----------
Parameters
----------
partner : Corr
Time symmetry partner of the Corr
partity : int
@ -309,8 +306,8 @@ class Corr:
def deriv(self, symmetric=True):
"""Return the first derivative of the correlator with respect to x0.
Attributes:
-----------
Parameters
----------
symmetric : bool
decides whether symmertic of simple finite differences are used. Default: True
"""
@ -414,8 +411,8 @@ class Corr:
def fit(self, function, fitrange=None, silent=False, **kwargs):
"""Fits function to the data
Attributes:
-----------
Parameters
----------
function : obj
function to fit to the data. See fits.least_squares for details.
fitrange : list
@ -447,8 +444,8 @@ class Corr:
def plateau(self, plateau_range=None, method="fit"):
""" Extract a plateu value from a Corr object
Attributes:
-----------
Parameters
----------
plateau_range : list
list with two entries, indicating the first and the last timeslice
of the plateau region.
@ -578,8 +575,8 @@ class Corr:
def dump(self, filename):
"""Dumps the Corr into a pickel file
Attributes:
-----------
Parameters
----------
filename : str
Name of the file
"""

View file

@ -59,7 +59,7 @@ class Fit_result(Sequence):
def least_squares(x, y, func, priors=None, silent=False, **kwargs):
"""Performs a non-linear fit to y = func(x).
Arguments:
Parameters
----------
x : list
list of floats.
@ -87,22 +87,23 @@ def least_squares(x, y, func, priors=None, silent=False, **kwargs):
enough.
silent : bool, optional
If true all output to the console is omitted (default False).
Keyword arguments
-----------------
initial_guess -- can provide an initial guess for the input parameters. Relevant for
initial_guess : list
can provide an initial guess for the input parameters. Relevant for
non-linear fits with many parameters.
method -- can be used to choose an alternative method for the minimization of chisquare.
The possible methods are the ones which can be used for scipy.optimize.minimize and
migrad of iminuit. If no method is specified, Levenberg-Marquard is used.
Reliable alternatives are migrad, Powell and Nelder-Mead.
resplot -- If true, a plot which displays fit, data and residuals is generated (default False).
qqplot -- If true, a quantile-quantile plot of the fit result is generated (default False).
expected_chisquare -- If true prints the expected chisquare which is
corrected by effects caused by correlated input data.
This can take a while as the full correlation matrix
has to be calculated (default False).
method : str
can be used to choose an alternative method for the minimization of chisquare.
The possible methods are the ones which can be used for scipy.optimize.minimize and
migrad of iminuit. If no method is specified, Levenberg-Marquard is used.
Reliable alternatives are migrad, Powell and Nelder-Mead.
resplot : bool
If true, a plot which displays fit, data and residuals is generated (default False).
qqplot : bool
If true, a quantile-quantile plot of the fit result is generated (default False).
expected_chisquare : bool
If true prints the expected chisquare which is
corrected by effects caused by correlated input data.
This can take a while as the full correlation matrix
has to be calculated (default False).
"""
if priors is not None:
return _prior_fit(x, y, func, priors, silent=silent, **kwargs)
@ -254,6 +255,8 @@ def odr_fit(x, y, func, silent=False, **kwargs):
def total_least_squares(x, y, func, silent=False, **kwargs):
"""Performs a non-linear fit to y = func(x) and returns a list of Obs corresponding to the fit parameters.
Parameters
----------
x : list
list of Obs, or a tuple of lists of Obs
y : list
@ -276,15 +279,14 @@ def total_least_squares(x, y, func, silent=False, **kwargs):
silent : bool, optional
If true all output to the console is omitted (default False).
Based on the orthogonal distance regression module of scipy
Keyword arguments
-----------------
initial_guess -- can provide an initial guess for the input parameters. Relevant for non-linear
fits with many parameters.
expected_chisquare -- If true prints the expected chisquare which is
corrected by effects caused by correlated input data.
This can take a while as the full correlation matrix
has to be calculated (default False).
initial_guess : list
can provide an initial guess for the input parameters. Relevant for non-linear
fits with many parameters.
expected_chisquare : bool
If true prints the expected chisquare which is
corrected by effects caused by correlated input data.
This can take a while as the full correlation matrix
has to be calculated (default False).
"""
output = Fit_result()

View file

@ -13,15 +13,16 @@ from ..fits import fit_lin
def read_rwms(path, prefix, version='2.0', names=None, **kwargs):
"""Read rwms format from given folder structure. Returns a list of length nrw
Attributes
-----------------
version -- version of openQCD, default 2.0
Keyword arguments
-----------------
r_start -- list which contains the first config to be read for each replicum
r_stop -- list which contains the last config to be read for each replicum
postfix -- postfix of the file to read, e.g. '.ms1' for openQCD-files
Parameters
----------
version : str
version of openQCD, default 2.0
r_start : list
list which contains the first config to be read for each replicum
r_stop : list
list which contains the last config to be read for each replicum
postfix : str
postfix of the file to read, e.g. '.ms1' for openQCD-files
"""
known_oqcd_versions = ['1.4', '1.6', '2.0']
if not (version in known_oqcd_versions):
@ -174,19 +175,25 @@ def extract_t0(path, prefix, dtr_read, xmin, spatial_extent, fit_range=5, **kwar
Parameters
----------
path -- Path to .ms.dat files
prefix -- Ensemble prefix
dtr_read -- Determines how many trajectories should be skipped when reading the ms.dat files.
Corresponds to dtr_cnfg / dtr_ms in the openQCD input file.
xmin -- First timeslice where the boundary effects have sufficiently decayed.
spatial_extent -- spatial extent of the lattice, required for normalization.
fit_range -- Number of data points left and right of the zero crossing to be included in the linear fit. (Default: 5)
Keyword arguments
-----------------
r_start -- list which contains the first config to be read for each replicum.
r_stop -- list which contains the last config to be read for each replicum.
plaquette -- If true extract the plaquette estimate of t0 instead.
path : str
Path to .ms.dat files
prefix : str
Ensemble prefix
dtr_read : int
Determines how many trajectories should be skipped when reading the ms.dat files.
Corresponds to dtr_cnfg / dtr_ms in the openQCD input file.
xmin : int
First timeslice where the boundary effects have sufficiently decayed.
spatial_extent : int
spatial extent of the lattice, required for normalization.
fit_range : int
Number of data points left and right of the zero crossing to be included in the linear fit. (Default: 5)
r_start : list
list which contains the first config to be read for each replicum.
r_stop: list
list which contains the last config to be read for each replicum.
plaquette : bool
If true extract the plaquette estimate of t0 instead.
"""
ls = []

View file

@ -11,8 +11,8 @@ from ..obs import Obs
def read_sfcf(path, prefix, name, **kwargs):
"""Read sfcf C format from given folder structure.
Keyword arguments
-----------------
Parameters
----------
im -- if True, read imaginary instead of real part of the correlation function.
single -- if True, read a boundary-to-boundary correlation function with a single value
b2b -- if True, read a time-dependent boundary-to-boundary correlation function
@ -113,15 +113,12 @@ def read_sfcf(path, prefix, name, **kwargs):
def read_sfcf_c(path, prefix, name, quarks='.*', noffset=0, wf=0, wf2=0, **kwargs):
"""Read sfcf c format from given folder structure.
Arguments
-----------------
Parameters
----------
quarks -- Label of the quarks used in the sfcf input file
noffset -- Offset of the source (only relevant when wavefunctions are used)
wf -- ID of wave function
wf2 -- ID of the second wavefunction (only relevant for boundary-to-boundary correlation functions)
Keyword arguments
-----------------
im -- if True, read imaginary instead of real part of the correlation function.
b2b -- if True, read a time-dependent boundary-to-boundary correlation function
names -- Alternative labeling for replicas/ensembles. Has to have the appropriate length
@ -236,8 +233,8 @@ def read_sfcf_c(path, prefix, name, quarks='.*', noffset=0, wf=0, wf2=0, **kwarg
def read_qtop(path, prefix, **kwargs):
"""Read qtop format from given folder structure.
Keyword arguments
-----------------
Parameters
----------
target -- specifies the topological sector to be reweighted to (default 0)
full -- if true read the charge instead of the reweighting factor.
"""

View file

@ -20,9 +20,6 @@ def derived_array(func, data, **kwargs):
automatic differentiation to work, all numpy functions have to have
the autograd wrapper (use 'import autograd.numpy as anp').
data -- list of Obs, e.g. [obs1, obs2, obs3].
Keyword arguments
-----------------
man_grad -- manually supply a list or an array which contains the jacobian
of func. Use cautiously, supplying the wrong derivative will
not be intercepted.

View file

@ -8,8 +8,8 @@ from .obs import Obs
def gen_correlated_data(means, cov, name, tau=0.5, samples=1000):
""" Generate observables with given covariance and autocorrelation times.
Arguments
-----------------
Parameters
----------
means -- list containing the mean value of each observable.
cov -- covariance matrix for the data to be geneated.
name -- ensemble name for the data to be geneated.

View file

@ -73,9 +73,12 @@ def inv_propagator(prop):
def Zq(inv_prop, fermion='Wilson'):
""" Calculates the quark field renormalization constant Zq
Attributes:
inv_prop -- Inverted 12x12 quark propagator
fermion -- Fermion type for which the tree-level propagator is used
Parameters
----------
inv_prop : array
Inverted 12x12 quark propagator
fermion : str
Fermion type for which the tree-level propagator is used
in the calculation of Zq. Default Wilson.
"""
_check_geometry()

View file

@ -53,7 +53,7 @@ class Obs:
def __init__(self, samples, names, idl=None, means=None, **kwargs):
""" Initialize Obs object.
Attributes
Parameters
----------
samples : list
list of numpy arrays containing the Monte Carlo samples
@ -150,57 +150,11 @@ class Obs:
res[e_name].append(e_name)
return res
def expand_deltas(self, deltas, idx, shape):
"""Expand deltas defined on idx to a regular, contiguous range, where holes are filled by 0.
If idx is of type range, the deltas are not changed
Parameters
----------
deltas -- List of fluctuations
idx -- List or range of configs on which the deltas are defined.
shape -- Number of configs in idx.
"""
if type(idx) is range:
return deltas
else:
ret = np.zeros(idx[-1] - idx[0] + 1)
for i in range(shape):
ret[idx[i] - idx[0]] = deltas[i]
return ret
def calc_gamma(self, deltas, idx, shape, w_max, fft):
"""Calculate Gamma_{AA} from the deltas, which are defined on idx.
idx is assumed to be a contiguous range (possibly with a stepsize != 1)
Parameters
----------
deltas -- List of fluctuations
idx -- List or range of configs on which the deltas are defined.
shape -- Number of configs in idx.
w_max -- Upper bound for the summation window
fft -- boolean, which determines whether the fft algorithm is used for
the computation of the autocorrelation function
"""
gamma = np.zeros(w_max)
deltas = self.expand_deltas(deltas, idx, shape)
new_shape = len(deltas)
if fft:
max_gamma = min(new_shape, w_max)
# The padding for the fft has to be even
padding = new_shape + max_gamma + (new_shape + max_gamma) % 2
gamma[:max_gamma] += np.fft.irfft(np.abs(np.fft.rfft(deltas, padding)) ** 2)[:max_gamma]
else:
for n in range(w_max):
if new_shape - n >= 0:
gamma[n] += deltas[0:new_shape - n].dot(deltas[n:new_shape])
return gamma
def gamma_method(self, **kwargs):
"""Calculate the error and related properties of the Obs.
Keyword arguments
-----------------
Parameters
----------
S : float
specifies a custom value for the parameter S (default 2.0), can be
a float or an array of floats for different ensembles
@ -378,6 +332,52 @@ class Obs:
self.ddvalue = np.sqrt(self.ddvalue) / self.dvalue
return
def expand_deltas(self, deltas, idx, shape):
"""Expand deltas defined on idx to a regular, contiguous range, where holes are filled by 0.
If idx is of type range, the deltas are not changed
Parameters
----------
deltas -- List of fluctuations
idx -- List or range of configs on which the deltas are defined.
shape -- Number of configs in idx.
"""
if type(idx) is range:
return deltas
else:
ret = np.zeros(idx[-1] - idx[0] + 1)
for i in range(shape):
ret[idx[i] - idx[0]] = deltas[i]
return ret
def calc_gamma(self, deltas, idx, shape, w_max, fft):
"""Calculate Gamma_{AA} from the deltas, which are defined on idx.
idx is assumed to be a contiguous range (possibly with a stepsize != 1)
Parameters
----------
deltas -- List of fluctuations
idx -- List or range of configs on which the deltas are defined.
shape -- Number of configs in idx.
w_max -- Upper bound for the summation window
fft -- boolean, which determines whether the fft algorithm is used for
the computation of the autocorrelation function
"""
gamma = np.zeros(w_max)
deltas = self.expand_deltas(deltas, idx, shape)
new_shape = len(deltas)
if fft:
max_gamma = min(new_shape, w_max)
# The padding for the fft has to be even
padding = new_shape + max_gamma + (new_shape + max_gamma) % 2
gamma[:max_gamma] += np.fft.irfft(np.abs(np.fft.rfft(deltas, padding)) ** 2)[:max_gamma]
else:
for n in range(w_max):
if new_shape - n >= 0:
gamma[n] += deltas[0:new_shape - n].dot(deltas[n:new_shape])
return gamma
def print(self, level=1):
warnings.warn("Method 'print' renamed to 'details'", DeprecationWarning)
self.details(level > 1)
@ -539,9 +539,10 @@ class Obs:
def dump(self, name, **kwargs):
"""Dump the Obs to a pickle file 'name'.
Keyword arguments
-----------------
path -- specifies a custom path for the file (default '.')
Parameters
----------
path : str
specifies a custom path for the file (default '.')
"""
if 'path' in kwargs:
file_name = kwargs.get('path') + '/' + name + '.p'
@ -836,7 +837,8 @@ def merge_idx(idl):
Parameters
----------
idl -- List of lists or ranges.
idl : list
List of lists or ranges.
"""
# Use groupby to efficiently check whether all elements of idl are identical
@ -893,14 +895,15 @@ def filter_zeroes(names, deltas, idl, eps=Obs.filter_eps):
Parameters
----------
names -- List of names
deltas -- Dict lists of fluctuations
idx -- Dict of lists or ranges of configs on which the deltas are defined.
Has to be a subset of new_idx.
Optional parameters
----------
eps -- Prefactor that enters the filter criterion.
names : list
List of names
deltas : dict
Dict lists of fluctuations
idx : dict
Dict of lists or ranges of configs on which the deltas are defined.
Has to be a subset of new_idx.
eps : float
Prefactor that enters the filter criterion.
"""
new_names = []
new_deltas = {}
@ -931,9 +934,6 @@ def derived_observable(func, data, **kwargs):
the autograd wrapper (use 'import autograd.numpy as anp').
data : list
list of Obs, e.g. [obs1, obs2, obs3].
Keyword arguments
-----------------
num_grad : bool
if True, numerical derivatives are used instead of autograd
(default False). To control the numerical differentiation the
@ -1072,10 +1072,13 @@ def reduce_deltas(deltas, idx_old, idx_new):
Parameters
----------
deltas -- List of fluctuations
idx_old -- List or range of configs on which the deltas are defined
idx_new -- List of configs for which we want to extract the deltas.
Has to be a subset of idx_old.
deltas : list
List of fluctuations
idx_old : list
List or range of configs on which the deltas are defined
idx_new : list
List of configs for which we want to extract the deltas.
Has to be a subset of idx_old.
"""
if not len(deltas) == len(idx_old):
raise Exception('Lenght of deltas and idx_old have to be the same: %d != %d' % (len(deltas), len(idx_old)))
@ -1109,9 +1112,6 @@ def reweight(weight, obs, **kwargs):
configurations in obs[i].idl for all i.
obs : list
list of Obs, e.g. [obs1, obs2, obs3].
Keyword arguments
-----------------
all_configs : bool
if True, the reweighted observables are normalized by the average of
the reweighting factor on all configurations in weight.idl and not
@ -1146,8 +1146,8 @@ def reweight(weight, obs, **kwargs):
def correlate(obs_a, obs_b):
"""Correlate two observables.
Attributes:
-----------
Parameters
----------
obs_a : Obs
First observable
obs_b : Obs
@ -1193,10 +1193,11 @@ def covariance(obs1, obs2, correlation=False, **kwargs):
is constrained to the maximum value in order to make sure that covariance
matrices are positive semidefinite.
Keyword arguments
-----------------
correlation -- if true the correlation instead of the covariance is
returned (default False)
Parameters
----------
correlation : bool
if true the correlation instead of the covariance is
returned (default False)
"""
for name in sorted(set(obs1.names + obs2.names)):
@ -1450,9 +1451,14 @@ def pseudo_Obs(value, dvalue, name, samples=1000):
def dump_object(obj, name, **kwargs):
"""Dump object into pickle file.
Keyword arguments
-----------------
path -- specifies a custom path for the file (default '.')
Parameters
----------
obj : object
object to be saved in the pickle file
name : str
name of the file
path : str
specifies a custom path for the file (default '.')
"""
if 'path' in kwargs:
file_name = kwargs.get('path') + '/' + name + '.p'
@ -1471,6 +1477,11 @@ def load_object(path):
def merge_obs(list_of_obs):
"""Combine all observables in list_of_obs into one new observable
Parameters
----------
list_of_obs : list
list of the Obs object to be combined
It is not possible to combine obs which are based on the same replicum
"""
replist = [item for obs in list_of_obs for item in obs.names]