{'fraction_non_alphanumeric': 0.043716684026714046, 'fraction_numerical': 0.02758715764965382, 'mean_word_length': 3.65224827193045, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 19, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Weakly-bound Hadronic Molecule near a 3-body Threshold

The X(3872) seems to be a loosely-bound hadronic molecule whose constituents are two charm mesons. A novel feature of this molecule is that the mass difference of the constituents is close to the mass of a lighter meson that can be exchanged between them, namely the pion. We analyze this feature in a simple model with spin-0 mesons only. Various observables are calculated to next-to-leading order in the interaction strength of the exchanged meson. Renormalization requires summing a geometric series of next-to-leading order corrections. The dependence of observables on the ultraviolet cutoff can be removed by renormalizations of the mass of the heaviest meson, the coupling constant for the contact interaction between the heavy mesons, and short-distance coefficients in the operator product expansion. The next-to-leading order correction has an unphysical infrared divergence at the threshold of the two heavier mesons that can be eliminated by a further resummation that takes into account the nonzero width of the heaviest meson.


I. INTRODUCTION
The X(3872) is a new hadronic resonance discovered in 2003 by the Belle Collaboration [1] and subsequently confirmed by the CDF, Babar, and D0 Collaborations [2][3][4]. All the properties of the X(3872) that have been measured thus far are compatible with its identification as a weakly-bound molecule whose constituents are a superposition of the charm meson pairs D * 0D0 and D 0D * 0 .
Because it is so weakly bound, the X(3872) has many features in common with the deuteron, which is a weakly-bound baryonic molecule consisting of a proton and a neutron. Their binding energies are both small compared to the natural energy scale associated with pion exchange: m 2 π /(2M 12 ), where M 12 is the reduced mass of the two constituents. The binding energy 2.2 MeV of the deuteron is small compared to the natural scale of about 20 MeV. After taking into account a recent precision measurement of the D 0 mass by the CLEO Collaboration [26], the difference between the measured mass of the X(3872) and the D * 0D0 threshold is M X − (M 1 + M 2 ) = −0.6 ± 0.6 MeV.
This is small compared to the natural energy scale of about 10 MeV. The Belle Collaboration has set an upper bound on the width of the X [1]: This is also small compared to the natural energy scale. The small width is most easily understood if the X is below the D * 0D0 threshold. Because the binding energies of the deuteron and the X(3872) are small compared to the natural energy scales, they have universal properties that are determined by the large scattering length a of the constituents [10]. The universality of few-body systems with a large scattering length has many applications in atomic, nuclear, and particle physics [27]. The universal features of the X(3872) were first exploited by Voloshin to describe its decays into D 0D0 π 0 and D 0D0 γ, which can proceed through decay of the constituent D * 0 orD * 0 [8]. Universality has also been applied to the production process B → KX [14,17], to the line shape of the X [18], and to decays of X into J/ψ and pions [20]. These applications rely on factorization formulas that separate the length scale a from all the shorter distance scales of QCD [18]. The factorization formulas can be derived using the operator product expansion for a low-energy effective field theory [24].
There is one feature of the X(3872) that is very different from the deuteron. The mass difference between the constituents D * 0 andD 0 , M 2 − M 1 = 142.12 ± 0.07 MeV, is very close to the π 0 mass: m = 134.98 MeV. Thus the D * 0 mass is very close to the D 0 π 0 threshold: M 2 − M 1 − m = 7.14 ± 0.07 MeV. This is important because the pion is the lightest meson that can be exchanged between the D * 0 andD 0 . The small energy gap between the D * 0 mass and the D 0 π 0 threshold makes it easy for the D * 0 to make a transition to D 0 π 0 . An important consequence is that the X could have a substantial 3-body component D 0D0 π 0 . Suzuki has argued that the near equality of M 2 − M 1 and m implies that the pion-exchange interaction is too weak to bind the X [21], but this conclusion has been criticized [23].
An analysis of the effect of the near equality M 2 − M 1 ≈ m in the case of the X(3872) is complicated by other effects that may be comparable in importance. There are other nearby two-body and three-body thresholds. The thresholds for the charged mesons D * + D − and D + D * − are higher than the D * 0D0 threshold by 8.08±0. 12 MeV. The D 0 D − π + and D +D0 π − thresholds are higher than the D * 0D0 threshold by only 2.23 ± 0.12 MeV. There are several other features that complicate the analysis of the effects of pions on the X. The constituent D * 's are spin-1 particles. The charge conjugation C = + requires the molecule to be a superposition of charm mesons: D * 0D0 + D 0D * 0 . The interaction of the charm mesons with pions depends on the pion momentum. A quantitative analysis of the effects of pions on the X might require all these effects to be taken into account.
Since the new feature of the X(3872) is the near equality between the mass difference M 2 − M 1 of the constituents and the mass m of a meson that can be exchanged between them, it is worthwhile to analyze this feature in a simple model that avoids all the other complications of the X. The simplest such model is one in which the primary constituents D 1 and D 2 are spin-0 mesons with a large positive scattering length and a momentumindependent coupling to a lighter spin-0 meson φ. We will calculate the effects of the exchanged meson in this model to next-to-leading order in the D 2 − D 1 φ coupling constant g.
In Sec. II, we write down the lagrangian for the model and discuss its parameters. In Sec. III, we summarize the results of the model at 0 th order in g that were derived previously in Refs. [18,24]. The results include the short-distance production rate for X and the line shape of X in a short-distance decay channel. The ultraviolet divergences can be removed by renormalization of the coupling constant for the D 1 D 2 contact interaction and by renormalization of short-distance coefficients in the operator product expansion. In Sec. IV, we calculate the observables in Sec. III to next-to-leading order in g. We also calculate the short-distance production rate for D 1 D 1 φ and the decay rate of X into D 1 D 1 φ. There are new ultraviolet divergences that can be removed by renormalization of the D 2 mass. We show that after summing a geometric series of next-to-leading order corrections, all remaining ultraviolet divergences can again be removed by renormalization of the coupling constant for the D 1 D 2 contact interaction and by renormalization of short-distance coefficients in the operator product expansion. A further resummation is required to eliminate an infrared divergence at the D 1 D 2 threshold in the next-to-leading order corrections. In Sec. V, we summarize our results. The two-loop Feynman diagrams that arise at next-to-leading order in g are calculated in an Appendix.

A. Lagrangian
We consider the simplest possible model with a bound state X whose mass is close to a 3-body threshold. The primary constituents D 1 and D 2 of the bound state are scalar mesons whose masses M 1 and M 2 satisfy M 1 < M 2 . There is also a scalar meson φ whose mass m is close to the mass difference M 2 − M 1 . We will refer to this model as the scalar meson model. It is defined by a nonrelativistic quantum field theory with three complex scalar fields D 1 , D 2 , and φ. The free terms in the lagrangian 1 are 1 The superscript † on a field indicates its complex conjugate. The interaction terms for the scalar meson model are The interaction vertices are illustrated in Fig. 1. The coupling constants λ 0 and g have mass dimensions −2 and −1/2, respectively. We assume that there is a fine-tuning that makes the D 1 D 2 scattering length large, which implies that the D 1 D 2 contact interaction with coupling constant λ 0 must be treated nonperturbatively. The subscript on λ 0 emphasizes that it is a bare coupling constant that requires renormalization. We assume that the D 2 − D 1 φ interaction can be treated perturbatively and we calculate its effects through order g 2 . At this order, no renormalization of g is required. However, we will see that nonperturbative renormalization of λ 0 requires summing a geometric series to all orders in g 2 . There are also ultraviolet divergent corrections to the mass of the D 2 which are cancelled by the mass counterterm δM in Eq. (4).

B. Masses
It is convenient to introduce concise notations for various combinations of the masses. The sum and difference of the D 1 and D 2 masses are The reduced masses of D 1 D 2 and of D 1 φ are If m is much smaller than M 1 , the natural energy and 3-momentum scales associated with exchange of the meson φ between the mesons D 1 and D 2 are m 2 /M 12 and m, respectively. We assume that the energy gap between the D 2 mass and the D 1 φ threshold is small compared to the natural energy scale: The total number of heavy mesons D 1 and D 2 is conserved by the interactions in Eq. (4). We are interested in the sector of the theory with two heavy mesons and, more specifically, in the threshold region where the invariant mass M of all the particles is very close to M 1+2 : This restricts the possible scattering states to D 1 D 2 and D 1 D 1 φ.
For numerical illustrations, we will use masses that correspond to the D * 0D0 /D 0D0 MeV. The condition in Eq. (7) is only marginally satisfied: M 12 (M 2−1 − m)/m 2 ≈ 0.4. We take the binding energy of X to be E X = 0.4 MeV, which is near the central value of the measurement of M X in Eq. (1). If we express this binding energy as E X = κ 2 X /(2M 12 ), the binding momentum is κ X = 27.8 MeV.

C. Coupling constants
We impose an ultraviolet cutoff on the momenta p of particles in the center-of-momentum frame that restricts them to the region | p | ≪ m. This guarantees that all the particles remain nonrelativistic. It was shown in Ref. [24] that at leading order g 0 , all ultraviolet divergences can be absorbed into λ 0 and into short-distance coefficients in the operator product expansion. At order g 2 , there are new ultraviolet divergences that can be absorbed into the D 2 mass. We will find that the remaining ultraviolet divergences can again be absorbed into λ 0 and into short-distance coefficients in the operator product expansion.
The coupling constant λ 0 must be tuned to near the critical value at which the D 1 D 2 scattering length diverges in order for the bound state X to have a small binding energy E X ≪ m 2 /M 12 . In order to compare results at order g 0 and order g 2 , we must have a welldefined renormalization prescription for λ 0 or, equivalently, a clear definition of the binding energy E X . We will find that all the amplitudes for processes near the D 1 D 2 threshold have a resonant term with a common factor A(E) that depends on the energy E in the rest frame. For example, the line shape in a short-distance decay channel of X produced by a short-distance process has a long-distance factor |A(M)| 2 , where M is the invariant mass of the decay products. The square of the resonant factor can be expressed in the form As a function of the real energy E, this resonance factor has a peak when E is close to M 1+2 −E X . The peak arises because Re A −1 (E) vanishes near that value of E and Im A −1 (E) is much smaller than Re A −1 (E) except near that value of E. As our renormalization prescription for λ 0 , we demand that Re A −1 (E) vanishes at an energy determined by a real variable κ X : The dependence of an observable on the bare parameter λ 0 can be eliminated in favor of the renormalization parameter κ X . We assume that κ X is small compared to the natural momentum scale: κ X ≪ m.
The amplitude A(E) has a pole at a complex energy E pole near the real energy at which |A(E)| 2 has its maximum. We choose to define the mass M X and the width Γ X of X by expressing that complex energy in the form We also choose to define the binding energy E X by We assume that the binding energy and the width of X are both small compared to the natural energy scale: E X , Γ X ≪ m 2 /M 12 . In Ref. [24], the authors proposed a renormalization scheme in which λ 0 was tuned to give a prescribed value for the complex parameter E pole . The binding energy E X and the width Γ X were taken as the input parameters that determine E pole . This renormalization scheme requires λ 0 to be a complex coupling constant. The imaginary part of λ 0 takes into account inelastic scattering channels for D 1 D 2 that are outside the energy range described by the effective field theory. In this paper, we assume that any such inelastic channels do not give a significant contribution to the width of X. The coupling constant λ 0 is therefore real valued. We could use the binding energy E X as the input parameter and eliminate λ 0 in favor of E X , and then Γ X would be calculable. We find it more convenient to eliminate λ 0 in favor of the variable κ X introduced in Eq. (10). In this renormalization scheme, E X and Γ X are both calculable in terms of κ X , g, and the masses M 1 , M 2 , and m.
The most convenient input for determining the coupling constant g is the width Γ 2 of the heavy meson D 2 . Since M 2 > M 1 + m, the D 1 D 2 φ vertex in Fig. 1(b) allows D 2 to decay into D 1 φ. Calculating the width using nonrelativistic phase space, we obtain where M 12 and m 1 are the reduced masses defined in Eqs. (6). The value of the coupling constant g can be determined from Γ 2 and the masses. We assume that the width of D 2 is small compared to the natural energy scale: Γ 2 ≪ m 2 /M 12 . We also assume Γ 2 is small enough to allow perturbation theory in the coupling constant g.
For numerical illustrations, we will use values of the parameters that correspond to the D * 0D0 /D 0D0 π 0 system. We take the width of D 2 to be Γ 2 = 0.07 MeV, which is approximately the total width for D * 0 [18]. This is much smaller than the natural energy scale m 2 /M 12 and also rather small compared to M 2−1 − m: Γ 2 /(M 2−1 − m) ≈ 0.01. Using Eq. (13), we determine the numerical value of the D 1 D 2 φ coupling constant to be g = 0.0064 MeV −1/2 . When considering the dependence of observables on g, it will be more convenient to regard them as functions of Γ 2 .
We now list all the energies that are assumed to be small compared to the natural scale m 2 /M 12 : • the energy M − M 1+2 of a D 1 D 2 or D 1 D 1 φ system relative to the D 1 D 2 threshold, • the energy difference M 2−1 − m between the D 2 mass and the D 1 φ threshold, • the decay width Γ 2 of D 2 , which is given in Eq. (13), • the binding energy E X of the molecule, which is approximately equal to κ 2 X /(2M 12 ), where κ X is a renormalization parameter.
We will carry out our calculations under the assumption that these energies are all comparable and that the masses M 2 , M 1 , and m are also comparable. The only approximations we will make are those that can be justified by a nonrelativistic approximation, such as At order g 2 , renormalization of the D 2 mass is necessary. In the on-shell renormalization scheme, the parameter M 2 in the free lagrangian in Eq. (3) is chosen to be the physical D 2 mass. The physical mass M 2 and the width Γ 2 of D 2 can be defined by specifying that the exact D 2 propagator at 0 momentum has a pole in the energy at p 0 = M 2 − iΓ 2 /2. The pole in the energy p 0 at 0 momentum is the solution to the equation where Σ(p 0 , p) is the D 2 self-energy. The exact D 2 self-energy is the sum of the one-loop diagram in Fig. 2(a) and the mass counterterm in Fig. 2 There are no other diagrams for the self-energy at higher order in g. The one-loop diagram is calculated in Section A 2 of the Appendix. The analytic expression for the self-energy diagram using dimensional regularization is given in Eq. (A9). The result in a general regularization scheme is The term Σ 2a (M 1 + m, 0) is real valued and includes a linear ultraviolet divergence. The divergence is cancelled by a divergent term in the mass counterterm δM in Eq. (15). In the on-shell renormalization scheme, the finite terms in δM are chosen so that the solution to Eq. (14) satisfies Re p 0 = M 2 . The wavefunction normalization factor for a D 2 with momentum p is Since Eq. (14) cannot be solved analytically, the on-shell renormalization scheme is not the most convenient choice. A simpler renormalization scheme for the D 2 mass is to choose the mass counterterm so the self-energy Σ 2a (p 0 , 0) vanishes at the D 1 D 2 threshold: δM = −Σ 2a (M 1 + m, 0). The self-energy is then With this renormalization prescription, the parameter M 2 differs at order g 4 from the physical D 2 mass, which is given by the real part of the solution to Eq. (14). The width Γ 2 , which is given by the imaginary part of the solution to Eq. (14), reduces at leading order in g to Γ 2 = −2 Im Σ 2a (M 2 , 0). This agrees with the explicit expression for the decay rate of D 2 at order g 2 obtained in Eq. (13). To order g 2 , the wavefunction normalization factor in Eq. (17) is Since the momentum p is restricted to the nonrelativistic region, the term p 2 /(2(M 1 +m)M 2 ) can always be neglected compared to 1. The resulting wavefunction normalization factor is independent of p: The correction term is pure imaginary.

E. The Complex Mass Scheme
The D 2 − D 1 φ interaction allows a virtual D 2 whose energy is greater than M 1 + m to decay into D 1 φ. In the amplitude for such a process at leading order in g, the propagator of the virtual D 2 diverges. Near such a point in momentum space, the corrections to the D 2 propagator from higher orders in g are not suppressed. This problem can be solved by summing the D 2 propagator corrections to all orders. The effect is to add −Σ(p 0 , p) to the denominator p 0 − p 2 /(2M 2 ) + iε of the D 2 propagator. The imaginary part of the self-energy Σ(p 0 , p) eliminates any divergence in the D 2 propagator for real values of the energy. One drawback of this solution is that it leads to unnecessarily complicated expressions for the integrals in loop diagrams with D 2 propagators.
Summing the D 2 self-energy corrections to all orders is essential only near the pole in the propagator. In other regions of the energy p 0 , Σ(p 0 , p) is a perturbative correction to the denominator of the propagator that is suppressed by a power of g 2 . In those regions, summing the D 2 propagator corrections to all orders is optional. At the pole in p 0 , the self-energy Σ(p 0 , p) reduces to −iΓ 2 /2. An alternative resummation that eliminates the divergence in the propagator for real values of the energy without changing the order in g of the truncation error is to sum only the −iΓ 2 /2 term in the self-energy to all orders and treat the remainder of the self-energy as a perturbation. An advantage of this partial resummation is that it leads to simpler expressions for the integrals in loop diagrams with D 2 propagators.
A systematic method for implementing this partial resummation is the complex mass scheme [28]. This method has been used in Standard Model calculations at next-to-leading order [29][30][31]. The application of this method to the scalar meson model involves adding cancelling terms to the free and interaction terms in the Lagrangian: where Γ 2 is the D 2 width and Z 2 is the D 2 wavefunction normalization factor. To order g 2 , Γ 2 and Z 2 are given by Eqs. (13) and (20). Since iΓ 2 /2 is pure imaginary and Z 2 is complex, the free Lagrangian consisting of the sum of Eqs. (3) and (21a) corresponds to a nonhermitian Hamiltonian. The effects of the nonhermiticity are cancelled exactly if the interaction terms in Eq. (21b) are calculated to all orders. At any finite order in perturbation theory, the complex mass scheme gives amplitudes that are uniformly accurate to the appropriate order in g, even if there is a virtual D 2 that decays into D 1 φ. The effect of adding Eq. (21a) to the free Lagrangian in Eq. (3) is to change the D 2 propagator: .
The effect of adding Eq. (21b) to the interaction Lagrangian in Eq. (4) is to change the counterterm vertex in Fig. 1: It is convenient to also rescale the fields D 2 (x) and D † 2 (x) by the same complex factor Z 1/2 2 . This eliminates the factor of Z 2 from the numerator of the D 2 propagator in Eq. (22). It also multiplies all the interaction terms in Eqs. (4) and (21b) by Z 2 or Z 1/2 2 . Since we will use nonperturbative renormalization for the coupling constant Z 2 λ 0 for the D 1 D 2 contact interaction, the factor of Z 2 only affects the value of the unphysical bare coupling constant λ 0 . Since Z 2 = 1+O(g 2 ), the other factors of Z 2 or Z 1/2 2 in the interaction vertices contribute only at order g 4 and higher. We will be calculating only to order g 2 , so we will ignore these factors of Z 2 and Z 1/2 2 . The perturbation expansion in the complex mass scheme corresponds to ordinary perturbation theory in g followed by a resummation to all orders of a subset of terms that are higher order in g. We will refer to calculations to order g 0 in this perturbation expansion as leading order (LO) in the complex mass scheme. We will refer to calculations through order g 2 in this perturbation expansion as next-to-leading order (NLO) in the complex mass scheme.

F. Propagator for X
We can use the scalar meson model to describe processes with asymptotic states that include not only the particles D 1 and φ, but also the bound state X. The local composite operator D † 1 D † 2 (x) has a nonzero amplitude to create X from the vacuum. Thus can be used as an interpolating field for X. It is more convenient to use the operator λ 0 D 1 D 2 (x), because this is a renormalized operator whose matrix elements do not depend on the ultraviolet cutoff [24]. The corresponding propagator for X is The momentum 4-vector in the Fourier transform is P µ = (E, P ), where E and P are the energy and momentum of X. At P = 0, this propagator has a pole at the complex energy E pole . Near the pole in the energy, the behavior of the propagator at P = 0 is The energy E pole determines the binding energy E X and the width Γ X of X through Eq. (11). The residue of the pole defines the wavefunction normalization factor Z X , which is complex valued. Because the propagator ∆ X (E, P ) has mass dimension −2, Z X has mass dimension −1.
Strictly speaking, since X has a nonzero width Γ X , there are no T-matrix elements for processes involving X in the initial or final state, because the decay of the X prevents it from being a truly asymptotic state. However if Γ X is sufficiently small, the X can propagate over a long time interval before decaying, and it can therefore be treated as a quasi-asymptotic state. We can use the Lehmann-Symanzik-Zimmermann (LSZ) formalism to define T-matrix elements for the X that are closely related to T-matrix elements for decay products of the X whose total energy is tuned to the peak of the X resonance. To define the T-matrix element T for a process with X in the initial or final state, we start with a Green's functions for the operator The connected Green's function at 4-momentum (E, P ) is amputated by multiplying by the inverse propagator ∆ X (E, P ) −1 , it is normalized by multiplying by Z 1/2 X , and then it is evaluated at the real energy M X + P 2 /(2M X ). This gives the T-matrix element multiplied by i for a state X with the standard nonrelativistic normalization. The T-matrix element T for a state X with the standard relativistic normalization is obtained by multiplying by √ 2M X . To justify this prescription, we consider a process A → B + C, where A and B both represent one or more particles and C denotes a set of particles that can be decay products of X. The T-matrix for this process will have a resonant enhancement when the invariant mass M of the particles in C is near the mass M X . The resonant contribution to the T-matrix element for A → B + C in the rest frame of C can be approximated by Our prescription for T-matrix elements such as T [A → B + X] and T [X → C] involves evaluating the amputated connected Green's function at the real energy E = M X . This ensures that the expression in Eq. (26) is as good an approximation as possible to the T-matrix element T [A → B + C] near the peak of the resonance. As a simple illustration of the prescription for T-matrix elements for processes involving X, we consider the forward scattering process X → X. The relevant Green's function in the rest frame of X is the negative of the propagator ∆ X (E, 0). This must be multiplied by two factors of ∆ X (E, 0) −1 Z 1/2 X , one for the X in the inital state and one for the X in the final state. The T-matrix element iT [X → X] is then obtained by evaluating this at E = M X and multiplying it by 2M X : G. Short-distance production and decay processes The scalar meson model defined by the lagrangian in Eqs. (3) and (4) can be a low-energy approximation to a more fundamental quantum field theory. The fundamental theory may include high energy processes that can create D 1 D 2 or D 1 D 1 φ with invariant mass near M 1+2 . If long-distance effects involving momenta much smaller than m can be separated from short-distance effects involving momenta of order m and larger, we can use the scalar meson model to calculate the long-distance effects. The basic tool required to separate long-distance effects from short-distance effects is the operator product expansion [24]. Applications of the operator product expansion to the scalar meson model at 0 th order in g were described in Ref. [24]. Below, we give the leading terms in the operator product expansions of the T-matrix elements for some of the processes considered in Ref. [24].
The general production process for X has the form A → B + X, where A and B each represent one or more particles. There can be analogous production processes for D 1 D 1 φ. In the case of D 1 D 1 φ, we assume that the invariant mass M satisfies the condition in Eq. (8). This implies that the relative momenta of the mesons, which we denote generically by p , are small compared to the natural momentum scale m. We call the production process a short-distance process if all the particles in A and B have momenta in the rest frame of X or D 1 D 1 φ that are of order m or larger. The T-matrix element for a short-distance production process can be expanded in powers of the small energy differences M − M 1+2 and M 2−1 − m divided by m 2 /M 12 and higher energy scales and in powers of the small relative momenta p divided by m and higher momentum scales. The operator product expansion can be used to organize the expansions of the T-matrix elements into sums of products of short-distance coefficients and matrix elements of local operators between the vacuum state |∅ and the appropriate final state. The leading terms in the expansions are those with the lowest dimension operator, which is D † 1 D † 2 (x). The corresponding renormalized operator is is represented by a dot from which a D 1 line and a D 2 line emerge, as illustrated in Fig. 3. If we keep only the leading terms, the factorization formulas for the T-matrix elements reduce to 2 A is a function of high energy scales, such as m and M 12 , and of the total momentum of the X or D 1 D 1 φ. The only dependence on whether the final state contains X or D 1 D 1 φ is in the matrix elements of the local operator λ 0 D † 1 D † 2 (0). Shortdistance and long-distance effects are separated in Eqs. (28). The decay of X into D 1 D 1 φ can be described within the scalar meson model. The fundamental theory may also allow the decay of X into other final states C. We call such a process a short-distance decay if all the particles in C have momenta in the rest frame that are of order m or larger. The T-matrix element for such a process can be expanded in powers of the small energy differences M X − M 1+2 and M 2−1 − m divided by m 2 /M 12 and higher energy scales. The operator product expansion can be used to organize the expansion of the T-matrix element into sums of products of short-distance coefficients and matrix elements of local operators between the state |X and the vacuum state. The leading term in the expansion is the one with the lowest dimension operator, which is D 1 D 2 (x). The corresponding renormalized operator is λ 0 D 1 D 2 (x). If we keep only this term, the factorization formula for the T-matrix element reduces to 2 The Wilson coefficient C C 12 is a function of high energy scales, such as m and M 12 . Shortdistance and long-distance effects are separated in Eqs. (29).
If the fundamental theory includes processes that allow the production of X via A → B + X and the decay of X via X → C, it also allows the process A → B + C, where C represents the same particles but with a variable invariant mass M instead of M X . This process has a resonant enhancement when M is near the D 1 D 2 threshold as specified by Eq. (8). We call it a short-distance process if each of the particles in A and B has momentum large compared to m in the rest frame of C and if the relative momentum between each particle in C and each particle in A or B is large compared to m. The T-matrix element for such a process can be described within the scalar meson model by a double operator product expansion. The leading terms in this expansion are where the 4-vector is P µ = (M, 0 ). The Wilson coefficients C B,12 A and C C 12 are the same ones that appear in the operator product expansions in Eqs. (28) and (29). The first term C B,C A on the right side of Eq. (30) takes into account the direct production of C at short distances. This term can be expanded in powers of the small energy differences M − M 1+2 and M 2−1 − m divided by energy scales of order m 2 /M 12 or higher. The leading term in this expansion is simply a constant. According to Eq. (24), the Fourier transform in Eq. (30) is just the X propagator evaluated at E = M. Short-distance and long-distance effects are not yet separated in Eq. (30), because the X propagator requires an additive renormalization [24].

III. LEADING ORDER IN g
At 0 th order in g, the only interaction in the scalar meson model is the contact interaction between D 1 and D 2 with coupling constant λ 0 . The D 1 D 1 φ states remain noninteracting at this order in g. In this section, we summarize some of the results at 0 th order in g that were obtained in Ref. [24] and we generalize those results to the complex mass scheme.
At 0 th order in g, all the observables for processes near the D 1 D 2 threshold are related in a simple way to the connected Green's function for D 1 D 2 → D 1 D 2 . The basic building block for this Green's function is the amplitude iL 0 (E) for the propagation of the D 1 D 2 pair between successive contact interactions, which is given by the Feynman diagram in Fig. 4. The amplitude L 0 (E) is calculated in Section A 3 of the Appendix. The analytic expression using dimensional regularization is given in Eq. (A13). The result in a general regularization scheme is where κ(E) is the energy variable defined by The renormalization prescription for λ 0 specified by Eq. (10) is that Re A −1 0 (E) must vanish at M 1+2 − κ 2 X /(2M 12 ). This requires the inverse bare coupling constant to be 1 Using this expression for 1/λ 0 and the expression for L 0 (E) in Eq. (A13), the amplitude in Eq. (33) reduces to If we were to use the complex mass scheme, the amplitude L 0 (E) would be given by Eq. (31) with M 1+2 replaced by M 1+2 − iΓ 2 /2. If we again demand that Re A −1 0 (E) must vanish at M 1+2 − κ 2 X /(2M 12 ), the renormalized expression for the amplitude would be where γ(E) is a complex energy variable defined by and γ X is its value at E = M 1+2 − κ 2 X /(2M 12 ): The variable γ(E) vanishes at the complex threshold E = M 1+2 − iΓ/2 and it has a negative imaginary part for all real values of E.

B. Propagator for X
If the local composite operator λ 0 D 1 D 2 (x) is used as an interpolating field for X, the propagator for X is given in Eq. (24). The diagrams for the propagator of X at 0 th order in g are shown in Fig. 6. These diagrams form a geometric series whose sum in the rest frame P = 0 is Using the expression for A 0 (E) in Eq. (33), we can express the propagator for X as To obtain a renormalized propagator for X that does not depend on the ultraviolet cutoff, an additive renormalization is necessary. After adding the constant iλ 0 , we obtain the renormalized propagator −iA 0 (E). The term −iA 0 (E) in the propagator for X in Eq. (40) has a pole at the complex energy E pole given in Eq. (11). Using the renormalized expression for A 0 (E) in Eq. (35), we find that the binding energy and width of X at 0 th order in g are The wavefunction normalization factor for X is If we were to use the complex mass scheme, the amplitude A 0 (E) would be given by Eq. (36). This amplitude has a pole at a complex energy E pole that determines the binding energy and width of X. The complex parameter E pole satisfies where γ(E) and γ X are defined in Eqs. (37) and (38). Comparing with Eq. (11), we find that at LO in the complex mass scheme, the binding energy and width of X are If we treat Γ 2 as being of order g 2 , the binding energy of X in Eq. (44a) differs from κ 2 X /(2M 12 ) only at order g 4 . The wavefunction normalization constant for X at LO in the complex mass scheme is Re γ X .
C. Short-distance Production of X The operator product expansion of the T-matrix element for a short-distance production process A → B + X is given in Eq. (28a). The vacuum-to-X matrix element can be calculated via the LSZ formalism using λ 0 D 1 D 2 (x) as an interpolating field for X [24]: At 0 th order in g, the normalization constant Z X is given in Eq. (42). At LO in the complex mass scheme, Z X is given in Eq. (45). The T-matrix element in Eq. (28a) can be expressed as the product of a short-distance factor and a long-distance factor: The rate for producing X is obtained by squaring the T-matrix element in Eq. (47) and integrating over the appropriate phase space. If A consists of a single particle, its decay rate into B + X can be expressed in the factored form We have followed Ref. [24] in choosing the long-distance factor in Eq. (48) to be M 12 |Z X |, which is dimensionless. The short-distance factor Γ B A is where the 4-vector P µ satisfies P 2 = M 2 1+2 .
D. Short-distance Decay of X The operator product expansion of the T-matrix element for a short-distance decay process X → C is given in Eq. (29). The X-to-vacuum matrix element can be calculated via the LSZ formalism using λ 0 D † 1 D † 2 (x) as an interpolating field for X [24]: 3 At 0 th order in g, the normalization constant Z X is given in Eq. (42). At LO in the complex mass scheme, Z X is given by Eq. (45). The T-matrix element in Eq. (29) can be expressed as the product of a short-distance factor and a long-distance factor: The decay rate of X into the particles represented by C is obtained by squaring the T-matrix element in Eq. (51) and integrating over the phase space of those particles. It can be expressed in the factored form We have followed Ref. [24] in choosing the long-distance factor to be M 12 |Z X |, which is dimensionless. The short-distance factor Γ C is where the 4-vector P µ satisfies P 2 = M 2 1+2 .

E. Line shape of X in a Short-distance Decay Channel
If X can be produced via the short-distance process A → B + X and if it can decay via the short-distance process X → C, then there can be resonant enhancement of the process A → B + C, where C represents the same particles but with a variable invariant mass M instead of M X . The T-matrix element for this process can be expressed as the double operator product expansion in Eq. (30). The Fourier transform of the matrix element is just the X propagator evaluated at (E, P ) = (M, 0 ). The T-matrix element therefore reduces to The T-matrix element can be expressed in a form in which the short-distance effects and long-distance effects are separated: where A 0 (E) is given in Eq. (35) and c B,C A is a complex constant with dimensions of length that is completely determined by short-distance factors: (56) The nonresonant term c B,C  (52), the short-distance factors cancel. This combination of observables has only a long-distance factor given by |A 0 (M)| 2 /(2π|Z X | 2 ).
The line shape of X as a function of the invariant mass E of the particles in C is given by the factor |A 0 (E)| 2 in Eq. (57). In Fig. 7, we compare the line shapes at 0 th order in g and at LO in the complex mass scheme. At 0 th order in g, the line shape diverges at M 1+2 − κ 2 X /(2M 12 ). In the complex mass scheme, the line shape is a resonance whose maximum is near M 1+2 − κ 2 X /(2M 12 ) and whose width is determined by Γ 2 . If Γ 2 ≪ κ 2 X /(2M 12 ), the line shape for |E − M 1+2 + κ 2 X /(2M 12 )| ≪ κ 2 X /(2M 12 ) is a nonrelativistic Breit-Wigner resonance whose full width at half-maximum is Γ 2 : However, a Breit-Wigner resonance has tails that fall off like E −2 and it is integrable. In contrast, the line shape |A 0 (E)| 2 with A 0 (E) given by Eq. (36) has tails that fall off only like E −1 and it is therefore not integrable. Thus the area under the resonance is not as well-defined as for a Breit-Wigner resonance.

IV. SECOND ORDER IN g
In this section, we calculate the results in Section III to next-to-leading order in g. Since the D 2 − D 1 φ interaction allows transitions to D 1 D 1 φ states, we also calculate the rates for processes whose final state includes D 1 D 1 φ.
One of the basic building blocks for the amplitudes for processes near the D 1 D 2 threshold is the amplitude for the propagation of D 1 D 2 between contact interactions. The term in this amplitude of order g 0 , iL 0 (E), is given in Eq. (31). The term of order g 2 , iL 2 (E), is the sum of the three Feynman diagrams in Fig. 8. These diagrams are calculated in Section A 4 of the Appendix. The final expression for L 2 (E) is the sum of the amplitudes L 2a (E) and L 2b (E) given by Eqs. (A21) and (A29): where κ(E) is the energy variable defined in Eq. (32), the function F (κ) is and κ 1 is the value of κ(E) at the D 1 D 1 φ threshold 2M 1 + m: The function t(z) in the integrand in Eq. (60) is a rational function of z that increases from 0 to 1 as z increases from 0 to 1: The term L 2 (2M 1 + m) in Eq. (59) is real valued and includes logarithmic ultraviolet divergences. The function F (κ) vanishes at the D 1 D 1 φ threshold κ = κ 1 by definition. The first few terms in the Laurent expansion of F (κ) around the D 1 D 2 threshold at κ = 0 can be obtained from Eqs. (A24) and (A31). The leading term is At 0 th order in g, the ultraviolet divergence in L 0 (E) was eliminated by summing a geometric series in L 0 and then renormalizing the coupling constant λ 0 . The ultraviolet divergence in L 2 (E) can be eliminated in a similar way. The geometric series of diagrams that must be summed are those that can be obtained by replacing each one-loop subdiagram in Fig. 5 by the sum of the one-loop subdiagram and the three diagrams for L 2 (E) in Fig. 8. The sum of that geometric series gives an amplitude iA 2 (E) that can be obtained by making the substitution L 0 → L 0 + L 2 in Eq. (33): Our renormalization prescription for λ 0 in Eq. (10) requires Re A −1 2 (E) to vanish at M 1+2 − κ 2 X /(2M 12 ). Using this prescription, the amplitude in Eq. (64) can be expressed as The quantity |A 2 (E)| 2 with A 2 (E) given in Eq. (65) is the line shape of the X resonance in a short-distance decay channel to order g 2 in ordinary perturbation theory. In Fig. 9, we compare this line shape with that at LO in the complex mass scheme, which is |A 0 (E)| 2 with A 0 (E) given by Eq. (36). They are shown for two values of the coupling constant g that correspond to the width of D 2 being Γ 2 = 0.07 MeV and 0.28 MeV. The two line shapes have the same qualitative behavior below the resonance and near the peak of the resonance. They have qualitatively different behavior near the D 1 D 2 threshold E = M 1+2 . Unlike |A 0 (E)| 2 , the function |A 2 (E)| 2 vanishes at the D 1 D 2 threshold and has a sharp peak just above the threshold. The inset in Fig. 9 shows the behavior near the D 1 D 2 threshold, which is clearly unphysical.
The reason |A 2 (E)| 2 vanishes at the D 1 D 2 threshold can be seen from the Laurent expansion of F (κ) around κ = 0. The leading term, which is given in Eq. (63), diverges like 1/κ as E → M 1+2 . No matter how small the coupling constant g, the term proportional to g 2 /κ cannot be treated as a perturbation in the region near κ = 0. The resolution of the problem can be found by noting that the diverging term in Eq. (63) can be expressed as where Γ 2 is the width of the D 2 from its decay into D 1 φ, which is given in Eq. (13). Thus the divergence at κ = 0 implies that near the D 1 D 2 threshold, the imaginary part of the D 2 self-energy must be resummed to all orders. That resummation shifts the branch point in the leading order amplitude from the real threshold E = M 1+2 to the complex threshold E = M 1+2 − iΓ 2 /2. Making this change and then imposing our renormalization scheme, the amplitude in Eq. (65) becomes where γ(E) and γ X are defined in Eqs. (37) and (38). This amplitude is a smooth function of E at κ = 0. A convenient way to implement the resummation in Eq. (67) that can be extended staightforwardly to higher orders in g is to use the complex mass scheme. This automatically gives an amplitude A 2 (E) that has smooth behavior at the D 1 D 2 threshold and is correct through relative order g 2 . In the complex mass scheme, the contribution to L 2 (E) from the diagrams in Figs. 8(a) and 8(b) are obtained by replacing the energy variable κ(E) by γ(E) defined in Eq. (37) and by replacing κ 1 by γ 1 , which is the value of γ(E) at the D 1 D 1 φ threshold 2M 1 + m: In the complex mass scheme, there are additional contributions to the function L 2 (E) coming from the one-loop diagram in Fig. 8(c) with the counterterm vertex in Eq. (21b). The additional contributions are given in Eq. (A33). Adding those terms is equivalent to making the substitution Using the expressions to order g 2 for Γ 2 in Eq. (13) and Z −1 2 in Eq. (20), we can see that the additional terms in Eq. (69) cancel the κ −1 and κ terms coming from the Laurent expansions of L 2b (E) in Eq. (A31). Our renormalization prescription for λ 0 in Eq. (10) requires Re A −1 2 (E) to have a zero at M 1+2 − κ 2 X /(2M 12 ). Using this prescription, the amplitude can be expressed as where γ(E) and γ X are given by Eqs. (37) and (38) and the function F cms (γ) is The function F cms (γ) can be obtained from F (κ) in Eq. (60) by subtracting several terms in the Laurent expansion in κ and then replacing κ(E) and κ 1 by γ(E) and γ 1 . The terms that are subtracted are the imaginary term of order κ −1 , the real terms of order κ 0 , and the imaginary term of order κ in the second of the two terms on the right side of Eq. (60). The real terms of order κ 0 are subtracted because they would cancel in the denominator in Eq. (70) anyway. The leading terms in the Laurent expansion of F cms (γ) around the complex D 1 D 2 threshold at γ = 0 are For real values of E, there is an ambiguity in the function F cms (γ(E)) in Eq. (70) for E < 2M 1 + m because of the branch point at γ = γ 1 in Eq. (71). This ambiguity can be resolved by the substitution γ → γ − iǫ. In Fig. 10, we compare the line shapes of X in a short-distance decay channel at NLO in the complex mass scheme, which is |A 2 (E)| 2 with A 2 (E) given by Eq. (70), and that at LO in the complex mass scheme, which is |A 0 (E)| 2 with A 0 (E) given by Eq. (36). They are shown for two values of the coupling constant g that correspond to the width of D 2 being Γ 2 = 0.07 MeV and 0.28 MeV. The curves are normalized so their maximum values are 1. The line shapes at LO and NLO have the same qualitative behavior. As g decreases, the quantitative differences between the line shapes at LO and NLO decrease. Thus the complex mass scheme seems to provide a good solution to the infrared problem near the D 1 D 2 threshold that arises in ordinary perturbation theory. We will therefore use the complex mass scheme in most of the calculations in the remainder of this Section.

B. Propagator for X
If the local composite operator λ 0 D 1 D 2 (x) is used as the interpolating field for X, the propagator for X is defined in Eq. (24). The expression for the propagator at 0 th in g is given in Eq. (39). A renormalizable expression for the X propagator that is accurate to order g 2 can be obtained simply by making the substitution L 0 → L 0 + L 2 in Eq. (39): Comparing with the expression for the amplitude A 2 (E) in Eq. (64), we see that the propagator in Eq. (73) can be expressed as The expression for A 2 (E) at NLO in the complex mass scheme is given in Eq. (70). This amplitude has a pole at a complex energy E pole that determines the binding energy and the width of X. The energy E pole is the solution to where γ(E) and γ X are defined in Eqs. (37) and (38). The wavefunction normalization constant for X is .
If we treat Γ 2 as order g 2 , the solution to Eq. (75) for the complex variable E pole through order g 2 is Using the expression for F cms (γ) in Eq. (71), we find that the binding energy and width of the X to order g 2 are where the upper endpoint z 0 of the integral is For κ X = 27.8 MeV, the second term on the right side of Eq. (78b) is 0.011 Γ 2 . Thus the width of the constituent D 2 accounts for most of the width of X at order g 2 .
We can obtain an analytic expression for the width Γ X if κ X ≪ (m 1 /M 12 ) 1/2 κ 1 . In this limit, the upper limit z 0 on the integral approaches 1 and the width reduces to As m/M 1 varies from 0 to ∞, the factor in square brackets ranges from π/2 to 1− M 1+2 /M 2 . Note that the expression for Γ X in Eq. (80) can be smaller than Γ 2 for some values of the parameters.

C. Short-distance Production of X
The operator product expansion of the T-matrix element for the short-distance production process A → B + X is given in Eq. (28a). At NLO in the complex mass scheme, the operator matrix element is given by Eq. (46), where Z X is now the wavefunction normalization constant for X in Eq. (76). The factored expression for the T-matrix element is given by Eq. (47). If A consists of a single particle, its decay rate into B + X can be expressed in the factored form in Eq. (48), where Z X is given in Eq. (76).

D. Short-distance decay of X
The operator product expansion for the T-matrix element for the short-distance decay process X → C is given in Eq. (29). At NLO in the complex mass scheme, the operator matrix element is given in Eq. (50), where Z X is the wavefunction normalization constant for X in Eq. (76). The factored expression for the T-matrix element is given by Eq. (51). The decay rate of X into the particles represented by C can be expressed in the factored form in Eq. (52), where Z X is given in Eq. (76).

E. Line shape of X in a Short-distance Decay Channel
The T-matrix element for the short-distance process A → B + C, where C is a set of particles with invariant mass M close to the D 1 D 2 threshold, can be expressed as the double operator product expansion in Eq. (30). The T-matrix element at order g 0 is given in Eq. (54). An expression for this T-matrix element that is accurate through order g 2 can be obtained by substituting L 0 → L 0 + L 2 : .
The expression for the T-matrix element in which short-distance and long-distance factors are separated is given by Eq. (55) with A 0 (M) replaced by A 2 (M). The invariant mass distribution for the short-distance decay channel is where Γ B A and Γ C are the short-distance factors defined in Eqs. (49) and (53). The expression for A 2 (M) at NLO in the complex mass scheme is given in Eq. (70). The line shapes at LO and NLO in the complex mass scheme are compared in Fig. 10.
If the invariant mass distribution in Eq. (82) is divided by the product of the decay rates Γ[A → B + X] and Γ[X → C] in Eqs. (48) and (52), the short-distance factors cancel, leaving only the long-distance factors |A 2 (M)| 2 /(2π|Z X | 2 ), where A 2 (E) and Z X are given in Eqs. (70) and (76). The corresponding quantity at LO in the complex mass scheme is |A 0 (M)| 2 /(2π|Z X | 2 ), where A 0 (E) and Z X are given in Eqs. (36) and (45). The dependence of this long-distance factor on the invariant mass M is illustrated in Fig. 10, where the curves are all normalized to maximum value 1. The peak value of |A(M)| 2 /(2π|Z X | 2 ) is also completely determined by long distances. For κ X = 27.8 MeV and Γ 2 = 0.07 MeV, the The fundamental theory may have short-distance production processes A → B + D 1 D 1 φ, where A and B both represent one or more particles whose momenta in the D 1 D 1 φ rest frame are of order m or larger. The D 1 D 1 φ invariant mass M is assumed to be close to M 1+2 , as specified by the condition in Eq. (8). The operator product expansion of the T-matrix element for this process is given in Eq. (28b).
We first calculate the rate for A → B + D 1 D 1 φ to order g 2 using ordinary perturbation theory. At order g, the operator matrix element in Eq. (28b) is given by the series of Feynman diagrams in Fig. 11 summed over the two permutations of the external D 1 lines. The sum of the geometric series of diagrams is A 0 (M)/(−λ 0 ) multiplied by the first diagram in the series. According to the Optical Theorem, the square of the first diagram integrated over the 3-body phase space for D 1 D 1 φ is given by the sum of the 3-body cuts through the Feynman diagrams for L 2 (E) in Fig. 8. If M < M 1+2 , the sum of the 3-body cuts gives the imaginary part of L 2 (E). Thus the square of the T-matrix element in Eq. (28b) integrated over the 3-body phase space of D 1 D 1 φ is 1 2 where p 1 , p 2 , and k are 4-momenta for D 1 , D 1 , and φ, respectively, and P 2 = M 2 . If A consists of a single particle, the invariant mass distribution for A → B + D 1 D 1 φ at order g 2 can be expressed as where Γ B A is the short-distance factor defined in Eq. (49) and we have replaced a multiplicative factor of M by M 1+2 .
If there is a short-distance process A → B + D 1 D 1 φ, then the process A → B + D 1 D 2 is also allowed. The invariant mass distribution for A → B + D 1 D 2 at order g 0 is given by the expression on the right side of Eq. (84) with Im L 2 (M) replaced by Im L 0 (M). The threshold for this process is M > M 1+2 . Since D 2 ultimately decays into D 1 φ, a decay into D 1 D 2 can also be regarded as a contribution to the inclusive decay into D 1 D 1 φ.
The result in Eq. (84) is an order g 2 correction to the inclusive invariant mass distribution for D 1 D 1 φ. It is nonzero for all M > 2M 1 + m. For M > M 1+2 , there are additional corrections of order g 2 that can be obtained by replacing any of the one-loop subdiagrams in Fig. 11 by one of the two-loop diagrams for L 2 (E) in Fig. 8. Nonperturbative renormalization of the coupling constant λ 0 then requires that a geometric series of these corrections be summed to all orders. The net effect is to replace A 0 (M) in Eq. (84) by A 2 (M). Thus the complete result for the inclusive invariant mass distribution for D 1 D 1 φ through order g 2 is  85) is The inclusive invariant mass distribution in Eq. (85) when calculated using ordinary perturbation theory has unphysical behavior at the threshold M = M 1+2 . The factor |A 2 (M)| 2 has a zero at M = M 1+2 and a sharp peak just above that threshold as illustrated in Fig. 9. The term −2 ImL 2 (M) also has unphysical behavior at this threshold. The imaginary part of L 2 (E) in the region 2M 1 + m < E < M 1+2 is given by the sum of Eqs. (A22) and (A30). The resulting expression for −2 Im L 2 (E) is where the upper limit z 0 of the integral is The function in Eq. (87) diverges like 1/κ(E) as E → M 1+2 . As a consequence, near the D 1 D 2 threshold, corrections of higher order in g 2 that come from the D 2 self energy are not suppressed. This is the same problem that prompted us to introduce the complex mass scheme in Section IV A. We next consider the inclusive rate for A → B + D 1 D 1 φ in the complex mass scheme. At LO in the complex mass scheme, the invariant mass distribution for D 1 D 1 φ can be obtained from Eq. (86) by replacing A 2 (M) by A 0 (M) in Eq. (36). This can be written more explicitly as where γ(E) and γ X are defined in Eqs. (37) and (38). This expression has a nonzero imaginary part at all energies E, even below the D 1 D 1 φ theshold. For E < 2M 1 + m, the imaginary part is small and it is cancelled by higher orders in the complex mass scheme. At NLO in the complex mass scheme, the invariant mass distribution for D 1 D 1 φ is obtained from Eq. (86) by using the expression for A 2 (M) in Eq. (70).
In Fig. 12, we compare the line shape of X in the D 1 D 1 φ channel, which is proportional to Im A 2 (E), with the line shape of X in a short-distance decay channel, which is proportional to |A 2 (E)| 2 . They are shown for two values of the coupling constant g that correspond to the width of D 2 being Γ 2 = 0.07 MeV and 0.28 MeV. The curves are normalized so their maximum values are 1. The line shapes have the same qualitative behavior. The line shape for D 1 D 1 φ is shifted towards higher energy by an amount that decreases as g decreases.
The Belle Collaboration has observed an enhancement in the production of D 0D0 π 0 near threshold in the decay B → K + D 0D0 π 0 [32]. The enhancement peaks at the invariant mass 3875.2 +1.1 −1.9 MeV, where we have combined the errors in quadrature. If we use the recent precision measurement of the D 0 mass by the CLEO Collaboration [26], the peak is above the D * 0D0 threshold by 3.4 +1.2 −1.9 MeV. In contrast, the peak observed in the short-distance decay channel J/ψ π + π − is below the D * 0D0 threshold by 0.6 ± 0.6 MeV. The difference is 4.0 +1.3 −2.0 MeV, which differs from 0 by two standard deviations. Part of the discrepancy may be due to a shift in the invariant mass distribution analogous to the one illustrated in Fig. 12.
The width Γ X of the X resonance was defined in Eq. (11) in terms of the energy E pole of the pole in the X propagator: Γ X = −2 Im E pole . At LO in the complex mass scheme, Γ X is simply Γ 2 . At NLO in the complex mass scheme, Γ X is obtained by solving Eq. (75) for E pole . If the width of the X resonance is sufficiently small, its width can also be interpreted as its decay rate. A direct calculation of that decay rate will give a result that is closely related to but not identical to Γ X . We will denote the result of the direct calculation of the inclusive decay rate by Γ[X → D 1 D 1 φ].
We can use the Optical Theorem to calculate the inclusive decay rate for X → D 1 D 1 φ: The forward T-matrix element for X → X can be defined using the LSZ formalism and it is given in Eq. (27). At NLO in the complex mass scheme, the normalization constant Z X and the propagator ∆ X (E, 0) are given in Eqs. (76) and (73). The resulting expression for the decay rate of X in Eq. (90) is We have used Eq. (64) to express this in a form with a factor of A −1 2 (M X ), because A 2 (E) does not depend on the ultraviolet cutoff Λ. The decay rate in Eq. (91) depends on Λ through the bare coupling constant λ 0 and through the ultraviolet-divergent additive constants in L 0 (M X ) and L 2 (M X ). From the renormalized expression for A 2 (E) in Eq. (64), we can see that the divergent terms in L 0 (E) + L 2 (E) must be cancelled by 1/λ 0 as Λ → ∞. The leading divergence is a linear divergence in L 0 (M 1+2 ) that is linear in Λ. Thus 1/λ 0 must include a cancelling linear divergence. This implies that λ 0 [L 0 (E) + L 2 (E)] must approach 1 as Λ → ∞. Upon taking the limit Λ → ∞, Eq. (91) reduces to An explicit expression for the inclusive decay rate of X can be obtained from Eq. (92) by inserting the renormalized expression for A 2 (E) in Eq. (70): where Z X , γ(E), γ X , and F cms (γ) are given in Eqs. (76), (37), (38), and (71). If we treat Γ 2 as order g 2 , the decay rate in Eq. (93) reduces at order g 2 to the corresponding result for the width Γ X in Eq. (78b). In Fig. 13, we compare the width Γ X obtained by solving Eq. (75) and the decay rate Γ[X → D 1 D 1 φ] in Eq. (93) as functions of Γ 2 . They agree for small Γ 2 , because the difference between Γ X and Γ[X → D 1 D 1 φ] is of order g 4 . They differ by less than 1% if Γ 2 < 0. 16 MeV.

