Anomalous quartic gauge couplings and unitarization for the vector boson scattering process $pp\rightarrow W^+W^+jjX\rightarrow \ell^+\nu_\ell\ell^+\nu_\ell jjX$

Weak vector boson scattering (VBS) at the LHC provides an excellent source of information on the structure of quartic gauge couplings and possible effects of physics beyond the SM in electroweak symmetry breaking. Parameterizing deviations from the SM within an effective field theory at tree level, the dimension-8 operators, which are needed for sufficiently general modeling, lead to unphysical enhancements of cross sections within the accessible energy range of the LHC. Preservation of unitarity limits is needed for phenomenological studies of the $VVjj$ events which signify VBS. Here we develop a numerical unitarization scheme for the full off-shell VBS processes and apply it to same-sign $W$ scattering, i.e. processes like $qq\to qqW^+W^+$. The scheme is implemented within the Monte Carlo program VBFNLO, including leptonic decay of the weak bosons and NLO QCD corrections. Distributions differentiating between higher dimensional operators are discussed.


Introduction
Among the scattering processes which can be studied at the CERN Large Hadron Collider (LHC), weak vector boson scattering (VBS) is particularly interesting as a probe of electroweak symmetry breaking. Within the Standard Model (SM), intricate cancellations between Feynman amplitudes involving quartic gauge boson interactions, trilinear gauge boson couplings, and Higgs exchange lead to scattering amplitudes for longitudinally polarized weak bosons which do not grow with energy and which, for a light Higgs boson, respect bounds derived from unitarity. Modifications of the weak boson couplings, among themselves or to the Higgs boson, spoil these cancellations and can lead to sizable cross section increases. For example, reduced weak boson couplings to the light, m h = 125 GeV Higgs boson and compensation by an additional heavy Higgs in a two-Higgs-doublet model would lead to a cross section increase at high energy, as would a change only in the quartic gauge couplings. Absent clear hints for a particular theory beyond the Standard Model (BSM), a bottom up approach is conveniently formulated within an effective field theory (EFT) approach [1,2]. Given the observation of a light Higgs boson at the LHC [3,4], we opt for a linear representation of the light fields in order to construct dimension-six and -eight operators for the EFT. To give just one example, a deviation in the Higgs sector could manifest itself via the dimension-8 term in the effective Lagrangian, which is present in the linearÉboli-basis [5]. Here, the covariant derivative of the Higgs-doublet field, D µ φ, contains W and Z fields, Λ, is the energy scale of new physics, and the coupling coefficient f S 1 is used later to allow different strengths for independent dimension-8 operators. This operator will induce an anomalous contribution to the four-W , four-Z, and W W ZZ vertices, which alter the scattering (predominantly) of the longitudinal degrees of freedom of the weak vector bosons. The impact of anomalous couplings on V V → V V scattering can be studied at the LHC via the full process pp → V V jjX as illustrated in Fig. 1, where the two final state vector bosons V can decay either leptonically or hadronically. The current, observed limits for f S 1 /Λ 4 , derived from same-sign W scattering by CMS, are [−21.6, 21.8] TeV −4 for 35.9 fb −1 of √ s = 13 TeV data [6]. We are not aware of new results for Run-II published by ATLAS. However, comparing old limits for f S 1 /Λ −4 of [−118, 120] TeV −4 from 19.4 fb −1 of √ s = 8 TeV data from CMS [7] and [−960, 960] TeV −4 from 20.3 fb −1 of √ s = 8 TeV data from ATLAS [8], one observes a substantial difference in precision. 1 This difference is mainly due to the different high-energy extrapolation of the EFT ansatz in the generation of BSM Monte-Carlo events. The EFT is only valid up to a certain energy scale Λ valid < Λ, where the operator product expansion breaks down. However, the experiment is only sensitive to the ratio f S 1 /Λ 4 and the scale Λ is a priori not known. Using just the EFT as input for the generation of Monte Carlo data will usually overshoot any result allowed by perturbative unitarity in the high energy region. Naturally this will result in more stringent limits for the EFT-coefficients. CMS is using this approach in presenting their limits. ATLAS on the other hand is using the T-matrix [12,13] unitarization scheme to provide a theoretically consistent description of the high energy region, where unitarity would otherwise be violated, with a proper interpolation to the low energy EFT. T-matrix unitarization leads to lower generated event rates for a given f S 1 /Λ 4 , which leads to weaker limits for this Wilson coefficient.
The effect is demonstrated in Fig. 2, where we compare the W + W + invariant mass distribution expected with the current CMS limit F S 1 = f S 1 /Λ 4 = 21.8 TeV −4 for a naive EFT description (dashed pink), with a T-matrix unitarized [12] prescription, with larger coupling F S 1 = 60 TeV −4 , which will give the same fiducial cross section (solid orange histogram). Also shown are the SM expectation (solid black) and the naive EFT expectation for the larger coupling of 60 TeV −4 (dashed brown), which agrees with the T-matrix unitarized expectation only at small invariant masses. In the energy range above M W W ≈ 1500 GeV, an L S 1 induced excess above the orange F S 1 = 60 TeV −4 curve violates unitarity, i.e. it is unphysical, and should therefore not be considered to estimate EFT coefficients. This is quite general: a pure effective Lagrangian/anomalous coupling analysis of LHC observables, with a finite set of terms in the effective Lagrangian, is insufficient in practice because the unbounded growth of amplitudes with energy typically corresponds to unitarity violation within the energy reach of the LHC. We thus need a general and versatile unitarization procedure for the naive EFT amplitudes at high momentum transfers, which smoothly interpolates to the pure EFT description well below Λ valid . In order to analyze V V jj production data, any unitarity considerations for V V → V V scattering must be extrapolated from on-shell bosons to the space-like incoming and time-like outgoing virtualities of the vector bosons which is implicit in the kinematics of Fig. 1. As we shall see, this extrapolation will require a few additional assumptions and will induce some model dependence. To obtain predictions which are compatible with unitarity, the T-matrix unitarization prescription can be used. So far, however, an implementation of this scheme is only available for a small number of effective Lagrangian operators for VBS due to the difficulty to handle VBS with arbitrarily polarized off-shell vector bosons in the full V V jj production process [12,13]. In this paper, we introduce a variant of the T/K-matrix unitarization scheme [12,14], called T u unitarization below, for general combinations of operators within VBS, for arbitrary space-like virtualities of the incoming vector bosons, and describe its implementation in the Monte Carlo generator VBFNLO [15,16,17].
The paper is organized as follows. We introduce the full set of bosonic dimension-8 operators with a list of current experimental limits in Section 2. In Section 3, we first consider how unitarity relations can be extended to off-shell VBS processes. Beyond the definition of off-shell polarization vectors, this entails partial wave decomposition for off shell sub-amplitudes for V V → V V scattering and its fast numerical implementation. The T u unitarization model, which we have implemented in VBFNLO for same-sign W scattering, is introduced in Section 3.3. Section 4 is devoted to numerical results for same-sign W -scattering, i.e the process pp → W + W + jjX → + ν + ν jjX at NLO QCD precision, which is now implemented in VBFNLO including T u unitarization for any combination of the dimension-8 operators listed in Section 2. We compare our T u unitarized model with naive EFT descriptions for different dimension-8 operators. Furthermore, we will also give examples of observables helping to distinguish experimentally between different subclasses of dimension-8 operators. Final conclusions are drawn in Section 5.

