pyerrors.misc

  0import pickle
  1import numpy as np
  2from .obs import Obs
  3
  4
  5def dump_object(obj, name, **kwargs):
  6    """Dump object into pickle file.
  7
  8    Parameters
  9    ----------
 10    obj : object
 11        object to be saved in the pickle file
 12    name : str
 13        name of the file
 14    path : str
 15        specifies a custom path for the file (default '.')
 16    """
 17    if 'path' in kwargs:
 18        file_name = kwargs.get('path') + '/' + name + '.p'
 19    else:
 20        file_name = name + '.p'
 21    with open(file_name, 'wb') as fb:
 22        pickle.dump(obj, fb)
 23
 24
 25def load_object(path):
 26    """Load object from pickle file.
 27
 28    Parameters
 29    ----------
 30    path : str
 31        path to the file
 32    """
 33    with open(path, 'rb') as file:
 34        return pickle.load(file)
 35
 36
 37def pseudo_Obs(value, dvalue, name, samples=1000):
 38    """Generate an Obs object with given value, dvalue and name for test purposes
 39
 40    Parameters
 41    ----------
 42    value : float
 43        central value of the Obs to be generated.
 44    dvalue : float
 45        error of the Obs to be generated.
 46    name : str
 47        name of the ensemble for which the Obs is to be generated.
 48    samples: int
 49        number of samples for the Obs (default 1000).
 50    """
 51    if dvalue <= 0.0:
 52        return Obs([np.zeros(samples) + value], [name])
 53    else:
 54        for _ in range(100):
 55            deltas = [np.random.normal(0.0, dvalue * np.sqrt(samples), samples)]
 56            deltas -= np.mean(deltas)
 57            deltas *= dvalue / np.sqrt((np.var(deltas) / samples)) / np.sqrt(1 + 3 / samples)
 58            deltas += value
 59            res = Obs(deltas, [name])
 60            res.gamma_method(S=2, tau_exp=0)
 61            if abs(res.dvalue - dvalue) < 1e-10 * dvalue:
 62                break
 63
 64        res._value = float(value)
 65
 66        return res
 67
 68
 69def gen_correlated_data(means, cov, name, tau=0.5, samples=1000):
 70    """ Generate observables with given covariance and autocorrelation times.
 71
 72    Parameters
 73    ----------
 74    means : list
 75        list containing the mean value of each observable.
 76    cov : numpy.ndarray
 77        covariance matrix for the data to be generated.
 78    name : str
 79        ensemble name for the data to be geneated.
 80    tau : float or list
 81        can either be a real number or a list with an entry for
 82        every dataset.
 83    samples : int
 84        number of samples to be generated for each observable.
 85    """
 86
 87    assert len(means) == cov.shape[-1]
 88    tau = np.asarray(tau)
 89    if np.min(tau) < 0.5:
 90        raise Exception('All integrated autocorrelations have to be >= 0.5.')
 91
 92    a = (2 * tau - 1) / (2 * tau + 1)
 93    rand = np.random.multivariate_normal(np.zeros_like(means), cov * samples, samples)
 94
 95    # Normalize samples such that sample variance matches input
 96    norm = np.array([np.var(o, ddof=1) / samples for o in rand.T])
 97    rand = rand @ np.diag(np.sqrt(np.diag(cov))) @ np.diag(1 / np.sqrt(norm))
 98
 99    data = [rand[0]]
