HYM-flation: Yang-Mills cosmology with Horndeski coupling

We propose new mechanism for inflation using classical SU(2) Yang-Mills (YM) homogeneous and isotropic field non-minimally coupled to gravity via Horndeski prescription. This is the unique generally and gauge covariant ghost-free YM theory with the curvature-dependent action leading to second-order gravity and Yang-Mills field equations. We show that its solution space contains de Sitter boundary to which the trajectories are attracted for some finite time, ensuring the robust inflation with a graceful exit. The theory can be generalized to include the Higgs field leading to two-steps inflationary scenario, in which the Planck-scale YM-generated inflation naturally prepares the desired initial conditions for the GUT-scale Higgs inflation.


Introduction
The 2015 Planck's collaboration release [1] confirmed that inflation had occurred at relatively low energy scale, 10 15 ÷ 10 16 GeV, with almost absent non-Gaussianity, and small tensor fluctuations being. A variety of phenomenological models with a scalar inflaton slowly rolling down in a flat potential describes well the data. These include Starobinsky R 2 model [2], Higgs inflation [3] and some other traditional models demanding tuning of parameters at the classical level. The tuning, however, is not protected from large quantum corrections, so various attempts were undertaken to find symmetries underlying the desired flatness of the potential making inflation "natural". One possibility, based on scalar fields only, invokes models with hidden conformal or shift symmetry [4,5,6,7]. The simplest such model [4] contains in the Jordan frame the following combination of the Einstein term and the conformally coupled scalar field term vanishing at the boundary φ 2 = 6. In the Einstein frame this leads to flattening of the scalar potential V (φ) → V ( √ 6 tanh(ϕ/ √ 6), which results in de Sitter solution as ϕ → ∞. Higgs inflation is particularly attractive since it identifies the inflaton with some known field. Note that this mechanism invokes the non-minimal coupling of Higgs to gravity via the curvature scalar.
The GUT-scale conformal inflation, however, raises the problem of initial conditions, discussed already in the early days of inflation scenario [8], and recently reconsidered again in the modern setting [9,10]. One of its solutions is a stage of preliminary inflation starting at Planck scale and driving the inflaton to the plateau of the GUT-scale observed inflation. The second inflaton is usually taken as another scalar field. Here we suggest to use for this purpose the vector YM field which is generically present in gauge/supergravity theories, so one does not need to introduce the second scalar by hand. To realize such a scenario one has to assume, likewise in the Higgs case, the non-minimal coupling of YM to gravity using the Horndeski prescription [11].
Recall that the SU(2) Yang-Mills field has an isotropic and homogeneous mode and satisfies (in the case of the standard YM lagrangian) the conformal equation of state p = ǫ/3, thus mimicking the hot Universe [12]. It was studied in 90-ies both in classical and quantum minisuperspace settings [13] and Euclidean quantum gravity [14]. Recently this idea was revived in anticipation of the future precise measurements of the primordial gravitational waves imprint on the cosmic microwave background (CMB) [15]. It was found that the perturbations of YM and tensor gravitational perturbations mix together, leading to difference in the evolution of right-and left-polarized gravitational waves (parity violation) which can become testable soon.
During past decade various modifications of the standard YM action breaking the conformal symmetry in a way consistent with the Standard Model and its extensions were introduced. Replacing the YM lagrangian by the Born-Infeld string motivated lagrangian, e. g., leads to an equation of state interpolating between that of the string gas and the photon gas −ǫ/3 < p < ǫ/3, but this is insufficient for inflation [16]. Phenomenologically, it was noticed that the lagrangian generically depending on two invariants L(f, g), f = F a µν F aµν , g =F a µν F aµν , leads to desired inflationary equation of state if the dependence on g is non-linear, i.e. ∂ 2 L/∂g 2 = 0 [17]. A particular model of this type with the quadratic g 2 −term in the lagrangian was called "gauge-flation" and received a lot of attention [18]. But physical origin of such a term, which, moreover, must enter with a large coefficient, is somewhat obscure.
Conformal symmetry is also broken once interaction of the YM field with the scalar fields is introduced. This is what happens in the full gauge theories. It was shown that the YM-Higgs model with the complex doublet Higgs leads to a kind of hybrid inflation scenario, [19] in which dynamics of the Higgs field is modulated by the YM component (this model was recently revived in [20]). Other typical interactions of YM fields include dilaton [21] and axion. The latter option attracted much attention as implementing the idea of "naturalness" [22] and became known as "chromo-natural" inflation [23].
The next class of models, closer to the present one, consists in exploiting the non-minimal coupling of the inflaton to gravity, which was applied to Higgs field under the name of Higgs inflation [3,24]. The idea was also applied to the YM field possibly together with modifying gravity lagrangian [25,26,27,28]. The potential danger of curvature-modified gravity is the emergence of higher derivative field equations plagued with the Ostrogradski ghosts. In attempts to avoid ghosts, new ideas associated with massive gravity and/or galileons were invoked [29,30]. General classes of couplings of vector fields to gravity involving curvature tensor couplings whose equations of motion does not contain ghosts were found by Horndeski [11]. Initially the non-minimal vector coupling to gravity was introduced as the extension of the Maxwell theory in curved space which preserves the second order equations of motion, admits the energy-momentum and charge conservation laws, and reduces to Maxwell theory in the flat space limit. Later it was revealed that Lovelock gravity [31], galileon models [29] and Horndeski theory [32,11] are strongly interrelated [30]. Using the Abelian vector fields in cosmology [33,34] leads either to isotropy or gauge invariance problems. So here we consider the unique case free from these complications which was not discussed before: the lowest order Hordeski coupling of the SU(2) YM field to the dual Riemann tensor. This model contains only one extra parameter of the dimension of mass which turns out to be the Hubble constant of the de Sitter stage in this model. We demonstrate that de Sitter solution is the boundary of the solution space which attracts a large set of trajectories, keeps them for some finite time and then relaxes to the hot universe state.
2 Non-minimal coupling of vector field to gravity General gauge-invariant curvature-dependent action quadratic in the vector field strength F µν and linear in the curvature can be written in the form where the susceptibility tensor R αβµν has the same index permutation symmetries as the Riemann tensor. It can be presented as the linear combination where R αβµν is the Riemann tensor, R αµ is the Ricci tensor, R is the scalar curvature, and the brackets [ ] mean an alternation over indices with the factor 1/2. Such a structure is typical for the one-loop corrections to the Maxwell action in curved space QED [35], where the coefficients q 1 , q 2 , q 3 have certain particular values. Here we consider this action as phenomenological, but subject to some theoretical restrictions. The field F µν in (2) can be either Abelian, or non-Abelian, in which case we will use the matrix notation, assuming the SU(2) gauge group and adding the trace operatorR before the lagrangian. For generic coefficients q 1 , q 2 , q 3 the resulting theory contains higher derivatives generating extra degrees of freedom which are plagued with Ostrogradski ghosts. The unique curvaturedependent coupling leading to the ghost-free theory was found by Horndeski [11]. It corresponds to q 1 = q 2 = q 3 , in which case the susceptibility tensor reduces to the double-dual Riemann tensor: where the Levi-Civita tensors contain suitable √ −g factors. This tensor satisfies the Bianchi which is crucial for making the theory ghost-free. Note that the Horndeski action can be written in two equivalent forms: using the dual field tensorsF αβ ≡ 1 2 ǫ αβµν F µν . This structure is reminiscent of the Gauss-Bonnet lagrangian from which it can be obtained replacing the Riemann tensor by the product of two field tensors. This is not accidental: the Horndeski action can be derived from the higher-dimensional Gauss-Bonnet theory by dimensional reduction [36]. It is worth noting that the vector Horndeski lagrangian in four dimensions is much simpler than the scalar Horndeski one [32] which was widely used recently in attempts to improve the simplest non-minimal Higgs inflation model [37].
Using the variation of the Riemann tensor and the Bianchi identity for the YM field where the gauge covariant derivative is introduced, one can write the variation of the action (8) over the metric in the form in which the absence of higher derivatives is manifest. On the contrary, one can notice the presence of the cubic term in F . We conclude this section with brief review of the earlier proposals to use non-minimally coupled vector fields in cosmology. The slow-roll inflation model with (generic) non-minimally coupled Maxwell field was suggested in [38], yet suffering the issue of anisotropy. The latter has been evaded in the non-Abelian case [39,28], however these earlier models either were loosing gauge invariance or contained ghosts. The Abelian vector model with Horndeski coupling both gauge invariant and without ghosts was investigated in [33], but this model was unable to provide de Sitter solutions unless the cosmological constant was added by hand [34]. The Yang-Mills-Higgs cosmology with general non-minimal coupling (3) was examined in [25], where de Sitter solutions were found, but no detailed investigation of inflation was undertaken.
In what follows we will study inflationary solutions in the ghost-free SU(2) Horndeski Yang-Mills model showing that robust inflation emerges because of general property of the Horndseki coupling, which closely resembles the coupling of the scalar field used in the models of conformal attractors. It turns out that the YM non-linearity is crucial for possibility of this scenario: it is impossible in the Maxwell case.

