From 109d758ff98ce4f475063f5be1712c7794d554a7 Mon Sep 17 00:00:00 2001 From: Peter Stangl Date: Wed, 3 Feb 2021 15:01:03 +0100 Subject: [PATCH] implement exact relations between (v,Mh2) and (m2,Lambda) (#59) * implement exact relations between (v,Mh2) and (m2,Lambda) * remove unnecessary parentheses * add comment on derivation of expressions * require positive arguments for `m2Lambda_to_vMh2` and `vMh2_to_m2Lambda` * use LF line endings --- wilson/run/smeft/smpar.py | 74 +++--- wilson/run/smeft/tests/test_smpar.py | 360 ++++++++++++++------------- 2 files changed, 224 insertions(+), 210 deletions(-) diff --git a/wilson/run/smeft/smpar.py b/wilson/run/smeft/smpar.py index 1fff99e..08114fa 100644 --- a/wilson/run/smeft/smpar.py +++ b/wilson/run/smeft/smpar.py @@ -37,45 +37,44 @@ def m2Lambda_to_vMh2(m2, Lambda, C): """Function to numerically determine the physical Higgs VEV and mass - given the parameters of the Higgs potential.""" - try: - v = (sqrt(2 * m2 / Lambda) + 3 * m2**(3 / 2) / - (sqrt(2) * Lambda**(5 / 2)) * C['phi'].real) - except ValueError: - v = 0 - Mh2 = 2 * m2 * (1 - m2 / Lambda * (3 * C['phi'].real - 4 * Lambda * C['phiBox'].real + - Lambda * C['phiD'].real)) - return {'v': v, 'Mh2': Mh2} + given the parameters of the Higgs potential. -def _vMh2_to_m2Lambda_SM(v, Mh2): - m2 = Mh2/2 - Lambda = 2 * m2 / v**2 - return {'m2': m2, 'Lambda': Lambda} + The relations used by this function have been derived from the Lagrangian + in eq. (3.1) of arXiv:1704.03888. This paper presents results expanded up to + linear order in the Wilson coefficients while the all-order expressions are + implemented here. + """ + if not m2 > 0 or not Lambda > 0: + raise ValueError('`m2` and `Lambda` are expected to be positive.') + Cphi = C['phi'].real + Ckin = C['phiBox'].real - C['phiD'].real / 4 + if abs(Cphi) < 1e-16: + v2 = 2 * m2 / Lambda + else: + sqrt_arg = Lambda**2 - 12 * Cphi * m2 + if not sqrt_arg >= 0: + raise ValueError("'Lambda**2 - 12 * Cphi * m2' must be positive.") + v2 = ( Lambda - sqrt(sqrt_arg) )/( 3 * Cphi ) + Mh2 = v2 * ( 1 + Ckin * v2 )**2 * ( Lambda - 3 * Cphi * v2 ) + return {'v': sqrt(v2), 'Mh2': Mh2} def vMh2_to_m2Lambda(v, Mh2, C): """Function to numerically determine the parameters of the Higgs potential - given the physical Higgs VEV and mass.""" - if C['phi'] == 0 and C['phiBox'] == 0 and C['phiD'] == 0: - return _vMh2_to_m2Lambda_SM(v, Mh2) - else: - def f0(x): # we want the root of this function - m2, Lambda = x - d = m2Lambda_to_vMh2(m2=m2.real, Lambda=Lambda.real, - C=C) - return np.array([d['v'] - v, d['Mh2'] - Mh2]) - dSM = _vMh2_to_m2Lambda_SM(v, Mh2) - x0 = np.array([dSM['m2'], dSM['Lambda']]) - try: - xres = scipy.optimize.newton_krylov(f0, x0) - except (scipy.optimize.nonlin.NoConvergence, ValueError) as e: - warnings.warn('Standard optimization method did not converge. The GMRES method is used instead.', Warning) - try: - xres = scipy.optimize.newton_krylov(f0, x0, method='gmres', - f_tol=1e-7) - except (scipy.optimize.nonlin.NoConvergence, ValueError) as e: - raise ValueError("No solution for m^2 and Lambda found. This problem can be caused by very large values for one or several Wilson coefficients.") - return {'m2': xres[0], 'Lambda': xres[1]} + given the physical Higgs VEV and mass. + The relations used by this function have been derived from the Lagrangian + in eq. (3.1) of arXiv:1704.03888. This paper presents results expanded up to + linear order in the Wilson coefficients while the all-order expressions are + implemented here.""" + if not v > 0 or not Mh2 > 0: + raise ValueError('`v` and `Mh2` are expected to be positive.') + v2 = v**2 + Cphi = C['phi'].real + Ckin = C['phiBox'].real - C['phiD'].real / 4 + Ckin_factor = ( 1 + Ckin * v2 )**2 + Lambda = ( Mh2 + 3 * Cphi * v2**2 * Ckin_factor ) / ( v2 * Ckin_factor ) + m2 = Lambda * v2 / 2 - 3/4 * Cphi * v2**2 + return {'m2': m2, 'Lambda': Lambda} def get_gpbar(ebar, gbar, v, C): r"""Function to numerically determine the hypercharge gauge coupling @@ -141,11 +140,10 @@ def smpar(C): """Get the running effective SM parameters.""" m2 = C['m2'].real Lambda = C['Lambda'].real - v = (sqrt(2 * m2 / Lambda) + 3 * m2**(3 / 2) / - (sqrt(2) * Lambda**(5 / 2)) * C['phi'].real) + vMh2 = m2Lambda_to_vMh2(m2, Lambda, C) + v = vMh2['v'] + Mh2 = vMh2['Mh2'] GF = 1 / (sqrt(2) * v**2) # TODO - Mh2 = 2 * m2 * (1 - m2 / Lambda * (3 * C['phi'].real - 4 * Lambda * C['phiBox'].real + - Lambda * C['phiD'].real)) eps = C['phiWB'].real * (v**2) gb = (C['g'] / (1 - C['phiW'].real * (v**2))).real gpb = (C['gp'] / (1 - C['phiB'].real * (v**2))).real diff --git a/wilson/run/smeft/tests/test_smpar.py b/wilson/run/smeft/tests/test_smpar.py index 9aab1d6..4badd59 100644 --- a/wilson/run/smeft/tests/test_smpar.py +++ b/wilson/run/smeft/tests/test_smpar.py @@ -1,172 +1,188 @@ -import unittest -from wilson import wcxf -import numpy as np -from wilson.run.smeft import smpar, SMEFT -from math import sqrt, pi -import ckmutil - - -np.random.seed(110) - - -def get_random_wc(eft, basis, scale=160, cmax=1e-2): - """Generate a random Wilson coefficient instance for a given basis.""" - basis_obj = wcxf.Basis[eft, basis] - _wc = {} - for s in basis_obj.sectors.values(): - for name, d in s.items(): - _wc[name] = cmax * np.random.rand() - if 'real' not in d or not d['real']: - _wc[name] += 1j * cmax * np.random.rand() - return wcxf.WC(eft, basis, scale, wcxf.WC.dict2values(_wc)) - - - -class Testgpbar(unittest.TestCase): - def test_sm(self): - gp = smpar.get_gpbar(0.3, 0.6, 246, {'phiWB': 0, 'phiB': 0}) - self.assertEqual(gp, 0.3*0.6/sqrt(-0.3**2+0.6**2)) - - def test_general(self): - C = {'phiWB': np.random.rand() / 200**2, 'phiB': np.random.rand() / 200**2} - v = 246 - gp = smpar.get_gpbar(0.3, 0.6, v, C) - gpb = gp / (1 - C['phiB'] * (v**2)) - gb = 0.6 - eps = C['phiWB'] * (v**2) - eb = (gb * gpb / sqrt(gb**2 + gpb**2) * - (1 - eps * gb * gpb / (gb**2 + gpb**2))) - self.assertAlmostEqual(eb, 0.3) - - -class TestMh2v(unittest.TestCase): - def test_sm(self): - v = 246 - Mh2 = 125**2 - d = smpar._vMh2_to_m2Lambda_SM(v, Mh2) - self.assertAlmostEqual(d['m2'], Mh2/2) - self.assertAlmostEqual(d['Lambda'], Mh2/v**2) - C = {k: 0 for k in ['phi', 'phiBox', 'phiD']} - d = smpar.vMh2_to_m2Lambda(v, Mh2, C) - self.assertAlmostEqual(d['m2'], Mh2/2) - self.assertAlmostEqual(d['Lambda'], Mh2/v**2) - - def test_Cphi0(self): - v = 246 - Mh2 = 125**2 - C = {k: np.random.rand() / 500**2 for k in ['phiBox', 'phiD']} - C['phi'] = 0 - d = smpar.vMh2_to_m2Lambda(v, Mh2, C) - d2 = smpar.m2Lambda_to_vMh2(d['m2'], d['Lambda'], C) - self.assertAlmostEqual(d2['v'], v) - self.assertAlmostEqual(d2['Mh2'], Mh2) - - def test_general(self): - v = 246 - Mh2 = 125**2 - C = {k: np.random.rand() / 500**2 for k in ['phi', 'phiBox', 'phiD']} - d = smpar.vMh2_to_m2Lambda(v, Mh2, C) - d2 = smpar.m2Lambda_to_vMh2(d['m2'], d['Lambda'], C) - self.assertAlmostEqual(d2['v'], v, places=6) - self.assertAlmostEqual(d2['Mh2'], Mh2, places=6) - - def test_noconvergence(self): - v = 246 - Mh2 = 125**2 - C = {k: 0 for k in ['phi', 'phiBox']} - C['phiD'] = 9.678e-06 - self.assertWarns(Warning, smpar.vMh2_to_m2Lambda, v, Mh2, C) - -class TestSMpar(unittest.TestCase): - def test_smeftpar_small(self): - wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-24) - smeft = SMEFT(wc=None) - smeft.scale_in = 160 - smeft._set_initial_wcxf(wc, get_smpar=False) - with self.assertRaises(ValueError): - smpar.smeftpar(smeft.scale_in, smeft.C_in, 'flavio') - CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') - p = smpar.p - self.assertAlmostEqual(CSM['m2'], p['m_h']**2/2) - v = sqrt(1 / sqrt(2) / p['GF']) - self.assertAlmostEqual(CSM['Lambda'], p['m_h']**2/v**2) - self.assertAlmostEqual(CSM['g'], 2*p['m_W']/v) - # self.assertAlmostEqual(CSM['gp']) - self.assertAlmostEqual(CSM['gs'], sqrt(4*pi*p['alpha_s'])) - self.assertAlmostEqual(CSM['Gd'][0, 0], p['m_d']/(v/sqrt(2))) - self.assertAlmostEqual(CSM['Gd'][1, 1], p['m_s']/(v/sqrt(2))) - self.assertAlmostEqual(CSM['Gd'][2, 2], p['m_b']/(v/sqrt(2))) - self.assertAlmostEqual(CSM['Ge'][0, 0], p['m_e']/(v/sqrt(2))) - self.assertAlmostEqual(CSM['Ge'][1, 1], p['m_mu']/(v/sqrt(2))) - self.assertAlmostEqual(CSM['Ge'][2, 2], p['m_tau']/(v/sqrt(2))) - UL, S, UR = ckmutil.diag.msvd(CSM['Gu']) - V = UL.conj().T - self.assertAlmostEqual(S[0], p['m_u']/(v/sqrt(2))) - self.assertAlmostEqual(S[1], p['m_c']/(v/sqrt(2))) - self.assertAlmostEqual(S[2], p['m_t']/(v/sqrt(2))) - self.assertAlmostEqual(abs(V[0, 2]), p['Vub']) - self.assertAlmostEqual(abs(V[1, 2]), p['Vcb']) - self.assertAlmostEqual(abs(V[0, 1]), p['Vus']) - - def test_smpar_small(self): - wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-24) - smeft = SMEFT(wc=None) - smeft.scale_in = 160 - smeft._set_initial_wcxf(wc, get_smpar=False) - CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') - Cboth = CSM.copy() - Cboth.update(smeft.C_in) - Cback = smpar.smpar(Cboth) - for k in smpar.p: - if k not in ['m_Z', 'delta']: - self.assertAlmostEqual(smpar.p[k], Cback[k], - msg="Failed for {}".format(k)) - - - def test_smpar_roundtrip(self): - wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-6) - smeft = SMEFT(wc=None) - smeft.scale_in = 160 - smeft._set_initial_wcxf(wc, get_smpar=False) - CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') - Cboth = CSM.copy() - Cboth.update(smeft.C_in) - Cback = smpar.smpar(Cboth) - for k in smpar.p: - if k in ['m_Z']: - self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, - msg="Failed for {}".format(k), - delta=0.05) - elif k in ['delta']: - self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, - msg="Failed for {}".format(k), - delta=1e-3) - else: - self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, - msg="Failed for {}".format(k), - delta=1e-5) - -class TestGetSMpar(unittest.TestCase): - def test_wcxf_smpar(self): - wc = get_random_wc('SMEFT', 'Warsaw', 1e5, 1e-11) - smeft = SMEFT(wc) - C_out = smeft._rgevolve(91.1876) - p_out = smpar.smpar(C_out) - for k in p_out: - self.assertAlmostEqual(p_out[k] / smpar.p[k], 1, - delta=0.05, - msg="Failed for {}".format(k)) - - def test_wcxf_smpar_incomplete(self): - wc = wcxf.WC('SMEFT', 'Warsaw', 160, {'qd1_1111': {'Im': 1e-6}}) - smeft = SMEFT(wc) - smeft.run(91.1876) - - def test_method(self): - wc = wcxf.WC('SMEFT', 'Warsaw', 1000, {'uG_33': 1e-6}) - smeft = SMEFT(wc) - p_out = smeft.get_smpar() - p_def = smpar.p - for k, v in p_out.items(): - self.assertAlmostEqual(v / p_def[k], 1, places=1, - msg="Failed for {}".format(k)) +import unittest +from wilson import wcxf +import numpy as np +from wilson.run.smeft import smpar, SMEFT +from math import sqrt, pi +import ckmutil + + +np.random.seed(110) + + +def get_random_wc(eft, basis, scale=160, cmax=1e-2): + """Generate a random Wilson coefficient instance for a given basis.""" + basis_obj = wcxf.Basis[eft, basis] + _wc = {} + for s in basis_obj.sectors.values(): + for name, d in s.items(): + _wc[name] = cmax * np.random.rand() + if 'real' not in d or not d['real']: + _wc[name] += 1j * cmax * np.random.rand() + return wcxf.WC(eft, basis, scale, wcxf.WC.dict2values(_wc)) + + + +class Testgpbar(unittest.TestCase): + def test_sm(self): + gp = smpar.get_gpbar(0.3, 0.6, 246, {'phiWB': 0, 'phiB': 0}) + self.assertEqual(gp, 0.3*0.6/sqrt(-0.3**2+0.6**2)) + + def test_general(self): + C = {'phiWB': np.random.rand() / 200**2, 'phiB': np.random.rand() / 200**2} + v = 246 + gp = smpar.get_gpbar(0.3, 0.6, v, C) + gpb = gp / (1 - C['phiB'] * (v**2)) + gb = 0.6 + eps = C['phiWB'] * (v**2) + eb = (gb * gpb / sqrt(gb**2 + gpb**2) * + (1 - eps * gb * gpb / (gb**2 + gpb**2))) + self.assertAlmostEqual(eb, 0.3) + + +class TestMh2v(unittest.TestCase): + def test_sm(self): + v = 246 + Mh2 = 125**2 + C = {k: 0 for k in ['phi', 'phiBox', 'phiD']} + d = smpar.vMh2_to_m2Lambda(v, Mh2, C) + d2 = smpar.m2Lambda_to_vMh2(d['m2'], d['Lambda'], C) + self.assertAlmostEqual(d['m2'], Mh2/2) + self.assertAlmostEqual(d['Lambda'], Mh2/v**2) + self.assertAlmostEqual(d2['v'], v) + self.assertAlmostEqual(d2['Mh2'], Mh2) + + def test_Cphi0(self): + v = 246 + Mh2 = 125**2 + C = {k: np.random.rand() / 500**2 for k in ['phiBox', 'phiD']} + C['phi'] = 0 + d = smpar.vMh2_to_m2Lambda(v, Mh2, C) + d2 = smpar.m2Lambda_to_vMh2(d['m2'], d['Lambda'], C) + self.assertAlmostEqual(d2['v'], v) + self.assertAlmostEqual(d2['Mh2'], Mh2) + + def test_general(self): + v = 246 + Mh2 = 125**2 + C = {k: np.random.rand() / 500**2 for k in ['phi', 'phiBox', 'phiD']} + d = smpar.vMh2_to_m2Lambda(v, Mh2, C) + d2 = smpar.m2Lambda_to_vMh2(d['m2'], d['Lambda'], C) + self.assertAlmostEqual(d2['v'], v, places=6) + self.assertAlmostEqual(d2['Mh2'], Mh2, places=6) + + def test_vreal(self): + m2 = 7812 + Lambda = 0.258 + C = {k: 0 for k in ['phiD', 'phiBox']} + C['phi'] = 1e-6 + with self.assertRaises(ValueError): + smpar.m2Lambda_to_vMh2(m2, Lambda, C) + + def test_args_positive(self): + v = 246 + Mh2 = 125**2 + m2 = 7812 + Lambda = 0.258 + C = {k: 0 for k in ['phi', 'phiBox', 'phiD']} + with self.assertRaises(ValueError): + smpar.vMh2_to_m2Lambda(-v, Mh2, C) + with self.assertRaises(ValueError): + smpar.vMh2_to_m2Lambda(v, -Mh2, C) + with self.assertRaises(ValueError): + smpar.m2Lambda_to_vMh2(-m2, Lambda, C) + with self.assertRaises(ValueError): + smpar.m2Lambda_to_vMh2(m2, -Lambda, C) + +class TestSMpar(unittest.TestCase): + def test_smeftpar_small(self): + wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-24) + smeft = SMEFT(wc=None) + smeft.scale_in = 160 + smeft._set_initial_wcxf(wc, get_smpar=False) + with self.assertRaises(ValueError): + smpar.smeftpar(smeft.scale_in, smeft.C_in, 'flavio') + CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') + p = smpar.p + self.assertAlmostEqual(CSM['m2'], p['m_h']**2/2) + v = sqrt(1 / sqrt(2) / p['GF']) + self.assertAlmostEqual(CSM['Lambda'], p['m_h']**2/v**2) + self.assertAlmostEqual(CSM['g'], 2*p['m_W']/v) + # self.assertAlmostEqual(CSM['gp']) + self.assertAlmostEqual(CSM['gs'], sqrt(4*pi*p['alpha_s'])) + self.assertAlmostEqual(CSM['Gd'][0, 0], p['m_d']/(v/sqrt(2))) + self.assertAlmostEqual(CSM['Gd'][1, 1], p['m_s']/(v/sqrt(2))) + self.assertAlmostEqual(CSM['Gd'][2, 2], p['m_b']/(v/sqrt(2))) + self.assertAlmostEqual(CSM['Ge'][0, 0], p['m_e']/(v/sqrt(2))) + self.assertAlmostEqual(CSM['Ge'][1, 1], p['m_mu']/(v/sqrt(2))) + self.assertAlmostEqual(CSM['Ge'][2, 2], p['m_tau']/(v/sqrt(2))) + UL, S, UR = ckmutil.diag.msvd(CSM['Gu']) + V = UL.conj().T + self.assertAlmostEqual(S[0], p['m_u']/(v/sqrt(2))) + self.assertAlmostEqual(S[1], p['m_c']/(v/sqrt(2))) + self.assertAlmostEqual(S[2], p['m_t']/(v/sqrt(2))) + self.assertAlmostEqual(abs(V[0, 2]), p['Vub']) + self.assertAlmostEqual(abs(V[1, 2]), p['Vcb']) + self.assertAlmostEqual(abs(V[0, 1]), p['Vus']) + + def test_smpar_small(self): + wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-24) + smeft = SMEFT(wc=None) + smeft.scale_in = 160 + smeft._set_initial_wcxf(wc, get_smpar=False) + CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') + Cboth = CSM.copy() + Cboth.update(smeft.C_in) + Cback = smpar.smpar(Cboth) + for k in smpar.p: + if k not in ['m_Z', 'delta']: + self.assertAlmostEqual(smpar.p[k], Cback[k], + msg="Failed for {}".format(k)) + + + def test_smpar_roundtrip(self): + wc = get_random_wc('SMEFT', 'Warsaw', cmax=1e-6) + smeft = SMEFT(wc=None) + smeft.scale_in = 160 + smeft._set_initial_wcxf(wc, get_smpar=False) + CSM = smpar.smeftpar(smeft.scale_in, smeft.C_in, 'Warsaw') + Cboth = CSM.copy() + Cboth.update(smeft.C_in) + Cback = smpar.smpar(Cboth) + for k in smpar.p: + if k in ['m_Z']: + self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, + msg="Failed for {}".format(k), + delta=0.05) + elif k in ['delta']: + self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, + msg="Failed for {}".format(k), + delta=1e-3) + else: + self.assertAlmostEqual(smpar.p[k]/Cback[k], 1, + msg="Failed for {}".format(k), + delta=1e-5) + +class TestGetSMpar(unittest.TestCase): + def test_wcxf_smpar(self): + wc = get_random_wc('SMEFT', 'Warsaw', 1e5, 1e-11) + smeft = SMEFT(wc) + C_out = smeft._rgevolve(91.1876) + p_out = smpar.smpar(C_out) + for k in p_out: + self.assertAlmostEqual(p_out[k] / smpar.p[k], 1, + delta=0.05, + msg="Failed for {}".format(k)) + + def test_wcxf_smpar_incomplete(self): + wc = wcxf.WC('SMEFT', 'Warsaw', 160, {'qd1_1111': {'Im': 1e-6}}) + smeft = SMEFT(wc) + smeft.run(91.1876) + + def test_method(self): + wc = wcxf.WC('SMEFT', 'Warsaw', 1000, {'uG_33': 1e-6}) + smeft = SMEFT(wc) + p_out = smeft.get_smpar() + p_def = smpar.p + for k, v in p_out.items(): + self.assertAlmostEqual(v / p_def[k], 1, places=1, + msg="Failed for {}".format(k))