Realistic compactification in spatially flat vacuum cosmological models in cubic Lovelock gravity: High-dimensional case

We investigate possible regimes in spatially flat vacuum cosmological models in cubic Lovelock gravity. The spatial section is a product of three- and extra-dimensional isotropic subspaces. This is the second paper of the series and we consider D=5 and general D>=6 cases here. For each D case we found critical values for $\alpha$ (Gauss-Bonnet coupling) and $\beta$ (cubic Lovelock coupling) which separate different dynamical cases and study the dynamics in each region to find all regimes for all initial conditions and for arbitrary values of $\alpha$ and $\beta$. The results suggest that for D>=3 there are regimes with realistic compactification originating from `generalized Taub' solution. The endpoint of the compactification regimes is either anisotropic exponential solution (for $\alpha>0$, $\mu \equiv \beta/\alpha^2<\mu_1$ (including entire $\beta<0$)) or standard Kasner regime (for $\alpha>0$, $\mu>\mu_1$). For D>=8 there is additional regime which originates from high-energy (cubic Lovelock) Kasner regime and ends as anisotropic exponential solution. It exists in two domains: $\alpha>0$, $\beta<0$, $\mu \leqslant \mu_4$ and entire $\alpha>0$, $\beta>0$. Let us note that for D>=8 and $\alpha>0$, $\beta<0$, $\mu<\mu_4$ there are two realistic compactification regimes which exist at the same time and have two different anisotropic exponential solutions as a future asymptotes. For D>=8 and $\alpha>0$, $\beta>0$, $\mu<\mu_2$ there are two realistic compactification regimes but they lead to the same anisotropic exponential solution. This behavior is quite different from the Einstein-Gauss-Bonnet case. There are two more unexpected observations among the results -- all realistic compactification regimes exist only for $\alpha>0$ and there is no smooth transition from high-energy Kasner regime to low-energy one with realistic compactification.


I. INTRODUCTION
It is not widely known but the idea of extra dimensions is older then General Relativity (GR) itself. Indeed, the first known extra-dimensional model was introduced by Nordström in 1914 [1] -it unified Nordström's second gravity theory [2] with Maxwell's electromagnetism. Soon after that Einstein proposed GR [3], but it took years before it was accepted: during the solar eclipse Einstein equations and Maxwell's electromagnetism. But for such decomposition to exist, the extra dimensions should be "curled" or compactified into a circle and "cylindrical conditions" should be imposed. The work by Kaluza was followed by Klein who proposed [5,6] a nice quantum mechanical interpretation of this extra dimension and so the theory, called Kaluza-Klein after its founders, was finalized. It is interesting to note that their theory unified all known interactions at their time.
With flow of time, more interactions were discovered and it became clear that to unify all of them, more extra dimensions are needed. At present, one of the promising theories to unify all interactions is M/string theory.
One of the distinguishing features of M/string theories is the presence in the Lagrangian of the corrections which are quadratic in curvature. Scherk and Schwarz [7] demonstrated that R 2 and extra-dimensional scale factors were studied for the first time in [43], and exponentially expanding (3+1)-dimensional part and an exponentially shrinking extra-dimensional scale factor were described. Power-law solutions have been considered in [18,44] and more recently in [20,[45][46][47][48] so that by now there is more-or-less complete description of the solutions of this kind (see also [49] for comments regarding different physical branches of the power-law solutions). Solutions with exponential scale factors [50] have also been studied in detail, namely, the models with both variable [51] and constant [52] volume; the general scheme for finding anisotropic exponential solutions in EGB gravity was developed and generalized for general Lovelock gravity of any order and in any dimensions [53]. The stability of the exponential solutions was addressed in [54] (see also [55] for stability of general exponential solutions in EGB gravity), and it was demonstrated that only a handful of the solutions found and described in [53] could be called "stable", while the most of them are either unstable or have neutral/marginal stability.
In order to find all possible cosmological regimes in EGB gravity, one needs to go beyond an exponential or power-law Ansatz and keep the scale factor generic. We are especially interested in models that allow dynamical compactification, so that we consider the metric to be the product of a spatially three-dimensional and extra-dimensional parts. In that case the three-dimensional part is "our Universe" and we expect for this part to expand while the extra-dimensional part should be suppressed in size with respect to the three-dimensional one. In [56] we demonstrated the there exist the phenomenologically sensible regime when the curvature of the extra dimensions is negative and the EGB theory does not admit a maximally-symmetric solution. In this case both the threedimensional Hubble parameter and the extra-dimensional scale factor asymptotically tend to the constant values. In [57] we continued investigation of this case and performed a detailed analysis of the cosmological dynamics in this model with generic couplings. Recent analysis of this model [58] revealed that, with an additional constraint on couplings, Friedmann-type late-time behavior could be restored.
With the exponential and power-law solutions described in the mentioned above papers, another natural question arise -could these solutions describe realistic compactification (i.e., with proper both past and future asymptotes) or are they just solutions with no connection to the reality? To answer this question, we have considered the cosmological model in EGB gravity with the spatial part being the product of three-and extra dimensional parts with both subspaces being spatially flat. As both subspaces are spatially flat, we can rewrite the equations of motion in terms of Hubble parameters, so that they become first order differential equations and could be analytically analyzed to find all possible regimes, asymptotes, exponential and power-law solutions. For vacuum EGB model it was done in [59] and reanalyzed in [60]. The results suggest that the vacuum model has two physically viable regimes -first of them is the smooth transition from high-energy GB Kasner (or generalized Taub solution) to low-energy GR Kasner. This regime exists for α > 0 (Gauss-Bonnet coupling) at D = 1, 2 (the number of extra dimensions) and for α < 0 at D 2 (so that at D = 2 it appears for both signs of α). Another viable regime is the smooth transition from high-energy GB Kasner to anisotropic exponential solution with expanding three-dimensional section ("our Universe") and contracting extra dimensions; this regime occurs only for α > 0 and at D 2.
The same analysis but for EGB model with Λ-term was performed in [61,62] and reanalyzed in [60]. The results suggest that the only realistic regime is the transition from high-energy GB Kasner to anisotropic exponential solution, it requires D 2, see [60][61][62] for exact limits on (α, Λ) in each particular D. The low-energy GR Kasner is forbidden in the presence of the Λ-term so the corresponding transition does not occur.
In these studies we have made two important assumptions -we considered both subspaces being isotropic and spatially flat. But what will happens in we lift these conditions? Indeed, the spatial section being a product of two isotropic spatially-flat subspaces could hardly be called "natural", so that we considered the effects of anisotropy and spatial curvature in [63]. The initial anisotropy affects the results greatly -indeed, say, in vacuum (4+1)-dimensional EGB gravity with Bianchi-I-type metric (all the directions are independent) the only future asymptote is nonstandard singularity [46]. Our analysis [63] suggest that the transition from Gauss-Bonnet Kasner regime to anisotropic exponential expansion (with expanding three and contracting extra dimensions) is stable with respect to breaking the symmetry within both three-and extra-dimensional subspaces.
However, the details of the dynamics in D = 2 and D 3 are different -in the latter there exist anisotropic exponential solutions with "wrong" spatial splitting and all of them are accessible from generic initial conditions. For instance, in (6 + 1)-dimensional space-time there are anisotropic exponential solutions with [3 + 3] and [4 + 2] spatial splittings, and some of the initial conditions in the vicinity of E 3+3 actually end up in E 4+2 -the exponential solution with four and two isotropic subspaces. In other words, generic initial conditions could easily end up with "wrong" compactification, giving "wrong" number of expanding spatial dimensions (see [63] for details).
The effect of the spatial curvature on the cosmological dynamics also could be dramatic -say, positive curvature changes the inflationary asymptotic [64,65]. In the case of EGB gravity the influence of the spatial curvature reveal itself only if the curvature of the extra dimensions is negative and D 3 -in that case there exist "geometric frustration" regime, described in [56,57] and further investigated in [58]. Full investigation of the spatial curvature effects on the existing regimes could be found in [63].
The current manuscript is a direct continuation of [66], where the same analysis was performed for low-D (D = 3, 4) cases. It also could be called a spiritual successor of [59][60][61][62] -now we are performing the same analysis but for cubic Lovelock gravity. In this paper we consider only vacuum case, Λ-term case, as well as possible influence of anisotropy, spatial curvature and different kinds of matter source are to be considered in the papers to follow.  [59,60]. At last, we draw conclusions and formulate perspective directions for further investigations.

II. EQUATIONS OF MOTION
Lovelock gravity [16] has the following structure: its Lagrangian is constructed from terms where δ i 1 i 2 ...i 2n j 1 j 2 ...j 2n is the generalized Kronecker delta of the order 2n. One can verify that L n is Euler invariant in D < 2n spatial dimensions and so it does not give nontrivial contribution into the equations of motion. Then the Lagrangian density for any given D spatial dimensions is sum of all Lovelock invariants (1) upto n = D 2 which give nontrivial contributions into equations of motion: where g is the determinant of metric tensor, c n are coupling constants of the order of Planck length in D dimensions and we assume summation over all n in consideration. The metric ansatz has the form g µν = diag{−1, a 2 1 (t), a 2 2 (t), . . . , a 2 n (t)}.
As we mentioned, we are interested in the dynamics in cubic Lovelock gravity, so we consider n up to three (n = 0 is the boundary term while n = 1 is Einstein-Hilbert, n = 2 is Gauss-Bonnet and n = 3 is cubic Lovelock contributions). Substituting metric (3) into the Lagrangian and following the standard procedure gives us the equations of motion: as the ith dynamical equation. The first Lovelock term-the Einstein-Hilbert contribution-is in the first set of brackets, the second term-Gauss-Bonnet-is in the second set and the third -cubic Lovelock term-is in the third set; α is the coupling constant for the Gauss-Bonnet contribution while β is the coupling constant for cubic Lovelock; we put the corresponding constant for Einstein-Hilbert contribution to unity 1 . Since we consider spatially flat cosmological models, scale factors do not hold much physical sense and the equations are rewritten in terms of the Hubble parameters . Apart from the dynamical equations, we write down the constraint equation As we mentioned in the Introduction, we want to investigate the particular case with the scale factors split into two parts -separately three dimensions (three-dimensional isotropic subspace), which are supposed to represent our Universe, and the remaining represent the extra dimensions (D designs the number of additional dimensions) and the equations take the following form: the dynamical equation that corresponds to H, the dynamical equation that corresponds to h, and the constraint equation, Looking at (6)-(8) one can notice that the structure of the equations depends on the number of extra dimensions D (terms with (D − 4) multiplier nullifies in D = 4 and so on). In previous papers, dedicated to study cosmological dynamics in EGB gravity, we performed analysis in all dimensions, sensitive to EGB case [59][60][61][62]. In the cubic Lovelock, the structure of the equations of motion is different in D = 3, 4, 5 and in the general D 6 cases. In the directly previous paper [66] we performed the analysis for D = 3 and D = 4 cases, leaving D = 5 and the general D 6 cases for this paper. Also, similar to [66], since the these papers are dedicated to the vacuum case, we have Λ ≡ 0.

III. GENERAL SCHEME
The procedure of the analysis is exactly the same as described in our previous papers [59][60][61][62].
In the directly previous paper [66] this scheme is described in great detail for D = 3 case. In the current paper we just list the scheme and use it with almost no details. So the scheme is as follows: • we solve (8)  If we take the discriminant of (8) with respect to H, and then its discriminant with respect to h, we obtain critical values for (α, β) which separate qualitatively different cases; • we find analytically isotropic exponential solutions: we substituteḢ =ḣ ≡ 0 as well as h = H into (6)-(8); the system simplifies into a single equation, we solve it and find not only roots but also the ranges of (α, β) where they exist; • we find analytically anisotropic exponential solutions: we substitutė H =ḣ ≡ 0 into (6) • first three steps provides us with a set of critical values for (α, β) which separate domains with different dynamics; • we solve (6)-(7) with respect toḢ andḣ; • we substitute obtained H i intoḢ andḣ and obtain the latter as a single-variable functions: • the obtainedḢ(h) andḣ(h) expressions and graphs are analyzed for all possible domains in (α, β) space to obtain all possible regimes; • obtained exponential regimes are compared with exact isotropic and anisotropic solutions (see [53]) to find the nature of the exponential regimes in question; • power-law regimes are analyzed in terms of Kasner exponents (p i = −H 2 i /Ḣ i ) to verify that low-energy power-law regimes are standard Kasner regimes with p i = p 2 i = 1 or "generalized Taub" regimes (see below) while high-energy power-law regimes are Lovelock Kasner regimes with The above scheme allows us to completely describe all existing regimes for a given set of the parameters (α, β). In the previous paper [66] we followed the scheme in great detail for D = 3 case, with full description of all steps, and in the current papers we just briefly mentioned the details and mainly focus on the results.
Before proceeding with the particular cases, it is useful to introduce the notations we are going to use through the paper. We denote Kasner regime as K i where i is the total expansion rate in terms of the Kasner exponents p i = (2n − 1) where n is the corresponding order of the Lovelock contribution (see, e.g., [20]). So that for Einstein-Hilbert contribution n = 1 and p i = 1 (see [67]) and the corresponding regime is K 1 , which is usual low-energy regime in vacuum EGB case (see [59,60]) and we expect it to remain here. For Gauss-Bonnet n = 2 and so p i = 3 and the regime K 3 is typical high-energy regime for EGB case (again, see [59][60][61][62]). Finally, for cubic Lovelock n = 3 and so p i = 5 and the regime K 5 is typical high-energy regime for this case in low D case [66].
Another power-law regime is what is called "generalized Taub" (see [68] for the original solution).
We mistakenly taken it for K 3 in [59], but then in [60] corrected ourselves and explained the details (they both have p i = 3 which causes misinterpreting). It is a situation when for one of the subspaces the Kasner exponent p is equal to zero and for another -to unity. So we denote P 1,0 the case with p H = 1, p h = 0 and P 0,1 the case with p H = 0, p h = 1.
The exponential solutions are denoted as E with subindex indicating its details -E iso is isotropic exponential solution and E 3+D is anisotropic -with different Hubble parameters corresponding to three-and extra-dimensional subspaces. But in practice, in each particular case there are several different anisotropic exponential solutions, so that instead of using E 3+D we use E i where i counts the number of the exponential solution (E 1 , E 2 etc). In case if there are several isotropic exponential solutions, we count them with upper index: E 1 iso , E 2 iso etc. And last but not least regime is what we call "nonstandard singularity" and denote it is as nS.
It is the situation which arise in Lovelock gravity due to its nonlinear nature. Since the equations (6)-(7) are nonlinear with respect to the highest derivative (Ḣ andḣ in our case), when we solve them, the resulting expressions are ratios with polynomials in both numerator and denominator.
So there exist a situation when the denominator is equal to zero for finite values of H andor h.
This situation is singular, as the curvature invariants diverge, but it happening for finite values of H andor h. Tipler [69] call this kind of singularity as "weak" while Kitaura and Wheeler [70,71] -as "type II". Our previous research demonstrate that this kind of singularity is widely spread in EGB cosmology -in particular, in totally anisotropic (Bianchi-I-type) (4 + 1)-dimensional vacuum cosmological model it is the only future asymptote [46].
Following the described above procedure, first we find the critical values for µ for H(h). We find the discriminant of (11) with respect to H and then its discriminant with respect to h -it is 16th-order polynomial with respect to µ = β/α 2 and it has the following roots: double root The next step is the abundance of the isotropic exponential solutions; we substituteḢ =ḣ ≡ 0 and h ≡ H into (9)-(11) to obtain single equation which governs the isotropic solutions: The analysis of its nontrivial solutions suggests that we have one isotropic exponential solution iff β < 0 (for any sign of α) and two isotropic exponential solutions for α < 0, β > 0, µ 5/8.
The final step before the regime analysis is the abundance of the anisotropic exponential solutions. We substituteḢ =ḣ ≡ 0 into (9)-(11) and solve the resulting system with respect to H and h. The resulting equation on h is bi-9-power and its discriminant is 24th-order polynomial on µ with the following roots: two roots µ 3 ≈ 0.0517, µ 4 ≈ 0.1328 belong to a certain six-order polynomial and cannot be expressed through elementary functions; double pair of roots µ 5,6 = 55/18 ± 10 √ 7/9 ≈ 0.1158, 5.9953; single roots µ 7 = 1681/8242 ≈ 0.1995 and µ = 5/8 ≡ µ 2 ; plus there is a triple set of imaginary roots from a fourth-order polynomial. The further analysis suggests that for α < 0, β < 0 there is one anisotropic exponential solution, for α < 0, β > 0 there are two anisotropic exponential solutions if µ > 5/8. For α > 0, β < 0 there are also two anisotropic exponential solutions while for α > 0, β > 0 the situation is more complicated, similar to the previous cases. So for µ < µ 3 there are three, for µ = µ 3 four and for µ 4 > µ > µ 3 five.
With further growth of µ we have decrease of the number of solutions: four for µ = µ 4 , three for µ 7 > µ > µ 4 , two for µ = µ 7 and only one for µ > µ 7 . One can see that, similar to the previous cases, there is a fine structure of the anisotropic exponential solutions at α > 0, β > 0 and we describe it separately.
With all preliminaries done, it is time to describe all individual regimes. We solve (11) with respect to H and plot the resulting H(h) curves in Fig. 1. There red curve corresponds to H 1 , blue to H 2 and green to H 3 . The panels layout is as follows: α < 0, β < 0 on (a) panel, α < 0, β > 0, We analyze the individualḢ(h),ḣ(h), p H , p h curves and plot the derived regimes directly on H(h) evolution curves (see Fig. 1). The resulting regimes for α < 0, β < 0, presented in Fig. 1(a), are: P 1,0 ↔ E iso and P 0,1 ↔ E iso on blue hyperbola-like H 2 curve, K 5 ↔ K 1 on green H 3 curve and K 5 ↔ E 1 , nS ↔ E 1 , nS ↔ nS and K 1 ↔ nS on H 1 . One can see that neither of the regimes have compactification feature - The next case to consider is α < 0, β > 0, µ < 5/8, presented in Fig. 1(b). The regimes there are: P 0,1 ↔ E 2 iso and nS ↔ E 2 iso on hyperbola-like part of H 1 , nS ↔ K 1 on the remaining part of H 1 , P 1,0 ↔ E 1 iso and K 5 ↔ E 1 iso on hyperbola-like H 2 , and K 5 ↔ K 1 on H 3 . For the same reasons as in the previous case, K 1 cannot be called realistic compactifications, and since there are no other candidates, there are no compactifications in this case either.
Previous case naturally followed by α < 0, β > 0, µ > 5/8, presented in Fig. 1(c). On the boundary value, µ = 5/8, the isotropic exponential solutions from Fig. 1(b) coincide so that we have one pair instead of two, but the other regimes are the same so we skipped µ = 5/8 from separate consideration. The µ > 5/8 case has the regimes: K 5 ↔ K 1 on H 3 and nS ↔ K 1 like in previous case, E 1 ↔ nS and E 1 ↔ P 1,0 on one banana-like curve and E 2 ↔ P 0,1 and E 2 ↔ K 5 on another banana-like curve. For the same reasons as in the previous cases, K 1 cannot give us realistic compactification, and E 1 and E 2 are located either in first, so that having H > 0, h > 0, or in third (H < 0, h < 0) quadrants, so that also cannot serve as realistic compactification either.
For α < 0, similar to the previous D cases (see [66]), there are no realistic compactification regimes.
We proceed with α > 0, and the first case to consider is β < 0, presented in Fig. 1(d). There on a hyperbola-like part of H 2 branch we have P 0,1 ↔ E iso and nS ↔ E iso and P 1,0 ↔ E 2 and nS ↔ E 2 on edge-shaped part of H 3 branch. On the latter, P 1,0 → E 2 is viable compactification regime - The final case to consider is α > 0, β > 0. The limiting cases -µ < µ 3 and µ > µ 7 are presented in Fig. 1 -Fig. 1(e) for µ < µ 3 and Fig. 1(f) for µ > µ 7 ; the fine structure between them is presented in Fig. 2: µ < µ 3 (the same as in Fig. 1 Let us now describe the regimes which appear in these cases. At µ < µ 3 (see Figs. 1(e) and 2(a)) we have K 1 → K 5 along the H 1 branch (we focus on the second quadrant; the regimes in the fourth quadrant are time-reversal of the described), Of these regimes only P 1,0 → E 3 from the H 3 branch has realistic compactification. For the next several cases the regimes along H 1 and H 3 do not change (and so the realistic compactification P 1,0 → E 3 from the H 3 branch is presented in them as well), the difference is only within H 2 branch and we describe only these changes. This way, for µ = µ 3 (Fig. 2(b)) along H 2 we have (the outmost nS turned into anisotropic exponential solution), but no new viable compactifications appear; for µ 4 > µ > µ 3 (Fig. 2(c)) we for for µ = µ 4 (Fig. 2(d)  to the different branches (red -to H 1 , blue -to H 2 and green -to H 3 , in accordance with the designation in Fig. 1) (see the text for more details).
we have K 1 ← nS → nS ← E 4 → nS ← E 1 → K 5 ; for µ 1 > µ > µ 4 (Fig. 2(e)) we have One can see that neither of the regimes which appear on H 2 have realistic compactification. The next case is µ = µ 1 , presented in Fig. 2(f) and this is the situation when the physical branches "touch" each other and "reconnecting", making different "routes", so that many of the existing regimes terminating and other regimes are formed instead.
New regimes could be found on Fig. 2(g), where we plot µ 7 > µ > µ 1 case. One can note that the abundance of the exponential solutions and nonstandard singularities is exactly the same as for µ 1 µ > µ 4 (Figs. 2(e, f)) but the regimes are completely different, since H(h) forming completely different physical branches. So the new regimes for µ 7 > µ > µ 1 (Fig. 2(g)) are: K 5 ← E 1 → K 5 on the outmost branch, K 5 ← E 2 → E 3 ← nS → nS ← K 1 on the middle branch and P 1,0 → K 1 on the innermost branch. One can see that only the latter gives us realistic compactification -in all other cases K 5 is future asymptote (or past asymptote but for the regimes in the fourth quadrant, which have H < 0). Further increase of µ changing the regimes along middle branch, leaving the outmost and innermost unchanged and so keeping P 1,0 → K 1 as a viable compactification. The change of the regimes along the middle branch is as follows: Fig. 2(h)) and K 5 ← nS → nS ← K 1 for µ > µ 7 (see Fig. 2(i)); one can see that there are no realistic compactification regimes among those within H 2 .

V. GENERAL D 6 CASE
In this case we use the core (6)-(8) equations -for all D 6 the structure of the equations of motion is unchanged. Following the procedure, we find the discriminant of (8) with respect to H and then its discriminant with respect to h -it is 16th-order polynomial with respect to µ = β/α 2 and it has the following roots: and µ 2 < µ 3 -the solutions of certain six-order polynomial with the coefficients made up to D 16 .
The expressions for µ 2,3 cannot be obtained in terms of elementary functions in the general case but for each D they could be derived, at least numerically. Similar to the previous cases, µ 1 separates two different H(h) regimes in α < 0, β > 0 while µ 2,3 -for α > 0, β > 0, Unlike previous D = 3, 4 cases (see [66]), where there is only one separation µ in the α > 0, β > 0 case, now we have two.

Analysis of its nontrivial solution
suggests that there is one isotropic exponential solution for β < 0 (for both signs of α) and two for α < 0, β > 0, µ < µ 1 . Let us note that this is the same scheme we had for all previous cases, making it true for all D.
The situation with anisotropic exponential solutions is more complicated. Following usual procedure, we obtain equation for h as bi-nine-power polynomial with the discriminant being 24th-order polynomial in µ with coefficients made up to D 130 . Interesting enough, in the general case there are only two roots of this discriminant which affect the structure of the anisotropic exponential solutions: µ 1 -the same as in H(h) and isotropic exponential solutions -and So that we can say that in the general D 6 case there is no "fine structure" of the anisotropic exponential solutions -at least not in the sense we have it for D = 3, 4, 5. One can also note that We perform the same procedure as in the previous cases and obtain the resulting regimes on the H(h) curves; we present them in Fig. 3. As always, different colors correspond to different branches: H 1 is red, H 2 is blue and H 3 is green. The panels layout is the following: α < 0, β < 0 presented on (a) panel, α < 0, β > 0, µ < µ 1 on (b) panel, α < 0, β > 0, µ > µ 2 on (c) panel, α > 0, β < 0 on (d) panel, α > 0, β > 0, µ < µ 2 on (e) panel, and α > 0, β > 0, µ > µ 3 on (f) panel. Let us now describe the regimes.
The first case to consider, α < 0, β < 0, is presented in Fig. 3(a). There one can see K 5 ↔ E iso and P 1,0 ↔ E iso on hyperbola-like H 2 branch, K 5 ↔ K 1 on H 3 and a number of regimes along H 1 (according to the second quadrant; in fourth quadrant they are time-reversed): One can see that neither of the listed regimes have viable compactification: E 1 is either unstable or has H < 0, as for K 1 , it is either past attractor or has H < 0, or has nS as a past attractor.
The next case is α < 0, β > 0, µ < µ 1 and it is presented in Fig. 3(b). There we have K 5 ↔ E 1 iso and nS ↔ E 2 iso on hyperbola-like part of H 1 branch and nS ↔ K 1 on the remaining part, K 5 ↔ E 2 iso and P 1,0 ↔ E 2 iso on H 2 and K 5 ↔ K 1 on H 3 . Similar to the previous case (and for the similar reasons) there are no realistic compactification schemes in this case as well.
Similarly to previously considered D cases, at µ = µ 1 the isotropic exponential solutions E 1 iso and E 2 iso "touch" each other and coincide; the regimes remain the same so we skip this case from consideration.
Again, similar to the previous D cases, isotropic exponential solutions switched into anisotropic ones; hyperbola-like branches, which "touch" each other at µ = µ 1 , "detouch" and form new physical "banana-shaped" branches with the exponential solutions are located on them -one at a branch. So on one of them we have K 5 ↔ E 1 with two different K 5 , on another nS ↔ E 2 and P 1,0 ↔ E 2 ; there are also nS ↔ K 1 and K 1 ↔ K 5 on H 3 -the sane as in previous cases. One can see that E 1 and E 2 in the third quadrant are unstable while those in the first are stable but have H > 0 and h > 0 so they cannot describe realistic compactification 2 . The comments about K 1 from the previous cases are true here as well, so there are no viable compactification in this case. To conclude, similarly to the previous D cases, α < 0 domain does not have realistic compactification regimes. Now let us turn to α > 0. The first of these cases is β < 0, presented in Fig. 3(d)  on the edge-shaped H 3 part; among them P 1,0 → E 2 is a realistic compactification regime. Other regimes include K 5 → K 1 on one and K 5 ← E 1 → nS ← nS → K 1 on another physical branch, listed according to the second quadrant. Let us note that the listed regimes along the second branch are exactly the same as in the described above α < 0, β < 0 case (see Fig. 3(a)). So that in this case we have one viable compactification regime -P 1,0 → E 2 .
What remains is the description of the cases for α > 0, β > 0. The first of the α > 0, β > 0 cases, µ < µ 2 , is presented in Fig. 3(e), the following cases -in Fig. 4 and the last case, µ > µ 3 , in Fig. 3(f). In this regard, we called the subcases in Fig. 4 as a "fine-structure" of the α > 0, β > 0 case, but it is not in the same sense as previous D cases. In the previous D cases we have fine-structure with respect to the rapid change in the number of exponential solutions -within a small µ region, the exponential solutions appear, disappear, merge, turn in nS -they have rich and rapid-changing dynamics. In this case we does not have all of this, but we have change of the physical branches. So the manner is the same and so we keep the name "fine-structure". The subcases within it are presented in Fig. 4 and the panel layout is as follows: µ = µ 2 is on (a) panel, Let us start with µ < µ 2 case, presented in Fig. 3(e). The regimes for this case are (we list them according to the second quadrant): Of all the regimes, only P 1,0 → E 3 have realistic compactification. The next case is µ = µ 2 , presented in Fig. 4(a). Here one can immediately see the difference between D = 4, 5 and the general D 6 cases -one can see that there are two touching points for the branches and in D = 4, 5 they touch and decouple at the same µ (see e.g. Fig. 2(f)). But in the general D 6 they touch and decouple at two different µ -one at a time -and this is the reason for "fine-structured" consideration of the α > 0, β > 0 case. The regimes in this case are the same as for µ < µ 2 . Once µ is increased to µ 4 > µ > µ 2 , the regimes changes and are presented in Fig. 4(b).
The regimes are: only P 1,0 → K 1 has realistic compactification. The next case is µ = µ 4 , presented in Fig. 4(c), and at this point the number of the exponential solution changes. So that instead of a pair of anisotropic exponential solutions E 2 and E 3 on H 3 , we have singe meta-stable ( H = 0 -constant-volume solution, see [52]) E 2 . The P 1,0 → K 1 and K 1 → K 5 regimes remain unchanged and the regimes list with anisotropic exponential solution is:  µ 3 > µ > µ 4 , presented in Fig. 4(d), even E 2 disappears, leaving K 5 ← nS → nS ← E 1 → K 5 in addition to P 1,0 → K 1 and K 1 → K 5 to the list of regimes; of them only P 1,0 → K 1 has realistic compactification. With further increase to µ = µ 3 , the second "touch" happening -see Fig. 4(e).
And finally the situation for µ > µ 3 is presented in Fig. 3(f). The regimes there are: P 1,0 → K 1 , One can see that of them only P 1,0 → K 1 has viable compactification.
This finalize our study of D = 6 vacuum case. We see that the dynamics and its features are different from the previous cases, but the regimes with realistic compactification and their abundance is the same: all of the regimes require α > 0 and there are two distinct regimes -P 1,0 → E 3+6 for µ < µ 2 (including β < 0) and P 1,0 → K 1 for µ > µ 2 .
The remaining cases for D = 7 are shown in Fig. 5. As always, different colors correspond to three different branches: H 1 is red, H 2 is blue and H 3 is green, and the panels layout is the following: Fig. 5(a), µ = µ 3 in Fig. 5(b), µ 4 > µ > µ 3 in Fig. 5(c) and µ = µ 4 in Fig. 5(d); for µ > µ 4 the situation and the regimes are the same as in µ > µ 3 for D = 6 and so are presented in Fig. 3(f).
Let us have a closer look on new regimes in this case. The first case, µ 3 > µ > µ 2 , is presented in Fig. 5(a) and is exactly the same as µ 4 > µ > µ 2 from D = 6, presented in Fig. 4(b) (one can compare and verify that the regimes are the same). Then, the regime with realistic compactification is also the same -it is P 1,0 → K 1 . The next case is µ = µ 3 , presented in Fig. 5(b) -the second "detouch" -and the regimes are the same. With further increase of µ to the µ 4 > µ > µ 3 range, the situation is presented in Fig. 5(c) and the regimes are: P 1,0 → K 1 , K 5 ← E 1 → K 5 (with two different K 5 ) and K 5 ← E 2 → E 3 ← nS → nS ← K 1 . One can see that the regimes along the last mentioned branch are the same as they are in the appropriate µ range for D = 6. The realistic compactification in this case is only P 1,0 → K 1 . The next case is µ = µ 4 , presented in Fig. 5(d), and the changes are the same as they are in D = 6 -a pair of anisotropic exponential solutions E 2 and E 3 collapsed to a single E 2 with H = 0; the regimes are K 5 ← E 2 ← nS → nS ← K 1 and the only realistic compactification regime is P 1,0 → K 1 . Finally, for µ > µ 4 E 2 also disappears and the situation and the regimes are identical to those in µ > µ 3 D = 6 case, presented in Fig. 3(f).
In Fig. 6(a) we have presented α > 0, β < 0, µ < µ 4 case. One can see that the behavior along hyperbola-like H 2 , edge-shaped part of H 3 and H 1 − H 3 physical branch is the same as in D = 6 case (compare with Fig. 3(d)), the changes are introduced to H 1 − H 2 physical branch. One of the changes regards the number of the anisotropic exponential solutions and another -the time direction of the evolution. One can see that ever since this branch (or its predecessor) appear in D = 4 (see [66]) and in following cases D = 5 (see Fig. 1(d)) and D = 6, 7 (see Fig. 3(d)), K 5 always was future attractor -the regime was E 1 → K 5 . But in D 8 it is reversed -K 5 now is past attractor and can serve as a cosmological singularity, making K 5 → E 1 successful compactification.
Apart from it, another realistic compactification regime is P 1,0 → E 2 from edge-shaped part of solution; the K 5 → E 1 regime remains, as well as P 1,0 → E 2 . With further growth of µ, to µ > µ 4 , the situation is presented in Fig. 6(c). One can see that E 1 disappeared and so we have K 5 → nS instead and only P 1,0 → E 1 remains as a realistic compactification regime. So that, for D 8 at α > 0, β < 0, µ µ 4 we have new realistic compactification regime - The change of the time direction of one of the branches, described just above, affect α > 0, β > 0 as well -K 5 on H 3 branch, which is future asymptote in D = 4 (see [66]), D = 5 (see Figs. 1(e, f)), D = 6 (see Figs. 3(e, f)) and D = 7 (see Figs. 5), now in D 8 it is a past asymptote, in analogue with α > 0, β < 0 case described above. So that, similar to the previous case, we have K 5 → E 2 realistic compactification, and it is present for all µ > 0 (see Figs. 6(d)-(f)). Another realistic compactification regime is P 1,0 → E 2 for µ < µ 2 (see Fig. 6(d)) and P 1,0 → K 1 for µ > µ 2 (see Figs. 6(e), (f)). The details and other regimes (without compactifications) are quite similar to the analogues in previous case so we skip their consideration.

VI. DISCUSSIONS
In the current paper we have analyzed the cosmological dynamics of the cubic Lovelock gravity, with Einstein-Hilbert and Gauss-Bonnet terms present as well. We have chosen a setup with a topology being a product of two isotropic subspaces -three-dimensional, representing our Universe, and D-dimensional, representing extra dimensions. Both subspaces are flat, which simplifies our equations of motion and makes it possible to analyze them analytically. Current paper is a direct continuation of [66] where we considered low-D (D = 3, 4) cases and in current paper we extend the analysis on the remaining high-D cases. In a sense, it is also a logical continuation of [59][60][61][62], where we considered the same problem but in EGB gravity -vacuum case in [59,60] and Λ-term case in [60][61][62]. In [60] we reviewed all the results for EGB from [59,61,62] and changed the visualization of the regimes -in the original papers [59,61,62] we use tables to list of all the regimes, and this way sometimes is not easy to read. On the contrary, in [60] we put all the regimes on H(h) curves and added arrows to demonstrate t → ∞ directed evolution. In the current paper we decided to keep visualization from [60].
First of all, let us summarize the results, as they are scattered over mini-conclusions in each particular sections here and in [66]. The fist case is D = 3, presented in [66] and it has interesting feature -since the equations of motion are cubic in both H and h, there could be up to three branches of the solutions. On the other hand, it is the lowest possible dimension for cubic Lovelock gravity, so there are no Kasner solutions (see [49]). Then the only possibility is what we call "generalized Taub" solution -the situation when the expansion in each direction is characterized by Kasner exponent p i = −H 2 i /Ḣ i equal to either 1 or 0; so that for our topology it is either P 1,0 (p H = 1, p h = 0 -expansion of the three-dimensional subspace and "static" extra dimensions) or P 0,1 (p H = 0, p h = 1 -expansion of the extra-dimensional subspace and "static" three dimensions).
The D = 4 case has one cubic Lovelock Kasner solution K 5 but it is still not enough for all branches, so we still nonstandard singularities at (α < 0, β > 0) and (α > 0, β < 0) while for (α < 0, β < 0) and (α > 0, β > 0) the evolution curves have complicated shapes. In D = 4 we still have P 1,0 regime, but not P 0,1 , and some of the nonstandard singularities have power-law behavior and so designated as nS/P . Unlike D = 3, where the regimes within the fine structure existed on an isolated H(h) curve, in D = 4 they are located on one of the physical branches connected with K 1 (see [66] for details). The realistic compactification regimes are P 1,0 → E 3+4 for α > 0, µ < µ 1 (including entire β < 0) and P 1,0 → K 1 , for α > 0, µ > µ 1 -exactly the same as in D = 3 case, and again both of the regimes exist only for α > 0.
In D = 5 (see Figs. 1 and 2) we have two different K 5 and a "return" of P 0,1 , and the general evolution curves and the regimes resemble previous cases; the same is true for the realistic compactification regimes -the same P 1,0 → E 3+5 for α > 0, µ < µ 1 (including entire β < 0) and P 1,0 → K 1 , for α > 0, µ > µ 1 . So that, we can conclude that in D = 3, 4, 5 the dynamics, with the exception of some features, is generally the same, so as the realistic compactification regimes.
Formally the equations of motion keep the functional form in all D 6 cases (meaning that no new terms appear, how it was in the lower dimensions), but our analysis suggests that the details of the dynamics -similar to the differences between D = 3, 4, 5 -are different in D = 6, 7 and D 8 cases, so that we considered them separately. With a minor details, the dynamics of D = 6 (see Fig. 3) and D = 7 (see Fig. 5) do not differ much from the previous cases, and the realistic compactification regimes are the same. The only difference is absence of the "finestructure" of anisotropic exponential solutions, but instead we received the "fine-structure" of the H(h) curves. On the other hand, general D 8 case brought us a whole new regime with realistic compactification -K 5 → E 3+D -a regime which originates from "normal" Kasner regime, instead of P 1,0 . This regime exist in two domains: α > 0, β < 0, µ µ 4 and entire α > 0, β > 0. The regimes P 1,0 → E 1 and P 1,0 → K 1 are also present in the general D 8 case.
The above-mentioned "generalized Taub" solution deserves additional comments. Formally it fits the description of the "generalized Milne" solution -the second branch of the power-law solutions in Lovelock gravity (see [20] for details), but only formally -it fits only because it is degenerative. As we demonstrated in [49], strict "generalized Milne" cannot exist in pure highestorder Lovelock gravity, as it leads to degeneracy in the equations of motion. But if additional (lower-order) Lovelock contributions are involved, it this branch of solutions could be restored, but it was never demonstrated before. So that on the particular example of [3 + D] spatial splitting we demonstrated this possibility. Still, a little is known about this regime and it deserves additional investigation in the separate papers.
When we consider this "generalized Taub" solution as a past asymptote -and this is the case for all possible realistic compactification models in D = 3 ÷ 7 -it feels unnatural. Indeed -the P 1,0 regime imply H → ∞ and h → 0 as t → 0 (by "0" we mean here initial cosmological singularity), so that we initially have "burst"-like expansion of three-dimensional subspace while the extradimensional subspace is almost static. In addition to the feeling of unnaturalness, it is a question if this regime could be reached from totally anisotropic space, in a manner it was done in [63] for Kasner regimes in EGB case. So that it gives additional reason to deeply investigate this regime and we are going to do it in the nearest time.
Similar to the results of [66], the results of this paper suggest that the variety and abundance of the regimes is closer to Λ-term EGB, rather then to the vacuum EGB models. The reasons for that are not exactly clear, but we suspect that number of the free parameters plays a role here.
Indeed, for vacuum EGB model there is only one parameter -α, Gauss-Bonnet coupling, while for Λ-term EGB and vacuum cubic Lovelock there are two -α and Λ for the former and α and β (cubic Lovelock coupling) for the latter. So that we expect that the dynamics of the Λ-term cubic Lovelock gravity to be even more interesting and we are going to consider this case shortly.

VII. CONCLUSIONS
This concludes our study of the cosmological models in vacuum cubic Lovelock gravity. We have found that in all D 3 there are compactification regimes of two kinds, the first of them originate from "generalized Taub" solution; for the future asymptote we have either Kasner regime or anisotropic exponential solution. In D 8 there appears another compactification scheme which originates from high-energy Kasner regime and has anisotropic exponential solution as future asymptote. So that for D 8 and some parameter the two of them coexist on different branchesthe situation we never had in EGB gravity.
In addition to the regimes with successful compactification, we described and plotted on H(h) curves all possible transitions for all initial conditions and all structurally different cases. The variety and abundance of the regimes exceed even Λ-term EGB case, featuring transition between two anisotropic exponential solutions and transition between two different "generalized Taub" solutions.
There are two interesting observations which require additional investigation, as both are quite unexpected. First of them is that all of the realistic compactification regimes have α > 0 requirement. This is unexpected, as in both vacuum and Λ-term EGB cases we have viable compactifications for both signs of α. We can note that for the Λ-term case the joint analysis of our cosmological bounds and those coming from AdS/CFT and other considerations allows us to conclude α > 0 (see [60,61]), but for that we involved external (to our analysis) results. On the contrary, in the current case without any external bounds we already have realistic compactification only for α > 0.
The second observation is that there is no K 5 → K 1 transition with realistic compactification.
In EGB vacuum case [59,60] we have the transitions of this kind, so we expected that in the higher-order Lovelock gravity they would also appear, but our investigation reveals that they do not. There is K 5 → K 1 transition, but with contracting three and expanding extra dimensions, so it