Spectroscopy of Masses and Couplings during Inflation

In this work, we extend the idea of Quasi Single Field inflation to the case of multiple isocurvaton fields with masses of order of Hubble, which are coupled kinetically to the inflaton field and have some interactions among themselves. We consider the effects of these massive modes in both the size and the shape of the bispectrum. We show that the shape of the bispectrum in the squeezed limit is dominated by the lightest field and is the same as in Quasi Single Field inflation. This is a generic feature of multiple isocurvaton fields and is independent of the details of the interactions among the massive fields. When the isocurvaton fields have similar masses, we can potentially distinguish two different shapes in the squeezed limit so that the shape of the bispectrum can act as a particle detector. However, in the presence of hierarchy among the massive fields, the dominant effect is due to the lightest field.


I. INTRODUCTION
Inflation is the leading paradigm for early universe and structure formation which is well consistent with recent observations from Planck satellite [2,3]. In simple models of inflation a scalar field rolls down towards the minimum of its potential supporting a long enough period of inflation. In order to solve the flatness and the horizon problem one requires 60 or so number of e-foldings. The near dS background of the inflationary mechanism usually guarantees that the primordial anisotropies imprinted on cosmic microwave background (CMB) map are nearly scale-invariant, nearly adiabatic and nearly Gaussian. These basics predictions of inflation are well consistent with the cosmological observations.
There have been much interests on non-Gaussianities during past few years, for a review see [4,5]. There was no detection of non-Gaussianity by Planck [6] and there are strong upper bound on the amplitude of local-type non-Gaussianities, f loc N L = 2.7 ± 5.8 (68 % CL) [6]. Having this said, the constraints on other shapes of non-Gaussianities, such as equilateral-type shape, are not as tight and there is still room for large observable non-Gaussianities in these shapes. On the other hand, there are also some constraints coming from the Large Scale Structure (LSS) analysis such as galaxy power spectrum as well as the scale dependence corrections to galaxy bias. As it is shown in [7] the galaxy halo bias and CMB analysis are complementray to each other. In the sense that, the halo bias is relatively more sensitive to the isocurvaton mass, since the scale dependent corrections depend on the behavior of curvature bispectrum in the squeezed limit. The sensitivity of CMB on the other hand, is more on the values of f N L . So combining both, depending on the amplitude of f N L , in principle we can distinguish the Quasi Single Field inflation (QSF) model, [1], as well is its extended version, in a significant fraction of mass scales.
While the simplest inflationary models are based on single field, one can easily conceive situations in which more then one fields are responsible for inflation dynamics and generating curvature perturbations. There have been much works on models of multiple field inflation, for a review see [11,12]. In multiple fields scenarios, there is one adiabatic direction responsible for curvature perturbations and several isocurvature directions which can generate entropy perturbations [13]. In these models, it is usually assumed that all fields are light, i.e. lighter than H, the Hubble expansion rate during inflation. This is because only the light fields with mass H are expected to play important roles. Very heavy fields with mass m H are expected to rapidly rolls down towards their minimum and play no important roles during inflation. Having this said, it is possible that one encounters situations in which two or more light or semi-heavy fields with some hierarchies of masses are present during inflation. Therefore, it is an important question as how the hierarchy of masses and couplings show themselves in inflationary predictions such as the power spectrum and the bispectrum.
The effects of an isocurvaton field, a semi-heavy field with the mass m 3H/2, on inflation bispectrum and trispectrum were studied by Chen and Wang in [1,14]. Phenomenologically, this model is very interesting. Because by adjusting the magnitude of the ratio m/H, the shape of the bispectrum in this model effectively interpolates between the shape of the bispectrum in single field and the shape of bispectrum in multi field inflationary models. On the other hand, this model is well-motivated from the theoretical point of view too. Indeed, many models of high energy physics such as string theory and supersymmetry contain many scalar degrees of freedom, for a review see [15][16][17][18][19][20][21][22]. In embedding inflation in these models of high energy physics one can naturally encounter several fields with different masses. Also in this setup it is quite natural to have massive fields, with masses at the order of Hubble parameter H or so. The authors in [1] considered the situation in which the isocurvaton field does not affect the power spectrum but can have significant contributions on bispectrum via turning trajectory. To see how this can happen suppose we denote the inflaton field by θ which is an angle in a two-dimensional field space trajectory while σ denotes the additional semi-heavy field with the self-coupling potential V (σ). At the background level, we can stabilize the semi-heavy field on its own minimum without any effect on the background trajectory. In [1] it is shown that the amplitude of f N L scales like θ H 3 V,σσσ H P ζ − 1 2 in which P ζ ∼ 10 −9 is the observed curvature perturbation power spectrum. With large enough self-coupling V ,σσσ one can easily saturate the current observational bound on the amplitude of non-Gaussianity. In addition, it is shown in [1] that depending on the σ field's mass, a new shape of non-Gaussianity is generated which continuously interpolate between the local shape and the equilateral shape. For other works with similar ideas see [23][24][25][26][27][28][29][30][31][32][33].
In this work we extend the analysis of [1] to the case in which two semi-heavy fields are present in the model, denoted by σ 1 and σ 2 . We show that while the correction to the power spectrum is quite negligible, these semi-heavy fields can have significant contributions in the bispectrum. It is interesting to study the effects of these massive fields in both the size and the shape of the bispectrum in the squeezed limit. It turns out that, while these fields have non trivial effects in the size of the bispectrum, they will not produce any new shapes in the squeezed limit. This means that in the squeezed limit, the leading shape is the same as in the quasi single field inflation, which scales as k −( 3 2 +νi) 3 with (k 3 k 1 k 2 ) and where we have defined ν i as, ν i ≡ 9 4 − m 2 i H 2 , (i = 1, 2). This is a generic feature of multiple isocurvaton fields and is independent of the details of the interactions among the massive modes.
The rest of this paper is organized as follows. In Section II we present our model containing the motivation as well as the phenomenological Lagrangian. We study the perturbation at the quadratic level, considering the free field action as well as the exchange vertex interactions which are necessary in order to convert the contribution from the isocurvaton fields into the curvature perturbation. In addition we investigate the interaction among the isocurvaton fields at the cubic order in perturbations. In Section III we present our power spectrum analysis. We show that the correction into the power-spectrum, which is controlled by the ratioθ 0 H , is small and can be neglected. In Section IV we present the bispectrum analysis employing the in-in formalism. As in [1] we use two different methods in order to calculate the in-in integrals. In Section V we consider the squeezed limit of the bispectrum followed by our conclusions and discussions in Section VI. We relegate many technical analysis of the in-in integrals into Appendices.