V. SUMMARY
We have examined the effects of a weakly-bound hadronic molecule on a nearby 3-body threshold, where the three bodies consist of a primary constituent of the molecule and the decay products of the other constituent. We studied these effects in the simplest possible model: the scalar meson model, with spin-0 constituents D 1 , D 2 , and φ and momentumindependent interactions. The primary constitutents of the molecule X are D 1 D 2 and the 3-body threshold is that for D 1 D 1 φ. The D 1 D 2 contact interaction with coupling constant λ 0 was treated nonperturbatively. The D 2 − D 1 φ interaction with coupling constant g was treated perturbatively. Several observables were calculated to all orders in λ 0 and to nextto-leading order in g. Both ultraviolet problems and infrared problems were encountered in these calculations.
We found that at next-to-leading order in g, the ultraviolet divergences can be removed by a perturbative renormalization of the D 2 mass, a nonperturbative renormalization of λ 0 , and a nonperturbative renormalization of short-distance coefficients in the operator product expansion. The nonperturbative renormalizations required the summing of a geometric series of order-g 2 subdiagrams to all orders. The subdiagrams represent order-g 2 contributions to the amplitude for the propagation of D 1 D 2 between contact interactions. Thus the renormalized amplitudes at next-to-leading order include all terms of order g 2 together with subsets of terms from all higher orders in g.
We found that the next-to-leading order amplitudes suffer from an infrared problem at the D 1 D 2 threshold. The problem is illustrated dramatically in Fig. 9, which shows the line shape |A(E)| 2 in a short-distance decay channel as a function of the energy E. The pathological behavior at E = M 1+2 arises because the D 2 − D 1 φ interaction shifts the D 1 D 2 threshold from the real value M 1+2 to the complex value M 1+2 − iΓ 2 /2. Since Γ 2 is of order g 2 , a strict expansion in powers of g leads to singularities at E = M 1+2 . One possible solution to this problem is to sum the geometric series of self-energy corrections to the D 2 propagator to all orders, but this leads to very complicated loop integrals. A simpler solution is to use the complex mass scheme, in which the width Γ 2 is included in the Feynman rule for the D 2 propagator and its effects are systematically compensated at higher orders in g through counterterms. This partial resummation of the D 2 propagator corrections leads to smooth behavior of the line shape |A(E)| 2 at the D 1 D 2 threshold as illustrated in Fig. 10. Having solved the infrared problems, we calculated several observables to NLO in the complex mass scheme. They include the line shapes of the X resonance in both a short-distance decay channel and in the D 1 D 1 φ channel.
The X(3872) is a weakly-bound hadronic molecule whose primary constitutents are a superposition of charm mesons: D * 0D0 + D 0D * 0 . Its mass is near the 3-body thresholds for D 0D0 π 0 , D 0 D − π + , and D +D0 π − . The methods we applied to the scalar meson model can be extended straightforwardly to this multichannel problem. A conventional perturbation expansion in the D * − Dπ coupling constant will have infrared problems at the D * 0D0 threshold which are associated with the decay D * 0 → D 0 π 0 . These infrared problems can be avoided by using the complex mass scheme. The ultraviolet problems can be expected to be more severe in this system, because the D * − Dπ interaction is proportional to the 3-momentum of the pion. Thus renormalization may not be as simple as in the scalar meson model. If the renormalization problem can be solved, it will be straightforward to calculate observables to NLO in the complex mass scheme. One of the most interesting applications will be to predict the difference between the line shapes of the X(3872) resonance in a short-distance decay channel, such as J/ψ π + π − , and in the D 0D0 π 0 channel.
where µ is the renormalization scale and (z) a is the Pochhammer symbol: (z) a = Γ(z + a) Γ(z) . (A2) The function of ǫ in the prefactor is designed to simplify the analytic expressions for dimensionally regularized integrals by cancelling the effects of the analytic continuation of angular integrals. Thus the integral of a scalar function of k is Loop integrals are evaluated by first inserting the appropriate nonrelativistic propagators: where i = 1, 2, φ and M φ = m. Integrals over the loop energy are then evaluated using the residue theorem. After using the Feynman parameter trick to combine denominators, the dimensionally regularized integrals over the 3-momenta can be evaluated. The final steps are the evaluation of the Feynman parameter integrals and the analytic continuation to ǫ = 0. It will be convenient to express the results in terms of the variables 2. Self-energy of D 2 The one-loop diagram for the self-energy Σ(p 0 , p) of D 2 is shown in Fig. 2(a). The expression for the diagram is Σ 2a (p 0 , p) = −ig 2 d 4 q (2π) 4 ∆ 1 (q 0 , q) ∆ φ (p 0 − q 0 , | p − q |).
The residue theorem can be used to evaluate the energy integral: If we implement dimensional regularization using the prescription in Eq. (A1), the analytic result is The self-energy in dimensional regularization is given by the analytic continuation of Eq. (A8) to ǫ = 0: The momentum integral in Eq. (A7) has a linear ultraviolet divergence that is set to zero by dimensional regularization. The self-energy in a general regularization scheme is given in Eq. (16). The linear ultraviolet divergence is contained in the extra term Σ 2a (M 1 + m, 0), which is real valued.

