Next-to-leading order QCD corrections to the production of $B_c$ and $B_c^*$ through $W^+$-boson decays

In this paper, we calculate the total decay widths for the $W^+$-boson decays, $W^+ \to B_c+b+\bar{s}+X$ and $W^+ \to B^*_c+b+\bar{s}+X$, up to next-to-leading order (NLO) accuracy within the framework of the nonrelativistic QCD theory. Both the fixed-order and the fragmentation approaches are adopted to do the calculation. Differential decay widths $d\Gamma/dz$ and $d\Gamma/ds_1$ are also given. We find that the NLO corrections are significant in those two $W^+$ decay channels. Our numerical results show that at the LHC, there are about $7.03\times 10^4$ $B_c$-meson events and $5.10\times 10^4$ $B^*_c$-meson events to be produced via the $W^+$-boson decays per operation year.

Besides the direct production mechanisms, the B c meson can also be indirectly produced through the topquark [27][28][29], the Z 0 -boson [30][31][32][33][34], the Higgs-boson [35,36] and the W -boson [37][38][39] decays.These indirect production channels can also generate abundant B c mesons at the LHC or the future high-energy colliders.The Wboson is the propagating media for the weak interaction, and the study on it is important for testing the SM.The LHC is a fruitful W -boson factory, there are about 3.07 × 10 10 W -bosons to be produced at the LHC per operation year [37].In the paper, we shall concentrate on the production of the B c meson and its first excited state B * c meson through the W -boson decays.The heavy constituent quarks move non-relativistically in the B c meson, and the processes involving the B c meson can be calculated within the nonrelativistic QCD (NRQCD) factorization formalism [40].Generally, the production cross-section or the decay width can be factorized into the product of the short-distance coefficients and the long-distance matrix elements.The short-tainty.This ambiguity introduces an important systematic error to pQCD predictions.It has been pointed out that one can use the higher-order β-terms to achieve an effective value of the strong running coupling α s , and the resultant conformal series is independent to the choice of renormalization scale and thus the conventional renormalization ambiguity is eliminated [41][42][43].The principle of maximum conformality (PMC) has been designed for such purpose [44][45][46][47][48], which provides a systemic way to eliminate the renormalization scheme-and-scale ambiguities simultaneously.The key idea of PMC is to set the correct momentum flow of the process by absorbing the non-conformal β-terms that govern the behavior of α s through the renormalization group equation (RGE).The β 0 -terms in the NLO coefficients can be adopted to set the α s value, thus in the paper, in addition to the conventional treatment, we shall also adopt the PMC to deal with the W + -boson decays, W + → B ( * ) c + b + s+ X.
In the decays, W + → B ( * ) c + b + s + X, the involved hard scales satisfy, m W ≫ m b , m c , so it is expected that the fragmentation mechanism dominates those decays.The NLO fragmentation functions for a heavy quark to a B c or B * c have recently been given by Ref. [23].It is interesting to apply those NLO fragmentation functions to the present processes, and compare the results from the fragmentation approach with those from the fixedorder approach.And in the present paper, besides the fixed-order approach, we shall also adopt the fragmentation approach to do the calculation.In usual cases, because the fragmentation probability of b → B c is about two orders of magnitude larger than that of c → B c , the b fragmentation is generally more important than the c fragmentation for the B c meson production.However, in cases with the W + -boson decays, the decay width for the b fragmentation shall be depressed due to the small values of the Cabibbo-Kobayashi-Maskawa (CKM) matrix elements |V cb | and |V ub |; Thus, it provides a good platform to test the fragmentation function for c → B c .
The paper is organized as follows.In Sec.II, we briefly present useful formulas at the LO accuracy under the fixed-order approach.In Sec.III, we present the formulas to calculate the NLO QCD corrections for the B ( * ) c meson production through the W -boson decays under the fixed-order approach.In Sec.IV, we present the useful formulas to calculate the decay widths under the fragmentation approach up to NLO accuracy.In Sec.V, numerical results and discussions are presented.Sec.VI is reserved as a summary.