II. THE SETUP
Here we present our setup. This is an extension of the work by Chen and Wang [1] known as Quasi Single Field Inflation (QSF). Our model contains two iso-curvaton fields (semi-heavy fields) σ 1 and σ 2 . Following the logic in [1] we want to study the effects of these iso-curvaton fields on power spectrum and bispectrum. First we start with the motivation.
In original QSF idea, it is imagined that the light field θ is moving along a constant radius circular trajectory with the fixed radius R determined by the vacuum expectation value of the heavy field σ, R = σ 0 . Here we extend this view to higher dimensional field space geometry. Consider the model containing a light field moving on the surface of a two-dimensional sphere. To be specific, let us denote the light field by the azimuthal direction φ, while the heavy fields are denoted by the radius r and the other angular direction θ. In this picture the heavy fields (r, θ) are stabilized around the background value (r 0 , θ 0 ), while the light field φ moves nearly freely on the surface of the sphere. The Lagrangian for the kinetic energy is Now consider the fluctuations r(t) = r 0 + δr(t), θ(t) = θ 0 + δθ(t) and φ(t) = φ 0 (t) + δφ(t). Plugging these back into the kinetic Lagrangian yields the quadratic Lagrangian plus the exchange vertex interactions ∆L 3 =φ 0 r 0 sin 2 θ 0 δr + r 2 0 sin 2θ 0 δθ δφ .
This is a multiple field extension of [1] with the exchange vertex interactions δrδφ and δθδφ. Note that from Eq. (1) we also obtain other sub-leading interactions, e.g. the exchange vertex between r and θ, which we discard. Having presented our motivation for the multiple massive fields kinetically coupled to the inflaton field, we proceed with our phenomenological model now. From now on we follow the notation of [1] with θ being the inflaton field while the heavy fields are denoted by σ 1 and σ 2 . The action is in which V sr (θ) is the slow-roll potential and V (σ 1 , σ 2 ) represents the potential between the fields σ 1 and σ 2 . Also c i , i = 1, 2 are both constants of the order R. At the background level, the Hubble equation and the continuity equation are where we have absorbed the background fields' values σ 10 and σ 20 into the net turning radius by defining Also the equations of motion for σ i0 (t), (i = 0, 1) and θ 0 (t) are In this picture σ i 's have been stabilized around their background values while the inflaton field θ is slowly rolling along the potential V sr . Now we can specify the form of the V (σ 1 , σ 1 ) up to the third orders in fields perturbations which will be used to calculate the bispectrum Furthermore, we chose the mass basis such that the mass matrix is diagonal so we can neglect the term m 12 . In what follows we assume 0 ≤ m i /H ≤ 3/2 so the iso-curvaton fields can be light or semi-heavy. This may have a natural interpretation for the super-symmetric completion of the theory. We choose the spatially flat gauge in which the scale factor of the metric is homogeneous, We expand only the matter part of the lagrangian, given by Eq. (4), while ignoring the perturbations in the gravitational sector. This is justified, since from Maldacena's analysis [34] it is expected that the gravitational back-reactions in bispectrum are slow-roll suppressed. This was specifically demonstrated in [1] and this is expected to be the case in our setup. For the future references, the slow-roll parameters are where We quantize the Fourier components of the free fields in the interaction picture denoted by δθ I k , δσ I 1k and δσ I 2k such that Here a k (a † k ) , b k (b † k ) and c k (c † k ) are the usual annihilation (creation) operators defined for δθ I k , δσ I 1k and δσ I 2k respectively. Since these fields are treated as independent random fields we assume that the annihilation and creation operators for different fields are independent of each other and they commute with each other. Furthermore, they satisfy the usual commutation relations The mode functions u k , v k and w k satisfy the following linear equations of motion obtained from the free Hamiltonian H 0 given in Eq. (15) The cubic Hamiltonian density in Fourier space is As in the power spectrum case, the transfer vertexes, the first two terms in (16), convert the above interactions to curvature perturbations. As it is mentioned in [1] the three point function for δθ is given by in which the symbols T andT represent the time ordering and anti-time ordering respectively. We can expand the above equation in two equivalent forms, which are called in [1] as the "factorized" and "commutator" forms. Each of these forms has its own advantage and disadvantages. Specially, as it is mentioned in [1], in the calculation of the integrals there are some unphysical divergences in both IR limit, τ → 0, and UV limit, τ → −∞ which can be neglected in commutator and factorized forms, respectively. In the following we briefly mention each of them and leave the details to Appendix B.
Let us start with the factorized form. The expansion of Eq. (37) is In each term, the self interactions between the δσ i 's can appear in one of four H I . For example, for H 3I 3 interaction we need one transfer vertex of the form H 2 1 ≡ −c 1 a(t) 3θ 0 δθδσ 1 and two of the form H 2 2 ≡ −c 2 a(t) 3θ 0 δθδσ 2 . After contractions we get, There are altogether 80 different terms in this format, which are given in details in Appendix B 1. Now let us look into the commutator form. The expansion of Eq. (37) in terms of the nested commutators is given by [35] Again we should consider each of H 3 i separately. There are 24 terms. The details are given in Appendix B 2. As mentioned in [1], when we calculate the integrals we enter unphysical singularities. To eliminate this spurious singularities we have to perform the integral in the so-called "mixed form" which is a combination of the factorized form and the commutator form. In principle one can find the full shape of the bispectrum in the mixed form in terms of different couplings as well as different masses of σ i fields. However, this involves very complicated analysis which is beyond the scope of this work. In the following we only consider the squeezed limit of the bispectrum which contains interesting physical information.