HYM cosmology
It is convenient to rescale coordinates and the YM potential as follows: GeV is the Planck mass, and g is a gauge coupling constant. Then we choose the units gM pl = 1. Actually, in most gauge theories the coupling constant is of the order of unity, so in what follows we assume gM pl ∼ M pl .
Adding the Einstein term and the standard YM term to the Horndeski action we obtain the total action where the "induction" tensor is introduced with the coupling µ −2 (dimensionless in the rescaled quantities): Even before passing to the Friedmann metrics, one can notice the following fundamental property of the HYM action: its matter part vanishes in de Sitter space with some curvature radius. Indeed, adjusting the Hubble parameter of de Sitter to be in which case Φ µν = 0. Obviously this has to be the boundary of the physical domain, beyond which we would get the phantom YM field (such an option is not considered here). We will see that the inflationary solutions are attracted to the phantom boundary, but do not cross it.
It is worth noting that the existence of the de Sitter boundary is encountered in the more general non-minimal theory (2) too (see [25]). In this case we get so in de Sitter space with the Hubble parameter, satisfying (6(q 1 − q 2 ) + q 3 )H 2 = µ 2 instead of (16), the induction tensor also vanishes. The Einstein equations of the HYM theory can be written in the usual form where the effective energy-momentum tensor reads The equations of motion for YM field are simply Now we pass to the homogeneous and isotropic cosmology, restricting for simplicity by the spatially flat metric: The YM matrix-valued one-form can be written in certain gauge in terms of a single function ψ(t) (for more general gauges and any spatial curvature see [40,16]): showing that the direction of A in the color space coincides with the space direction. Choosing the proper time gauge N = 1 and introducing the 'electric' and 'magnetic' effective fields: E ≡ψ + Hψ, H ≡ ψ 2 we can present the standard YM Lagrangian in the Maxwell form: The effective one-dimensional Lagrangian of the SU (2) HYM model (14) then reads: where one can notice the factor 1 − H 2 /µ 2 in the gauge sector, indicating on the de Sitter boundary described above.
It is easy to check that the stress-energy tensor is diagonal and isotropic, with the energy density and pressure In the limit of vanishing coupling, µ → ∞, the system represents radiation with the equation of state p g = ρ g /3.
The dynamics of the system is governed by Friedmann equations, and the gauge field equation, An important characteristic of the system (28)(29) is the determinant of the matrix of coefficients before the derivativesḢ,ψ : When this quantity vanishes, the solution meets the singularity. One can show that the boundary H 2 = µ 2 separates the domain of non-singular solutions from that of singular ones. Consider a solution crossing the boundary at the moment t = t 1 , so that H 2 (t 1 ) = µ 2 ,Ḣ(t 1 ) = 0. Then the gauge field equation (29) implies that the expression in the round brackets in the second term vanishes: Hψ + H 2 ψ + ψ 3 = 0 . This implies vanishing of the determinant D, indicating the singularity. Thus the non-singular trajectories should not cross the boundary H 2 = µ 2 . The physical trajectory must be non-singular and reach the flat space asymptotic: H = ψ =ψ = 0, which implies D → 1. Therefore, physical initial states should reside in the following domain of the phase space: to which the flat space asymptotic belongs. Note that with H < µ the sign of the kinetic term of the Lagrangian (25) remains positive. Thus the boundary H 2 = µ 2 also preserves the system from falling into the phantom state.