Propagation between contact interactions at leading order
The amplitude for the propagation of D 1 D 2 between contact interactions is given at order g 0 by the one-loop diagram in Fig. 4. In the center-of-momentum frame where the 4-momentum is (E, 0), this amplitude is The residue theorem can be used to evaluate the energy integral: If we implement dimensional regularization using the prescription in Eq. (A1), the analytic result is where κ is defined in Eq. (A5a). The amplitude in dimensional regularization is given by analytic continuation to ǫ = 0: The momentum integral in Eq. (A11) has a linear ultraviolet divergence that is set to zero by dimensional regularization. The amplitude in a general regularization scheme is given in Eq. (31). The linear ultraviolet divergence is contained in the extra term L 0 (M 1+2 ), which is real valued.

Propagation between contact interactions at next-to-leading order
The amplitude for the propagation of D 1 D 2 between contact interactions has contribution of order g 2 from the two-loop diagrams in Fig. 8. In the center-of-momentum frame where the 4-momentum is (E, 0), the order-g 2 amplitude iL 2 (E) is a function of E only. We proceed to calculate the contributions from each of the three diagrams.

a. Correction from φ exchange
The contribution to the amplitude iL 2 (E) from the two-loop diagram in Fig. 8(a) in which φ is exchanged between the D 1 and D 2 is The integrals over the loop energies k 0 and q 0 can be evaluated by closing the contours in the upper half-planes. After introducing Feynman parameters and evaluating the integrals over the loop momenta, the amplitude reduces to where κ and κ 1 are defined in Eqs. (A5). One of the Feynman parameter integrals can be evaluated analytically by changing variables to z = x + y and w = (x − y)/(x + y). The integral over w gives a hypergeometric function: where t(z) is the rational function of z given in Eq. (62). The pole at ǫ = 0 in Eq. (A16) is a logarithmic ultraviolet divergence. The amplitude L 2a (E) in dimensional regularization is defined by the Laurent expansion of the expression in Eq. (A16) to order ǫ 0 . The expansion in powers of ǫ is facilitated by using a transformation formula to replace the hypergeometric function in Eq. (A16) by one of the form 2 F 1 (1, ǫ, 3 2 ; t): The expansion of 2 F 1 (1, ǫ, 3 2 ; t) to order ǫ is simple: where κ and κ 1 are defined in Eqs. (A5) and M 11 is given by the one-loop diagram in Fig. 8(c) is iL 2c (E) = −i d 4 k (2π) 4 ∆ 1 (E − k 0 , k)∆ 2 (k 0 , k) 2 (δM + iΓ 2 /2) In the term with the factor Z −1 2 − 1, the integral is proportional to L 0 (E). This amplitude in a general regularization scheme is given in Eq. (31). In the term with the factor δM + iΓ 2 /2, the integral is convergent. The complete result for the diagram is