V. SQUEEZED LIMIT OF BISPECTRUM
The goal of this section is the calculation of the bispectrum in the squeezed limit, k 3 k 1 k 2 . The importance of doing analysis in this limit is that, we can analytically estimate the effect of the massive fields in both the size and the shape of the bispectrum. As we will show explicitly, while these massive fields have non-trivial effect in the amplitude of the bispectrum, they will not generate any new shapes as compared with the quasi single field inflation, [1]. Before going through the analysis, let us review the case of quasi single field inflation. In quasi single field inflation, we can understand the squeezed limit behavior in the following way, [19]. First, remember that the squeezed limit corresponds to the correlation between a long mode and two short modes. Let us assume that this long mode crosses the horizon at τ 1 while short modes cross the horizon some time after, say τ 2 . Since the massive modes decay after exiting the horizon, at τ 2 amplitude has been reduced as, Which as we see is the leading shape of the bispectrum in the squeezed limit. We argue that this result is generic and is due to a couple of points. First, although there is a mixing between the massive modes with inflaton field, the approximate scale invariance of the model is still alive, the Goldstone boson will not be massive. Second, there is not any mixing between massive modes in the free field action, which means that at the zeroth order in the quadratic action, these massive modes do not communicate with each other and freely decay after the horizon crossing.
In the following we are going to show the above results in more details.
In this limit we can use the small argument of the Hankel functions which is given in the following, In what follows, we only consider the commutator forms, since the results of the factorized form is exactly the same as this form, [1]. we separate the contributions from three types of terms denoted by A, B and C. First we look at the A term, (B43). As we see from this equation, there are four different contributions originating from H 3 i , i = 1..4. But as we mentioned in the Appendix B 2 only two of them are independent, namely H 3 1 and H 3 3 while the contribution from the other terms can be easily calculated. So we can skip them and only consider the independent terms. We just mention one term here and leave the details in the appendix, B 3.
This part is very similar to [1], so we summarize the calculations and express the final result. With the definition Now as in [1], the term H x 2 ) in the 3rd line can be approximated in the small −k 3 /k 1 x 2 limit. However, the term H in the 3rd line can not be approximated. So by redefining y 4 ≡ k 3 /k 1 x 4 , we get The scaling behavior of the above term is We next look at the term with the permutation k 1 ↔ k 3 . With the same redefinition of x i (i = 1, 2, 3, 4), we get Again as in [1], in the 3rd line the first three functions can be approximated in the small −(k 3 /k 1 )x i , (i = 1, 2) limit, and in order to get a non-zero result, one of the Hankel functions should be expanded to O(x ν1 i ) . By this approximation, the scaling behavior of the above term is obtained to be Then we look at the term with the permutation k 2 ↔ k 3 Again in the 4th line, we can approximate the three functions in the small −x i k 3 /k 1 , (i = 2, 3) limit. For ν 1 > 1 2 , we use the leading term for the two Hankel functions and sub-leading term for the exponential function. Then the result is The scaling behavior of this term is P 1 ≡ k 5−2ν1 it is possible to choose some parts of the parameter space such that P 1 becomes larger than N 2 ≡ k This happens when in the ratio P1 N2 = k3 k1 3/2+ν2−2ν1 the expression 3/2 + ν 2 − 2ν 1 becomes negative. From Fig.4 we see that this actually happens for some values of ν i in the parameter space. So we should consider this term as well as N 1 and N 2 . On the other hand, for ν 1 < 1 2 one of the Hankel functions should be expanded to O(x ν1 i ) and by using the leading term for the exponential function, the scaling behavior is Finally, the other permutation k 1 ↔ k 2 gives each term a factor of 2. The final result can be obtained by summing up all of the above contributions from A, B and C terms.
This is the main result of this Section. Although the full analysis has been given in Appendix (B 3), it is worth to elaborate on how different terms are obtained, especially it is important to connect them to the Feynman diagrams in Fig.3. As it can be seen from Appendix (B 3), there are different contributions from the H 3I i s in the above shapes. Since S 1 (ν 1 , ν 2 ) is related to the leading shape, we express its form here, while leave the details of S 3 (ν 1 , ν 2 ) in the appendix, B 4.
Eq. (54) is one of our main results in this paper. There are different shapes in the bispectrum, which depending on the value of ν 1 and ν 2 i.e. different values of the isocurvatons mass and the coupling constants for the interaction among σ i fields, can be dominant or subdominant. It is an interesting question to study the full dependence of the bispectrum shape on the model parameters. This requires extensive analysis which is beyond the scope of this work.
Having this said, we still can obtain important information from our current results as follows.
First of all, as it was mentioned in [1], the origins of non-Gaussianity in this model are the self-interactions between δσ i 's which scales like k −2νi 3 , [1]. However, it has to convert to the curvature perturbation via the exchange vertex. While converting to the curvature bispectrum, for a typical value of ν i , the leading shape will be changed from its original contribution,k −2νi , which is more like the local shape, k −3 3 . However, as sub-leading term compared to this leading contribution, there is still k −2νi 3 shape which means that in the conversion procedure there is not any change in the shape of the bispectrum. Therefore, we conclude that the projection effect from δσ i bispectrum to that of curvature bispectrum can change both of the amplitude and the shape of bispectrum. In every case the amplitude will be changed. However, while for the leading contribution we have a change in the shape of three point function, there is still another sub-leading term which does not contain any change in its shape, as compared to the original contribution of σ i bispectrum.
The second point is that, as it is mentioned in [1], these shapes are between the equilateral-type shape, k −1 3 , and the local-type shape, k −3 3 . For small values of ν i they are more like equilateral shape while for the larger values of ν i they go to the local shape. It is physically understandable, since for the small values of ν i i.e. for higher values of the isocurvatons' masses, the isocurvaton wave function decays rapidly after horizon crossing which means that effectively we are in the single field regime. On the other hand, for the higher values of the ν i , or in another word for small values for the isocurvatons' masses, the wave functions decay slowly after horizon crossing and we are in multiple field regime with local shape.
Third, for the spectroscopy of different masses, it is worth to consider two different cases in the following, A) Suppose that we are in the regime in which both isocurvaton fields have masses of the same order. In this case, the first two terms in Eq. 54 will dominate and effectively we will have, B) In this case, we assume that there is a hierarchy between different mass scales, e.g. m σ2 m σ1 which means that ν 1 ν 2 . So as it is shown in Fig. 4, the leading scaling dependence of ν 2 becomes smaller than the subleading shape coming from σ 1 . As a result the shape is dominated with the lightest field and we can neglect the contribution from the heavy field, in which we have kept the contribution from the subleading term to specify that potentially one can detect this shape in the squeezed limit, although it is subleading.
Fourth, let us look at the order of magnitude of S i (ν 1 , ν 2 ), i = 1, .., 4. In our case, since S 3 (ν 1 , ν 2 ) and S 4 (ν 1 , ν 2 ) contain four layer integrals the detail calculation of these terms are tremendous and are left to the future work. However, compared to [1] we can estimate their order of magnitudes. We have checked that as we vary ν i from zero to 3 2 , the values of S 3 (ν 1 , ν 2 ) and S 4 (ν 1 , ν 2 ) change from zero to around 30, in the limit in which all of the coupling constants are of the same order.
Fifth, although these massive modes have non-trivial effects in the amplitude of the bispectrum, they will not produce any new shapes in this limit. In the sense that they are completely separated at the level of the shape and there will not be any mixing between their masses at this level. As we showed, it is a generic feature of this kind of models and is due to the scale invariance of this model and the fact that at the zeroth order in show roll these massive fields do not communicate with each other and decay after horizon crossing.
The last point is the comparison of this model with the PLANCK data. Because there was no detection of non-Gaussianity by PLANCK, there is not any preference for the values of ν i , see Fig. 5 but more strong constraints on the amplitude of f N L . However, it seems that there are still some rooms for the equilateral shapes to work after PLANCK. Naively speaking this means that our model has a better fit to the data for the smaller values of ν i i.e. for larger values of the isocurvatons' masses. In addition, there are also some constraints from the LSS analysis [7][8][9][10]. Interestingly the halo bias is more sensitive to the scale dependence of the bispectrum in the squeezed limit while the CMB analysis are more based on the Bispectrum amplitude. So these different measurements can act as a complementary to each other in order to distinguish QSF models among other models.