HYM-flation temporary attractor
Contrary to earlier negative verdict concerning the Abelian Horndeski-Maxwell cosmology [33,41] we would like to show that non-linearity of the YM theory makes the proposed HYM model much more promising. Our claim is that de Sitter solution H 2 = µ 2 is an inflationary attractor in non-Abelian Horndeski model; robust inflation emerges without fine-tuning of parameters or initial conditions.
To show this we first observe that the YM Eq. (29) is satisfied if H 2 = µ 2 ,Ḣ = 0. Then one can solve the first Friedmann equation in (28) as a quadratic equation inψ: Since only two of the three equations (28,29) are independent, the remaining second Friedmann equation will also be satisfied with H = µ andψ given by (32). Obviously, the square root in (32) is real only if the YM function is above the critical value ψ > ψ cr , satisfying ψ 6 cr + (3/2)ψ 4 cr µ 2 − µ 4 = 0. This solution is possible due to the YM non-linearity, which manifests itself in presence of the terms ψ 6 , ψ 4 under the square root. In the Abelian case the first Friedmann equation would imply µ 2 = −E 2 for the ansatz H 2 = µ 2 ,Ḣ = 0.
The critical value, ψ cr , is proportional to µ 2/3 for µ ≪ 1, and to µ 1/2 for µ ≫ 1. The large field limit, ψ ≫ max(µ, µ 2/3 ), corresponds to the dominance of ψ 6 under the square root in (32) and always satisfies ψ ≫ ψ cr , as well. Then the two branches of the solution (32) simplify and can be easily integrated:ψ The first branch corresponds to dominance of the kinetic term, E ≫ H, while in second case the YM potential prevails, H ≫ E. Consider now small deviations (δH, δψ, δψ) from these solutions and compute the eigenvalues of the corresponding linearized systems. These values can be viewed as local Lyapunov exponents [42] which describe the growth rate of the deviations in a given mode (note that the solutions (33) are not the stationary point of the system). The result reads: From this we deduce that ψ + is an unstable singular solution, while ψ − describes the mode which is stable for some period of time. For this reasons we call the solution ψ − an 'inflationary attractor', though strictly speaking it is not an attractor in the sense of the theory of dynamical system. The universe filled with the supercritical YM condensate at high density, ψ ≫ max(µ, µ 2/3 ), experiences inflation with H = µ and ψ ∝ e −µt/4 . As expansion is going on, the gauge condensate monotonously decays. Eventually, when ψ drops below ψ cr , the de Sitter stage ends and transition takes place to the universe filled with radiation.