II. DECAY WIDTHS AT THE LO LEVEL
According to the NRQCD factorization, the differential decay width for the B c -meson production from the W + -boson decays can be written as At the LO level, there are two Feynman diagrams for the (c b)[n]-pair production from the W + -boson decays, which are shown in Fig. 1.The LO amplitude for the (c b)[n] production can be written as the sum of two terms (M LO = M 1 + M 2 ) corresponding to two Feynman diagrams in Fig. 1, and we have where p 11 and p 12 are momenta of the c and b quarks in (c b)[n]-pair, where r c = m c /(m b + m c ). Π denotes the spin projector, for 1 S 0 state, and for 3 S 1 state, Λ 1 is color-singlet projector, and where 1 denotes the unit matrix of the color SU (3) group.
Using those amplitudes at the LO level, the LO decay width for (c b)-pair production can be calculated through where denotes the sum over the color and spin states of the initial and final particles, 1/3 comes from the spin average of the initial W + -boson.dΦ 3 denotes the differential phase space at the LO level, where d stands for the dimension of the space-time.With these formulas, the LO decay width for W + → (c b)[n] + b + s + X can be calculated directly.

III. THE NLO QCD CORRECTIONS
The NLO QCD corrections to the decay widths include virtual and real corrections.There are ultraviolet (UV) and infrared (IR) divergences in virtual correction, and IR divergence in real correction.We adopt the conventional dimensional regularization approach with d = 4 − 2ǫ to regulate these divergences.Then the UV and IR divergences appear as pole terms in 1/ǫ.We shall sketch the calculations for the virtual and real corrections in the following subsections.