Effective Field Theory Description of Anomalous Quartic Gauge Couplings
The bottom-up EFT framework is useful to quantify deviations from the SM in a model independent way and, once experimental evidence for such deviations is discovered, it gives hints, from which BSM effect a possible anomaly might originate. Short of such a desirable situation, experimental limits on the Wilson coefficients serve as a measure of the experimental precision. Two EFT representations are mainly used to describe BSM contributions for anomalous quartic gauge couplings, the linear and non-linear representation. They can be distinguished due to the different ordering of the EFT expansion which is written in terms of operators O i of energy dimension d i , corresponding Wilson coefficients f i (which allow for variations in importance of the individual operators) and energy scale of new physics, Λ. In the non-linear representation, the Higgs couplings are treated as additional free parameters and deviations in the Higgs sector can already be introduced at lowest order [18]. This was well motivated before the experimental Higgs discovery in case of a heavy or strongly interacting Higgs [19,20]. However, since no deviations from the light SM Higgs predictions have been observed so far, we choose the linear Higgs representation, where deviations from the SM predictions for Higgs couplings and trilinear or quartic gauge boson first appear at energy dimension d i = 6 [21,22,23]. Anomalous quartic gauge couplings (aQGC) are induced at the dimension-6 level already. However, they are not independent of changes in the Higgs couplings or of anomalous trilinear gauge couplings. These three-boson couplings are most easily measured in Higgs production or decay or in vector boson pair production (qq → V 1 V 2 ), at the LHC, and little additional information is to be expected from the measurement of the significantly smaller VBS cross sections. Also, the tensor structure of dimension-6 operators is not general enough to allow for sufficiently uncorrelated variations of the 81 helicity amplitudes which, in principle, can be probed in a VBS process, V 1 V 2 → V 3 V 4 , with massive vector bosons. In this paper, we study aQGC which enter the EFT at lowest order at dimension-8 without contributing to anomalous trilinear gauge interactions or to HV V couplings.
The contributing CP conserving operators can be assembled from three SM building blocks. One building block is the covariant derivative acting on the Higgs doublet field, which affects the coupling of longitudinal modes of the gauge bosons. Here, the Higgs, H is embedded in the Higgs doublet field in the unitary gauge: The other building blocks are the field strength tensors which are normalized such that [D µ , D ν ] = W µν + B µν for the covariant derivative in Eq. (3). The abelian parts of these field strength tensors lead to couplings of the transverse degrees of freedom of the gauge fields. The dimension-8 operators are separated into longitudinal, transverse, and mixed contributions, corresponding to the occurrence of the building blocks above. A revised list of dimension-8 operators from [5] and [24] is given in Eqs. (6,7,8). In comparison to the operators defined in [5], we choose a different normalization for the field strength in Eq. (5), which is accompanied by an additional factor of ig or ig /2. These normalization choices are labeled as "Éboli" for [5] and the normalization in Eq. (5) as "VBFNLO", in the following.
For the longitudinal operators the two normalization choices coincide: and has to be included to complete the operator set. Finally, the purely transverse operators are Same-sign W boson scattering is the VBS process which can be measured with the highest precision, due to a sizable signal cross section and a particularly low QCD background [6,7,8]. In Section 4 we will concentrate on this process, to which only operators with exactly four W ± fields can contribute at tree level. This eliminates all operators with a hypercharge field strength, B µν . The remaining ones have been probed by ATLAS [8] and CMS [6,7] in same-sign W scattering, and the results are summarized in Table 1 In comparing the different normalization ofÉboli and VBFNLO, one finds that the limits for VBFNLO differ by about one order of magnitude only, whereas the limits in theÉboli normalization vary by up to two orders of magnitude. The difference is simply due to consistently factorizing the small electroweak couplings, which are expected for any model explaining the EFT, into the definition of the operators for the VBFNLO normalization convention. For the 8 TeV data, the ATLAS result incorporates a unitarization model to prevent the generation of unphysical events at high energy, which violate unitarity constraints. The corresponding bound on f S 1 /Λ 4 for the same-sign W scattering process observed by ATLAS [8] is approximately one order of magnitude weaker than the CMS √ s = 8 TeV limit, which indicates the impact that unitarization can have on quoted experimental results.