VI. SUMMARY AND DISCUSSIONS
In this work we have extended the analysis of [1] to the model containing two semi-heavy isocurvaton fields. These kind of models are well motivated in the UV completion of the string theory and supersymmetry, when naturally we can have heavy fields with masses of the order of the Hubble constant. In this model, due to the turning trajectory, the perturbations from the isocurvaton fields are converted to the curvature perturbation which can have non-trivial effects in both the power spectrum as well as the bispectrum. The correction in the power spectrum is negligible. However, because of their strong self-interactions, their contributions in the amplitude and the shape of the bispectrum of the curvature perturbation are significant.
The origin of the non-Gaussianity comes from isocurvatons' bispectrum. However they have to convert to the curvature bispectrum via the turning trajectory. While converting to the curvature bispectrum, the leading shape also goes more to local-type shape. The amplitude of f N L is controlled with the strength of these interactions as well as the velocity of the turning trajectory. On the other hand, adding new massive field does not produce any new shape in the squeezed limit which is due to the scale invariance of the model and the fact that free propagating isocurvaton modes do not communicate with each other. However, squeezed limit can act as a particle detector to recognize the particles with masses of the same order.
It is worth to mention that, while the squeezed limit of the bispectrum is dominated by the lightest isocurvaton, it does not mean that the entire bispectrum signal should be dominated by the lightest isocurvaton field. It would be interesting to consider the analysis of the full shape of bispectrum. This is however beyond the scope of this paper and is left for the future work.
From the observational point of view, there are strong upper bound on the amplitude of the local non-Gaussianity from PLANCK. This means that the strength of the self-interactions between isocurvaton fields as well as the velocity of the turning trajectory can not be large. However, still there is room for the equilateral shape to be detected in the future observations such as in large scale structure surveys which are more sensitive tot he scale dependence of the Bispectrum curvature perturbation.
this work under the "Sandwich Training Educational Program" (STEP) fellowship.

