lattice/pyerrors
1
0
Fork 0
mirror of https://github.com/fjosw/pyerrors.git synced 2025-03-15 14:50:25 +01:00
Code Issues Projects Releases Packages Wiki Activity

Vectorized gamma_method (#207)

Browse source 
* feat: vectorized gamma_method function added.

* feat: vectorized gamma method generalized to also work on other objects
like Corr or Fit_result.

* feat: alias gamma_method for vectorized gamma_method added.

* docs: example 5 updated to include vectorized gamma_method.

* docs: output of example 5 updated.
This commit is contained in:
Fabian Joswig 2023-07-19 12:13:20 +01:00 • committed by GitHub
parent 6b7846486d
commit 1e438356fd
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
3 changed files with 76 additions and 46 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

73
examples/05_matrix_operations.ipynb
View file

@ -7,8 +7,7 @@
"outputs": [], "outputs": [],
"source": [ "source": [
"import pyerrors as pe\n", "import pyerrors as pe\n",
"import numpy as np\n", "import numpy as np"
"import scipy"
] ]
}, },
{ {
@ -27,8 +26,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[4.10(20)] Obs[-1.00(10)]]\n", "[[4.10(20) -1.00(10)]\n",
" [Obs[-1.00(10)] Obs[1.000(10)]]]\n" " [-1.00(10) 1.000(10)]]\n"
] ]
} }
], ],
@ -63,8 +62,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[17.81] Obs[-5.1]]\n", "[[17.81 -5.1]\n",
" [Obs[-5.1] Obs[2.0]]]\n" " [-5.1 2.0]]\n"
] ]
} }
], ],
@ -88,8 +87,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[4.1] Obs[-1.0]]\n", "[[4.1 -1.0]\n",
" [Obs[-1.0] Obs[1.0]]]\n" " [-1.0 1.0]]\n"
] ]
} }
], ],
@ -113,8 +112,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[78.12099999999998] Obs[-22.909999999999997]]\n", "[[78.12099999999998 -22.909999999999997]\n",
" [Obs[-22.909999999999997] Obs[7.1]]]\n" " [-22.909999999999997 7.1]]\n"
] ]
} }
], ],
@ -138,8 +137,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[30.161857460980094] Obs[-1.1752011936438014]]\n", "[[30.16185746098009 -1.1752011936438014]\n",
" [Obs[-1.1752011936438014] Obs[1.1752011936438014]]]\n" " [-1.1752011936438014 1.1752011936438014]]\n"
] ]
} }
], ],
@ -163,7 +162,7 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[Obs[2.00(40)] Obs[1.00(10)]]\n" "[2.00(40) 1.00(10)]\n"
] ]
} }
], ],
@ -171,8 +170,7 @@
"vec1 = pe.pseudo_Obs(2, 0.4, 'e1')\n", "vec1 = pe.pseudo_Obs(2, 0.4, 'e1')\n",
"vec2 = pe.pseudo_Obs(1, 0.1, 'e1')\n", "vec2 = pe.pseudo_Obs(1, 0.1, 'e1')\n",
"vector = np.asarray([vec1, vec2])\n", "vector = np.asarray([vec1, vec2])\n",
"for (i), entry in np.ndenumerate(vector):\n", "pe.gm(vector)\n",
" entry.gamma_method()\n",
"print(vector)" "print(vector)"
] ]
}, },
@ -192,14 +190,13 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[Obs[7.2(1.7)] Obs[-1.00(46)]]\n" "[7.2(1.7) -1.00(45)]\n"
] ]
} }
], ],
"source": [ "source": [
"product = matrix @ vector\n", "product = matrix @ vector\n",
"for (i), entry in np.ndenumerate(product):\n", "pe.gm(product)\n",
" entry.gamma_method()\n",
"print(product)" "print(product)"
] ]
}, },
@ -245,15 +242,14 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[2.025(49)] Obs[0.0]]\n", "[[2.025(49) 0.0]\n",
" [Obs[-0.494(50)] Obs[0.870(29)]]]\n" " [-0.494(50) 0.870(29)]]\n"
] ]
} }
], ],
"source": [ "source": [
"cholesky = pe.linalg.cholesky(matrix)\n", "cholesky = pe.linalg.cholesky(matrix)\n",
"for (i, j), entry in np.ndenumerate(cholesky):\n", "pe.gm(cholesky)\n",
" entry.gamma_method()\n",
"print(cholesky)" "print(cholesky)"
] ]
}, },
@ -273,8 +269,8 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[-8.881784197001252e-16] Obs[0.0]]\n", "[[-8.881784197001252e-16 0.0]\n",
" [Obs[0.0] Obs[0.0]]]\n" " [0.0 0.0]]\n"
] ]
} }
], ],
@ -299,18 +295,17 @@
"name": "stdout", "name": "stdout",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"[[Obs[0.494(12)] Obs[0.0]]\n", "[[0.494(12) 0.0]\n",
" [Obs[0.280(40)] Obs[1.150(39)]]]\n", " [0.280(40) 1.150(38)]]\n",
"Check:\n", "Check:\n",
"[[Obs[1.0] Obs[0.0]]\n", "[[1.0 0.0]\n",
" [Obs[0.0] Obs[1.0]]]\n" " [0.0 1.0]]\n"
] ]
} }
], ],
"source": [ "source": [
"inv = pe.linalg.inv(cholesky)\n", "inv = pe.linalg.inv(cholesky)\n",
"for (i, j), entry in np.ndenumerate(inv):\n", "pe.gm(inv)\n",
" entry.gamma_method()\n",
"print(inv)\n", "print(inv)\n",
"print('Check:')\n", "print('Check:')\n",
"check_inv = cholesky @ inv\n", "check_inv = cholesky @ inv\n",
@ -335,21 +330,19 @@
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"Eigenvalues:\n", "Eigenvalues:\n",
"[Obs[0.705(57)] Obs[4.39(19)]]\n", "[0.705(56) 4.39(20)]\n",
"Eigenvectors:\n", "Eigenvectors:\n",
"[[Obs[-0.283(26)] Obs[-0.9592(75)]]\n", "[[-0.283(26) -0.9592(75)]\n",
" [Obs[-0.9592(75)] Obs[0.283(26)]]]\n" " [-0.9592(75) 0.283(26)]]\n"
] ]
} }
], ],
"source": [ "source": [
"e, v = pe.linalg.eigh(matrix)\n", "e, v = pe.linalg.eigh(matrix)\n",
"for (i), entry in np.ndenumerate(e):\n", "pe.gm(e)\n",
" entry.gamma_method()\n",
"print('Eigenvalues:')\n", "print('Eigenvalues:')\n",
"print(e)\n", "print(e)\n",
"for (i, j), entry in np.ndenumerate(v):\n", "pe.gm(v)\n",
" entry.gamma_method()\n",
"print('Eigenvectors:')\n", "print('Eigenvectors:')\n",
"print(v)" "print(v)"
] ]
@ -371,9 +364,9 @@
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"Check eigenvector 1\n", "Check eigenvector 1\n",
"[Obs[-5.551115123125783e-17] Obs[0.0]]\n", "[-5.551115123125783e-17 0.0]\n",
"Check eigenvector 2\n", "Check eigenvector 2\n",
"[Obs[0.0] Obs[-2.220446049250313e-16]]\n" "[0.0 -2.220446049250313e-16]\n"
] ]
} }
], ],
@ -407,7 +400,7 @@
"name": "python", "name": "python",
"nbconvert_exporter": "python", "nbconvert_exporter": "python",
"pygments_lexer": "ipython3", "pygments_lexer": "ipython3",
"version": "3.8.10" "version": "3.11.0"
} }
}, },
"nbformat": 4, "nbformat": 4,