100    for i in range(1, samples):
101        data.append(np.sqrt(1 - a ** 2) * rand[i] + a * data[-1])
102    corr_data = np.array(data) - np.mean(data, axis=0) + means
103    return [Obs([dat], [name]) for dat in corr_data.T]
104
105
106def _assert_equal_properties(ol, otype=Obs):
107    if not isinstance(ol[0], otype):
108        raise Exception("Wrong data type in list.")
109    for o in ol[1:]:
110        if not isinstance(o, otype):
111            raise Exception("Wrong data type in list.")
112        if not ol[0].is_merged == o.is_merged:
113            raise Exception("All Obs in list have to have the same state 'is_merged'.")
114        if not ol[0].reweighted == o.reweighted:
115            raise Exception("All Obs in list have to have the same property 'reweighted'.")
116        if not ol[0].e_content == o.e_content:
117            raise Exception("All Obs in list have to be defined on the same set of configs.")
118        if not ol[0].idl == o.idl:
119            raise Exception("All Obs in list have to be defined on the same set of configurations.")

 6def dump_object(obj, name, **kwargs):
 7    """Dump object into pickle file.
 8
 9    Parameters
10    ----------
11    obj : object
12        object to be saved in the pickle file
13    name : str
14        name of the file
15    path : str
16        specifies a custom path for the file (default '.')
17    """
18    if 'path' in kwargs:
19        file_name = kwargs.get('path') + '/' + name + '.p'
20    else:
21        file_name = name + '.p'
22    with open(file_name, 'wb') as fb:
23        pickle.dump(obj, fb)

26def load_object(path):
27    """Load object from pickle file.
28
29    Parameters
30    ----------
31    path : str
32        path to the file
33    """
34    with open(path, 'rb') as file:
35        return pickle.load(file)

38def pseudo_Obs(value, dvalue, name, samples=1000):
39    """Generate an Obs object with given value, dvalue and name for test purposes
40
41    Parameters
42    ----------
43    value : float
44        central value of the Obs to be generated.
45    dvalue : float
46        error of the Obs to be generated.
47    name : str
48        name of the ensemble for which the Obs is to be generated.
49    samples: int
50        number of samples for the Obs (default 1000).
51    """
52    if dvalue <= 0.0:
53        return Obs([np.zeros(samples) + value], [name])
54    else:
55        for _ in range(100):
56            deltas = [np.random.normal(0.0, dvalue * np.sqrt(samples), samples)]
57            deltas -= np.mean(deltas)
58            deltas *= dvalue / np.sqrt((np.var(deltas) / samples)) / np.sqrt(1 + 3 / samples)
59            deltas += value
60            res = Obs(deltas, [name])
61            res.gamma_method(S=2, tau_exp=0)
62            if abs(res.dvalue - dvalue) < 1e-10 * dvalue:
63                break
64
65        res._value = float(value)
66
67        return res

 70def gen_correlated_data(means, cov, name, tau=0.5, samples=1000):
 71    """ Generate observables with given covariance and autocorrelation times.
 72
 73    Parameters
 74    ----------
 75    means : list
 76        list containing the mean value of each observable.
 77    cov : numpy.ndarray
 78        covariance matrix for the data to be generated.
 79    name : str
 80        ensemble name for the data to be geneated.
 81    tau : float or list
 82        can either be a real number or a list with an entry for
 83        every dataset.
 84    samples : int
 85        number of samples to be generated for each observable.
 86    """
 87
 88    assert len(means) == cov.shape[-1]
 89    tau = np.asarray(tau)
 90    if np.min(tau) < 0.5:
 91        raise Exception('All integrated autocorrelations have to be >= 0.5.')
 92
 93    a = (2 * tau - 1) / (2 * tau + 1)
 94    rand = np.random.multivariate_normal(np.zeros_like(means), cov * samples, samples)
 95
 96    # Normalize samples such that sample variance matches input
 97    norm = np.array([np.var(o, ddof=1) / samples for o in rand.T])
 98    rand = rand @ np.diag(np.sqrt(np.diag(cov))) @ np.diag(1 / np.sqrt(norm))
 99
100    data = [rand[0]]
101    for i in range(1, samples):
102        data.append(np.sqrt(1 - a ** 2) * rand[i] + a * data[-1])
103    corr_data = np.array(data) - np.mean(data, axis=0) + means
104    return [Obs([dat], [name]) for dat in corr_data.T]