A. The virtual correction
Four typical one-loop diagrams are shown in Fig. 2. The NLO virtual corrections come from the interference of those one-loop Feynman diagrams with the LO Feynman diagrams shown in Fig. 1.
The virtual corrections can be calculated through where M Virtual denotes the amplitude for the virtual corrections, dΦ 3 is the LO differential phase space.As a subtle point, there are Coulomb divergences in the hard part of the NLO amplitudes by using the traditional matching procedure.One may observe hat they appear in both the virtual corrections to the (c b)[n] production and the virtual corrections to the long-distance (2) [n] (n) , which shall be canceled by each other.As a result, no Coulomb divergence appears in the resultant pQCD series.In dimensional regularization, there is a simpler way to extract the NRQCD shortdistance coefficients using the method of regions [49].In this method, one can calculate the hard region contributions directly by expanding the relative momentum of the (c b)[n] pair before carrying out the loop integration (More explicitly, under the present lowest-order nonrelativistic approximation, one just needs to set q = 0 before the loop integration).Then the Coulomb divergences, which are power IR divergences, vanish in dimensional regularization.We adopt this new treatment, and the Coulomb divergences shall not appear in our present calculation.
There are UV and IR divergences in the loop-diagram contributions.The IR divergences in the virtual correction shall be canceled by the IR divergences in the real correction.The UV divergences should be removed through renormalization.We carry out the renormalization using counterterm approach, where the decay widths are calculated in terms of the renormalized quark mass m, the renormalized quark field Ψ r , the renormalized gluon field A µ r , and the renormalized coupling constant g s .The relations between the renormalized quantities and their corresponding bare quantities are where Z i = 1 + δZ i with i = m, 2, 3, g are renormalization constants, and they are fixed by the renormalization scheme.The renormalization scheme is adopted as follows: The renormalization of the heavy quark mass, the heavy quark field, the light quark field and the gluon field are performed in the on-shell scheme, whereas the renormalization of the strong coupling constant is performed in the MS scheme.The quantities δZ i can be calculated and they are where Q(= c, b) in the subscripts denotes a heavy quark, and q(= s) denotes a light quark.µ R is the renormalization scale, γ E is the Euler constant.For QCD, (3) (4) The real corrections come from the decay process ).The Feynman diagrams for the real corrections can be obtained through the LO Feynman diagrams by adding an additional gluon in the final state.Typical real correction Feynman diagrams are shown in Fig. 3.The real correc-tion can be calculated through where dΦ 4 denotes the differential phase space for the real corrections, and There are IR divergences in the real correction, which come from the phase-space integration.However, integrating the squared amplitude directly over the phase space in d dimensions is too difficult to be practical.In order to isolate the divergent and finite terms, we adopt the two-cutoff phase-space slicing method [50] to calculate the real correction.Following this method, the phase space for the real correction is decomposed into three regions by introducing two small cutoffs, δ s and δ c , which should satisfy the requirement δ c ≪ δ s [50].The three regions are soft region (S) with , where E 4 is defined in the rest frame of the initial W + boson.The IR finite hard-non-collinear part can be calculated numerically in four dimensions.
Applying the soft approximation to the soft part (e.g.all the terms of order δ s are neglected), we obtain where , and E 2 is defined in the rest frame of the W + boson.θ is the angle between p 2 and p 3 in the rest frame of the W + boson.
Applying the collinear approximation (e.g.terms of order δ c are neglected), we obtain where . Summing up three parts from the soft, hard-collinear and hard-non-collinear regions, we obtain the required real correction.Separate contributions from three regions depend on one or both of the two cutoff parameters δ s and δ c .However, the sum of those three contributions should be independent to the choices of δ s and δ c .The verification of this cut-off independence provides an important check for the correctness of the numerical calculation.We have checked this independence and have indeed found that the results are independent of the δ s and δ c by varying δ s from 10 −3 to 10 −7 with δ c = δ s /50.For clarity, we fix δ s = 10 −5 and δ c = 2 × 10 −7 to do the following numerical calculation.
Total NLO corrections can be obtained by summing up virtual and real corrections.The UV and IR divergences are exactly canceled after summing the real and virtual corrections, and the finite NLO corrections are obtained.
Then the decay widths , where n = 1 S 0 for B c and n = 3 S 1 for B * c , respectively.R S (0) denotes the radial wave function at the origin of the B c (B * c ) meson.In the calculation, we adopt the FeynArts package [51] to generate the Feynman diagrams and the corresponding amplitudes, and the FeynCalc package [52,53] to carry out the Dirac and color traces.Then we use the $Apart package [54] and the FIRE package [55] to do partial fraction and integration-by-parts (IBP) reduction of the loop integrals.After the IBP reduction, there are only few master integrals (e.g.A 0 , B 0 , C 0 , and D 0 functions) need to be calculated, which shall be dealt with by using the LoopTools package [56].Numerical phase-space integrations are carried out by the VEGAS program [57].

IV. DECAY WIDTHS UNDER THE FRAGMENTATION APPROACH
We take the process W + → B c + X as an example to illustrate the calculation under the fragmentation approach.The formulas for the B * c production are similar to the B c case.
The differential decay width for W + → B c + X under the fragmentation approach can be written as where d ΓW + →i+X (y, µ F ) denotes the decay width (coefficient function) for a W + to a parton i 1 , D i→Bc (z/y, µ F ) 1 Due to the coefficient function d ΓW + →i+X (y, µ F ) is IR safe, the denotes the fragmentation function for a parton i into a B c , and µ F is the factorization scale which separates the energy scales of two parts.
For comparison, we adopt several strategies to obtain the fragmentation predictions.More details about those strategies can be found in Refs.[23,58].For convenience, we denote them as "Frag, LO", "Frag, NLO" and "Frag, NLO+NLL", respectively.For the case of "Frag, LO", where Γ LO W + →c+s is the LO decay width for W + → c + s and D LO c→Bc (z, µ F ) is the LO fragmentation function.In the calculation, the factorization and renormalization scales are set as For the case of "Frag, NLO", where the NLO fragmentation function D NLO c→Bc (z, µ F ) can be found in Ref. [23].In the calculation, the factorization and renormalization scales are set as µ F = 2m b + m c and µ R = 2m b , and the factorization scheme is chosen as the MS scheme.
For the case of "Frag, NLO+NLL", dΓ Frag,NLO+NLL Q /m 2 W ). In the following fragmentation calculations, we shall adopt this approximation for simplicity.The neglected higher-power terms will be included in the results by combining the fixed-order and fragmentation approaches.
for the Mellin transformation method can be found in Ref. [65] 2 .

