Study of the weak annihilation contributions in charmless $B_s\to VV$ decays

In this paper, in order to probe the spectator-scattering and weak annihilation contributions in charmless $B_s\to VV$ (where $V$ stands for a light vector meson) decays, we perform the $\chi^2$-analyses for the end-point parameters within the QCD factorization framework, under the constraints from the measured $\bar B_{s}\to$$\rho^0\phi$, $\phi K^{*0}$, $\phi \phi$ and $K^{*0}\bar K^{*0}$ decays. The fitted results indicate that the end-point parameters in the factorizable and nonfactorizable annihilation topologies are non-universal, which is also favored by the charmless $B\to PP$ and $PV$ (where $P$ stands for a light pseudo-scalar meson) decays observed in the previous work. Moreover, the abnormal polarization fractions $f_{L,\bot}(\bar B_{s}\to K^{*0}\bar K^{*0})=(20.1\pm7.0)\%\,,(58.4\pm8.5)\%$ measured by the LHCb collaboration can be reconciled through the weak annihilation corrections. However, the branching ratio of $\bar B_{s}\to\phi K^{*0}$ decay exhibits a tension between the data and theoretical result, which dominates the contributions to $\chi_{\rm min}^2$ in the fits. Using the fitted end-point parameters, we update the theoretical results for the charmless $B_s\to VV$ decays, which will be further tested by the LHCb and Belle-II experiments in the near future.