Constant-roll regime
One can use numerics to explore behavior of solutions with different initial data demonstrating explicitly that their choice within the substantial region of the physical domain D phys (31) ensures qualitatively similar behavior. Namely, after some time the solution gets attracted to the de Sitter stage of finite duration. To show robustness of this process we choose zero initial value for the YM field ψ i = 0. Then one can easily see that two parameters defining D phys , The trajectories are temporarily attracted to the boundary H 2 = µ 2 . Then transition to radiation dominated universe, H = 1/2t, occurs.
satisfy the desired conditions 0 ≤ H 2 i < µ 2 /3, D i ≥ 1 for anyψ i . On the Fig. (1, 2) we present evolution of the YM function and the Hubble parameter starting from the initial states ψ i = 0,ψ i /µ 2 = 0..100. Such initial states correspond to domination of the kinetic energy, so the Hubble parameter rapidly decreases. Therefore, initially the system resides in the phase space region distant from the inflationary attractor and even moves further away from it. Nonetheless, soon after, the trajectories starting withψ i ≫ µ 2 become attracted to the boundary H = µ, signalling the onset of inflation. The finite duration of this stage is ensured by a subsequent exponential decay of the YM function along with the exponential growth of the scale factor. Once the field value eventually falls below the critical value ψ cr , the effect of the non-minimal coupling becomes negligible. Then the system undergoes transition to the radiation dominated stage H → 1/(2t), corresponding to an oscillating YM function.
Despite the fact that dynamical system (28,29) looks rather complicated, the physics behind it resembles much the very first model of chaotic inflation [43] with λϕ 4 potential, which is just an oscillator with the Hubble friction. In our case the self-interaction of non-Abelian gauge fields plays crucial role, generating the ψ 4 potential. The cosmology in the non-Abelian Horndeski theory can be viewed as an oscillator with the Hubble friction amplified by non-minimal coupling to gravity, provided the initial state belongs to physical domain, D phys . The precise initial conditions are totally irrelevant, only the energy of oscillations matters. If this energy allows the field to climb high enough above the non-minimal coupling scale, ψ ≫ µ, the system at maximal deviation will be inevitably attracted by the solution H = µ, and the amplified Hubble friction will impose the subsequent steady downhill motion. Therefore the domain of inflationary initial states can be presented as Comparing our model with scalar inflation with a monomial potential V (ϕ) ∝ ϕ n one observes the following important difference. In the scalar slow-roll case the Hubble parameter rapidly grows with increasing field value, H 2 ∼ V (ϕ). In the HYM model, the H 2 dependence is flattened near the boundary of the phase space: qualitatively H 2 ∼ f (tanh ψ/µ). This resembles the conformal attractors, in which case the effective potential is flattened near the phase space boundary, so that H 2 ∼ V (tanh ϕ/ √ 6). In the present model the potential ψ 4 is quite steep, while the Hubble friction is constant during the inflationary stage. This is why one observes the constant-roll motion [44], whenψ/Hψ = n = const. In our case n = −1/4, while the slow-roll motion corresponds to n ≃ 0 (and n = −3 is 'ultra slow-roll'). The peculiar feature of the constant-roll inflation is the absence of eternal inflation. Indeed, quantum fluctuations should be of the order δψ ≈ H/2π ≈ µ/6. But according to classical motion, for the time interval t = H −1 = µ −1 , the YM field decreases by the value ∆ψ ≈ ψ (1−e −1/4 ) ≈ ψ/5. During most of the inflation stage, one has ψ ≫ µ, therefore quantum fluctuations provide just a minor correction to the classical motion: ∆ψ ≫ δψ. Contrary to the standard chaotic inflation scenario [45], here the Hubble parameter (and therefore the amplitude of quantum fluctuations) remains constant, while the potential ψ 4 becomes steeper as the field value grows. Hence the field can not climb up due to quantum fluctuations. Only at the exit from the inflation stage one has ψ ∼ µ, and the quantum fluctuations may become comparable with the classical values. This can probably generate the large scale structure of the universe, but it can not lead to an eternal self-reproduction. It is a matter of opinion whether the eternal inflation is objectionable [46], or not. But anyway, our inflationary scenario seems to be free of it.
To summarize the results, on the Fig. (3) we schematically compare the chaotic HYM-flation emerging in the Horndeski vector-tensor theory with the ψ 4 self-coupling potential of the YM field, and the chaotic inflation in the minimally coupled scalar field theory with the potential λϕ 4 .