V. NUMERICAL RESULTS AND DISCUSSIONS
To do the numerical calculation, the input parameters are taken as follows: where G F is the Fermi coupling constant.The input value for |R S (0)| 2 is taken from the potential-model calculation [66].For the strong coupling constant, we use the two-loop formula , The decay widths for W + → B c +b+ s+X and W + → B * c + b + s + X under the fixed-order approach are given 2 In Ref. [65], the authors adopted an approximation that the sigularity of the NLO splitting function P Q→Q (z) at z = 1 is regularized by an overall "+" prescription.Then, 1 0 Pc→c(z) dz = 0.This approximation is good, since the probability to produce an additional heavy quark pair in the evolution process is very low at the energy of m W .For simplicity, we shall adopt this approximation in the present paper.
in Tables I and II, where Γ NLO = Γ LO + Γ Cor.
NLO denotes the NLO corrections, Γ Cor.NLO = Γ virtual + Γ Real .In the calculation, we take two typical energy scales (2m b and m W ) as the renormalization scale, and we have α s (2m b ) = 0.180 and α s (m W ) = 0.121.Tables I and II show that the NLO corrections are significant.After including the NLO corrections, the total decay width for W + → B c (B * c ) + b + s + X is increased by 69% (43%) for µ R = 2m b and 109% (92%) for µ R = m W .The momentum of the produced b-quark jet can be measured using vertex tagging technology, so the invariant mass of the B ( * ) c and b quark in the final state can be determined experimentally.We present the differential decay widths dΓ/ds 1 for W + → B c (B * c ) + b + s + X in Figs. 4 and 5, where µ R = 2m b and µ R = m W , respectively.Here s 1 ≡ (p 1 + p 2 ) 2 .Figs. 4 and 5 show that there is a peak near the minimum value of s 1 , indicating the dominant contributions of the decay processes come from the phase-space region near the threshold of producing the B c (B * c ) and b quark.This property may be helpful to distinguish the B c or B * c mesons produced through the W + -boson decays from the those produced through other production mechanisms at the high-energy colliders.
Tables I and II show that after including the NLO corrections, the renormalization scale dependence is softened.However, such scale dependence is still very large, e.g.W + → B c (B * c ) + b + s + X, the NLO total decay width decreases by 44% (39%), when µ R changes from 2m b to m W .
As mentioned in the Introduction, the PMC scalesetting approach provides a way to eliminate the renormalization scale ambiguity [44][45][46][47][48].As an attempt of showing how the PMC affects the decay width, we present the PMC predictions in the following.
To apply the PMC, we first schematically rewritten the NLO decay width as Using the RGE, the non-conformal term (− 3b 2 β 0 ) can be adopted to fix the strong running coupling.A PMC scale µ PMC is then determined, which corresponds to the (correct) typical momentum flow of the process.Then, following the standard PMC procedures, the NLO decay width changes to where µ PMC = µ R e 3b/2 .It is interesting to find that the PMC scale µ PMC is independent to any choice of renormalization scale µ R , e.g.µ PMC ≡ 6.67 GeV for B c and µ PMC ≡ 7.17 GeV for B * c , thus the conventional renormalization scale ambiguity is really eliminated.The PMC scales are closer to µ R = 2m b than µ R = m W , the conditions for the total decay widths are similar, thus the usual guessing choice of µ R = 2m b is more reasonable for conventional prediction.Thus in the following analysis, we fix µ R = 2m b for predictions when use conventional pQCD series.Numerical results for the total decay widths of W + → B c (B * c ) + b + s + X up to NLO accuracy under the PMC are shown in Table III.After applying the PMC scalesetting approach, the convergence is slightly better than conventional series, e.g. after including the NLO QCD corrections, the decay width for W + → B c (B * c )+b+s+X is increased by 58% (32%).It is interesting to know the differential distributions of those decay processes.We define the energy fraction , where E 1 and E max 1 are the energy and the maximum energy of the B c (B * c ) meson in the rest frame of the initial W + boson.The differential decay widths dΓ/dz for W + → B * c + b + s + X are presented in Figs. 6 and 7.In addition to the fixed-order results, we also present the results from the fragmentation approach up to NLO level.Here, in order to know whether the fragmentation mechanism dominates the decay processes, we do not resum the leading logarithms of m 2 Q /m 2 W , i.e., the coefficient function and the fragmentation functions are both calculated at the LO or NLO level without the DGLAP evolution.The factorization and renormalization scales are set as 2m b + m c and 2m b respectively in the fragmentation calculation, and the renormalization scale is set as 2m b in the fixed-order calculation.FIG. 8.The NLO differential decay width dΓ/dz for W + → Bc+b+s+X under the fixed-order approach (FO), fragmentation approach (Frag) and the combination of two approaches (Combined), respectively.FIG. 9.The NLO differential decay width dΓ/dz for W + → B * c + b + s+ X under the fixed-order approach (FO), fragentation approach (Frag) and the combination of two approaches (Combined), respectively.
Figs. 6 and 7 show that the fragmentation mechanism dominates the decay processes W + → B c (B * c )+b+ s+X, since the fixed-order and fragmentation shapes are close at the LO and NLO levels.Differences appear in small z region, indicating in this z-region, the non-fragmentation  terms become important.The neglected logarithms of m 2 Q /m 2 W may give sizable contributions, which can be resumed in fragmentation approach by using the DGLAP evolution equation.More explicitly, the NLO fragmentation results without or with resummation labelled as "Frag, NLO" and "Frag, NLO+NLL" are presented in Figs. 8 and 9, which are calculated by using Eqs.(20,21), respectively.The choices of the factorization and renormalization scales for the "Frag, NLO" and "Frag, NLO+NLL" calculations have been given below Eq.( 20) and Eq.( 21) accordingly.Here the NLO fixed-order results labelled as "FO, NLO" are presented as a comparison.The decay widths for "Frag, NLO" and "Frag, NLO+NLL" are presented in Table IV.One may observe that by resuming the next-to-leading logarithms of m 2 Q /m 2 W , more accurate behavior in large z region can be achieved, and the total decay widths shall be reduced by about 2% for both B c and B * c productions.In the fixed-order prediction, the large logarithms of m 2 Q /m 2 W may appear in specific kinematic region, and the fragmentation approach provides us a way to give a reasonable contribution in this region by resuming all the large logarithms.Thus a combination of those two approaches may be helpful.As an attempt, we combine the NLO results from the fixed-order and fragmentation approaches in the following way, Here, the choices of the factorization and renormalization scales for the "Frag, NLO" and "Frag, NLO+NLL" calculations have been given below Eq.( 20) and Eq.( 21) accordingly, and the renormalization for the "FO" calculation is set as µ R = 2m b .The differential decay widths dΓ/dz for "Combined, NLO+NLL" are also presented in Figs. 8 and 9.The total decay widths are Γ Combined,NLO+NLL Γ Combined,NLO+NLL The combined results are the definitive results in the presented paper because the large logarithms of m 2 Q /m 2 W have been resummed and the theoretical uncertainties are reduced compared to the fixed-order results.