Appendix A: The full Lagrangian up to third order
In this Appendix we derive the full third order action in the spatially flat gauge. We also show that (14) is the leading order interaction. We set M P = 1 in this Appendix.
The full action is where is the gravitational action and represents the matter part os the action which is given by Eq. (4). Using the ADM formalism the action becomes Here the index of N i is lowered by the three-dimensional metric h ij and R (3) is the three-dimensional Rici scalar constructed from h ij . The definitions of the extrinsic curvature E ij and its trace E are We choose the following spatially flat gauge: For the constant turn case σ i0 (t) = constant. The constraint equations for the Lagrangian multipliers N and N i are Now we expand N i and N up to first order Plugging them into (A8) and (A9), the solutions with proper boundary conditions are Now we plug these solutions into the action and expand up to third order in perturbations. The first order terms give the equation of motion for θ 0 (t) and σ 0 (t). The quadratic action is The third order action is As we see the terms coming from the gravity back-reaction are smaller than the matter effects.

Appendix B: Details of integrals in various forms
In this Appendix we give the details of the in-in integrals in the factorized and commutator forms.

All terms in the factorized form
Here we calculate the independent contributions to the bispectrum in the factorized form. Since H 3 2 can be obtained by changing (c 1 , λ 1 , v ki ) to (c 2 , λ 2 , w ki ) we skip it. Also the terms H 3 4 is not independent of H 3 3 so we skip it too. Different terms related to Eq.(38) and H 3 1 are: Different terms related to Eq.(38) and H 3 3 are: Different terms related to Eq.(39) and H 3 1 are: Different terms related to Eq.(39) and H 3 3 are: Different terms related to Eq.(40) and H 3 1 are: Different terms related to Eq.(40) and H 3 3 are: The three point function is the sum of the above 40 terms as,