HYM-flation as pre-inflation
During the constant-roll motion ψ exponentially decreases with the rate −µ/4, while the scale factor grows with the rate µ. So the number of e-folds gained by the scale factor at the moment t e is just the number of e-folds lost by ψ with a factor of four: Interpreting HYM-flation as the observed inflation, we will be faced with the problem of perturbations, however. If one considers µ as the scale parameter of the observed inflation, so that µ ∼ 10 −6 M pl , then ψ e , being slightly larger than µ for the inflationary solution, must be of the order 10 −6 ÷ 10 −5 M pl . Naturally choosing the initial conditions at the Planckian scale, ψ i ∼ M pl , one has N e−folds ≈ 46 ÷ 55 . Though technically this can be regarded satisfactory as the minimally required number of e-folds, one can barely obtain the observed perturbation spectrum in such a model. An accurate analysis of perturbation spectra in YM cosmology is rather involved (see calculations in the context of the "gauge-flation" [18]), but qualitative estimates do not seem to be in favor of the HYM-flation as the model of the observed inflation. Indeed, the amplitude of fluctuations is inversely proportional toψ, which exponentially decreases with time. Hence, the power spectrum should be rather blue-tilted, what contradicts to the Planck's data. Adopting a view that the observed inflation should occur at the GUT scale and be described by the scalar field slow-rolling in the plateau potential, we can suggest the Planck scale HYM accelerated expansion as the pre-inflation. It could help to solve the issue of initial conditions for the observed inflation. Recall that this problem was raised long ago in the context of the so-called new inflation [8] scenario. That time it was alleviated with invention of the chaotic inflation [43] in which the de Sitter expansion starts immediately after Planck era. Nowadays, observational data give preference to the low scale inflation with the plateau potential, during which the energy density is nearly ρ ≃ V infl ∼ 10 −10 . Thus a pre-inflationary evolution of the universe from some initial state with Planck energy ρ ∼ 1 is claimed again.
Note that the problem of an initial excess of kinetic energy, K ∼ 1 ≫ V infl ∼ 10 −10 , is not a critical one. Kinetic energy density drops as fast as K ∝ a −6 , while the value of the inflaton field changes not significantly during the pre-inflation epoch [9]. So, after a period of the post-planckian expansion, the inflaton can likely be found on a plateau, where the potential energy dominates. Of course, the closed universe with the Planck-size volume filled with non-inflating matter will collapse long before the low scale inflation could start. But this objection can be evaded arguing that universe was born non-closed.
The problem of inhomogeneities is more thorny [10]. If the universe in the pre-inflationary epoch was dominated by the kinetic energy or radiation, the scale factor would grow as a ∝ t 1/3 or a ∝ t 1/2 . The observed slow-roll inflation with the plateau potential V infl ∼ 10 −10 have to start when the kinetic energy density decreases by nearly eleven or twelve orders of magnitude, which requires the scale factor to gain roughly two orders of magnitude. This would take quite a long time, t infl ∼ 10 4 ÷ 10 5 . But in the decelerating universe, a ∝ t n , n < 1, the cosmological horizon grows with time: If the universe was born at Planck time, t 1 = t P l = 1, with Planck energy density, H(t 1 ) = 1/ √ 3, the distance at which the initial inhomogeneities could spread during till the start of inflation t 2 = t infl then would be d hor (t P l , t infl ) ∼ 10 2 ÷ 10 3 .
So, either the Universe was created homogeneous in a domain containing millions or billions of the identical Planck volumes (implausible), or some special topology has to be assumed [47,48], or the preliminary chaotic-type inflation [49,50] must be introduced which started immediately after the Planck era and ended before the stage of observed inflation. HYM-flation as pre-inflation looks as the most economic way to solve the problem of initial conditions for the GUT-scale inflation. This does not require the second scalar field, giving job to gauge fields already present in the GUT or supergravity models. Such a scenario thus may be viewed as an extension of the Higgs inflation within the full gauge theory involving vector fields. Let us take the natural value of coupling parameter µ = 1 and assume that the universe was born with the Planck energy density, ρ i ≃ 1. According to (36), for vanishing initial field value, ψ i = 0, the energy density ρ i = 1 requiresψ 2 i ≫ 1. Conversely, withψ 2 i ≪ 1, one has ρ i → 0, which is unlikely for the quantum creation process. Thus the universe with large probability emerges in a state belonging to HYM-flation domain, D infl , for whichψ 2 i ≫ 1. Of course, this is a simplified picture which does not take into account the non-vanishing initial gauge field value, the contribution of other fields into the energy density and so on. But, anyway, in the chaotic-like approach with randomly distributed initial conditions there should be non-small probability to find the universe (or some part of the universe) in a state belonging to HYM-flation domain, D infl .
As was argued above, the issue of initial conditions for the low scale inflation can be resolved if the scale factor grows by the factor of 10 2 ÷ 10 3 during the pre-inflation. This provides a large enough homogeneous patch, the core of which remains homogeneous until the time when the low scale inflation can start. In the non-Abelian Horndeski theory, according to Eq. (38), it is enough to have ψ i = 4 ÷ 6 in order to get the robust inflation with the desired gain of the scale factor. Mention that the initial domination of the potential term is not required: large kinetic energy of oscillating motion in the quartic potential due to the YM self-action will be transformed into the large potential energy, and then inflation starts. This extended HYM-Higgs two-stage inflationary model will be considered in more details in a separate publication.

