From 9eb1dc439049435deb3b1f4996b8ca0a5568f606 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 15:47:11 +0200 Subject: [PATCH 01/12] Interface QED and N3LO --- .../unpolarized/space_like/__init__.py | 9 +++ .../unpolarized/space_like/as4/__init__.py | 81 +++++++++++++++++++ 2 files changed, 90 insertions(+) diff --git a/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py b/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py index 6a3c6cf06..0fbceb805 100644 --- a/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py +++ b/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py @@ -155,6 +155,11 @@ def gamma_ns_qed(order, mode, n, nf): gamma_ns[3, 0] = as3.gamma_nsp(n, nf, cache) elif mode in [10202, 10203]: gamma_ns[3, 0] = as3.gamma_nsm(n, nf, cache) + if order[0] >= 4: + if mode in [10102, 10103]: + gamma_ns[4, 0] = as4.gamma_nsp(n, nf, cache) + elif mode in [10202, 10203]: + gamma_ns[4, 0] = as4.gamma_nsm(n, nf, cache) return gamma_ns @@ -280,6 +285,8 @@ def gamma_singlet_qed(order, n, nf): gamma_s[0, 2] = aem2.gamma_singlet(n, nf, cache) if order[0] >= 3: gamma_s[3, 0] = as3.gamma_singlet_qed(n, nf, cache) + if order[0] >= 4: + gamma_s[4, 0] = as4.gamma_singlet_qed(n, nf, cache) return gamma_s @@ -314,4 +321,6 @@ def gamma_valence_qed(order, n, nf): gamma_v[0, 2] = aem2.gamma_valence(n, nf, cache) if order[0] >= 3: gamma_v[3, 0] = as3.gamma_valence_qed(n, nf, cache) + if order[0] >= 4: + gamma_v[4, 0] = as4.gamma_valence_qed(n, nf, cache) return gamma_v diff --git a/src/ekore/anomalous_dimensions/unpolarized/space_like/as4/__init__.py b/src/ekore/anomalous_dimensions/unpolarized/space_like/as4/__init__.py index 10c4cb57c..56c2c16ff 100644 --- a/src/ekore/anomalous_dimensions/unpolarized/space_like/as4/__init__.py +++ b/src/ekore/anomalous_dimensions/unpolarized/space_like/as4/__init__.py @@ -55,3 +55,84 @@ def gamma_singlet(N, nf, cache, variation): np.complex_, ) return gamma_S_0 + + +@nb.njit(cache=True) +def gamma_singlet_qed(N, nf, cache, variation): + r"""Compute the leading-order singlet anomalous dimension matrix for the unified evolution basis. + + .. math:: + \\gamma_S^{(3,0)} = \\left(\begin{array}{cccc} + \\gamma_{gg}^{(3,0)} & 0 & \\gamma_{gq}^{(3,0)} & 0\\ + 0 & 0 & 0 & 0 \\ + \\gamma_{qg}^{(3,0)} & 0 & \\gamma_{qq}^{(3,0)} & 0 \\ + 0 & 0 & 0 & \\gamma_{qq}^{(3,0)} \\ + \\end{array}\right) + + Parameters + ---------- + N : complex + Mellin moment + nf : int + Number of active flavors + cache: numpy.ndarray + Harmonic sum cache + + Returns + ------- + numpy.ndarray + Leading-order singlet anomalous dimension matrix :math:`\\gamma_{S}^{(3,0)}(N)` + + """ + gamma_np_p = gamma_nsp(N, nf, cache) + gamma_qq = gamma_np_p + gamma_ps(N, nf, cache, variation[3]) + gamma_S = np.array( + [ + [ + gamma_gg(N, nf, cache, variation[0]), + 0.0 + 0.0j, + gamma_gq(N, nf, cache, variation[1]), + 0.0 + 0.0j, + ], + [0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j], + [gamma_qg(N, nf, cache, variation[2]), 0.0 + 0.0j, gamma_qq, 0.0 + 0.0j], + [0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j, gamma_np_p], + ], + np.complex_, + ) + return gamma_S + + +@nb.njit(cache=True) +def gamma_valence_qed(N, nf, cache): + r"""Compute the leading-order valence anomalous dimension matrix for the unified evolution basis. + + .. math:: + \\gamma_V^{(3,0)} = \\left(\begin{array}{cc} + \\gamma_{nsV}^{(3,0)} & 0\\ + 0 & \\gamma_{ns-}^{(3,0)} + \\end{array}\right) + + Parameters + ---------- + N : complex + Mellin moment + nf : int + Number of active flavors + cache: numpy.ndarray + Harmonic sum cache + + Returns + ------- + numpy.ndarray + Leading-order singlet anomalous dimension matrix :math:`\\gamma_{V}^{(3,0)}(N)` + + """ + gamma_V = np.array( + [ + [gamma_nsv(N, nf, cache), 0.0], + [0.0, gamma_nsm(N, nf, cache)], + ], + np.complex_, + ) + return gamma_V From 089391fa2e9a0f6c4c24736ab923f858bb634e34 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 15:57:03 +0200 Subject: [PATCH 02/12] Add n3lo_ad_variation key to missing functions --- src/eko/evolution_operator/__init__.py | 5 ++++- .../anomalous_dimensions/unpolarized/space_like/__init__.py | 6 ++++-- 2 files changed, 8 insertions(+), 3 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index eab3c64de..27b8933dc 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -450,6 +450,7 @@ def quad_ker_qed( ev_op_max_order, sv_mode, is_threshold, + n3lo_ad_variation, ): """Raw evolution kernel inside quad. @@ -489,6 +490,8 @@ def quad_ker_qed( scale variation mode, see `eko.scale_variations.Modes` is_threshold : boolean is this an itermediate threshold operator? + n3lo_ad_variation : tuple + |N3LO| anomalous dimension variation ``(gg_var, gq_var, qg_var, qq_var)`` Returns ------- @@ -497,7 +500,7 @@ def quad_ker_qed( """ # compute the actual evolution kernel for QEDxQCD if ker_base.is_QEDsinglet: - gamma_s = ad_us.gamma_singlet_qed(order, ker_base.n, nf) + gamma_s = ad_us.gamma_singlet_qed(order, ker_base.n, nf, n3lo_ad_variation) # scale var exponentiated is directly applied on gamma if sv_mode == sv.Modes.exponentiated: gamma_s = sv.exponentiated.gamma_variation_qed( diff --git a/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py b/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py index 0fbceb805..6358a5af8 100644 --- a/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py +++ b/src/ekore/anomalous_dimensions/unpolarized/space_like/__init__.py @@ -255,7 +255,7 @@ def choose_ns_ad_aem2(mode, n, nf, cache): @nb.njit(cache=True) -def gamma_singlet_qed(order, n, nf): +def gamma_singlet_qed(order, n, nf, n3lo_ad_variation): r""" Compute the grid of the QED singlet anomalous dimensions matrices. @@ -267,6 +267,8 @@ def gamma_singlet_qed(order, n, nf): Mellin variable nf : int Number of active flavors + n3lo_ad_variation : tuple + |N3LO| anomalous dimension variation ``(gg_var, gq_var, qg_var, qq_var)`` Returns ------- @@ -286,7 +288,7 @@ def gamma_singlet_qed(order, n, nf): if order[0] >= 3: gamma_s[3, 0] = as3.gamma_singlet_qed(n, nf, cache) if order[0] >= 4: - gamma_s[4, 0] = as4.gamma_singlet_qed(n, nf, cache) + gamma_s[4, 0] = as4.gamma_singlet_qed(n, nf, cache, n3lo_ad_variation) return gamma_s From 9db7ac84a8345db686bdd5dbabeb9f9b0caa5c4b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 16:05:28 +0200 Subject: [PATCH 03/12] Add n3lo_ad_variation key to missing functions. again --- src/eko/evolution_operator/__init__.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index 27b8933dc..122bfbc69 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -307,6 +307,7 @@ def quad_ker( ev_op_max_order, sv_mode, is_threshold, + n3lo_ad_variation, ) # recombine everything From 24bbfd30b5a5d0024abcda90b77bb7119591078a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 16:14:18 +0200 Subject: [PATCH 04/12] Fix test_kernels_QEDsinglet.py --- tests/eko/kernels/test_kernels_QEDsinglet.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/eko/kernels/test_kernels_QEDsinglet.py b/tests/eko/kernels/test_kernels_QEDsinglet.py index 302bb941b..f5803750c 100644 --- a/tests/eko/kernels/test_kernels_QEDsinglet.py +++ b/tests/eko/kernels/test_kernels_QEDsinglet.py @@ -79,7 +79,7 @@ def test_zero_true_gamma(monkeypatch): for qed in range(1, 2 + 1): order = (qcd, qed) n = np.random.rand() - gamma_s = ad.gamma_singlet_qed(order, n, nf) + gamma_s = ad.gamma_singlet_qed(order, n, nf, (0, 0, 0, 0)) # monkeypatch.setattr( # ad, # "exp_matrix_2D", From 350871307f88cfebda1ac63446b6a00a23be918f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 16:30:38 +0200 Subject: [PATCH 05/12] Add order=4 in fixed_alphaem_exact --- src/eko/kernels/non_singlet_qed.py | 25 +++++++++++++++++++++++++ 1 file changed, 25 insertions(+) diff --git a/src/eko/kernels/non_singlet_qed.py b/src/eko/kernels/non_singlet_qed.py index 86ec1a9ff..36d4d3c35 100644 --- a/src/eko/kernels/non_singlet_qed.py +++ b/src/eko/kernels/non_singlet_qed.py @@ -120,6 +120,29 @@ def as3_exact(gamma_ns, a1, a0, beta): return ns.nnlo_exact(gamma_ns, a1, a0, beta) +@nb.njit(cache=True) +def as4_exact(gamma_ns, a1, a0, beta): + """O(as3aem1) non-singlet exact EKO. + + Parameters + ---------- + gamma_ns : numpy.ndarray + non-singlet anomalous dimensions + a1 : float + target coupling value + a0 : float + initial coupling value + beta : list + list of the values of the beta functions + + Returns + ------- + e_ns^3 : complex + O(as4aem1) non-singlet exact EKO + """ + return ns.n3lo_exact(gamma_ns, a1, a0, beta) + + @nb.njit(cache=True) def dispatcher( order, @@ -216,6 +239,8 @@ def fixed_alphaem_exact(order, gamma_ns, a1, a0, aem, nf, mu2_from, mu2_to): qcd_only = as2_exact(gamma_ns_list[1:], a1, a0, betalist) elif order[0] == 3: qcd_only = as3_exact(gamma_ns_list[1:], a1, a0, betalist) + elif order[0] == 4: + qcd_only = as4_exact(gamma_ns_list[1:], a1, a0, betalist) else: raise NotImplementedError("Selected order is not implemented") return qcd_only * apply_qed(gamma_ns_list[0], mu2_from, mu2_to) From 0118d11807948113b59bc7c8149eb43990d52c04 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 16:45:39 +0200 Subject: [PATCH 06/12] Fix test_kernels_QEDns.py --- tests/eko/kernels/test_kernels_QEDns.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/eko/kernels/test_kernels_QEDns.py b/tests/eko/kernels/test_kernels_QEDns.py index e84956084..5574b1a12 100644 --- a/tests/eko/kernels/test_kernels_QEDns.py +++ b/tests/eko/kernels/test_kernels_QEDns.py @@ -304,7 +304,7 @@ def test_error(): for running in [True, False]: with pytest.raises(NotImplementedError): ns.dispatcher( - (4, 2), + (5, 2), "iterate-exact", np.random.rand(4, 3), [0.1, 0.2], From 361afa8d2b4cef4843e459374cb6e45437991fba Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 16:56:16 +0200 Subject: [PATCH 07/12] Fix again test --- tests/eko/kernels/test_kernels_QEDns.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/eko/kernels/test_kernels_QEDns.py b/tests/eko/kernels/test_kernels_QEDns.py index 5574b1a12..f2bd45731 100644 --- a/tests/eko/kernels/test_kernels_QEDns.py +++ b/tests/eko/kernels/test_kernels_QEDns.py @@ -302,7 +302,7 @@ def test_ode_true_gamma(): def test_error(): for running in [True, False]: - with pytest.raises(NotImplementedError): + with pytest.raises(ValueError): ns.dispatcher( (5, 2), "iterate-exact", From 17404d359beb475c67ea2a29a9f52be9ebf795c1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Fri, 6 Oct 2023 17:09:22 +0200 Subject: [PATCH 08/12] Fix another test --- .../anomalous_dimensions/unpolarized/space_like/test_init.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_init.py b/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_init.py index 2175eddcd..c7c955960 100644 --- a/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_init.py +++ b/tests/ekore/anomalous_dimensions/unpolarized/space_like/test_init.py @@ -208,7 +208,7 @@ def test_dim_singlet(): nf = 3 N = 2 cache = h.cache.reset() - gamma_singlet = ad_us.gamma_singlet_qed((3, 2), N, nf) + gamma_singlet = ad_us.gamma_singlet_qed((3, 2), N, nf, (0, 0, 0, 0)) assert gamma_singlet.shape == (4, 3, 4, 4) gamma_singlet_as1 = ad_us.as1.gamma_singlet_qed(N, cache, nf) assert gamma_singlet_as1.shape == (4, 4) From 263f6063cca6e23535209feb2cc525f3bb14d40e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Sat, 7 Oct 2023 19:44:38 +0200 Subject: [PATCH 09/12] Change nf->beta into ns kernels at n3lo --- src/eko/kernels/non_singlet.py | 24 ++++++++---------------- 1 file changed, 8 insertions(+), 16 deletions(-) diff --git a/src/eko/kernels/non_singlet.py b/src/eko/kernels/non_singlet.py index 24706a794..23b2ab47a 100644 --- a/src/eko/kernels/non_singlet.py +++ b/src/eko/kernels/non_singlet.py @@ -186,7 +186,7 @@ def nnlo_expanded(gamma_ns, a1, a0, beta): @nb.njit(cache=True) -def n3lo_expanded(gamma_ns, a1, a0, nf): +def n3lo_expanded(gamma_ns, a1, a0, beta): """|N3LO| non-singlet expanded EKO. Parameters @@ -206,12 +206,8 @@ def n3lo_expanded(gamma_ns, a1, a0, nf): |N3LO| non-singlet expanded EKO """ - beta0 = beta.beta_qcd((2, 0), nf) - b_list = [ - beta.b_qcd((3, 0), nf), - beta.b_qcd((4, 0), nf), - beta.b_qcd((5, 0), nf), - ] + beta0 = beta[0] + b_list = [betas / beta0 for betas in beta[1:]] j12 = ei.j12(a1, a0, beta0) j13 = as4_ei.j13_expanded(a1, a0, beta0, b_list) j23 = as4_ei.j23_expanded(a1, a0, beta0, b_list) @@ -225,7 +221,7 @@ def n3lo_expanded(gamma_ns, a1, a0, nf): @nb.njit(cache=True) -def n3lo_exact(gamma_ns, a1, a0, nf): +def n3lo_exact(gamma_ns, a1, a0, beta): """|N3LO| non-singlet exact EKO. Parameters @@ -245,12 +241,8 @@ def n3lo_exact(gamma_ns, a1, a0, nf): |N3LO| non-singlet exact EKO """ - beta0 = beta.beta_qcd((2, 0), nf) - b_list = [ - beta.b_qcd((3, 0), nf), - beta.b_qcd((4, 0), nf), - beta.b_qcd((5, 0), nf), - ] + beta0 = beta[0] + b_list = [betas / beta0 for betas in beta[1:]] roots = as4_ei.roots(b_list) j12 = ei.j12(a1, a0, beta0) j13 = as4_ei.j13_exact(a1, a0, beta0, b_list, roots) @@ -420,6 +412,6 @@ def dispatcher( "decompose-expanded", "perturbative-expanded", ]: - return n3lo_expanded(gamma_ns, a1, a0, nf) - return n3lo_exact(gamma_ns, a1, a0, nf) + return n3lo_expanded(gamma_ns, a1, a0, betalist) + return n3lo_exact(gamma_ns, a1, a0, betalist) raise NotImplementedError("Selected order is not implemented") From 9393ea2c125cab80f970e62b9328f74059118458 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Sun, 8 Oct 2023 15:43:31 +0200 Subject: [PATCH 10/12] Add pto = 3 in some qed tests --- tests/eko/evolution_operator/test_init.py | 2 +- tests/eko/kernels/test_kernels_QEDns.py | 14 +++++++++----- tests/eko/kernels/test_kernels_QEDsinglet.py | 4 ++-- tests/eko/kernels/test_kernels_QEDvalence.py | 4 ++-- 4 files changed, 14 insertions(+), 10 deletions(-) diff --git a/tests/eko/evolution_operator/test_init.py b/tests/eko/evolution_operator/test_init.py index 642de5bb7..47d5e7001 100644 --- a/tests/eko/evolution_operator/test_init.py +++ b/tests/eko/evolution_operator/test_init.py @@ -382,7 +382,7 @@ def test_compute(self, monkeypatch, theory_ffns, operator_card, tmp_path): err_msg=k, ) - for n in range(1, 3 + 1): + for n in range(1, 4 + 1): for qed in range(1, 2 + 1): g.config["order"] = (n, qed) o1 = Operator(g.config, g.managers, Segment(2.0, 2.0, 3)) diff --git a/tests/eko/kernels/test_kernels_QEDns.py b/tests/eko/kernels/test_kernels_QEDns.py index f2bd45731..e6179ba99 100644 --- a/tests/eko/kernels/test_kernels_QEDns.py +++ b/tests/eko/kernels/test_kernels_QEDns.py @@ -25,7 +25,7 @@ def test_zero(): nf = 3 ev_op_iterations = 2 running_alpha = [True, False] - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) gamma_ns = ( @@ -73,7 +73,7 @@ def test_zero_true_gamma(): for mode in br.non_singlet_pids_map.values(): if mode in [10201, 10101, 10200]: continue - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) n = np.random.rand() @@ -136,7 +136,8 @@ def test_ode(): for mode in br.non_singlet_pids_map.values(): if mode in [10201, 10101, 10200]: continue - for qcd in range(1, 3 + 1): + max_qcd = 4 + for qcd in range(1, max_qcd + 1): for qed in range(1, 2 + 1): order = (qcd, qed) sc = Couplings( @@ -154,14 +155,17 @@ def test_ode(): ) a1 = sc.a_s(mu2_to) gamma_ns = ( - np.random.rand(3 + 1, 2 + 1) + np.random.rand(3 + 1, 2 + 1) * 1j + np.random.rand(max_qcd + 1, 2 + 1) + + np.random.rand(max_qcd + 1, 2 + 1) * 1j ) gamma_ns[0, 0] = 0.0 gamma_ns[2, 1] = 0.0 gamma_ns[3, 1] = 0.0 + gamma_ns[4, 1] = 0.0 gamma_ns[1, 2] = 0.0 gamma_ns[2, 2] = 0.0 gamma_ns[3, 2] = 0.0 + gamma_ns[4, 2] = 0.0 gammatot = 0.0 for i in range(0, order[0] + 1): for j in range(0, order[1] + 1): @@ -226,7 +230,7 @@ def test_ode_true_gamma(): for mode in br.non_singlet_pids_map.values(): if mode in [10201, 10101, 10200]: continue - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) sc = Couplings( diff --git a/tests/eko/kernels/test_kernels_QEDsinglet.py b/tests/eko/kernels/test_kernels_QEDsinglet.py index f5803750c..b953c053c 100644 --- a/tests/eko/kernels/test_kernels_QEDsinglet.py +++ b/tests/eko/kernels/test_kernels_QEDsinglet.py @@ -21,7 +21,7 @@ def test_zero(monkeypatch): nf = 3 ev_op_iterations = 2 ev_op_max_order = (3, 0) - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) gamma_s = ( @@ -75,7 +75,7 @@ def test_zero_true_gamma(monkeypatch): nf = 3 ev_op_iterations = 2 ev_op_max_order = (3, 0) - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) n = np.random.rand() diff --git a/tests/eko/kernels/test_kernels_QEDvalence.py b/tests/eko/kernels/test_kernels_QEDvalence.py index 95d0cbc4c..a27bc6ef6 100644 --- a/tests/eko/kernels/test_kernels_QEDvalence.py +++ b/tests/eko/kernels/test_kernels_QEDvalence.py @@ -21,7 +21,7 @@ def test_zero(monkeypatch): nf = 3 ev_op_iterations = 2 ev_op_max_order = (3, 0) - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) gamma_v = ( @@ -62,7 +62,7 @@ def test_zero_true_gamma(monkeypatch): nf = 3 ev_op_iterations = 2 ev_op_max_order = (3, 0) - for qcd in range(1, 3 + 1): + for qcd in range(1, 4 + 1): for qed in range(1, 2 + 1): order = (qcd, qed) n = np.random.rand() From 569a0fe8f972eaf0b71549c0f8969542b9e38336 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Tue, 10 Oct 2023 16:42:26 +0200 Subject: [PATCH 11/12] restructure self.compute_aem_list --- src/eko/evolution_operator/__init__.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index 122bfbc69..29e67c19c 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -626,7 +626,7 @@ def __init__( self.alphaem_running = self.managers["couplings"].alphaem_running if self.log_label == "Evolution": self.a = self.compute_a() - self.compute_aem_list() + self.as_list, self.a_half_list = self.compute_aem_list() @property def n_pools(self): @@ -701,8 +701,8 @@ def compute_aem_list(self): """ ev_op_iterations = self.config["ev_op_iterations"] if self.order[1] == 0: - self.as_list = np.array([self.a_s[0], self.a_s[1]]) - self.a_half_list = np.zeros((ev_op_iterations, 2)) + as_list = np.array([self.a_s[0], self.a_s[1]]) + a_half = np.zeros((ev_op_iterations, 2)) else: as0 = self.a_s[0] as1 = self.a_s[1] @@ -722,7 +722,7 @@ def compute_aem_list(self): couplings = self.managers["couplings"] mu2_steps = utils.geomspace(self.q2_from, self.q2_to, 1 + ev_op_iterations) mu2_l = mu2_steps[0] - self.as_list = np.array( + as_list = np.array( [ couplings.compute( a_ref=a_start, nf=self.nf, scale_from=mu2_start, scale_to=mu2 @@ -738,7 +738,7 @@ def compute_aem_list(self): ) a_half[step] = [a_s, aem] mu2_l = mu2_h - self.a_half_list = a_half + return as_list, a_half @property def labels(self): From 2f65ee975f869bad50e79c17fea3a3a9a13427e8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Niccol=C3=B2=20Laurenti?= Date: Tue, 10 Oct 2023 17:41:23 +0200 Subject: [PATCH 12/12] Remove unused default parameters and remove qed argument from apply_pdf --- src/eko/basis_rotation.py | 4 ++-- src/eko/evolution_operator/flavors.py | 8 ++++---- src/eko/evolution_operator/matching_condition.py | 2 +- src/eko/evolution_operator/physical.py | 2 +- src/ekobox/apply.py | 11 ++++------- src/ekomark/benchmark/runner.py | 1 - tests/eko/evolution_operator/test_flavors.py | 16 +++++++++------- .../test_matching_condition.py | 10 +++++----- tests/eko/test_basis_rotation.py | 8 ++++---- 9 files changed, 30 insertions(+), 32 deletions(-) diff --git a/src/eko/basis_rotation.py b/src/eko/basis_rotation.py index 4b972592c..e31a7beb6 100644 --- a/src/eko/basis_rotation.py +++ b/src/eko/basis_rotation.py @@ -275,7 +275,7 @@ } -def ad_projector(ad_lab, nf, qed=False): +def ad_projector(ad_lab, nf, qed): """ Build a projector (as a numpy array) for the given anomalous dimension sector. @@ -355,7 +355,7 @@ def select_light_flavors_uni_ev(ad_lab, nf): return map_ad_to_evolution[ad_lab] -def ad_projectors(nf, qed=False): +def ad_projectors(nf, qed): """ Build projectors tensor (as a numpy array), collecting all the individual sector projectors. diff --git a/src/eko/evolution_operator/flavors.py b/src/eko/evolution_operator/flavors.py index 7ad7b5801..fac94a608 100644 --- a/src/eko/evolution_operator/flavors.py +++ b/src/eko/evolution_operator/flavors.py @@ -52,7 +52,7 @@ def pids_from_intrinsic_evol(label, nf, normalize): return weights -def get_range(evol_labels, qed=False): +def get_range(evol_labels, qed): """Determine the number of light and heavy flavors participating in the input and output. Here, we assume that the T distributions (e.g. T15) appears *always* @@ -73,7 +73,7 @@ def get_range(evol_labels, qed=False): nf_in = 3 nf_out = 3 - def update(label, qed=False): + def update(label, qed): nf = 3 if label[0] == "T": if not qed: @@ -129,7 +129,7 @@ def rotate_pm_to_flavor(label): return l -def rotate_matching(nf, qed=False, inverse=False): +def rotate_matching(nf, qed, inverse=False): """Rotation between matching basis (with e.g. S,g,...V8 and c+,c-) and new true evolution basis (with S,g,...V8,T15,V15). Parameters @@ -206,7 +206,7 @@ def rotate_matching(nf, qed=False, inverse=False): return l -def rotate_matching_inverse(nf, qed=False): +def rotate_matching_inverse(nf, qed): """Inverse rotation between matching basis (with e.g. S,g,...V8 and c+,c-) and new true evolution basis (with S,g,...V8,T15,V15). Parameters diff --git a/src/eko/evolution_operator/matching_condition.py b/src/eko/evolution_operator/matching_condition.py index 059617fbf..8df4775ff 100644 --- a/src/eko/evolution_operator/matching_condition.py +++ b/src/eko/evolution_operator/matching_condition.py @@ -21,7 +21,7 @@ def split_ad_to_evol_map( nf, q2_thr, intrinsic_range, - qed=False, + qed, ): """ Create the instance from the |OME|. diff --git a/src/eko/evolution_operator/physical.py b/src/eko/evolution_operator/physical.py index c963e2a26..1ee501a2a 100644 --- a/src/eko/evolution_operator/physical.py +++ b/src/eko/evolution_operator/physical.py @@ -22,7 +22,7 @@ class PhysicalOperator(member.OperatorBase): """ @classmethod - def ad_to_evol_map(cls, op_members, nf, q2_final, intrinsic_range, qed=False): + def ad_to_evol_map(cls, op_members, nf, q2_final, intrinsic_range, qed): """ Obtain map between the 3-dimensional anomalous dimension basis and the 4-dimensional evolution basis. diff --git a/src/ekobox/apply.py b/src/ekobox/apply.py index c392ffe78..f59928186 100644 --- a/src/ekobox/apply.py +++ b/src/ekobox/apply.py @@ -12,7 +12,6 @@ def apply_pdf( lhapdf_like, targetgrid=None, rotate_to_evolution_basis=False, - qed=False, ): """ Apply all available operators to the input PDFs. @@ -34,13 +33,14 @@ def apply_pdf( out_grid : dict output PDFs and their associated errors for the computed mu2grid """ + qed = eko.theory_card.order[1] > 0 if rotate_to_evolution_basis: if not qed: rotate_flavor_to_evolution = br.rotate_flavor_to_evolution else: rotate_flavor_to_evolution = br.rotate_flavor_to_unified_evolution return apply_pdf_flavor( - eko, lhapdf_like, targetgrid, rotate_flavor_to_evolution, qed + eko, lhapdf_like, targetgrid, rotate_flavor_to_evolution ) return apply_pdf_flavor(eko, lhapdf_like, targetgrid) @@ -48,9 +48,7 @@ def apply_pdf( CONTRACTION = "ajbk,bk" -def apply_pdf_flavor( - eko: EKO, lhapdf_like, targetgrid=None, flavor_rotation=None, qed=False -): +def apply_pdf_flavor(eko: EKO, lhapdf_like, targetgrid=None, flavor_rotation=None): """ Apply all available operators to the input PDFs. @@ -65,8 +63,6 @@ def apply_pdf_flavor( if given, interpolates to the pdfs given at targetgrid (instead of xgrid) flavor_rotation : np.ndarray Rotation matrix in flavor space - qed : bool - activate qed Returns ------- @@ -97,6 +93,7 @@ def apply_pdf_flavor( if error_final is not None: out_grid[ep]["errors"] = dict(zip(eko.bases.targetpids, error_final)) + qed = eko.theory_card.order[1] > 0 # rotate to evolution basis if flavor_rotation is not None: for q2, op in out_grid.items(): diff --git a/src/ekomark/benchmark/runner.py b/src/ekomark/benchmark/runner.py index 8c07e7153..e9ae298ee 100644 --- a/src/ekomark/benchmark/runner.py +++ b/src/ekomark/benchmark/runner.py @@ -243,7 +243,6 @@ def log(self, theory, _, pdf, me, ext): pdf, xgrid, flavor_rotation=rotate_to_evolution, - qed=qed, ) for q2, ref_pdfs in ext["values"].items(): log_tab = dfdict.DFdict() diff --git a/tests/eko/evolution_operator/test_flavors.py b/tests/eko/evolution_operator/test_flavors.py index 81bb35fc4..2f3761404 100644 --- a/tests/eko/evolution_operator/test_flavors.py +++ b/tests/eko/evolution_operator/test_flavors.py @@ -61,13 +61,15 @@ def get(d): def test_get_range(): - assert (3, 3) == flavors.get_range([]) - assert (3, 3) == flavors.get_range([member.MemberName(n) for n in ["S.S", "V3.V3"]]) + assert (3, 3) == flavors.get_range([], False) + assert (3, 3) == flavors.get_range( + [member.MemberName(n) for n in ["S.S", "V3.V3"]], False + ) assert (3, 4) == flavors.get_range( - [member.MemberName(n) for n in ["S.S", "V3.V3", "T15.S"]] + [member.MemberName(n) for n in ["S.S", "V3.V3", "T15.S"]], False ) assert (4, 4) == flavors.get_range( - [member.MemberName(n) for n in ["S.S", "V3.V3", "T15.T15"]] + [member.MemberName(n) for n in ["S.S", "V3.V3", "T15.T15"]], False ) assert (3, 3) == flavors.get_range( [member.MemberName(n) for n in ["S.S", "Td3.Td3"]], True @@ -113,7 +115,7 @@ def test_rotate_pm_to_flavor(): def test_rotate_matching(): - m = flavors.rotate_matching(4) + m = flavors.rotate_matching(4, False) assert len(list(filter(lambda e: "c+" in e, m.keys()))) == 2 assert len(list(filter(lambda e: "b-" in e, m.keys()))) == 1 @@ -156,8 +158,8 @@ def load(m): return mm for nf in range(4, 6 + 1): - m = load(flavors.rotate_matching(nf)) - minv = load(flavors.rotate_matching_inverse(nf)) + m = load(flavors.rotate_matching(nf, False)) + minv = load(flavors.rotate_matching_inverse(nf, False)) np.testing.assert_allclose(m @ minv, np.eye(len(br.evol_basis)), atol=1e-10) diff --git a/tests/eko/evolution_operator/test_matching_condition.py b/tests/eko/evolution_operator/test_matching_condition.py index 8cbd6b537..c854fb90c 100644 --- a/tests/eko/evolution_operator/test_matching_condition.py +++ b/tests/eko/evolution_operator/test_matching_condition.py @@ -38,7 +38,7 @@ def update_intrinsic_OME(self, ome): def test_split_ad_to_evol_map(self): ome = self.mkOME() - a = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, []) + a = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [], False) triv_keys = [ "V.V", "T3.T3", @@ -64,7 +64,7 @@ def test_split_ad_to_evol_map(self): ome[(200, 200)].value, ) # # if alpha is zero, nothing non-trivial should happen - b = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, []) + b = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [], False) assert sorted(str(k) for k in b.op_members.keys()) == sorted( [*triv_keys, *keys3] ) @@ -74,7 +74,7 @@ def test_split_ad_to_evol_map(self): # ) # nf=3 + IC self.update_intrinsic_OME(ome) - c = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [4]) + c = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [4], False) assert sorted(str(k) for k in c.op_members.keys()) == sorted( [*triv_keys, *keys3, "S.c+", "g.c+", "c+.c+", "c-.c-"] ) @@ -83,7 +83,7 @@ def test_split_ad_to_evol_map(self): b.op_members[member.MemberName("V.V")].value, ) # nf=3 + IB - d = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [5]) + d = MatchingCondition.split_ad_to_evol_map(ome, 3, 1, [5], False) assert sorted(str(k) for k in d.op_members.keys()) == sorted( [*triv_keys, *keys3, "b+.b+", "b-.b-"] ) @@ -92,7 +92,7 @@ def test_split_ad_to_evol_map(self): np.eye(self.shape[0]), ) # nf=4 + IB - d = MatchingCondition.split_ad_to_evol_map(ome, 4, 1, [5]) + d = MatchingCondition.split_ad_to_evol_map(ome, 4, 1, [5], False) assert sorted(str(k) for k in d.op_members.keys()) == sorted( [ *triv_keys, diff --git a/tests/eko/test_basis_rotation.py b/tests/eko/test_basis_rotation.py index 680ef065d..194669d46 100644 --- a/tests/eko/test_basis_rotation.py +++ b/tests/eko/test_basis_rotation.py @@ -9,19 +9,19 @@ def test_ad_projector(): g = br.rotate_flavor_to_evolution[2] v3 = br.rotate_flavor_to_evolution[br.evol_basis.index("V3")] - s_to_s = br.ad_projector((100, 100), nf=6) + s_to_s = br.ad_projector((100, 100), nf=6, qed=False) np.testing.assert_allclose(s @ s_to_s, s) np.testing.assert_allclose(g @ s_to_s, 0.0) np.testing.assert_allclose(v3 @ s_to_s, 0.0) - g_to_s = br.ad_projector((21, 100), nf=6) + g_to_s = br.ad_projector((21, 100), nf=6, qed=False) np.testing.assert_allclose(s @ g_to_s, 0.0) np.testing.assert_allclose(g @ g_to_s, s) np.testing.assert_allclose(v3 @ g_to_s, 0.0) - ns_m = br.ad_projector((br.non_singlet_pids_map["ns-"], 0), nf=6) + ns_m = br.ad_projector((br.non_singlet_pids_map["ns-"], 0), nf=6, qed=False) np.testing.assert_allclose(s @ ns_m, 0.0, atol=1e-15) np.testing.assert_allclose(g @ ns_m, 0.0) @@ -74,7 +74,7 @@ def test_ad_projectors(): for nf in range(3, 6 + 1): diag = np.array([0] * (1 + 6 - nf) + [1] * (1 + 2 * nf) + [0] * (6 - nf)) identity = np.diag(diag) - projs = br.ad_projectors(nf) + projs = br.ad_projectors(nf, qed=False) # sum over diagonal projectors form an identity np.testing.assert_allclose(