Unitarity for VBS: going off-shell
We need to apply unitarity considerations to electroweak processes of the type pp →ψ 1 ψ 2ψ3 ψ 4 jj at O(α 6 ) (LO) and at O(α 6 α s ) (NLO), i.e. including QCD corrections. At the parton level, the ψ i represent decay leptons of two vector bosons, the initial pp state represents the scattering partons (quarks or anti-quarks in the LO case) and jj stands for the final state partons yielding two tagging jets. Representative Feynman graphs for the 8-fermion processes at LO are given in Fig. 3 and include vector boson emissions off quark lines as in Fig. 3a as well as VBS contributions as in Fig. 3b. The BSM physics, which we consider via the introduction of bosonic operators, will only contribute to the VBS subprocess V V → V V . The SM contributions to the complete process are gauge invariant by themselves, they are "small" and they respect perturbative unitarity. Splitting the full amplitude into the SM and a BSM piece, it is, therefore, sufficient to unitarize the BSM piece only, via the VBS subprocess, which means that we neglect the interference of SM and BSM amplitudes for unitarization.

Identification of the V V → V V Subamplitude
Within the VBFNLO approach, the entire V V → V V subprocess is contained inside a leptonic tensor, which then is contracted with quark currents, J µ p→jV . These represent the emission of a virtual vector boson, V , off an initial parton in Fig. 3b. This structure is well suited to implement the BSM amplitude as where M BSM V 1 V 2 →4f µν denotes the BSM contribution to the leptonic tensor. The quark currents, J µ p→jV , and the decay currents, J µ V →f f , are conserved, since we are neglecting fermion masses, and this allows for a simple expansion of the off-shell vector boson propagators in terms of polarization vectors of fixed helicity. When writing