16
pyerrors/obs.py
View file

@ -1032,6 +1032,22 @@ class CObs:
return f"({self.real:{format_type}}{self.imag:+{significance}}j)" return f"({self.real:{format_type}}{self.imag:+{significance}}j)"
def gamma_method(x, **kwargs):
"""Vectorized version of the gamma_method applicable to lists or arrays of Obs.
See docstring of pe.Obs.gamma_method for details.
"""
return np.vectorize(lambda o: o.gm(**kwargs))(x)
def gm(x, **kwargs):
"""Short version of the vectorized gamma_method.
See docstring of pe.Obs.gamma_method for details
"""
return gamma_method(x, **kwargs)
def _format_uncertainty(value, dvalue, significance=2): def _format_uncertainty(value, dvalue, significance=2):
"""Creates a string of a value and its error in paranthesis notation, e.g., 13.02(45)""" """Creates a string of a value and its error in paranthesis notation, e.g., 13.02(45)"""
if dvalue == 0.0 or (not np.isfinite(dvalue)): if dvalue == 0.0 or (not np.isfinite(dvalue)):

33
tests/obs_test.py
View file

@ -1,4 +1,5 @@
import autograd.numpy as np import numpy as np
import autograd.numpy as anp
import os import os
import copy import copy
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
@ -126,9 +127,9 @@ def test_function_overloading():
fs = [lambda x: x[0] + x[1], lambda x: x[1] + x[0], lambda x: x[0] - x[1], lambda x: x[1] - x[0], fs = [lambda x: x[0] + x[1], lambda x: x[1] + x[0], lambda x: x[0] - x[1], lambda x: x[1] - x[0],
lambda x: x[0] * x[1], lambda x: x[1] * x[0], lambda x: x[0] / x[1], lambda x: x[1] / x[0], lambda x: x[0] * x[1], lambda x: x[1] * x[0], lambda x: x[0] / x[1], lambda x: x[1] / x[0],
lambda x: np.exp(x[0]), lambda x: np.sin(x[0]), lambda x: np.cos(x[0]), lambda x: np.tan(x[0]), lambda x: anp.exp(x[0]), lambda x: anp.sin(x[0]), lambda x: anp.cos(x[0]), lambda x: anp.tan(x[0]),
lambda x: np.log(x[0]), lambda x: np.sqrt(np.abs(x[0])), lambda x: anp.log(x[0]), lambda x: anp.sqrt(anp.abs(x[0])),
lambda x: np.sinh(x[0]), lambda x: np.cosh(x[0]), lambda x: np.tanh(x[0])] lambda x: anp.sinh(x[0]), lambda x: anp.cosh(x[0]), lambda x: anp.tanh(x[0])]
for i, f in enumerate(fs): for i, f in enumerate(fs):
t1 = f([a, b]) t1 = f([a, b])
@ -306,9 +307,9 @@ def test_derived_observables():
test_obs = pe.pseudo_Obs(2, 0.1 * (1 + np.random.rand()), 't', int(1000 * (1 + np.random.rand()))) test_obs = pe.pseudo_Obs(2, 0.1 * (1 + np.random.rand()), 't', int(1000 * (1 + np.random.rand())))
# Check if autograd and numgrad give the same result # Check if autograd and numgrad give the same result
d_Obs_ad = pe.derived_observable(lambda x, **kwargs: x[0] * x[1] * np.sin(x[0] * x[1]), [test_obs, test_obs]) d_Obs_ad = pe.derived_observable(lambda x, **kwargs: x[0] * x[1] * anp.sin(x[0] * x[1]), [test_obs, test_obs])
d_Obs_ad.gamma_method() d_Obs_ad.gamma_method()
d_Obs_fd = pe.derived_observable(lambda x, **kwargs: x[0] * x[1] * np.sin(x[0] * x[1]), [test_obs, test_obs], num_grad=True) d_Obs_fd = pe.derived_observable(lambda x, **kwargs: x[0] * x[1] * anp.sin(x[0] * x[1]), [test_obs, test_obs], num_grad=True)
d_Obs_fd.gamma_method() d_Obs_fd.gamma_method()
assert d_Obs_ad == d_Obs_fd assert d_Obs_ad == d_Obs_fd
@ -1283,12 +1284,14 @@ def test_format_uncertainty():
pe.obs._format_uncertainty(1, np.NaN) pe.obs._format_uncertainty(1, np.NaN)
pe.obs._format_uncertainty(np.NaN, np.inf) pe.obs._format_uncertainty(np.NaN, np.inf)
def test_format(): def test_format():
o1 = pe.pseudo_Obs(0.348, 0.0123, "test") o1 = pe.pseudo_Obs(0.348, 0.0123, "test")
assert o1.__format__("+3") == '+0.3480(123)' assert o1.__format__("+3") == '+0.3480(123)'
assert o1.__format__("+2") == '+0.348(12)' assert o1.__format__("+2") == '+0.348(12)'
assert o1.__format__(" 2") == ' 0.348(12)' assert o1.__format__(" 2") == ' 0.348(12)'
def test_f_string_obs(): def test_f_string_obs():
o1 = pe.pseudo_Obs(0.348, 0.0123, "test") o1 = pe.pseudo_Obs(0.348, 0.0123, "test")
print(f"{o1}") print(f"{o1}")
@ -1297,6 +1300,7 @@ def test_f_string_obs():
print(f"{o1:-1}") print(f"{o1:-1}")
print(f"{o1: 8}") print(f"{o1: 8}")
def test_f_string_cobs(): def test_f_string_cobs():
o_real = pe.pseudo_Obs(0.348, 0.0123, "test") o_real = pe.pseudo_Obs(0.348, 0.0123, "test")
o_imag = pe.pseudo_Obs(0.348, 0.0123, "test") o_imag = pe.pseudo_Obs(0.348, 0.0123, "test")
@ -1307,7 +1311,24 @@ def test_f_string_cobs():
print(f"{o1:-1}") print(f"{o1:-1}")
print(f"{o1: 8}") print(f"{o1: 8}")
def test_compute_drho_fails(): def test_compute_drho_fails():
obs = pe.input.json.load_json("tests/data/compute_drho_fails.json.gz") obs = pe.input.json.load_json("tests/data/compute_drho_fails.json.gz")
obs.gm() obs.gm()
assert np.isclose(obs.dvalue, 0.0022150779611891094) assert np.isclose(obs.dvalue, 0.0022150779611891094)
def test_vec_gm():
obs = pe.misc.gen_correlated_data(np.arange(3), np.array([[0.0364 , 0.03627262, 0.03615699],
[0.03627262, 0.03688438, 0.03674798],
[0.03615699, 0.03674798, 0.03732882]]), "qq", 3.8, 1000)
pe.gm(obs[0], S=0)
assert obs[0].S["qq"] == 0
pe.gm(obs, S=1.3)
assert np.all(np.vectorize(lambda x: x.S["qq"])(obs) == 1.3)
aa = np.array([obs, obs, obs])
pe.gamma_method(aa, S=2.2)
assert np.all(np.vectorize(lambda x: x.S["qq"])(aa) == 2.20)
cc = pe.Corr(obs)
pe.gm(cc, S=4.12)
assert np.all(np.vectorize(lambda x: x.S["qq"])(cc.content) == 4.12)