All of terms in the commutator form
Here we present the form of the integrals in commutator form. Again we only calculate the independent parts of bispectrum in this form.
Different terms related to Eq. (42) are: in which we have defined As in the factorized form, we only calculate the independent parts of bispectra in the commutator method. Therefore, in the following, we specify the contribution of different terms of H 3 i , i = 1, 3 in A ,B and C: As in the factorized case, we have three terms for A H 3 3 , related to the different possible positions for H 2 Also for B we have: Finally; for C we have: As we mentioned above the summation of these A, B and C terms give us the final result for the bispestrum in the commutator form.

Different terms in the squeezed limit
In the following, we are going to present the whole independent contributions in the squeezed limit in the commutator form.
(2 A )A H 3 2 : Since this term can be easily obtained by replacing (ν 1 , c 1 , λ 1 ) with (ν 2 , c 2 , λ 2 ), we skip it now and insert it in the final result.
This term is related to our new interaction and is divided into three different parts as follows.
Again by approximating the small limit of − k3 k1 x 2 in the third line we have, We next look at the term with the permutation k 1 ↔ k 3 . Again, the scaling behavior of this term is Which, as we argued before, is negligible. Also for the term with the permutation k 2 ↔ k 3 and ν 2 > 1 2 we have, (3 A2 )A H 3 32 : For this term we can use the results of the above case, A H 3 31 : As we can see the scaling behavior of this term is as N 2 , which means that this term is not negligible. We next look at the term with the permutation k 1 ↔ k 3 . Again the scaling behavior of this term is which, as we argued before, is negligible. Also for the term with the permutation k 2 ↔ k 3 and ν 1 > 1 2 , we have As we can see the scaling behavior of this term is as N 1 which means that, depending on the parameter space, this term can be significant. We next look at the term with the permutation k 1 ↔ k 3 . Again, the scaling behavior of this term is which, as we argued before, is negligible. Also for the term with the permutation k 2 ↔ k 3 and ν 2 > 1 2 we have Now as in [1], the term H in the 4th line can be approximated in the small −k 3 /k 1 x 3 limit. However, the term H (2) ν1 (− k3 k1 x 4 ) in the 4th line can not be can approximated. Then by redefining y 4 ≡ k 3 /k 1 x 4 , we get The scaling behavior of the above term is as N 1 . We next look at the term with the permutation k 1 ↔ k 3 : Again as in [1], the terms containing −k 3 /k 1 x i (i = 1, 3) can be approximated in the small argument limit. The reason for i = 1 is the following. In the integrand we have the factor H x 2 ) so if |x 2 | 1, this term becomes fast-oscillating and hence suppresses the integration. On the other hand, since the upper bound of the integral of x 2 is x 1 , |x 1 | < |x 2 |. So the terms containing −k 3 /k 1 x 1 is small. In addition, due tot he term H (1) ν1 (−x 3 ), the smallness of the term containing −k 3 /k 1 x 3 is somewhat more clear. Therefore, we have Then we look at the term with the permutation k 2 ↔ k 3 : Again as in [1], the terms containing −k 3 /k 1 x i , (i = 2, 3) can be approximated in the small argument limit. The smallness of the term containing −k 3 /k 1 x 3 is more clear due to the term H On the other hand, the reason for the smallness for i = 2 is that the upper bound of the integral of x 3 is x 2 . This means that |x 2 | < |x 3 |. So the terms containing −k 3 /k 1 x 2 is small too. Therefore, we have As in the above cases, for the small value of k3 k1 x 3 , we have The scaling behavior of the above term is as N 2 .
Again, the terms containing −k 3 /k 1 x i , (i = 1, 3) can be approximated in the small argument limit. The proof is the same as in the above cases: Then we look at the term with the permutation k 2 ↔ k 3 : Again, the terms containing −k 3 /k 1 x i , (i = 2, 3) can be approximated in the small argument limit: The scaling behavior of the above term is as N 2 . We next look at the term with the permutation k 1 ↔ k 3 : Then we look at the term with the permutation k 2 ↔ k 3 : (3 B3 )B H 3 33 : Since the analysis is very similar to the previous terms, we just mention the final result, Then we look at the term with the permutation k 1 ↔ k 3 : Finally we look at the term with the permutation k 2 ↔ k 3 : (4 B )B H 3 4 : Again, since this term can be obtained from B H 3 3 , we skip it now and insert it in the final result. Now we move on to the final case, C.
(1 C )C H 3 1 : Unlike the case [1], this term can be relevant in our case. The reason is that we have two different ν i and we can choose our parameter space such that the related term can be large enough. So we should consider all of the terms related to this part carefully.
One can easily see that the terms containing −k 3 /k 1 x i (i = 3, 4) can be approximated in the small argument limit. The result is We next look at the term with the permutation k 1 ↔ k 3 :  Then we look at the term with the permutation k 2 ↔ k 3 : (2 C )C H 3 2 : Since this term is very similar to the above case, C H 3 1 , we skip it at the moment and insert it in the final result.
(3 C )C H 3 3 : This part is related to our new interaction and is divided in three terms.
Then we look at the term with the permutation k 2 ↔ k 3 : (3 C2 )C H 3 32 : Again we just mention the final result: We next look at the term with the permutation k 1 ↔ k 3 :  Then we look at the term with the permutation k 2 ↔ k 3 : (3 C3 )C H 3 33 : Again we just mention the final result: