Moduli trapping mechanism in modular flavor symmetric models

We discuss how the moduli in modular flavor symmetric models dynamically select enhanced symmetry points at which the residual modular symmetry renders extra matter fields massless. The moduli dynamics non-perturbatively produces the extra matter particles, which gives (time-dependent) effective potential that traps the moduli to enhanced symmetry points. We show analytic estimates of particle production rate consistent with numerical results, and the dynamics of moduli based on the analytic estimates.


Introduction
In effective theory of superstring, moduli fields, light fields associated e.g. with higherdimensional gravitational degrees of freedom, play crucial roles in constructing realistic models of particle physics as well as cosmology.
In the modular flavor symmetric models, Yukawa couplings correspond to modular forms of complex structure moduli τ , and the structure of Yukawa couplings is determined by the vacuum expectation values (VEVs) of complex structure moduli, which are complex scalar fields in 4D effective theory.Therefore, moduli stabilization plays a crucial role in such flavor symmetry models.Indeed, moduli stabilization has been studied in modular flavor symmetric models [52,57,[80][81][82][83][84][85].
One of the ways to stabilize the moduli is the use of fluxes of p-form fields along compactified dimensions that yield effective potential to complex structure moduli.If somehow the chosen flux potential yields preferable VEV to the complex structure moduli, it would explain the flavor structure in the standard model by chance.Another possibility is due to non-perturbative effects.On the other hand, there might be some dynamical origin why the complex structure moduli may take some special value in the moduli space: In [95], it was pointed out that particle production takes place when moduli dynamically crosses the special point called the enhanced symmetry point (ESP) at which particles coupled to the moduli become massless.The produced particle gives effective potential to the moduli such that the moduli are attracted to the ESP and trapped around it.Such mechanism may explain why the VEV of moduli chooses some particular value from the landscape of moduli space.
In this work, we consider the dynamical moduli trapping mechanism in modular flavor symmetric models.First, with simple toy models we briefly review and develop a numerical approach to discuss the moduli dynamics back-reacted by the production of spectator fields.Despite its straightforwardness, we find a common difficulty due to the computational costs in the simulations.Furthermore, the modular flavor symmetric models contain extra complications due to the structure of Yukawa couplings as well as moduli being coupled to matters via Planck suppressed operators.Therefore, we will develop a semi-analytic approach where only the first particle production event is taken into account, and contributions from all the momentum modes are included.Using the analytic result, we show that the moduli trapping mechanism works even for (complicated) modular flavor symmetric models.
This paper is organized as follows.In Sec. 2, we review and discuss a fully numerical approach, which is in principle applicable to any dynamical systems.In Sec. 3, we review the modular flavor symmetry , and then specify the system we consider in this work and derive the equations of motion that we need to solve.We then analytically evaluate the amount of particle production when moduli crosses the ESP in Sec. 4. We also numerically examine the moduli dynamics using the analytic expressions we derived.Finally, we summarize our results and discuss implications to more realistic models in Sec. 5.In Appendix A, we give a brief review on particle productions in time-dependent background.Modular forms, which we use, are shown in Appendix B. In Appendix C, we review field dynamics in expanding Universe.
Throughout this paper, we take a natural unit convention = 1, c = 1 but we write the Planck scale M pl ∼ 2.4 × 10 18 GeV explicitly to clarify the hierarchy between various scales.

A toy model
We briefly review some ingredients that we will use in modular flavor symmetric models.For simplicity, we discuss the following system consists of two massless real scalar fields, We assume that φ behaves as a classical homogeneous field φ = φ(t), the background to be the Friedman-Robertson-Walker (FRW) spacetime ds 2 = −dt 2 + a 2 (t)dx 2 and χ to be a quantum field, and then the coupling term λφ(t) 2 χ 2 behaves as an effective time-dependent mass term of χ.We quantize χ as follows: First, the action of χ is which, by introducing χ = a 3 2 χ, can be rewritten as Then, the quantum field χ is written by where the mode function f k (t) satisfies and The mode function f k (t) satisfies the normalization condition which implies the canonical commutation relation for all t.In general, it is impossible to solve (2.5) analytically, but we are able to introduce a formal adiabatic solution (2.9) The auxiliary functions α k (t), β k (t) satisfy and 1 for all t, which is equivalent to (2.7).Note that the choice of the adiabatic solution is not unique, and we have chosen the zeroth order adiabatic solution with V k (t) = 0 in [96][97][98].We choose the initial condition α k (t) → 1 and β k (t) → 0 as t → −∞ which realizes the adiabatic past vacuum. 2In general, β k (t) becomes non-vanishing due to the time-dependence of background fields, which can be physically understood as "particle production from vacuum".We briefly review it in Appendix A. The production of particles back-reacts to the dynamics of the background field.Throughout this work, we assume that the backreaction affects only the modulus φ, and the background spacetime is intact.
Let us consider the dynamics of φ back-reacted by χ-particles, whose E.O.M. is given by φ where χ2 (x) ren is an expectation value of a renormalized χ2 operator with the adiabatic vacuum state.The vacuum expectation value of χ2 without renormalization is explicitly given by (2.13) We note that β k (t) is vanishing unless particle production occurs, which cannot be removed by local counter terms.Therefore, we may identify the β k -independent term to be a "vacuum" contribution.We first evaluate the "vacuum" contribution as where we have used the dimensional regularization H 2 is the effective mass, 3 and µ is renormalization scale parameter.We here applied dimensional regularization such that the invariance under trivial scale transformation a → ca, k → ck (x → c −1 x) with a constant c holds.Thus, with appropriate counter terms, we can formally write the E.O.M. of φ as where δV (φ) denotes quantum corrected potential terms of φ and where n k (t) ≡ |β k (t)| 2 and we have dropped fast-oscillatory terms since it would be averaged to be zero. 4We note that |β k (t)| typically decays faster than any powers of k, and therefore, there would be no UV divergence associated with it. 5 We would like to give some comments about solving the dynamics of the system under consideration: 1.In the following, we assume that the quantum-corrected effective potential δV (φ) = 0 throughout this work.Such a situation is effectively realized e.g. if the theory is supersymmetric.Strictly speaking, the time-dependence of backgrounds break supersymmetry spontaneously and there would be some potential.We will leave the effect of the quantum corrected potential for future work.
2. We use Emarkov-Milne equation, which directly yields the time-dependent particle number density: We use f k (t) = ξ k (t)e −iλ k (t) as an Ansatz for the mode function.
From the normalization condition (2.7) the function λ k (t) satisfies and the mode equation reads 3 The effective mass M 2 eff may become tachyonic, which leads to imaginary part to the effective potential.Such instability simply shows the existence of tachyonic modes that cannot be naively integrated out.Therefore, one has to treat such modes separately if exist.In the following, we consider only the cases where such instability shows up. 4 The neglected terms contain memory of the past, which makes the E.O.M. of φ an integro-differential equation. 5For adiabatic vacuum states β(−∞) = 0, β k (t) becomes non-zero due to Stokes phenomena/particle production, which is not a local but a global property with respect to time.Therefore, it is reasonable that the terms associated with such global (non-local) effects do not lead to UV divergences which are local effects due to short wave length modes.
with the initial values In terms of ξ k (t), the particle number density n k (t) is simply given by which is much easier to evaluate, since it can be solved simply as a set of differential equations.It is crucial to neglect the oscillatory terms in F (t) to avoid an integrodifferential equation.
3. As done in [95] we discretize the momentum integration appearing in F (t), which allows us to perform numerical simulations.However, even with such approximation, it is still a hard problem unless we reduce the number of k-modes in numerical simulations.To do so, it would be useful to rescale the momentum k by some reference scale.Assuming that the effective mass M 2 eff is dominated by the modulus coupling λφ 2 , the most relevant scale turns out to be v = λ Noting that where ξk = √ vξ k , we rewrite F (t) as Let us discuss how to perform the discretization of the integral effectively.We note that the particle number density after the first crossing of m 2 (t) = 0 is approximately given by This is why we have taken v as the reference scale.The mode π k2 ∼ O (10) does not contribute to the integration and here we take k = N to be the effective cut-off of the momentum integration where N = O(1) is a positive integer.Thus, F (t) is approximately given by where n is a positive integer characterizing the lattice spacing and kj = N j n . (2.26) We may recover the continuous case by taking N, n → ∞.Thus, we have reduced the problem to be a set of n differential equations. 7sing above approximation, we find a set of equations: ) As quoted in [95], in this model, parametric resonance occurs after the second and subsequent zero crossings, which extremely enhances the number density of low k modes and strengthens the trapping effect.
In numerical simulations, we have observed that the numerical solutions with different numbers of the mode number n are qualitatively similar to each other but quantitatively different.As another illustration, we consider the background where β is a positive constant, and we have normalized the scale factor such that a( t0 ) = 1 at the initial time t = t0 .Then, In Figs 2,3 we show the numerical solution for the matter dominated universe β = 2 3 and the radiation dominated universe β = 1/2, respectively.These figures show that the trapping mechanism works even in the expanding background.Note however that, if the initial Hubble parameter is large enough ( t0 ∼ O(1)), it is possible to slow the trapping by the dilution of particles as well as Hubble friction.Let us consider a more involved case that instead of λφ 2 the mass term is given by where A is a real parameter.In this case, there is a true symmetric point φ = n ∈ Z where χ becomes massless, but there is a fake symmetric point as shown in Fig. 4. The E.O.M. of φ turns to φ + 2πµ 2 (sin πφ + A sin 3πφ)(cos πφ + 3A cos 3πφ) χ 2 = 0. (2.32) The behavior of the φ-dependent mass (2.31) for µ 2 = 1, A = 0.7.φ = n ∈ Z are the true symmetric points and there are fake points, where the mass becomes small but non-zero.
We show two numerical solutions with different parameters in Minkowski spacetime a(t) = 1 in Figs 5, 6.We note that the false vacua are φ = 2n+1 2 , (n ∈ Z).In the former case, the modulus is immediately trapped to the true symmetric vacuum φ = 1, whereas the latter shows that the modulus is trapped both at the false and the true symmetric vacua, but finally reaches the true symmetric vacuum.We expect that modulus trapping at the true vacuum is not a generic property since the effective potential vanishes even at the false vacua, which can be seen from the E.O.M. of φ.Nevertheless it is true that the particle production takes place more efficiently near the true vacua since χ becomes much lighter at the point.Therefore, if the true and false vacua are sufficiently separated, we expect the modulus to be trapped near the enhanced symmetry point.
Before closing this section, we finally comment on some technical issues remaining in our numerical approach shown here.We have tried our numerical simulations with different set of parameters or initial conditions several times.Although the behavior qualitatively converges, the results are quantitatively unstable under the change of parameters such as mode numbers.Such a behavior may be understood from the fact that parametric resonance become important when the modulus φ is almost trapped to the ESP, and the resonance is sensitive to the parameters such as momentum.As we will discuss later, the modular flavor symmetric models contain more issues in performing fully numerical simulations.Therefore, we will propose a semi-analytic approach to capture the moduli dynamics.Nevertheless, the numerical approach we have taken in this section would still be useful for some class of models.