C. Uncertainty analysis
In this subsection, we shall estimate the theoretical uncertainties for these decay widths.The main uncertainty sources for these decay widths include the heavy quark masses (m c and m b ), the renormalization scale and the radial wave function at the origin |R S (0)|.|R S (0)| 2 is an overall factor in the calculation, the uncertainty due to it can be figured out easily.So we shall not consider the uncertainty of |R S (0)| 2 and concentrate our attention on the uncertainties from the heavy quark masses and the renormalization scale.For clarify, when considering the uncertainty caused by one input parameter, the other input parameters are fixed to their center values.We first consider the uncertainties caused by the heavy quark masses.We estimate the uncertainties by varying the heavy quark masses with m c = 1.5 ± 0.1 GeV and m b = 4.9 ± 0.2 GeV.The uncertainties caused by m c and m b are presented in Tables V and VI respectively.For comparison, the fixed-order results and the combined NLO+NLL results are presented explicitly.From the two tables, we can see that the decay widths decrease with the increment of m ).For the fragmentation contribution, there are lower factorization and renormalization scales, and upper factorization and renormalization scales involved in the calculation.The lower scales should be O(m Q ) and the upper scales should be O(m W ). Since we adopt the approximation used in Ref. [65] where  The results in Tables VII and VIII show that the combined NLO+NLL results are more sensitive to the lower renormalization scale than the upper renormalization scale.Comparing the results in Tables VII and VIII with those in Tables I and II, we can see that after resumming the large logarithms of m 2 Q /m 2 W , the uncertainties caused by the choice of the renormalization scale are reduced explicitly.