Outlook
Previous attempts to construct models of inflation using YM fields were based on some ad hoc assumptions about mechanism of the conformal symmetry breaking [18,17], and, as a consequence, an introduction into the lagrangian of the new terms whose theoretical origin remained obscure. Our present mechanism looks more natural, appealing to now very popular non-minimal gravity couplings. Moreover, even within the realm of the Horndeski-inspired models, the present one is especially attractive by it simplicity and uniqueness. It has an intrinsic de Sitter attractor which is manifest already at the level of the lagrangian, the property which can hardly be underestimated. Its second crucial feature is the presence of the quartic self-interaction term which was absent in the previous Abelian non-minimal models, in which case stable inflation could not be achieved [25,41]. We have demonstrated that in the HYM theory the corresponding homogeneous and isotropic cosmology has a robust inflationary stage starting from a large variety of initial data. During inflation, the YM field decays down to some limiting value where the exponential expansion stops, ensuring a natural graceful exit. This HYM-flationary solution shares the features of both the plateau inflation and the chaotic inflation with quartic potential: the Hubble parameter is nearly constant (as for the plateau potential), but the potential remains steep. This results in a specific constant-roll motion, instead of the slow-roll, and leads to absence of eternal inflation. We propose this model on a role of the Planck-scale pre-inflation preparing the initial conditions for the GUT-scale observed inflation which can be both incorporated into the Yang-Mills-Higgs model with Horndeski coupling of the YM field to gravity. Such model would combine advantages of the Higgs conformal inflation and the preliminary chaotic-like inflation in a very natural and economic way. Similar construction seems to be possible in the context of the supergravity models where the corresponding pair of the YM field and a scalar also can be found.