Modular flavor symmetric models
We show the modular flavor symmetric model considered in the following sections.In this section, we first briefly review modular flavor symmetry and then derive the E.O.M. of complex structure moduli as well as a spectator scalar χ while we will study in detail the χ-particle production in the next section.

Modular flavor symmetry
The SL(2, Z) group, where a, b, c, d are integer and ad − bc = 1, is generated by two generators, S, and T , They satisfy the following algebraic relations: This group is referred to the homogeneous modular group, Γ = SL(2, Z).
Under the modular symmetry, the modulus τ transforms as Note that τ is invariant under S 2 .That is, the generators satisfy on the modulus τ .This group is referred to the inhomogeneous modular symmetry, Γ = P SL(2, Z) = SL(2, Z)/Z 2 .Also we define the congruence subgroup Γ(N ), which includes T N .Similarly, we can define Γ(N ).
The modular forms f (τ ) i are the holomorphic functions of τ , which transform under the modular symmetry, where k is the modular weight and ρ(γ) ij is a unitary matrix.Suppose that for γ ∈ Γ(N ).Then, the matrix ρ(γ) ij represents the quotient, Γ N = Γ/ Γ(N ), where T N = 1.Interestingly, these quotients Γ N with N = 2, 3, 4, 5 are isomorphic to S 3 , A 4 , S 4 , A 5 , respectively.There are fixed points on τ , i.e. τ = i, ω = e 2πi/3 , i∞, where residual symmetries remain.That is, Z S 2 , Z ST 3 , Z T N symmetries remain at τ = i, ω = e 2πi/3 , i∞, respectively.The modular form f (τ ) has definite charges under these residual symmetries and their behaviors are determined by their charges.
As illustrating examples, we study modular A 4 symmetric models in this paper.Modular flavor symmetric superpotential can be written by where Φ i are the chiral superfields with modular weights k i and they have some representations under A 4 .The Yukawa couplings Y ijk (τ ) are modular forms.The superpotential must be invariant under the modular symmetry including the A 4 symmetry.The Kähler potential of chiral matter fields can be written by The fundamental modular forms of A 4 have the modular weight 2 and they correspond to the [27].Their explicit forms are shown in Appendix B. The modular forms of higher modular weights can be obtained by their tensor products.For example, we use the modular forms of weight 8, because all of three A 4 singlets, 1, 1 , and 1 appear as modular forms when the weight is 8, i.e.Y As a simple illustrating model, in the following sections we consider a single complex scalar spectator field that has the mass term where M χ (τ, τ )| 2 is given by modular forms.Such a complex scalar field can appear in supersymmetric models with where Φ is a chiral superfield whose scalar component is χ, k is the modular weight of Φ, m a mass parameter, and Y (τ ) a holomorphic function of τ given by modular forms.S is an additional chiral superfield which is a singlet under modular symmetry. 8Although there are fermions in such a model, which are as light as the scalar χ, we only discuss the complex scalar χ only.The fermionic particle production can be similarly discussed along the line e.g. of [99][100][101] but would be more involved than the bosonic case we consider here because of the chiral structure.
We would like to comment on the relation to the standard model.In the supersymmetric standard model, mass terms of scalars and fermions through Yukawa couplings are not generated until the electroweak symmetry is spontaneously broken.Therefore, we do not expect standard model particles to contribute the moduli trapping. 9Instead, we can identify Φ as a superfield that obtains its mass term through the GUT symmetry breaking 10 .Therefore, the mass parameter m can be sufficiently large (but smaller than M pl ).