VI. SUMMARY
In the present paper, we have calculated the W + -boson decays, W + → B c (B * c ) + b + s + X, up to NLO QCD corrections under the NRQCD framework.Both the fixed-order and fragmentation approaches are adopted for the calculations.The theoretical uncertainties are analyzed through varying the heavy quark masses and the renormalization scale.Our results show the NLO corrections are significant.Under conventional scale-setting approach, the decay widths for W + → B c (B * c )+b+ s+X shall be increased by 69% (43%) for the case of µ R = 2m b after including the NLO corrections.The scale dependence can be suppressed after including the NLO corrections, even though it is still large.By using the PMC, we show that the renormalization scale ambiguity can be eliminated for those two decay processes; thus they provide another successful applications of the PMC.
The differential distributions dΓ/dz for the W + -boson decays W + → B c (B * c ) + b + s + X have been given under the fixed-order and the fragmentation approaches, respectively.Our results show that both decay processes are dominated by the fragmentation mechanism, and the differences exist in small z region.By combining the fixed-order prediction with the fragmentation approach to resum the leading and next-to-leading logarithms of m 2 Q /m 2 W , a more accurate distribution to the fixed-order prediction can be achieved.
The total W -boson decay width is Γ W = 2.09GeV [67].Using the combined results (25,26) from the fixed-order and fragmentation approaches, we obtain the branching fractions for the considered channels, e.g.
If the LHC runs with a luminosity of 10 34 cm −2 s −1 [37], the expected W + -boson events is about 3.07×10 10 events per operation year.Then, there are about 7.03 × 10 4 B c mesons and 5.10×10 4 B * c mesons to be produced through the W + -boson decays per operation year.Two 2S-level excited states B c (2 1 S 0 ) and B c (2 3 S 1 ) have recently been observed by CMS and LHCb collaborations [68,69].Using |R 2S (0)| 2 = 0.983 GeV 3 [66], there are about 4.21 × 10 4 B c (2 1 S 0 ) and 3.05 × 10 4 B * c (2 3 S 1 ) to be produced through the W + -boson decays per operation year.Those excited states can also be important sources for the ground-state B c meson.Thus by carefully measuring the W + -boson decays, we may have a good chance to study the B c meson properties.
where d ΓW + →(c b)[n]+b+s+X denotes the decay width for the production of a perturbative state (c b)[n] with quantum numbers [n].The long-distance matrix element O Bc (n) is the transition probability for a (c b)[n] state to the B c meson.In the lowest-order nonrelativistic approximation, only color-singlet contributions need to be considered, and the long-distance matrix elements for the color-singlet contributions can be determined through potential models.Practically, we first calculate the decay width for an on-shell (c b)-pair, i.e., dΓ W + →(c b)[n]+b+s+X .Then the decay width for the B c meson, i.e., dΓ W + →Bc+b+s+X , can be obtained from dΓ W + →(c b)[n]+b+s+X by replacing