vector boson propagators may be taken as
since the q µ terms are contracted with a conserved current and, thus, vanish. The indices J and M on the polarization vectors distinguish between those that are contracted with the currents, J , and the polarization vectors contracted with the VBS matrix element, M . Furthermore, we generalize the definition of the polarization vectors for off-shell vector bosons with four-momentum q to where N M and N J have to fulfill N J N M = 1/q 2 as normalization factors for the longitudinal polarization vectors. One could choose the individual factors to be equal in magnitude. However, in order to match the proper normalization for on-shell weak bosons and thus to reproduce the correct normalization of V V → V V scattering amplitudes for longitudinal V , we set With these definitions, the BSM contribution to the leptonic tensor becomes The full anomalous VBS information is contained in the helicity amplitudes For on-shell vector boson momenta, they correspond to the normal V V → V V helicity amplitudes induced by the dimension-8 operators. For the full pp → 4f jj process, however, we are dealing with incoming space-like four-momenta q 1 and q 2 , and outgoing time-like four-momenta q 3 and q 4 , which means that the initial and final states of even the elastic W + W + → W + W + process do not properly match.
In the presence of dimension-8 operators, the tree level VBS amplitudes M V BS can rise with the fourth power of the center-of-mass energy, √ s. For example, the BSM part of the helicity amplitude M V BS +−;+− , for the operator O T 0 , is given by where the time-like momenta are approximated as on-shell, q 2 3 = q 2 4 ≈ m 2 W . This example also shows that unphysical, strong enhancements are possible for large s and for large q 2 i , independently.
In order to avoid unphysical behaviour, within the energy range probed by the LHC, the subprocess amplitudes M BSM pp→4f jj of Eq. (10) need to be replaced by unitarized versions, for Wilson coefficients of practical interest. Since we intend to describe BSM interactions of the known SM bosons, the unitarization has to act at the level of V V → V V scattering instead of the full pp → 4f jj subprocess: working at the latter level, unitarity alone would e.g. allow replacement of the well-known narrow Breit Wigner propagator for the W by a broad spectral function, keeping the leptons produced in W → ν on top of resonance for virtuality ranges of hundreds of GeV, which would also result in large cross section increases. Our choice of the physics which we want to describe, forces us to match the offshell VBS amplitudes M V BS to unitarized on-shell V V → V V scattering amplitudes, which can be defined from first principles. The guiding principles here are • In the on-shell limit, q 2 i → m 2 V i , the unitarized off-shell amplitude must reduce to the corresponding unitarized on-shell amplitude.
• For large virtualities, q 2 i , and modest s (which is allowed for the incoming space-like bosons) the unitarized off-shell amplitude should not exceed the corresponding on-shell unitarity bound.
• The unitarization procedure must reduce to the EFT limit when the absolute values of all invariants (q 2 i , and the Mandelstam variables, s, t and u) are small compared to Λ 2 , which sets the new physics scale.
These principles must now be applied to the unitarization of the off-shell VBS amplitudes.

Unitarity Relation for 2 → 2 Amplitudes
It is useful to briefly recall the derivation of unitarity relations as exposed, for example, in Ref. [25]. Starting point of any unitarization procedure is the unitarity of the scattering matrix S where, exploiting momentum conservation, the i → f matrix elements are given by Truncating the sum over intermediate states to the two-boson subspace, the elements of the 2 → 2 scattering matrix, T f i (q 3 , q 4 ← q 1 , q 2 ), have to fulfill the condition where S n is the statistical factor for identical particles, S n = 1 2 for the W + W + case to be concentrated on later. Exploiting angular momentum conservation, every 2 → 2 helicity amplitude M λ 3 λ 4 ←λ 1 λ 2 = T f i can be expanded in corresponding partial wave amplitudes A j where d j λ 12 λ 34 denotes a Wigner d-function, λ ij = λ i − λ j , and N f i = N f i (q 3 , q 4 ; q 1 , q 2 ) is a normalization factor. Note that for the dimension-8 operators described in Section 2, only partial waves up to j max = 2 contribute. 3 Performing the angular integral in Eq. (22), the partial wave amplitudes A j are found to satisfy the relation 3 Since at most three partial waves contribute, knowledge of the helicity amplitude M V BS λ3λ4←λ1λ2 (Θ) at three angles is sufficient to determine all partial wave amplitude A j λ3λ4←λ1λ2 for a given set of helicities. We have implemented this procedure in VBFNLO.
where we have separated the sum over intermediate states into an explicit helicity sum and a sum over n, which corresponds to a sum over possible boson flavor combinations. Choosing N ni = s λ 1/4 (s, q 2 n,3 , q 2 n,4 )λ 1/4 (s, q 2 1 , q 2 2 ) with Källén function and analogously for N nf and N f i , the phase-space factor in Eq. (24) is canceled, resulting in a form analogous to Eq. (19). Note that for the case at hand, W ± W ± → W ± W ± scattering, the statistical factors are all equal, S i = S f = S n = 1/2. However, the above description readily generalizes to more complex cases like W + W − → W + W − , ZZ, Zγ, HH etc.. Diagonalizing the partial wave helicity amplitudes, the eigenvalues a j (s) will lie on an Argand circle of radius unity, which implies |Re(a j (s))| ≤ 1 .
We will refer to this limit as the unitarity bound on the scattering amplitude. Alternatively one could use |a j (s)| ≤ 2, which is reached for a purely imaginary scattering amplitude. This comparison shows that the precise place at which a (real) tree level amplitude violates unitarity is somewhat ambiguous. However, a polynomial growth with energy, as implied by a truncated EFT, is clearly forbidden by the unitarity relation of Eq. (24). For on-shell W ± W ± → W ± W ± scattering the specification of virtualities in the above equations is superfluous, of course. However, we want to extend the formalism to the unitarization of the off-shell amplitudes M V BS of Eq. (15), with space-like momenta q 1 and q 2 and time-like momenta q 3 and q 4 which are somewhat off the W Breit-Wigner peak. 4 For this general case, Eq. (23) together with the normalization factor of Eq. (25) defines the partial wave amplitudes to be used below.
Allowing free virtualities of the external particles leads to a new problem, however: already at tree level the scattering amplitudes A j λ 3 λ 4 ←λ 1 λ 2 no longer form normal matrices, i.e. TT † = T † T when states with different virtualities are identified, i.e. when they are associated with a single on-shell state. While the mismatch becomes sub-dominant for virtualities much smaller than the center of mass energy, i.e. for q 2 i /s 1, we here need an interpolation which also works for modest center of mass energies and virtualities, reproducing the EFT results, and which allows us to take the exact, off-shell helicity amplitudes M V BS as input for the unitarization in all regions of phase space. The proposed generalization will be described in the next section.