Dynamics of moduli and spectator scalars
We consider quantization of a spectator complex scalar field χ that couples to the complex structure modulus τ , and the action of χ is given by (3.14) We introduce quasi-canonical real fields which will become a choice that makes φ canonical.We assume that the "bare mass" M 2 (τ, τ ) is given by modular forms as where m denotes some mass scale such as the VEV of a GUT Higgs field, and the background spacetime is ds 2 = −dt 2 + a 2 (t)dx 2 .Introducing χ = χ1 + i χ2 , we rewrite the action as which can be made canonical by introducing a new basis χi = αχ i where 2 kφ/M pl . (3.18) With the new basis, we find where we have assumed that the background field depends only on t, and used integration by parts in the second equality.Thus, the effective mass of the canonically normalized spectator field is (3.20) 9 Scalars would have mass terms by supersymmetry breaking and they generally depend on moduli.
Therefore, the moduli dynamics may generate scalars in the supersymmetric standard model. 10Another candidate would be right-handed (s)neutrinos, where the superpotential can be written by The time derivative terms of φ in the second line may be interpreted as the higherdimensional curvature induced mass as Imτ is related to the area of the extra dimensional torus.We will neglect such terms as well as Hubble induced terms in the following. 11We also notice that the effective mass is modular invariant as it should be.
Let us consider the effective action of the complex structure modulus τ given by Note that the appearance of non-covariant expression +a −3 M 2 eff 2 i=1 χ 2 i is the result of making χ i canonical. 12The equations of motions of φ, θ are given by It is straightforward to generalize our previous discussion to this case by using The vacuum expectation value of χ2 1 is given by Thus, the renormalized effective potential is found to be where δ CT denotes possible local counter-terms eliminating divergent pieces. 13As the previous sections, we assume the cancellation of all but the last term in (3.26).
In analysing the dynamics of the moduli fields φ, θ, there are several technical issues in this model in addition to the ones discussed within toy models: Since the coupling between χ and moduli are suppressed by Planck scale which causes hierarchically small or large numbers, and the moduli dependence is quite complicated, the numerical costs become more than that in the toy models.More specifically, the effective field theory description requires the field velocity and the Hubble parameter to be much smaller than the Planck scale.Furthermore, if the moduli trapping occurs after inflation, the Hubble scale needs to be less than about O(10 13 GeV) by the constraint on the tensor-to-scalar ratio.Therefore, it would be natural to assume the ratio between the initial field velocity and Planck scale to be below O(10 −5 ), which appear in couplings if we normalize all dimensionful parameters by the scale of initial field velocity.The complication of the mass term given by modular forms shown later further makes the numerical simulations difficult.Therefore, we will develop semi-analytic approaches to discuss the moduli dynamics within modular flavor symmetric models in the next section.

