Issue with the running of WCs in 2HDM+VLQs

[ajueid]

Dear authors of Wilson,

We are using wilson+flavio for a global fit of the 2HDM+VLQs. The WCs in the Warsaw basis were calculated using MatchMakerEFT and then translated to python format. However, for many points we got several issues with some WCs that lead to the abort of the runs. The issue is summarized below (an example of a problematic benchmark point).

Parameters for general potential in generic basis

3        1.00000000e+00   # tan(beta)

18 2.00000000e+02 # m_12
20 9.99000000e-01 # sin (beta - alpha)
34 0.95000000e+00 # CTLu indpe
35 0.95000000e+00 # CTLd
24 2 # Yukawas Type
Block MASS # Mass spectrum (kinematic masses)

PDG Mass

 24        8.03980000e+01   # MW
 25        1.25090000e+02   # Mh1, lightest CP-even higgs
 35        5.50000000e+02   # Mh2, heaviest CP-even higgs
 36        5.50000000e+02   # Mh3, CP-odd higgs
 37        4.00000000e+02   # Mhc

6000006 1.50000000e+03 # MU1
6000007 1.50000000e+03 # MD1

Some Wilson coef:

OH 2.5047355741326885e-05 (max |WC| value)
OHD -1.5580641729755483e-07
OHWB -6.314669514751517e-10
OHW -3.525591158856639e-10
OHB -1.1683500926778163e-09
OHbox -2.41631797138483e-06

[FLAVPREC] compute_obs failed: No solution for gp found. This problem can be caused by very large values for one or several Wilson coefficients.

Traceback (most recent call last):
File "/lustre/ific.uv.es/prj/gl/som.mosayco/rruiz/software/miniconda3/envs/2hdm/lib/python3.8/site-packages/wilson/run/smeft/smpar.py", line 99, in get_gpbar
gpbar = scipy.optimize.brentq(f0, 0, 3)
File "/lustre/ific.uv.es/prj/gl/som.mosayco/rruiz/software/miniconda3/envs/2hdm/lib/python3.8/site-packages/scipy/optimize/zeros.py", line 776, in brentq
r = _zeros._brentq(f, a, b, xtol, rtol, maxiter, args, full_output, disp)
ValueError: f(a) and f(b) must have different signs

Some numbers from the failing function:

phiWB = -6.314669514751517e-10
phiB = -1.1683500926778163e-09
v = 246.21965079413738 ebar = 0.31345100004952897 gbar = 0.651450846764905
f0(0) = -0.31345100004952897
f0(3) = 0.32316825751986145
phiWB = (15049.805902382412+0.00012567165595918193j)
phiB = (83536.16327705246+0.0006975590276445177j)
v = 246.21965079413738 ebar = 0.31345100004952897 gbar = 0.651450846764905
f0(0) = -0.31345100004952897
f0(3) = -120448774.1479608

I attach also the wcs.py (a module containing the expressions of the WCs in our model). Please let me know if you require any info that helps to identify the problem.

Thanks a lot!

All the best,
Adil

flavour.py
wcs.py

[peterstangl]

Hi @ajueid,

it doesn't look like there is any issue with wilson, which seems to do exactly what it should. The error message you get tells you what's the problem:

No solution for gp found. This problem can be caused by very large values for one or several Wilson coefficients.

The extraction of the hypercharge gauge coupling gp depends on the values of the Wilson coefficient phiWB and phiB. If your phiWB and phiB have values of 15000/GeV^2 and 80000/GeV^2, you are far outside the region of parameter space in which the SMEFT can be used in any meaningful way. You get an error because the extraction of gp fails for such large values of Wilson coefficients. But the problem is actually more severe than that since you should not use the SMEFT at all if you get such large values from the matching. This problem has nothing to do with wilson and is either caused by your choice of UV model parameters (which might be outside the SMEFT validity region) or by a potential bug in your matching implementation.

[ajueid]

Hi @peterstangl,

Thanks a lot for your reply. Apologies for my slow reply. I think we now have an idea about why these coefficients explode. Some of them depend on differences of the squared masses (Bare masses of VLQ singlets and the mass of the second doublet in the Higgs basis) in the denominator. If they come close to each other then the corresponding WCs explode. We will try to find a way to regularize them or probably apply a requirement on the spectrum prior to the calculation of WCs. Have you ever encountered issues with models containing multiple mass scales?

Thanks again!

All the best,
Adil

[peterstangl]

Hi @ajueid,

OK, this makes sense. If you do one-loop matching, the loop functions depending on e.g. two different masses m1 and m2 have denominators containing the difference of the masses squared. But the limit m1->m2 actually exists and is finite. In this limit, the loop function just reduces to a simpler loop function depending on a single mass. So the problem you have is entirely numerical since you take the ratio of two very small numbers instead of taking the actual limit. So if both masses are exactly equal, you can easily do your analysis by using the simpler one-mass loop functions. If the mass difference is small but non-zero, you can perform a Taylor expansion in the small mass difference.

So in general, if you implement loop functions depending on several different masses, you have to carefully treat the limit in which the masses are close to each other (or even equal) to avoid numerical issues from taking the ratio of two small numbers.