Implementation of Unitarization: the T u Model
Using a truncated EFT model at tree level for large energy scales will violate unitarity above a certain energy. For current experimental limits on the EFT coefficients, this unphysical behavior happens within the energy reach of the LHC, as demonstrated in Fig. 2. Therefore, an extended model must be used to ensure that generated differential cross sections are not becoming unphysically large. Several procedures, with different high energy behavior, are available to extrapolate the EFT beyond its validity range. One possibility, which has been used in VBFNLO in the past, is the introduction of (somewhat ad hoc) form-factors which multiply the full BSM amplitude M VBS of Eq. (15) to ensure the unitary bound of Eq. (27).
Theoretically more attractive is the substitution of the tree level amplitudes by versions, which, at least approximately, satisfy the unitarity condition of Eq. (24). One such procedure is the linear T-matrix projection for the intermediate 2 → 2 interaction matrix that is introduced in [12,13]. With this projection, the 2 → 2 scattering amplitudes will approach the perturbative unitarity bound at high energies and are matched to the naive EFT at low energies. Given the starting point of a normal 5 tree level interaction matrix T 0 , the procedure corresponds to the substitution of T 0 by The T-matrix unitarization model has been implemented for the O S operators [12,13,27], which enhance the scattering of longitudinal vector bosons. In these implementations, an analytical approach has been chosen to provide T-matrix unitarized results at high center of mass energies. The next step in this program is the expansion of the method for operator classes O M and O T , i.e. the implementation for additional helicity combinations of the vector bosons [28].
Contrary to the analytical ansatz chosen in [12,13,27], which requires approximations which become exact only in the limit of s m 2 V , |q 2 i |, we here opt for a numerical approach, which gives us greater versatility for the additional dimension-8 operator classes, allowing investigations of arbitrary regions of phase space. As mentioned in the introduction, we here limit ourselves to the doubly charged channels, i.e. to scattering of two same-sign W ± bosons.
In the case of on-shell scattering, the interaction matrix becomes hermitian, at tree level, and we can expand the denominator in Eq. (28), to improve numerical stability, as As mentioned in the last section, the interaction matrix of the 2 → 2 vector boson scattering subprocess, within the process pp → W ± W ± jj, is not normal, because the momenta of the incoming vector bosons q 1 , q 2 are space-like and the momentum of the outgoing vector bosons q 3 , q 4 are time-like and almost on-shell. Although an extended procedure for non-normal interaction matrices is provided in [12], it is not feasible for a numerical approach.
To generalize Eq. (29) for off-shell sub-amplitudes M VBS , we distinguish states with time-like and space-like bosons as separate classes, labeling the corresponding matrix elements with s for space-like and t for time-like momenta. This leads us to consider three cases for the partial wave amplitudes A j defined in Eq. (23), which correspond to the amplitudes of the actual physical subprocess, with time-like final momenta and space-like initial momenta, its hermitian adjoint, and an approximately on-shell amplitude, respectively. A s←s is omitted, because a purely space-like 4-point function does not appear as a sub-amplitude in a scattering process initiated by two particles only. The additionally introduced time-like momenta k 1 , k 2 and space-like momenta k 3 , k 4 in Eq. (30) point in the same direction in 3-space as the original q i , but with swapped virtualities. More precisely, the invariant mass of the scattering weak boson pair, √ s, is kept fixed and We can identify the matrices of the right hand side of Eq. (29) by following the guiding principles introduced in section 3.1. The matrix T 0 in the numerator has to be A t←s to guarantee a reduction to the EFT limit for low energy scales. Additionally, the virtualities of polarization vectors in the sum over intermediate states have to be the same in order to guarantee reduction to the correct vector boson propagator (see Eq. (12)), i.e. in matrix multiplication of the helicity amplitudes in Eq. (30) only the products A i←t A t←j or A i←s A s←j are allowed. The unitarized interaction matrix has to be of transition type t ← s and, thus, the matrix product in the numerator is determined to be A t←t A t←s . The denominator has to behave as t ← t, which leaves only open the possibility of a linear combination of A t←t A t←t and A t←s A s←t for the matrix product in the denominator. This linear combination has to suppress both the polynomial rise with the invariant mass of the scattering system, √ s, as well as the rise with the space-like virtualities q 2 1 and q 2 2 . Time-like virtualities are of no concern once the dependence on s is addressed, because s provides an upper limit for q 2 3 and q 2 4 . Contributions involving high virtuality space-like momenta, especially for the transverse operators, will eventually lead to a unphysical cross section growth at s |q 2 1 | , |q 2 2 |. An example is given in Eq. (17). Either the s q 2 1/2 or the q 4 1/2 m 2 W /s term could become dominant at low √ s. To ensure that the unitarized amplitude will not rise due to un-suppressed space-like virtualities and therefore become unphysical, the denominator has to contain at least as many space-like states as the numerator. Hence, the matrix product A t←t A t←t has to be omitted and we arrive at the unitarization formula In this off-shell extension of the linear T-matrix unitarization, the eigenvectors of denominator and numerator will only align exactly in the on-shell limit. In fact, since A t←s is not normal, (non-aligned) eigenvectors can only be defined for the hermitian and the anti-hermitian parts of A t←s separately. As a result, the suppression of large enhancements cannot be guaranteed for states which fall along eigenvectors of small eigenvalues of the denominator. For the case at hand, W + W + scattering, this is not problematic for the operators O S , where only one helicity combination, namely the purely longitudinal ones, will receive a leading contribution, proportional to s 2 . However, multiple helicity combinations will receive a strong enhancement if at least one coefficient of transverse or mixed dimension-8 operators is non-zero. Therefore, the formula in Eq. (32) is still not satisfactory. Using the maximal eigenvalue a 2 max of the matrix product A t←s A s←t instead, individually for each j = 0, 1, 2 partial wave, will ensure that the resulting amplitudes are always below the unitarity limit. Our final unitarization formula, which we call the T u model, reads and fulfills all the guiding principles listed at the end of Section 3.1. Note that other choices would be possible for the suppression factors 1/(1 + a 2 max /4). For example, a 2 max could be taken the same for the j = 0, 1, 2 partial waves. This would correspond to a common overall form-factor, i.e. the dynamical suppression of the EFT growth would set in at a unique scale of new physics for all helicity combinations and partial waves. Clearly, such changes correspond to different models of the BSM dynamics. Here, we use the T u model because it is closer to the previous T-matrix unitarization model of Ref. [12].