Introduction
Non-leptonic B-meson weak decays play an important role in testing the flavor dynamics of the Standard Model (SM) and exploring possible hints of New Physics beyond it.Theoretically, one of the main obstacles for a reliable prediction on these decays is how to evaluate precisely the hadronic matrix elements of local operators between the initial and final hadronic states, especially due to the nontrivial QCD dynamics involved.To this end, several attractive QCDinspired approaches, such as QCD factorization (QCDF) [1,2], perturbative QCD (pQCD) [3,4] and soft-collinear effective theory (SCET) [5][6][7][8], have been proposed in the last decades.
However, the convolution integrals of the parton-level hard kernels with the asymptotic forms of light-cone distribution amplitudes (LCDAs) generally suffer from the end-point divergence in the weak annihilation (WA) amplitudes.This divergence limits the predictive power and introduces large theoretical uncertainties.
In the QCDF approach, the end-point divergent integrals, signalling of infrared-sensitive contributions, are usually parameterized by two complex quantities X A and X L that are defined, respectively, by [9,10] where Λ h = 0.5 GeV, and the two phenomenological parameters ρ A and φ A account for the strength and possible strong phase of WA contributions near the end-point, respectively.In addition, the spectator-scattering amplitudes also involve the end-point divergence, which is dealt with the same manner by introducing the complex quantity X H (ρ The numerical values of ρ A ,H and φ A ,H are unknown and can only be inferred by fitting them to the experimental data so far.While weakening the predictive power of the QCDF approach, the parameterization scheme provides a feasible way to explore the WA effects from a phenomenological point of view. Theoretically, the WA contributions with possible strong phase have attracted a lot of attention in the past few years, for instance, in Refs.[11][12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28][29][30].Traditionally, both ρ A and φ A are treated as universal parameters for different kinds of annihilation topologies.However, a global fit for the end-point parameters indicates that, while a relatively large end-point parameter is needed for the decays related by isospin symmetry, there exist some tensions in B → φK * and B → πK decays [22], with the latter exhibiting the so-called "πK CP-asymmetry puzzle". 1 In Refs.[23,24], after studying carefully the flavor dependence of the end-point parameters in charmless B → P P (where P stands for a light pseudo-scalar meson) decays, the authors suggest that the end-point parameters should be topology-dependent.Such a topologydependent parameterization scheme is also favored by most of the charmless B → P P and P V (where V stands for a light vector meson) decays, as demonstrated in Refs.[26][27][28], and it could provide a possible solution to the well-known "πK CP-asymmetry puzzle" [25].In addition, using the recent measurements of the pure annihilation B s → π + π − and B d → K + K − decays by the CDF [36], Belle [37] and LHCb collaborations [38,39], the authors of Refs.[23,24] find significant flavor-symmetry breaking effects in the nonfactorizable annihilation contributions.
Experimentally, due to the rapid development of dedicated heavy-flavor experiments, more precise measurements of non-leptonic B decays will be available.As reported in Ref. [40], for instance, over 10 11 b b quark pairs are produced per fb −1 of data at the LHCb experiment.
Furthermore, after the high-luminosity upgrade, a dataset of 50 fb −1 will be accumulated [41][42][43][44].In addition, most recently, the SuperKEKB/Belle-II experiment has started test operations and succeeded in circulating and storing beams in the electron and positron rings.The annual integrated luminosity is expected to reach up to 13 ab −1 , and over 10 10 samples of b b quark pairs will be accumulated by the Belle-II experiment [45].Thus, the forthcoming measurements of not only B u,d but also B s decays are expected to reach a high accuracy, which could provide us with a clearer picture of the WA contributions in these decays.
In this paper, motivated by the recent theoretical studies and the bright experimental prospects, we will investigate the WA contributions in charmless B s → V V decays, which involve more observables than in charmless B → P P and P V decays and may provide much stronger constraints on the end-point parameters.In addition, using the obtained values of the end-point parameters, we will update the theoretical results for charmless B s → V V decays within the QCDF framework.
Our paper is organized as follows.In section 2, we review briefly the WA amplitudes within 1 The direct CP asymmetries to be roughly the same [31][32][33].However, the current experimental data show a significant difference between them, [34,35], deviating from zero at about 5.6 σ level.
the QCDF framework and observables for charmless B s → V V decays.Section 3 is devoted to the numerical results and discussions.Finally, we give our conclusion in section 4.
2 Brief review of the theoretical framework for charm-

Amplitudes in QCD factorization
In the SM, the effective weak Hamiltonian for non-leptonic B-meson decays is given by [46,47] where V qb V * qp , with q ∈ {u, c, t} and p ∈ {d, s}, are products of the Cabibbo-Kobayashi-Maskawa (CKM) matrix elements, and C i the Wilson coefficients of the effective operators O i .Starting with the effective Hamiltonian and following the strategy proposed in Ref. [48], Beneke et al.
The theoretical framework for charmless B → V V decays has also been fully developed (cf. Refs. [19,20,53] for details).In this paper, we follow the same conventions as in Refs.[53,56].
Within the QCDF framework, after performing the convolution integrals of the O(α s ) hard kernels with the asymptotic forms of the light-meson LCDAs, one gets the following basic building blocks of the WA amplitudes in charmless B → V V decays [53,56]: for the non-vanishing longitudinal contributions, and for the transverse ones, where the superscripts 0 , ± refer to the vector-meson helicities.Here V denoting the longitudinal and transverse vector-meson decay constants; m B and m 1,2 are the masses of the initial and final states, respectively.
The two complex quantities X A and X L are introduced in Eqs. ( 4)- (10) to parameterize the end-point divergence (cf.Eqs. ( 1) and ( 2)).In addition, we have distinguished the WA contributions with the gluon emitted either from the initial (marked by the superscript "i") or from the final state (marked by the superscript "f "), corresponding to the nonfactorizable and the factorizable annihilation topologies, respectively.For the factorizable annihilation topologies, as argued in Refs.[23,24], since all decay constants have been factorized out from the hadronic matrix elements, the building blocks A f,0 3 and A f,− 3 are independent of the initial states.However, for the nonfactorizable annihilation topologies, X i A,L are generally non-universal for B u,d and B s decays [23,24].Besides, an additional complex quantity, X H (ρ H , φ H ) = X A (ρ A , φ A )| A→H , is introduced to parameterize the end-point divergence in the hard spectator-scattering (HSS) amplitudes (cf.Refs.[9,53] for details).
With the above prescriptions for the WA amplitudes, we consider in this paper the following B s → V V decay modes:2 (i) ∆D = 1 transition: the color-allowed tree-dominated Bs → ρ − K * + , the color-suppressed tree-dominated Bs → ρ 0 K * 0 and ωK * 0 , as well as the penguin-dominated Bs → φK * 0 decay, the amplitudes of which are given, respectively, as [9, 53] (ii) ∆S = 1 transition: the penguin-dominated Bs → K * + K * − , K * 0 K * 0 , and φφ, ρ 0 φ, ωφ, as well as the pure annihilation Bs → ρ + ρ − , ρ 0 ρ 0 , ρω, ωω decays, the amplitudes of which are given, respectively, as [9, 53] In the above decay amplitudes, the vertex, penguin and spectator-scattering corrections are encoded in the effective coefficients α p,h i (cf.Ref. [9,53] for details), and the WA contributions are denoted by β p,h i (or b p,h i ), which are defined, respectively, by [9,53] Based on the previous studies [22,53] and the amplitudes given above, B s → V V decay modes can be classified as follows according to their sensitivities to the WA and/or HSS corrections: • The pure annihilation decay modes: Because both WA and HSS corrections involve the undetermined end-point contributions, the interference between them presents an obstacle for precisely probing the WA contributions from data.Fortunately, such problem can be avoided by using the pure annihilation decay modes, which can be easily seen from Eqs. ( 20), ( 21), ( 22) and (23).In this paper, the Bs → ρ 0 ρ 0 , ρ + ρ − , ρ 0 ω and ωω decays belong to this category, but unfortunately, they have not been measured for now.
• The color-suppressed tree-and electroweak or QCD flavor-singlet penguin-dominated decays: Only two decay modes fall into this category, namely, Bs → ρ 0 φ and ωφ.From Eqs. ( 18) and (19), one can find that these decays are characterized by an interplay of the color-and CKM-suppressed tree amplitude, α 2 , electroweak penguin amplitude, α c 3,EW , and flavor-singlet QCD penguin amplitude, α c 3 .For Bs → ωφ decay, due to a partial cancellation between the QCD and electroweak penguin contributions, the largest partial amplitude is α 2 [53], which is very sensitive to the HSS corrections.For Bs → ρ 0 φ decay, α 2 is also nontrivial even though it is numerically smaller than α c 3,EW when ρ H is small.More importantly, a remarkable feature of such two decays is that their amplitudes are irrelevant to the WA contributions and, therefore, very suitable for probing the HSS corrections.Recently, the branching ratio of Bs → ρ 0 φ decay has been measured by the LHCb collaboration with a statistical significance of about 4σ [57].
• The color-suppressed tree-dominated ∆D = 1 decays: This class includes Bs → ρ 0 K * 0 and ωK * 0 decays, whose amplitudes are given by Eqs. ( 12) and ( 13), respectively.The CKM-factors relevant to the effective coefficients in their amplitudes are at the same order, ∼ λ 3 , and therefore one can roughly find that α 2 dominates their amplitudes.Considering further the fact that the HSS contribution in α 2 is proportional to the largest Wilson coefficient C 1 , we can generally expect that these decays present strong constraints on the HSS end-point parameters even though they are not as "clean" as Bs → ρ 0 φ decay due to the interference induced by β c 3 .However, there is no available data for these decays for now.
• The QCD penguin-dominated decays: This class contains the residual decays, except for Bs → ρ − K * + , considered in this paper, among which Bs → φK * 0 , φφ and K * 0 K * 0 decays have been measured.From their amplitudes, Eqs. ( 14), ( 16) and ( 17), one can find that the effective color-allowed QCD penguin amplitude, αc 4 ≡ α c 4 + β c 3 , plays a dominant role [53]; and the penguin-annihilation amplitude, β c 4 , presents the first subdominant contribution besides αc 4 for the longitudinal amplitude of the last two decay modes [53].Considering further the fact that the HSS contribution in α c 4 is trivial compared to the LO contribution C 4 + C 3 /N C , we can generally conclude that such three QCD penguindominated decays are suitable for probing the WA contribution from data even through they are not as "clean" as the pure annihilation decay modes.
It should be noted that, such an expectation or such a conclusion is valid only when ρ H is not very large especially for the decays involving α c 3 , for instance, Bs → φK * 0 and φφ decays.The HSS correction in α c 3 with a large ρ H can bring about a significant correction.
• The color-allowed tree-dominated ∆D = 1 decay: In this paper, only Bs → ρ − K * + decay belongs to this class.For this decay mode, the effects of HSS and WA contributions are generally not very significant, because of the dominant role played by the color-allowed tree amplitude, α 1 .It also has not been measured.

Observables
Using the amplitudes given above, the observables for B s → V V decays can be defined as follows.The most important observables are the CP-averaged branching ratio and direct CP asymmetry, which are defined theoretically as [22] respectively.The decay rates should be summed over the polarization state (h = L, , ⊥) for evaluating the "whole" observables.Following the convention of Ref. [53], the polarization amplitudes can be easily obtained from the helicity amplitudes through the relations, ĀL = Ā0 and Ā ,⊥ = Ā− ± Ā+ √ 2 .Besides branching ratio and CP asymmetry given by Eqs. ( 30) and ( 31), the two-body B s → V V decays with cascading decays V → P P provide additional observables in the full angular analysis of the 4-body final state [58].There are polarization fractions and relative phases defined by for Bs decays.The same quantities for B s decays are obtained by the replacement Āh → A h .
The CP-averaged polarization fractions and CP asymmetries are given, respectively, by among which only two of the polarization fractions are independent due to the normalization condition f L + f + f ⊥ = 1; such a definition for A h CP is in fact the same as Eq. ( 31) for a given h .The CP-averaged and CP-violating observables for the relative phases can be constructed, respectively, as for h = , ⊥.This phase convention for the amplitudes implies φ h = ∆φ h = 0 at leading order, where all strong phases are zero [53], while it should be noted that the sign of A L relative to the transverse amplitudes differs from the experimental convention, which leads to an offset of π for φ ,⊥ [22,53].
It should be noted that the above "theoretical" definitions are in the flavor-eigenstate basis and at t = 0. where with the heavy and light mass-eigenstates, |B H,L s = p|B s ∓q| Bs ; is the CP asymmetry due to the width difference; y s is the parameter proportional to the width difference Then the relation between the experimentally measurable and theoretically calculated branch- , can be written as [22,59] Here the decay width parameter y s is universal for B s decays and has been well measured, y s = 0.063 ± 0.005 [34].However, the CP asymmetry H f h is generally non-universal, not only for various B s decay modes, but also for various polarization states.Moreover, its values in most of B s decays are not measured.Therefore, we take the SM prediction [22] H Accordingly, the experimentally measurable polarization fractions should also be modified as they still satisfy the normalization condition In addition, such a correction induced by the B s oscillation does not affect the definition for the three polarization-dependent CP asymmetries given by Eq. ( 33).

Numerical results and discussions
With the theoretical formulae given above, we now present our numerical results and discussions.The values of the input parameters used in our evaluation are summarized in Table 1.
So far, only some observables of Bs →ρ 0 φ, φK * 0 , φφ and K * 0 K * 0 decays, including the CPaveraged branching ratios, the polarization fractions and the relative phases between different Table 2: Experimental data for the measured observables of Bs → ρ 0 φ, K * 0 K * 0 , φK * 0 and φφ decays, as well as the deviations of theoretical results from data in cases I and II, i.e., the χ i σ i (i denotes a given observable) evaluated at the best-fit points of (ρ i,f A , φ i,f A ) with the other inputs given by Eq. ( 42) and Table 1.The bold value denotes the largest deviation in the fit.

Observable
Decay mode Exp.[ [34] are listed in the "Exp."column of Table 2, and will be used as constraints in the following χ 2 -fits.
In order to probe the HSS and WA contributions in charmless B s → V V decays, we perform χ 2 -analyses for the end-point parameters, adopting the statistical fitting approach illustrated in our previous work [25] (cf.Appendix C of Ref. [25] for detail).In the following fits and posterior predictions, we have to evaluate the theoretical uncertainties induced by the inputs listed in Table 1.The total theoretical errors are obtained by evaluating separately the uncertainties induced by each input parameter and then adding them in quadrature.
Our χ 2 -fits are based on the topology-dependent parametrization scheme for the end-point divergence [23,24].This implies that we need four free parameters (ρ i,f A , φ i,f A ) (where the superscripts i and f , as introduced in section 2, correspond to the nonfactorizable and factorizable annihilation topologies, respectively) to describe the WA contributions.Besides, we also need two free parameters (ρ H , φ H ) to describe the HSS contributions.
3.1 Constraints on (ρ H , φ H ) from Bs → ρ 0 φ decay As has been illustrated in Refs.[1,2,22], the factor ρ H e iφ H summarizes the remainder of the non-perturbative contribution including a possible strong phase; the numerical size of such a complex parameter is unknown.However, a too large ρ H will give rise to numerically enhanced subleading Λ QCD /m b contributions compared with the formally leading terms.Thus, the size of ρ H should be carefully coped with.
As analyzed in the last section, the Bs → ρ 0 φ decay is independent of WA contributions and sensitive to HSS corrections.Therefore, it provides an ideal channel for probing the end-point parameters in the HSS amplitudes.Recently, the branching ratio of Bs → ρ 0 φ decay has been measured by the LHCb collaboration [57], with a significance of about 4σ.
Taking ρ H = 0 , 0.5 , 1 and using the central values of input parameters in Table 1, the dependence of the theoretical result for B( Bs → ρ 0 φ) on φ H is shown in Fig. 1 It has been noted that, besides the Bs → ρ 0 φ decay, the large ρ H is also disfavored by the color-suppressed tree-dominated B d → ρ 0 ρ 0 decay [29] even though a large HSS correction with large ρ H is helpful for explaining the "πK and ππ puzzle" [25].In the following analysis and evaluation, we take a conservative choice that as inputs.Even though such a φ H has a large uncertainty, its effect on the following analysis for (ρ i,f A , φ i,f A ) would be not significant because the HSS contribution with ρ H = 0.5 is severely suppressed.

Case I: constraints on (ρ A , φ A ) from Bs → φK * 0 and φφ decays
As can be seen from Eqs. ( 14) and ( 17), the Bs → φK * 0 and φφ decays have similar amplitude structures.However, being a ∆D = 1 transition, the Bs → φK * 0 decay amplitude is suppressed by one power of the Wolfenstein parameter λ 0.23 compared with that of Bs → φφ.This explains why the branching ratio B( Bs → φK * 0 ) should be much smaller than B( Bs → φφ).
One can see from Tables 2 and 3 that the previous QCDF [19,53] and pQCD [65] predictions are in good agreement with each other for B( Bs → φφ) and B( Bs → φK * 0 ); in addition, their predictions are consistent with the data for the former, but are much smaller than the data for the latter.Such a deviation could possibly be moderated by the different WA contributions involved in these two decays (cf.Eqs. ( 14) and ( 17)).To this end, we firstly take the measured Bs → φK * 0 , φφ decays as constraints to fit the WA contributions, which is named as "case I" for convenience of discussion.
Under the constraints from the measured Bs → φK * 0 and φφ decays (there are totally eight  observables available, see Table 2), the allowed spaces of (ρ i A , φ i A ) and (ρ f A , φ f A ) are shown in Fig. 2. We find that: • As shown in Figs.2(a) and 2(b), the spaces of (ρ i A , φ i A ) and (ρ f A , φ f A ) are bounded into three and two separate regions, respectively, at 68% C.L., which are labeled as SI-1, 2, 3 and SF-1, 2 for convenience of discussion.We do not find any direct correspondence between SF-1, 2 for (ρ f A , φ f A ) and SI-1, 2, 3 for (ρ i A , φ i A ). 3   • For (ρ f A , φ f A ), as shown in Fig. 2(b), the space of SF-1 is strictly bounded at 68% C.L., while 3 For instance, if we pick out the allowed space of SF-1 at 68% C.L. for (ρ f A , φ f A ), all of the three solutions SI-1,2,3 for (ρ i A , φ i A ) are still allowed.
the constraint on the one of SF-2 is very loose.The best-fit point with χ 2 min /n dof = 8.6/4 falls in SF-1; numerically, It should be noted that we can also find a point in SF-2, (0.45, −40 • ) [SF-2], having a χ 2 value similar to the χ 2 min in SF-1.The situation for (ρ i A , φ i A ) is similar to the one for (ρ f A , φ f A ), but is much more complicated as Fig. 2(a) shows.The best-fit point in SI-1 is These allowed spaces in Figs.2(a) and 2(b) will be further confronted with the measured observables of Bs → K * 0 K * 0 decay in the next subsection.
• The correlation, ρ i A vs. ρ f A , is shown by Fig. 2(c).One can see again that ρ f A is significantly divided into two parts.The relation between ρ f A and ρ i A is not clear due to the large uncertainties except that they can not be equal to zero simultaneously.The correlation, φ i A vs. φ f A , shown in Fig. 2(d) is very interesting.One can clearly see that the φ i A can be well determined except when φ f A ∼ −30 • or −190 • and vice versa; the case for φ f A is similar.Hence, the phases φ i,f A are expected to be well determined when more constraints are considered.
As argued in Refs.[23,24], the parameters (ρ f A , φ f A ) are expected to be universal for B u,d and B s systems, while (ρ i A , φ i A ) are flavor-dependent on the initial states.Comparing with the fitted results in B u,d system [29], we find that the best-fit values, Eqs. ( 43) and (44), are very similar to the results, ( , solution C given by Eq. ( 14) in Ref. [29]) obtained by fitting to B u,d → ρK * and K * K * decays.However, it should be noted that the results in Ref. [29] are based on the assumption (ρ H , φ H ) = (ρ i A , φ i A ), which is not employed in this paper because ρ H is strictly constrained by B s → ρφ decay as analyzed in the last subsection; thus, the flavor dependence of (ρ i A , φ i A ) is indeterminable here.The goodness of the fit can be characterized by χ 2 min /n dof and p-value. 4Numerically, we obtain χ 2 min /n dof = 8.6/4 , p-value = 0.07 (45) at the best-fit point given by Eqs. ( 43) and (44).In order to find which observables lead to the large χ 2 min and small p-value, we summarize the deviations of theoretical results from data for the considered observables in the fourth column of Table 2.It can be clearly seen that the tension between the theoretical result and data for B( Bs → φK * 0 ), ∼ −2.57σ, dominates the contributions to χ 2 min .Numerically, one can find χ 2 B( Bs→φK * 0 ) /χ 2 min = 77%.Such a tension implies that the problem of large B( Bs → φK * 0 ) mentioned in the beginning of this subsection is hardly to be moderated by the WA contribution due to the constraints from the other measured observables.In addition, from Table 2, we also find some significant tensions in the Bs → K * 0 K * 0 decay, which is not considered in the fit in this case (case I).This implies that the best-fit points given by Eqs. ( 43) and ( 44) might be excluded when the constraints from Bs → K * 0 K * 0 decay are considered; and the other fitted spaces in this case may also suffer challenges from Bs → K * 0 K * 0 decay, which will be studied in detail in the next two subsections.

decay
In the fit of case I, we do not include the measured Bs → K * 0 K * 0 decay, because it is difficult to understand its polarizations measured by the LHCb collaboration [66].It is well-known that, due to the (V − A) nature of the SM weak interactions, the hierarchical pattern among the three helicity amplitudes, A 0 : A − : A + = 1 : [67].Even after the QCD corrections are taken into account, the charmless B → V V decay amplitudes are generally still dominated by the longitudinal polarization component.For the penguin-dominated B → V V decays, the longitudinal polarization fraction is generally predicted at the level of about 50%, for instance in the B → φK * decays [68][69][70][71][72][73][74][75].Consistent with the above expectation, the longitudinal (transverse) polarization fraction, f L(⊥) ( Bs → K * 0 K * 0 ) ∼ 50% (25%), is predicted both in the QCDF [53] and in the pQCD approach [65].However, the obviously different experimental results have been measured by the LHCb collaboration [66], -350 -300 -250 -200 -150 -100 -50 0 0 10 20 30 40   which imply f ⊥ ( Bs → K * 0 K * 0 ) = (58.4± 8.5)%.Furthermore, these measurements are also inconsistent with the previous theoretical expectation, is satisfied by most of the charmless B → V V decays [53]).As a consequence, these possible anomalies present a challenge to the current theoretical predictions.Therefore, we would like to check if the modifications of end-point parameters could reconcile these anomalies.
In Fig. 3, we plot the dependences of B( Bs ) and (0.45, −40 • ), which are the best-fit values in SF-1 and -2, respectively, in case I. From Fig. 3, one can see that: • The small ρ i A 1 is perhaps allowed by the measured B( Bs → K * 0 K * 0 ) with φ i A around 0 • (or −360 • ) as Fig. 3(a) shows, but is excluded by both f L ( Bs → K * 0 K * 0 ) and f ⊥ ( Bs → K * 0 K * 0 ) as Figs.3(b) and 3(c) show.
• With ρ i A ∼ 3, B( Bs → K * 0 K * 0 ) and f L,⊥ ( Bs → K * 0 K * 0 ) present different requirements for φ i A , which can be seen by comparing Fig. 3(a) with Figs.3(b) and 3(c).This implies that the choice ρ i A ∼ 3 is also excluded by the anomalies of f L,⊥ ( Bs → K * 0 K * 0 ).
• Only a large ρ i A ∼ 6 with the phase φ i A ∼ −240 • or −130 • could possibly account for the current LHCb measurements for f L,⊥ ( Bs → K * 0 K * 0 ).It is very interesting that such possible solutions are similar to SI-2 and -3 shown in Fig. 2(a).However, the SI-1 is possibly excluded by f L,⊥ ( Bs → K * 0 K * 0 ).In order to further check such possible solutions, in the next subsection we will perform a combined fit for the end-point parameters with the measured observables of Bs → φK * 0 , φφ and K * 0 K * 0 decays as constraints.

Case II: combined fit to the measured B s → V V decays
Under the combined constraints from Bs → φK * 0 , φφ and K * 0 K * 0 decays (i.e., the measured 11 observables of B s → V V decays are now included), our fitted results for the end-point parameters are shown in Fig. 4, which is named as "case II" for convenience of discussion.It can be seen that: • For (ρ i A , φ i A ), because of the constraints from f L,⊥ ( Bs → K * 0 K * 0 ), the space SI-1 favored in case I is excluded at 68% C.L. in case II, which can be clearly seen from Fig. 4(a) and easily understood from the analysis in the subsection 3.3, while both SI-2 and -3 survive, but they are further restricted; the large ρ i A is required to fit f L,⊥ ( Bs → K * 0 K * 0 ).The spaces SI-2 and -3 are located symmetrically at two sides of φ i A ∼ −180 • For (ρ f A , φ f A ), the two spaces, SF-1 and -2, are still allowed and further restricted in case II as Fig. 4  For SF-2, even though the space of (ρ f A , φ f A ) is further restricted compared with case I, the constraints are still very loose.In addition, it also can be found that the relation (ρ ) is required at 68% C.L., which implies that the end-point parameters are topology-dependent and, therefore, confirms the suggestion proposed in Refs.[23,24].
In this case, we obtain χ 2 min /n dof = 11.2/7, p-value = 0.13 (50) at the best-fit point.The deviations of the theoretical results from data are summarized in the fifth column of Table 2.We again find that the observable B( Bs → φK * 0 ) results in the large χ 2 min and small p-value.Numerically, one can find χ 2 B( Bs→φK * 0 ) /χ 2 min = 90%, which is similar to case I.If we disregard B( Bs → φK * 0 ), we can find χ 2 min /n dof = 1.65/6 , p-value = 0.95 .
Comparing case II with case I, we find from Table 2 that their main difference is the deviations for the observables of Bs → K * 0 K * 0 decay: 0.00 σ vs. +4.31σ for B, +0.13 σ vs. +5.99σ for f L , and −1.21 σ vs. −4.96σ for f ⊥ , which again indicates that the abnormal data for Bs → K * 0 K * 0 decay can be explained through the WA contributions.However, in both cases I and II, the deviations for B( Bs → φK * 0 ), −2.57σ and −3.17 σ, are very large; this implies that the measured large B( Bs → φK * 0 ) = (1.13 ± 0.30) × 10 −6 [34], which is much larger than all of the current predictions, ∼ 0.4 × 10 −6 , in pQCD [65] and QCDF [19,53], is hardly to be accommodated by the WA contributions due to the constraints from the other observables and decay modes.

Updated results for charmless B s → V V decays
Using the fitted results given by Eqs. ( 43), (44) (case I) and ( 48), (49) (case II) for (ρ i,f A , φ i,f A ), Eq. ( 42) for (ρ H , φ H ), and the other input parameters listed in Table 1, we now present in Tables 3, 4, 5 and 6 our updated theoretical results (posterior predictions5 ) for the branching ratios, CP asymmetries, polarization fractions and relative phases in Bs → ρK * , ωK * , φK * , K * K * , φφ, ρφ, ωφ, ρρ, ρω, ωω decays, where the previous predictions in the QCDF [19,53] and pQCD [65] approaches are also listed for comparison.The first uncertainty for the results of cases I and II in these tables is caused by the input parameters listed in Table 1 and φ H given by Eq. ( 42); the second uncertainty for the results of case II corresponds to the uncertainties of (ρ i,f A , φ i,f A ) given by Eqs. ( 48) and (49).One can find from these tables that most of our results for the observables are in agreement with the current experimental data including the abnormal f L,⊥ ( Bs → K * 0 K * 0 ); the only Table 4: Theoretical results for the b → d induced Bs → ρ − K * + , ρ 0 K * 0 and ωK * 0 decays.The other captions are the same as in Table 3.

--
Table 5: Theoretical results for the b → s induced Bs → K * − K * + , ρ 0 φ and ωφ decays.The other captions are the same as in Table 3.
Table 6: Theoretical results for the pure annihilation Bs → ρρ, ρω and ωω decays.Our results A 0 CP = 0, ∆φ ⊥ = 0, ∆φ = 0 are in agreement with previous ones and are, therefore, not listed here.The other captions are the same as in Table 3.
exception is for B( Bs → φK * 0 ).A detailed discussion has been presented in the last subsections.
More theoretical and experimental efforts are needed to confirm or refute this possible puzzle.
Our results are also generally in consistence with the previous theoretical predictions in QCDF [19,53] and pQCD [65] within the theoretical uncertainties.The most obvious differences are the results for the pure annihilation B s decays, which can be seen from Table 6.Our results for the branching ratios of Bs → ρρ, ρω and ωω decays are about one order larger than the previous predictions; moreover, our results, f L,⊥ ∼ 55% , 25%, are also obviously different from f L,⊥ ∼ 100% , 0% [53,65].These differences in fact can be easily understood from the following: (i) The best-fit value of ρ i A is very large in order to fit the abnormal f L,⊥ ( Bs → K * 0 K * 0 ) ∼ 20% , 58% in case II as discussed above; it results in sizable nonfactorizable annihilation contributions.(ii) The pure annihilation decays, Bs → ωφ, ρρ, ρω, ωω, are only relevant to the nonfactorizable annihilation amplitudes.Therefore, it can be briefly concluded that, if one requires the WA corrections to account for the abnormal f L,⊥ ( Bs → K * 0 K * 0 ) measured by the LHCb collaboration, the large branching ratios and transverse polarization fractions of Bs → ρρ, ρω and ωω decays will be expected accordingly.Interestingly, the large nonfactorizable annihilation contributions have been observed in the pure annihilation B s → π + π − decay [22][23][24][36][37][38].
Finally, we would like to point out that the allowed spaces for the end-point parameters are still very large, especially in case I; our results are only based on the best-fit points, and the other allowed spaces for (ρ i,f A , φ i,f A ) shown in Figs. 2 and 4 are not taken into account here; moreover, it is also not clear whether the annihilation corrections should account for the abnormal f L,⊥ ( Bs → K * 0 K * 0 ).More data on B s → V V decays are needed for a definite conclusion.The pure annihilation decays mentioned above, except for Bs → ρω, having branching ratios O(10 −7 ) are in the scope of the LHCb and Belle-II experiments.Hence, a much clearer picture of the WA contributions in charmless B s → V V decays is expected to be obtained from these dedicated heavy-flavor experiments in the near future.

Conclusion
In summary, we have studied the HSS and WA contributions in charmless B s → V V decays.
In order to probe their strength and possible strong phase, we have performed χ 2 -analyses for the end-point parameters under the constraints from the measured Bs →ρ 0 φ, φK * 0 , φφ and K * 0 K * 0 decays.It is found that the end-point parameters in the factorizable and nonfactorizable annihilation topologies are non-universal at 68% C.L. due to the constraint from Bs → K * 0 K * 0 decay; this further confirms the findings in the previous work.Moreover, the abnormal polarization fractions f L,⊥ ( Bs → K * 0 K * 0 ) = (20.1 ± 7.0)% , (58.4 ± 8.5)% measured by the LHCb collaboration can be reconciled through the weak annihilation corrections with a large ρ i A .However, the B( Bs → φK * 0 ) exhibits a significant tension between the data and theoretical results, which dominates the contributions to χ 2 min in the fits.Using the best-fit end-point parameters, we have also updated the theoretical results for the charmless B s → V V decays within the framework of QCDF.It is found that the large branching fractions and transverse polarization fractions for the pure annihilation decays are possible if we require the WA contributions to account for the abnormal f L,⊥ ( Bs → K * 0 K * 0 ).Our results and findings will be further tested by the LHCb and Belle-II experiments in the near future.
and the difference between B[B s → f ] and B[B s → f ] can, therefore, reach up to O(10%) for final states that are CP-eigenstates, as has been observed

Figure 1 :
Figure 1: Fig.(a) shows the dependence of B( Bs → ρ 0 φ) on φ H with different ρ H labeled in the figure; the gray band is the experimental data within 1σ error bars.Fig.(b) shows the allowed spaces of (ρ H , φ H ) at 68% C.L. and 95% C.L. under the constraint from the measured B( Bs → ρ 0 φ); the dashed line corresponds to ρ H = 0.5, and the best-fit point corresponds to χ 2 min = 0.19.

Figure 2 :
Figure 2: The allowed spaces of (ρ i,f A , φ i,f A ) ( Figs.(a) and (b) ) at 68% C.L. and 95% C.L. under the constraints from the measured Bs → φK * 0 and φφ decays, namely, case I.The fitted spaces are also shown in (ρ i A , ρ f A ) and (φ i A , φ f A ) planes (Figs.(c) and (d) ) in order to show the possible correlation.The best-fit point corresponds to χ 2 min /n dof = 8.6/4.

•
The correlations, ρ i A vs. ρ f A and φ i A vs. φ f A , are shown by Figs.4(c) and 4(d), respectively.The main difference between case I and case II is that the ranges ρ i A 5.5 and φ i A ≈ −140 • or −220 • are excluded by the Bs → K * 0 K * 0 decay at 68% C.L. in case II.
[59]0]ct complicating the concept of B s decay observables is caused by the significant effects of time-dependent oscillation between Bs and B s states.Concerning the decays of Bs and B s mesons into a common final state, f = f , the untagged decay rate is the sum of the two time-dependent components, Γ[ Bs (t) → f ] + Γ[B s (t) → f ] [59,60], which yields the averaged and time-integrated branching ratio[59]

Table 1 :
Values of the input parameters: Wolfenstein parameters, pole and running quark masses, decay constants, form factors, Gegenbauer moments and decay width parameter y s .
(a).It can be clearly seen that the measured B( Bs → ρ 0 φ) presents a very stringent constraint on ρ H ; the large ρ H should obviously be ruled out.The fitted space for (ρ H , φ H ) is shown in Fig.1(b).H cannot be well determined due to the lack of data for the other observables; however, if ρ H 0.5, values of φ H around −180 • are favored.