Particle production at ESPs in modular flavor symmetric models
In order to overcome some technical difficulties within our setup, we consider a semi-analytic approach, where we analytically estimate the number density of the spectator scalar particle produced at the first particle production time at which M 2 χ (t 0 ) ≈ 0. Using the estimate, we analytically evaluate the effective potential arising after the first particle production event, and numerically solve the E.O.M. of moduli fields with the estimated effective potential, which does not contain the difficulties mentioned before.Although this approach cannot capture the subsequent particle production events, we expect such events just strengthen the trapping potential, which does not change the dynamics of moduli qualitatively.
We first show analytic formulas for particle number density produced at the crossing of the ESP.To do so, we need to know the behavior of the effective mass near the ESPs.In modular symmetric models, the ESPs are τ = i, e 2iπ/3 (≡ ω), +i∞ where some modular forms vanish.In the following discussion, we focus on τ = ω as a representative case, but similar analysis can be done for any other critical points in a similar way.
In the following subsections, we first discuss the behavior of the effective mass of the spectator scalar field near the ESP τ ∼ ω, which enables us to analytically estimate the amount of particle number density produced after the crossing with the ESP.
With the aid of the analytic expression of the effective mass around τ ∼ ω, we discuss one-dimensional dynamics of θ or φ where one of them is fixed to a constant value.We give analytic formulas for the particle number density as well as the effective forces of the produced particle on the modulus.