Consequences for LHC Physics
Both the newly introduced T u -model and T-matrix unitarization modify the naive EFT description in slightly different ways, but we expect both to agree at asymptotically large energies, s m 2 W , |q 2 i |, when a single helicity configuration and, thus, a single large eigenvalue of the scattering matrix dominates the high energy behavior. In order to demonstrate these features, we start with a comparison of the three models, using Wilson coefficients near the present experimental limits for the dimension-8 operators. Next we discuss their impact on various observables, with an eye to distinctions between the different operator classes.
The Monte-Carlo generator VBFNLO is used to calculate distributions and fiducial cross sections for the vector boson scattering process pp → W + W + jjX → + ν + ν jjX at NLO QCD for √ s = 13 TeV. Here, + denotes a positron or muon in the final state. The jets are defined by anti-k t clustering [29] with radius R = 0.4. They are ordered by transverse momenta and the tagging jets at NLO are defined as the two hardest jets. As default, we use the CT10 PDF set [30], and electroweak parameters are determined within the G F -scheme with the measured values of G F , m W , m Z and m H as input. For the fiducial cross section we follow the recent  CMS analysis [6] and use the following cuts, dubbed VBF cuts: m > 20 GeV, m jj > 500 GeV, In Fig. 4, we compare the prediction of the naive EFT (dashed pink), the linear T-matrix unitarization (solid orange), the newly introduced T u model (solid blue) and the SM (solid black) for the longitudinal operator O S 1 as a function of the invariant mass of the vector boson pair. The coefficient for the naive EFT is chosen as the current experimental limit, TeV. For the unitarized models, T-matrix and T u , we choose the coupling such that the fiducial cross section of the naive EFT and the T u model coincide within the VBF cuts of Eq. (34). The number of produced events in the unitarized model are therefore nearly identical to the naive EFT expectation within the high energy region used by CMS to set the experimental bound, and thus the chosen value of F S 1 = 60 TeV −4 approximates the present bound on this coupling for the two unitarized models. Note that the coupling for the unitarized models is approximately a factor 3 larger than for the naive EFT description. The expected excess of events with invariant masses above 2 TeV in the naive EFT description will violate unitarity and is therefore unphysical. Fig. 4b shows that the events of this unphysical high energy tail need to be redistributed to energies between 800 GeV and 2 TeV for the unitarized models, leading to the weaker limit on the Wilson coefficient. Limits on the dimension-8 coefficient derived with the naive EFT model overestimate the sensitivity of experiments to the scale of high energy BSM effects. As displayed in Fig. 4a, the T u model reproduces the linear T-matrix unitarization prescription very well in the high energy range, with barely visible differences at intermediate energies, below M W W ∼ 2 TeV, which can be traced to subleading effects in q 2 i /s. The deviation of the unitarization models from the SM at high energies is greatly reduced as compared to the naive EFT, and is a valid description beyond the point of unitarity violation in the naive EFT model.
In Table 2, we list the estimated bounds on the full set of dimension-8 coefficients for the T u model, derived as above by matching the fiducial cross section to the one obtained in the naive EFT model. We stress that these numbers should be taken as rough estimates only, to be superseded by full experimental analyses. Their main purpose here is for use in subsequent figures. They illustrate deviations from the SM which are at the edge of what is presently allowed experimentally. The bounds on these Wilson coefficients in the T u model are about a factor of 3 weaker than the corresponding bounds derived within the naive EFT, for all three types of dimension-8 operators. Note, also, that the normalization conventions differ betweenÉboli and VBFNLO definitions for mixed and transverse operators. In the following we use the notation  Table 2: Experimental limits (in TeV −4 ) on dimension-8 operators from the observation of pp → W ± W ± jjX by CMS [6] (first column) and corresponding estimates for the bounds on the Wilson coefficients f i /Λ 4 in the T u model (second column). Columns three and four give the corresponding numbers for the VBFNLO normalization of operators. See text for further details.
The differential cross section as a function of the invariant mass of the W-pair in Fig. 4 cannot be accessed experimentally, because the 4-momentum of the neutrinos is not measurable. A more readily accessible observable is the differential distribution in the invariant mass of the two charged leptons, which is correlated to a sufficient degree to the invariant mass of the two vector bosons. Fig. 5 shows the corresponding distribution for one non-zero coefficient of each class, namely the longitudinal O S 1 (Fig. 5a), the transverse O T 0 (Fig. 5b) and the mixed O M 0 (Fig. 5c) operator. For all three coefficients, the T u model is suppressed at high energy scales where the EFT description violates unitarity. However, the high energy tails differ by approximately one order of magnitude between the longitudinal operator O S 1 and the transverse operator O T 0 . The differential cross section of the mixed operator in the the T u model lies between these. Below 500 GeV, the event production is mainly driven by SM contributions, as indicated by the fact that the SM curve coincides with the T u model curve for all three operators.
In an attempt to distinguish the different operator types, we study the transverse momentum of the leading tagging jet, p T,max (j), and the difference of the two lepton transverse momenta, To optimize the ratio of BSM to SM events for the following study, Fig. 5 suggests a cut, m + + > 500 GeV, on the charged lepton pair invariant mass. We show only the transverse and the longitudinal operators in the following and omit the mixed operators, which fall somewhat in between. In Fig. 6a and Fig. 7a, the differential cross sections as a function of p T,max (j) and ∆p T, , respectively, are plotted. On the right-hand-side, in Fig. 6b and Fig. 7b the same curves are shown as normalized distributions, which helps to better expose differences in shape. We compare the slope of the SM (solid black), the longitudinal operator O S 1 (T u : solid blue, naive EFT: dashed pink) and the transverse operator O T 0 (T u : solid brown, naive EFT: dashed purple).
Since incoming transversely polarized weak bosons lead to a harder jet p T distribution than longitudinally polarized bosons [31], and since the transverse operators enhance the transverse components, we expect more events at larger p T,max (j) for the transverse operators as compared to the longitudinal ones. This is clearly borne out in Fig. 6, which shows a considerably harder p T,max (j) spectrum for the O T 0 operator than for O S 1 . The cross section enhancement for the longitudinal O S 1 operator occurs at small p T , which is typical for incident longitudinal bosons. At large p T,max (j), where incident transversely polarized W s dominate, the SM and O S 1 curves coincide, indicating that the underlying anomalous quartic gauge coupling is mostly longitudinal.
Anomalous transverse operators produce cross section enhancements also at large p T,max (j). Here, an interesting difference can be observed between the naive EFT model and our T u model: T u unitarization considerably softens the p T,max (j) spectrum (dashed purple to solid brown curves). This effect is caused by the suppression of any large enhancement of the  V V → V V partial wave amplitudes, irrespective of its origin. For the transverse operators one finds unphysically large enhancements also at high virtualities of the incoming W s, while the 2 → 2 center of mass energy, √ s = m W W , remains small. Such an enhancement would not be corrected by a unitarization attempt which relies only on suppression at large s, such as the form-factor unitarization implemented previously in VBFNLO. 6 Thus, one needs to be cautious when devising observables for transversely polarized scattering based on a naive EFT approach: The large enhancement at high p T,max (j) for the purple O T 0 curve is an artifact of the missing unitarization. The properly unitarized distribution has a shape which is almost identical to the SM curve in Fig. 6b, which is also dominated by incoming transversely polarized W s. Rather, the distinction between incoming longitudinal and transverse weak bosons has to rely on the differences in the 0 < p T,max (j) < 200 GeV region, where, fortunately, also the bulk of the cross section is concentrated in all cases.
A transversely polarized W + with helicity λ = +1 tends to emit the charged anti-lepton in the forward direction relative to the W -momentum, which leads to a high lepton p T . This is in contrast to a negatively polarized W + , which produces a relatively soft + , and the nearly equally shared energy between + and neutrino for the decay of an energetic, longitudinally polarized W . Thus, ∆p T, promises to distinguish (λ 3 , λ 4 ) = (+, +) helicities (high ∆p T, ) from e.g. the (0, 0) or (+, −) helicity combinations at lower average ∆p T, .
The corresponding differences are clearly exhibited in Fig. 7. Also for the ∆p T, -distribu- 6 This problem for a form-factor implementation can easily be cured by generalizing the functional dependence of the form-factor, e.g. to F (s, q 2 1 , q 2 2 ) = (1 + z(s 2 + (q 2 1 + q 2 2 ) 2 )/Λ 4 F F ) n , with n ≤ −1 and a phase factor z = 1 or z = i. tions, the slopes of the dimension-8 operator enhancements are noticeably influenced by the T u model. The unphysical events at larger ∆p T, are suppressed because ∆p T, and √ s = m W W are highly correlated. The T u model prediction for the longitudinal operators and the SM have more events at ∆p T, below 1000 GeV and receive a large suppression for larger ∆p T, . As expected, the transverse operator produces a broader distribution, i.e. the enhancement due to (λ 3 , λ 4 ) = (+, +) polarization is clearly visible. For this observable, the unitarization model even increases the discrimination power between different operators types.

