Merge branch 'develop' into documentation

This commit is contained in:
fjosw 2021-11-18 10:52:22 +00:00
commit 95fa13aaac
2 changed files with 24 additions and 17 deletions

View file

@ -174,20 +174,19 @@ def matmul(*operands):
return derived_array(multi_dot, operands) return derived_array(multi_dot, operands)
def jack_matmul(a, b): def jack_matmul(*operands):
"""Matrix multiply both operands making use of the jackknife approximation. """Matrix multiply both operands making use of the jackknife approximation.
Parameters Parameters
---------- ----------
a : numpy.ndarray operands : numpy.ndarray
First matrix, can be real or complex Obs valued Arbitrary number of 2d-numpy arrays which can be real or complex
b : numpy.ndarray Obs valued.
Second matrix, can be real or complex Obs valued
For large matrices this is considerably faster compared to matmul. For large matrices this is considerably faster compared to matmul.
""" """
if any(isinstance(o[0, 0], CObs) for o in [a, b]): if any(isinstance(o[0, 0], CObs) for o in operands):
def _exp_to_jack(matrix): def _exp_to_jack(matrix):
base_matrix = np.empty_like(matrix) base_matrix = np.empty_like(matrix)
for (n, m), entry in np.ndenumerate(matrix): for (n, m), entry in np.ndenumerate(matrix):
@ -201,10 +200,10 @@ def jack_matmul(a, b):
import_jackknife(entry.imag, name)) import_jackknife(entry.imag, name))
return base_matrix return base_matrix
j_a = _exp_to_jack(a) r = _exp_to_jack(operands[0])
j_b = _exp_to_jack(b) for op in operands[1:]:
r = j_a @ j_b r = r @ _exp_to_jack(op)
return _imp_from_jack(r, a.ravel()[0].real.names[0]) return _imp_from_jack(r, op.ravel()[0].real.names[0])
else: else:
def _exp_to_jack(matrix): def _exp_to_jack(matrix):
base_matrix = np.empty_like(matrix) base_matrix = np.empty_like(matrix)
@ -218,10 +217,10 @@ def jack_matmul(a, b):
base_matrix[n, m] = import_jackknife(entry, name) base_matrix[n, m] = import_jackknife(entry, name)
return base_matrix return base_matrix
j_a = _exp_to_jack(a) r = _exp_to_jack(operands[0])
j_b = _exp_to_jack(b) for op in operands[1:]:
r = j_a @ j_b r = r @ _exp_to_jack(op)
return _imp_from_jack(r, a.ravel()[0].names[0]) return _imp_from_jack(r, op.ravel()[0].names[0])
def inv(x): def inv(x):

View file

@ -418,9 +418,17 @@ class Obs:
""" """
return self.is_zero() or np.abs(self.value) <= sigma * self.dvalue return self.is_zero() or np.abs(self.value) <= sigma * self.dvalue
def is_zero(self): def is_zero(self, rtol=1.e-5, atol=1.e-8):
"""Checks whether the observable is zero within machine precision.""" """Checks whether the observable is zero within a given tolerance.
return np.isclose(0.0, self.value) and all(np.allclose(0.0, delta) for delta in self.deltas.values())
Parameters
----------
rtol : float
Relative tolerance (for details see numpy documentation).
atol : float
Absolute tolerance (for details see numpy documentation).
"""
return np.isclose(0.0, self.value, rtol, atol) and all(np.allclose(0.0, delta, rtol, atol) for delta in self.deltas.values())
def plot_tauint(self, save=None): def plot_tauint(self, save=None):
"""Plot integrated autocorrelation time for each ensemble. """Plot integrated autocorrelation time for each ensemble.