The behavior of the effective mass near the critical point
The modular forms can be classified as representations of subgroup.The singlet modular forms of weight 8, Y (8) which can be equivalently written as where is a variable for the deviation from the symmetric point τ = ω.Therefore, we obtain where r (u).Expanding both sides with respect to ω yields where q r = 0, 1, 2 for 1, 1 , 1 , respectively.This relation implies that = 0, and we find ) Note also that The leading order terms of each singlet 1, 1 , 1 are as follows: Y ( 8) and the last one shows that 1 cannot produce the particles at ESP τ = ω, because χ is still massive. 14.Thus, the effective mass can be approximated as where C ≡ , and ellipses denote the terms higher order in (τ − ω).Before going to the details of particle production, we would like give a few comments on the behavior of the effective mass for each representation.We notice that depending on the representation r, the behavior of the mass term near the ESP τ = ω changes.Since the moduli dependence comes with powers of χ/M pl or θ/M pl , we can say the particle production in r = 1 case should be much smaller than that in r = 1 .Therefore, the moduli trapping effect becomes more significant for r = 1 .Such observation is based on the fact that moduli are gravitationally coupled to matter fields, which are completely independent of particle production dynamics discussed below.

1D dynamics: φ-fixed
In this subsection, we discuss the particle production due to the θ dynamics while φ is fixed at the ESP φ = φ 0 where e 2 , and we parametrize θ as where t 0 is the time at which θ crosses the enhanced symmetry point and v > 0 is the velocity of θ at t = t 0 .This parametrization is a good approximation if Hubble parameter is smaller than √ v. (See appendix C.) Then, near t ∼ t 0 the effective mass can be approximated by (for r = 1), (4.17) and the effective frequency for k-mode is 14 More precisely speaking, for r = 1 , particles are produced if τ = ω is the local minimum of Y (8) 1 , but as the leading term is non vanishing, the particle production would be less than that of r = 1 or 1 where we have approximated the scale factor by its value at t = t 0 .We will show the analytic estimate of the particle production taking place around t ∼ t 0 separately for both r = 1 and r = 1 respectively.r = 1 case In order to estimate the particle production rate, we need to find the turning points at which ω 2 k (t) = 0 in a complex t-plane.(See [97,98,[102][103][104][105][106] for review.)With (4.17) the turning points are found to be , (4.20) with which the effective frequency can be rewritten as where t * can be any of turning points and Notice that there are two pairs of turning points (t + i , t − i ) (i = 1, 2) which are complex conjugate to each other, and Stokes lines connecting them crosses the real t-axis.The amount of particle production can be approximately given by a simple formula (see e.g.[97,98,[102][103][104][105][106] for details 15 ) where we have assumed that t is sufficiently later than t 0 , and ) More explicitly we find Therefore, the particle number density after crossing two Stokes lines can be approximately estimated as where We have numerically checked the validity of this formula, and find a very good agreement of the numerical result and the analytic formula as shown in Fig. 7. Thus, we can estimate the total particle production as well as the effective potential on the basis of the approximate formula (4.28).
0.5 1.0 1.5 Log 10 n p Figure 7.Comparison between the analytic formula (4.28) (the blue line) and numerical results (red dots).The momentum p is in units of the velocity |v| at the particle production time t 0 .In this example, we have taken vt 0 = 50, vt ini = 15, β = 1 2 , |v|/M pl = 5 × 10 −4 , m/ |v| = 5 × 10 2 .We have also checked that the agreement is quite well for various sets of parameters.
With the approximate particle number formula for a k-mode (4.28), we are able to estimate the total number density to be 16 Even though we have an approximate formula of the produce particle number density, it is still difficult to evaluate the effective force given by Since the effective support of the integrand is localized to small k, we may approximate ω p ≈ M eff in the integrand, which allows us to approximate F θ (t) as where we have introduced the Heaviside theta function Θ(t − t 0 ). 16The oscillatory part e γ k cos γ k turns out to be vanishing in the momentum integration.r = 1 case Next, we consider the case r = 1 where the effective mass is given by (4.18), and find a pair of turning points Similarly to the previous case, the particle number produced at t = t 0 can be evaluated by We show the comparison of our analytic formula (4.34) and numerical results in Fig. 8. Again, we have found an excellent agreement between them.Log 10 n p Figure 8.Comparison between the analytic formula (4.34) (the blue line) and numerical results (red dots).The momentum p is in units of the velocity |v| at the particle production time t 0 .
Accordingly, the total particle density is found to be and the effective force is approximately given by

1D dynamics: θ-fixed
We consider the case that θ is fixed at θ = − 2 and parametrize φ as with which the effective mass can be rewritten as where we have expanded the terms higher order in v(t − t 0 )/M pl since it is generally small.We notice that these expressions are the same as the φ-fixed case (4.17), (4.18) by replacing |C| 2 → 9 16 |C| 2 and |D| 2 → 3 4 |D| 2 , respectively.Thus, the particle number density of the k-mode is given by Therefore, the effective force on φ due to particle production can be written as 17 or We note that if there are N fields that has the same Yukawa coupling, the effective force would be multiplied by N , which enhances the effect.We have checked the agreement between analytic formulas (4.40) and numerical results shown in Fig. 9.We show a numerical example of the moduli dynamics by using our analytic approximation in Fig 10 .In the example, we have used the 1 -model.We have used the formulas (C.3),(C.4) to choose appropriate initial conditions such that φ crosses the ESP at t 0 . 18As a cross-check of our numerical solution, we show the behavior of the field velocity 17 Recall that R ∝ |C| 2 and also that F φ does not have an extra factor (Imτ ) 2 while F θ does. 18Another technical note: We have modified the denominator of (4.42) as M eff → 10 −6 v + M 2 eff such that the singular behavior at the critical point is avoided.We have checked that the result is not affected by the IR cut-off.φ in Fig 11 .As is clear from this simulation, the moduli trapping works despite Planck suppressed couplings between the matter field χ and the moduli.Thus, the moduli fields seem to prefer the ESP if they cross such a point along their time-evolution.We emphasize that here only the first particle production event is taken into account, but as we see, there would be secondary and more particle production events when the modulus crosses the ESP, which further strengthen the trapping effect as quoted before.

Summary and discussion
In this work, we have studied the moduli trapping mechanism due to the particle production near ESPs within modular flavor symmetric models.We have reviewed and developed a general method for numerical simulations of models that the classical background fields and quantum fields interact with each other.Despite generality and simplicity, we have found that such an approach is not suitable for the application to modular flavor symmetric models The time dependence of the field velocity φ with (left) or without (right) particle production effect.We have taken v as the unit of the velocity.The parameters are the same as that in Fig. 10.We see that in both cases the velocity is because of the complexity of the couplings between background fields and quantum fields.Therefore, we have developed a semi-analytic approach where we analytically evaluate the particle production near the first crossing with the ESP, which yields the effective potential arising from produced particles.Although this approach looses the effect of subsequent particle production events, we have found the expected behavior, namely that the moduli in the modular flavor symmetric models can be trapped around the ESP at which there is a residual discrete symmetry.
Although we have studied moduli trapping effects due to the scalar field χ, we could discuss effects due to spinor and vector fields in a similar way.As illustrating models, we have used the A 4 modular forms of weight 8, which have suppressed values around τ = ω.Similar results would be obtained for modular forms of generic weights and other finite groups, if they have suppression behavior around τ = ω and τ = i.Furthermore, concrete modular flavor symmetric models include many fields, whose masses are determined by similar modular forms.Thus several modes would be produced around fixed points, which enlarges the effect of moduli trapping.
There are several issues that should be addressed in future work.One is to embed this mechanism into more realistic models of particle physics.In particular, it would be interesting to study the moduli dynamics within magnetized orbifold models where the standard model flavor structure including flavor mixings as well mass hierarchies has been realized [108][109][110][111][112][113][114].It is also important to notice that the moduli trapping due to particle production never completes the moduli stabilization since the effective potential disappears as the particle numbers are diluted by expansion of the Universe, which therefore requires the introduction of additional potential terms.In particular, we have not taken into account the 1-loop effective potential such as the Coleman-Weinberg type potential that we have shown.Although it becomes small in supersymmetric models, it would be important to investigate the effect of such potential and the fate of moduli in the late Universe.
Another interesting issue would be the cosmological implications of moduli dynamics.It has been known that the particle production during inflation may give imprints on the curvature perturbation spectrum that can be seen from observations of cosmic microwave background [115][116][117][118][119]. In particular, moduli in modular flavor symmetric models may become inflaton directions.(See e.g.[120,121].)Furthermore, the presence of the dynamically changing CP phase for the standard model matter fields may realize baryogenesis.In such a case, the method we have applied in this work would be useful to discuss the dynamics of matter and moduli simultaneously.We leave these interesting questions and model buildings for future work.

A Review of particle production in time-dependent backgrounds
We give a short review of particle production in time-dependent backgrounds.Let us consider a scalar field having a time-dependent mass, where f (t) satisfies a mode equation (2.5) with a(t) = 1 and the normalization condition (2.7).We have "formally" introduced creation and annihilation operators, which satisfies the canonical commutation relation (2.8).This is yet insufficient to give any meanings to the vacuum state that is annihilated by âk since we have not determined f k (t).More precisely speaking, unless the boundary condition of f k (t) is specified, the above expansion has no physical meaning.In Minkowski spacetime, the annihilation operators are introduced as coefficients of positive frequency modes ∼ 1 √ 2ω k e −iω k t .In time-dependent backgrounds, it is impossible to define a "global" positive frequency mode, but locally (in time) it is possible to find an approximate solution.For instance, the formal solution (2.9) becomes asymptotically a positive frequency mode if α k (t) → 1 and β k (t) → 0 as t → −∞. 19Indeed, under this condition, the scalar operator becomes and the vacuum state that satisfies âk |0 in = 0 is understood as a vacuum state in the asymptotic past.Thus we have defined a past (adiabatic) vacuum.However, the "past" vacuum state is not a "vacuum" for a future observer.In general, β k (t) becomes non-zero due to Stokes phenomena20 , t being asymptotic future time.In such a case, it is not appropriate to call f k (t) as a "positive frequency mode".Notice that the scalar operator can be expanded as The new set of creation and annihilation operators bk , b † k defines a future vacuum state |0 out .Now, we find that the future creation and annihilation operators are given by linear combinations of the past creation annihilation operators.Note that there is no momentum exchange and the linear combination is diagonal with respect to momenta k.
How can we find particle production from vacuum?We define the future particle number density as We would like to know how much "future" particles are contained in the past vacuum state |0 in , which can be explicitly evaluated as where we have used lim k→0 δ 3 (k) = V (2π) 3 and V being spatial volume.Thus, we have found that the number density for (future) k-mode is given by (A.6) It turns out that the past "vacuum" state is not a "vacuum" for a future observer viewpoint.
The "particle production from vacuum" would be more appropriate to be understood as the ambiguity of energy, which does not allow us to globally define positive and negative modes.
We emphasize that the definition of "particle" is quite ambiguous except for constant backgrounds.In particular, when time-dependence is not turned off, we have to introduce adiabatic solutions such as WKB solutions, but there are infinitely many choices of adiabatic solutions.There is no clear answer to the question "which solution should we take?" but it is known that the optimal definition of the "adiabatic particle number" is related to the Stokes phenomena [97,98]. 21

C Field dynamics without potential in expanding Universe
Here we discuss the free field dynamics in the expanding Universe.The E.O.M. of a massless free scalar field in the FRW background is given by φ + 3H φ = 0, (C.1) which can also be written as Now since H(t r ) ∼ 1/t r , we can neglect the third and higher order terms as long as Hubble parameter at t = t r is sufficiently smaller than | φ(t r )|.
The relation between the initial velocity v ini = φ(t ini ) and v = φ(t 0 ) at some reference time t 0 > t ini is This would be useful to relate the initial velocity with the field velocity at the crossing with enhanced symmetry points particularly in estimating the particle production rate.

Figure 9 .
Figure 9.Comparison between the analytic formula (4.40) (blue lines) and numerical results (red dots).The left panel shows the r = 1-case and the right the r = 1 -case.We see a good agreement for both cases.

Figure 10 .
Figure 10.A numerical solution of the E.O.M. of φ with (left) or without (right) particle production effect.Here we have taken the following parameters: The expansion parameter β = 1 2 , the initial time √ vt ini = 5, the time of particle production √ vt 0 = 50, the Yukawa coupling parameter m/ √ v = 100, and the species number to be N = 1.The field velocity at the particle production is √ v = 10 −4 M pl .The blue solid curve is the trajectory of φ and the red dashed line is the critical point Imτ = √ 3 2 .

Figure 11 .
Figure 11.The time dependence of the field velocity φ with (left) or without (right) particle production effect.We have taken v as the unit of the velocity.The parameters are the same as that in Fig.10.We see that in both cases the velocity is 1 at vt = vt 0 = 50.