2 F /m 2 W or µ 2 R /m 2 W
where the upper factorization and renormalization scales in the coefficient function d ΓNLO W + →c+X (y, µ F )/dy are set as µ F = µ R = m W , so as to avoid the large logarithms of µ appear in the coefficient function.The fragmentation functionD NLO+NLL c→Bc (z, µ F = m W ) isobtained through solving the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution equation [59-61] with NLO splitting function for c → c [62-64], where the NLO fragmentation function D NLO c→Bc (z, µ F 0 = 2m b + m c ) with the lower factorization and renormalization scales µ F 0 = 2m b + m c and µ R0 = 2m b is used as the boundary condition.To solve the DGLAP evolution equation, the Mellin transformation method is adopted, and the related formulas heavy-quark mass m Q in the coefficient function can be approximately set to 0, and this approximation brings only a small error of O(m 2
c and m b , but they are more sensitive to m b than m c .Adding these two uncertainties caused by m c and m b in quadrature, we obtain the uncertainties caused by the heavy quark masses Γ FO,NLO W + →Bc+X = 4.89 +0.82 −0.65 keV, Γ FO,NLO W + →B * c +X = 3.55 +0.59 −0.46 keV, (27) and Γ Combined,NLO+NLL W + →Bc+X = 4.78 +0.80 −0.63 keV, Γ Combined,NLO+NLL W + →B * c +X = 3.47 +0.58 −0.45 keV.(28)Then we consider the uncertainties caused by the renormalization scale and the factorization scale.For the fixed-order results under the conventional scale setting, the uncertainties caused by the renormalization scale can be estimated by varying the renormalization scale between the two typical energy scales 2m b and m W .The fixed-order results with µ R = 2m b and µ R = m W have been presented in Tables I and II, i.e., Γ FO,NLO W + →Bc+X ∈ [2.72, 4.89] keV and Γ FO,NLOW + →B * c +X ∈ [2.15, 3.55]  keV with a variation of µ R ∈ [2m b , m W ]. These uncertainties caused by the choice of the renormalization scale in the fixedorder results are big.For the combined NLO+NLL results, there are several renormalization and factorization scales involved in the calculation.We can decompose the combined NLO+NLL results into fragmentation contribution dΓFrag,NLO+NLL

1 0
dzP c→c (z) = 0, the fragmentation probabilities for c → B c (B * c ) do not change under the DGLAP evolution.Thus, for the total decay widths, we need not to consider the uncertainties caused by the lower and upper factorization scales, and only need to consider the uncertainties caused by the lower and upper renormalization scales.For the power-correction contribution, the renormalization scales in dΓ FO,NLO same.We take the renormalization scales in the power-correction contribution as the same as the lower renormalization scale in the fragmentation contribution for simplicity.In order to estimate the uncertainties for the combined NLO+NLL results, we vary the lower and upper renormalization scales by a factor of 2 from their center values.The decay widths with lower renormalization scale µ R0 = 2m b and µ R0 = 4m b are presented in

TABLE II .
[67]he two-loop coefficient of the QCD β-function.According to α s (m Z ) = 0.1185[67], we obtain Λ The total decay widths of W + → B * c + b + s + X up to NLO level under the fixed-order approach.Two typical renormalization scales are adopted.
TABLEI.The total decay widths of W + → Bc + b + s + X up to NLO level under the fixed-order approach.Two typical renormalization scales are adopted.

TABLE IV .
The total decay widths (in unit: keV) for W + → Bc(B * c ) + b + s + X under the fragmentation approach.

TABLE V .
The decay widths (in unit: keV) for W + → Bc(B * c ) + b + s + X with a variation of mc = 1.5 ± 0.1 GeV.

Table VII ,
and the decay widths with upper renormalization scale µ R = m W and µ R = m W /2 are presented in Table VIII.More explicitly, we obtain Γ Combined,NLO+NLL R ∈ [m W /2, m W ].

TABLE VIII .
The "combined, NLO+NLL" decay widths (in unit: keV) with a variation of the upper renormalization scale µR.