Discussion and Conclusions
The parameterization of new physics effects in vector boson scattering via anomalous quartic gauge couplings or an effective field theory, including operators up to dimension-8, is a useful tool for analyzing VBS at the LHC. However, because of the large energy reach of hadron colliders, which spans from low energies and momentum transfers where the pure EFT description is valid, to regions of phase space where the polynomial growth of partial wave amplitudes with energy exceeds unitarity limits, the naive EFT description must be generalized to a model which respects unitarity bounds. In this paper, we have developed the T u model, which is one such generalization and which closely mirrors a K-matrix or linear T-matrix unitarization of anomalous VBS amplitudes. The T u model has been implemented as a purely numerical procedure in the Monte Carlo program VBFNLO [15], which allows to analyze VBS at NLO QCD precision, for arbitrary dimension-8 operators [17]. The T u model is constructed such that it reduces to the naive EFT approximation in all phase space regions where this description is valid, and it smoothly interpolates to a unitarized description for VBS at high virtualities. These high virtualities may either correspond to high boson-pair invariant masses, m V 3 V 4 , signified by high energy and transverse momentum of the produced vector bosons in V 1 V 2 → V 3 V 4 , or to highly off-shell incoming V 1 or V 2 , i.e. large space-like q 2 i , which corresponds to pp → V 3 V 4 jjX events with tagging jets at very high transverse momentum. Unphysical growth of VBS cross sections at high tagging jet p T , (see Fig. 6) which is present in a naive EFT implementation even at small m V 3 V 4 , also needs to be suppressed, and the T u model does provide this regularization.
The purely numerical implementation grants great versatility and avoids analytical approximations, like neglecting m 2 W /s or q 2 i /s suppressed terms in a high energy approximation. It allows for arbitrary combinations of dimension-8 operators to be present in the effective Lagrangian and thus provides a general unitarized framework to analyze the effects of dimension-8 operators in VBS at the LHC. In addition, the numerical isolation of off-shell V V → V V helicity amplitudes at intermediate steps of the calculation, allows, with little additional effort, to generate events for selected center of mass helicities in the BSM V V → V V contribution, similar to a recent implementation in the PHANTOM Monte Carlo [32,33]. So far, the implementation of the T u model has been tested and is available for same-sign W -boson scattering, more precisely for pp → W ± W ± jjX → ± ν ± ν jjX. However, the generalization to single charged VBS (W Z-scattering) and neutral channels will become available soon [34].
For same-sign W scattering we have analyzed distributions which promise a differentiation between individual tensor structures of the operators in the EFT expansion, beyond the only theoretically accessible di-boson invariant mass distribution in Fig. 4 or the invariant mass distribution of the two same-sign charged leptons in Fig. 5. The transverse momentum distribution of the tagging jets, e.g. p T,max (j), which is shown in Fig. 6 is a good separator between longitudinal and transverse polarization of the incident weak bosons. The charged lepton transverse momentum difference, ∆p T, , which is shown in Fig. 7, can be used to distinguish different combinations of W polarizations in the final state. For the same-sign W case considered in this paper, we have shown in Fig. 4 that the T u model closely agrees with the T-matrix unitarization discussed by the WHIZARD group [12] for longitudinal W + W + scattering. However, the treatment of subleading, m 2 W /s or q 2 i /s suppressed terms is different and means that the two schemes provide different unitarization models. The numerical framework which is now set up in the VBFNLO program allows for easy implementation of variants of the T u model, such as taking into account more than just the largest eigenvalue of the tree level scattering matrix for the denominator when going from Eq. (32) to Eq. (33), or by exploring other mappings of these real eigenvalues onto the Argand circle. We leave such investigations to future work.