A 3D Phase Space Analysis of Scalar Field Potentials

In this study, we present the phase-space analysis of Quintessence models specified by the choice of two potentials, namely the Recliner potential and what we call the broken exponential-law potential, which is a new proposal. Using a dynamical system analysis we provide a systematic study of the cosmological evolution of the two models and their properties. We find new scaling solutions characterised by a constant ratio between the energy density of the scalar field and that of the matter component. These solutions are of high interest in light of the possibility to alleviate the coincidence problem. Additionally, the models also show attractor solutions. We finally construct concrete models built using a double potential according to which one potential realises the early-time scaling regime and the second one allows to exit this regime and to enter in the epoch of cosmic acceleration driven by a scalar-field dominated attractor point.


Introduction
Since the discovery of the late time cosmic acceleration from Supernovae Ia, a lot of proposals to model this phenomenon have been put forward and tested with the large amount of data we have collected so far. The most successful one [1,2], based on General Relativity (GR) with the further assumption of an homogeneous and isotropic geometry, is the standard cosmological model (the ΛCDM model), where the late time acceleration is sourced by the cosmological constant Λ. The latter is characterised by an equation of state w Λ ≡ P Λ /ρ Λ = −1. However, the disagreement between the observed value of Λ and its theoretical one, suggested by quantum field theory, has cast doubts on the viability of such scenario. This is known as the cosmological constant problem [3][4][5]. Alternatively, the source of the late-time acceleration can be modelled with a dark energy fluid characterised by a time-dependent equation of state w DE (t). This framework is useful as null test of the ΛCDM model and enables to identify the tendencies of data with few free parameters.
One easy way to realise a dark energy scenario is to consider a dynamical scalar field φ, and the simplest Lagrangian description is given by Quintessence [6][7][8][9][10][11][12][13]. Quintessence is free from theoretical problems such as ghosts and Laplacian instabilities and the late-time acceleration is realised by a slowly varying field along a potential V(φ). This model, besides evolving differently at the background level due to a time-dependent equation of state, also predicts a different growth of perturbations compared to ΛCDM and its phenomenology depends on the form of the potential. Cosmological data such as galaxy clustering, redshift space distortion, supernovae, baryonic acoustic oscillation, and cosmic microwave background radiation data can be very efficient in constraining the form of the potential [14][15][16][17][18]. arXiv:2204.06420v1 [gr-qc] 13 Apr 2022 In this paper we perform a phase-space analysis [10,[19][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34] of the background cosmology of Quintessence. We identify a set of dimensionless variables that allows us to write the background equations as an autonomous system of first-order differential equations. To make the system autonomous, one has to make a choice on the form of the potential. Phasespace analysis in the context of Quintessence has been widely performed, considering many functional forms for the scalar-field potential [35][36][37][38][39][40][41][42][43][44][45][46][47][48]. The simplest option is to consider an exponential potential, as it reduces the dimensionality of the system to be two and its analysis is rather simple [19,[49][50][51]. In this work, we are interested in exploring different forms of the potential. In particular, we will explore the case in which the system is three-dimensional, given the feature of the tracker parameter (∝ d 2 V/dφ 2 ) to be a function of the roll parameter (∝ dV/dφ). These parameters depend, respectively, on the second and first derivative of the potential with respect to the scalar field. We then apply our method to some specific potentials in order to elucidate it.
The paper is organised as follows. In Section 2 we provide an overview of Quintessence equations on a homogeneous, isotropic, and flat space-time. We also present two different forms for the potential we aim at studying in detail. In Section 3 we present the general framework of the phase-space analysis, which will be used to study the evolution of the system associated with the Quintessence model. We then specialise the discussion on two potentials. In Section 4, we show the evolution of the scalar field density (ρ φ ) for the two models characterised by an early-time scaling behaviour (ρ φ proportional to the matter density). We provide the conclusions in Section 5.

Quintessence Models
Quintessence models are described by the following action where g is the determinant of the metric g µν , R is the Ricci scalar, and L Q is the Quintessence Lagrangian defined as where φ is the scalar field, V(φ) its potential, and ∇ µ is the covariant derivative. Finally, L m is the Lagrangian of matter fields, χ i . If we consider the flat Friedmann-Lemaître-Robertson-Walker (FLRW) metric, given by: where a(t) is the scale factor, t the cosmic time, x i the spatial coordinates and δ ij the Kronecker delta, the Friedmann and Klein-Gordon equations read, respectively, In the equations above, M 2 pl = 1/(8πG N ) is the reduced Planck mass squared in natural units whereh = c = 1 and with G N being Newton's constant, dots stand for the derivative(s) with respect to cosmic time, V φ ≡ dV dφ , H ≡˙a a is the Hubble parameter, ρ m represents the sum of the matter densities associated to cold dark matter and baryons. Similarly, we can associate an energy density ρ φ and a pressure P φ to the scalar field, which read, respectively, From them we can define the equation of the state of the scalar field Note that for a very shallow potential, whereφ 0, w φ −1, it mimics the case of the cosmological constant.
In these models, the potential V(φ) has to be specified in order to study the dynamics of the scalar field. In this respect, one can follow two different approaches: the first and most common approach specifies the form of the potential and then solves the Klein-Gordon equation to obtain the evolution of the scalar field and finally of H [7] (see also [52] for an explanation of how to construct the initial conditions); the second approach, instead, assumes a functional form for w φ and then using the Friedmann equations one can reconstruct V(φ) [45,53,54]. The approaches are, obviously, equivalent.
A fundamental theory which, from the first principles, allows us to specify either the functional form of the potential or of the equation-of-state parameter does not exist. All the functional forms proposed for w φ and V(φ) in the literature are, albeit motivated, simply phenomenological. Here we follow the first approach and we consider two forms of the potential. We will then analyse the resulting physical system using the technique of the phase-space dynamics.
In the present analysis we will focus on the potentials defined as follows: • The Recliner potential [55]: where V 0 > 0 and α > 0 are constant. The Recliner model is characterised by tracker solutions [55] and as such it belongs to a class of models known as α-attractors [56,57]. They have been proposed to study inflation and have the interesting feature of linking inflation and the present dark-energy dominated era. • The broken exponential-law potential. We propose the following potential: where p, q, and n are positive constants and A arbitrary, and χ = exp (φ/M pl ). It represents a simple approximation of α-attractor models and is inspired by a general functional form of the potential defining this class of models [58,59]. To the best of our knowledge, this is the first time such a potential is considered. As we will see in Section 4, this model gives results very similar to those of the Recliner model. This can somehow be expected as this potential can recover many sub-classes of the α-attractor models, including the Recliner. Therefore, the advantage of this potential is to allow a general and systematic investigation, albeit only formally correct when φ/M pl 1, of several models at the same time.
In Figures 1 and 2 we show the functional form of the Recliner and the broken exponentiallaw potentials for different sets of the free parameters characterising each model. Different line styles and colours correspond to different parameter sets. For simplicity, we only show positive values of the scalar field φ, but we will also comment on the behaviour of the potentials when φ < 0.  In Figure 1 we show the Recliner potential in units of V 0 . We can individuate three different regimes for positive values of the scalar field: for values φ/M pl 1 (which are characteristic at early times) the potential shows a plateau due to the fact that the exponential function is constant and approximately one. In this limit we have V(φ)/V 0 2. Larger values of φ/M pl are reached at later times and assuming α 1, the contribution of the exponential function is negligible and the Recliner potential flattens again, leading to an equation of state w φ ≈ −1, as shown in [55]. At late times, therefore, V ≈ V 0 . In the intermediate regime, the potential undergoes a transition from V = 2V 0 to V = V 0 . For different values of the parameter α, the potentials are clearly distinguishable during the tracker regime, whilst behaving in the same way at both early and late times. This shows the attractor nature of this potential: different initial conditions (values of the parameter α) will lead to the same final result, regardless of their exact values. For negative values of the scalar field, the evolution of the potential has some peculiar features. For |φ/M pl | 1, the plateau at V 2V 0 still exists and for smaller values of the scalar field (higher absolute values) the potential diverges as the argument of the exponential function is positive.
In Figure 2 we show the evolution of the broken exponential-law potential with respect to the scalar field φ. As we said before, this is inspired by a very general form of the α-attractor potentials and we can expect a wide variety of behaviours when choosing different values for the four free parameters. For small, positive values of the scalar field, the potential is flat and later it decreases. This happens because the exponential function is of the order unity when the scalar field is small and does not contribute to the evolution of the potential. When the exponential dominates, the potential diminishes. Note that we do not reach any plateau, so in general we do not expect the late-time dynamics to be similar to the cosmological constant. It is rather straightforward to establish limiting behaviours for this potential. When φ/M pl 1, we , with χ defined above. This means that when p − 2nq > 0 (p − 2nq < 0) the potential grows (decreases). When we consider negative values of the scalar field, the situation is more interesting. For |φ/M pl | 1, our figure would just extend to the left as the potential remains flat. For intermediate values the potential grows, reaches a maximum (its exact position and amplitude depend on the chosen parameters), and then declines to zero as χ p .

Phase-Space Analysis for Quintessence Models
The phase-space analysis is often used in cosmology to study the evolution of a dynamical system around critical points, p c,i . These points satisfy the equilibrium condition dp/dN| c,i = 0. In this section, we derive the equations necessary to the study of the dynamical system, and as such we consider a system made up of matter (either baryons and cold dark matter or radiation) and the scalar field. The first step is to rewrite the equations for the background introduced in Section 2 into an autonomous system of first-order differential equations. To this purpose, we make use of the following dimensionless variables: From a physical point of view, x represents the evolution of the kinetic term, y of the scalar field potential, and σ describes the matter density. The quantities λ and Γ are the roll and tracker parameters, respectively [60], which determine how fast the scalar field rolls down its potential (in other words, what is the shape of the potential) and whether the tracking behaviour exists. Using these variables we can derive the expressions for the scalar field density parameter Ω φ , its equation of state w φ and the total (also called effective) equation of state of the system, w eff , which read, respectively: In the previous set of equations, 0 ≤ w m < 1 is the equation of state of the matter component. As we are considering only a single matter fluid, w m can be either radiation (w m = 1/3) or dark matter plus baryons (w m = 0). However the system can be easily extended to two or more matter components. Moreover, 0 ≤ Ω φ ≤ 1 and both x and y are finite. While in the two-dimensional analysis they both lie inside a semicircle of radius r = 1, in three dimensions the situation is more complex, as we will show for each potential separately. Using Friedmann's equations for Quintessence, we finally write the dynamical system as follows where N = ln a. The system (12) is non-linear and non-autonomous. For our purpose to make it autonomous, we can proceed following two options: the first is to assume Γ = 1, and in this case it immediately follows that the form of the potential V is fixed and it is found to be an exponential potential [10]. Additionally, λ is a constant and this makes the system bi-dimensional. The second assumes Γ − 1 to be a function of λ. In this case the system becomes three-dimensional and gives rise to an interesting phenomenology associated to different forms of the potential V. The latter is the procedure we will consider in this work. In the following we will show that the choice of the potentials in Section 2 will indeed satisfy the condition Γ(λ) and we will proceed by identifying the relation between Γ and λ. Then, to investigate the dynamics of the system (12), we will determine the evolution around the critical points, p c,i . These points, as said before, are characterised by the equilibrium condition dp/dN| c,i = 0. Therefore, we will solve the following system of algebraic equations dx dN = dy dN = dλ dN = 0. The values p c,i ∈ {x c , y c , λ c } satisfying this system will be the critical points we are looking for. Let us notice that we will consider only the points with y c > 0 as all the functions appearing in its definition are positively defined. However, mathematically the system is symmetric with respect to y, and also the points with y c < 0 will satisfy the equilibrium condition. The final step consists in evaluating the stability of the system around the critical points. To this extent, we compute the eigenvalues µ i of the Jacobian matrix M. Critical points can be classified in hyperbolic and non-hyperbolic points according to whether all eigenvalues of the Jacobian matrix have, respectively, Re(µ i ) = 0 and in the second case at least one eigenvalue for which Re(µ i ) = 0. The dynamical system linearised around hyperbolic critical points is a good approximation of the full nonlinear system and small perturbations do not change qualitatively the phase portrait near the equilibrium point. For a 3D system one has three eigenvalues for each critical point and the stability depends on the sign of the eigenvalues: if all the eigenvalues are real and negative, the critical point is stable/attractor; if they are all real and positive the point is unstable; if the eigenvalues are real and at least one is positive while the others negative, the point is a saddle; if at least one eigenvalue is real and there are pairs of complex eigenvalues, the stability depends on the sign of their real part: negative real parts will give a stable point, positive real parts an unstable one, and if one is positive and the others are negative the critical point is a saddle point. These are the points we are more interested in as the analysis around them is robust. However, we will face some cases in which one of the eigenvalues is null. In this case, linear stability is not sufficient to infer the stability of the critical points and one has to resort either to a numerical investigation or to other mathematical methods, such as the Lyapunov method or the centre manifold theory [49,61]. For an application of the latter method to dynamical systems in cosmology, we refer the reader to [42,62,63].
Finally, note that when discussing the phase-space analysis of Quintessence models, we do not need to take into account stability conditions coming from linear perturbations, as they propagate luminally (c 2 s = 1). In addition, in determining the dynamics of the system, the exact value of V 0 is irrelevant, as it simplifies in the calculations.
Let us note that when discussing the critical points, we will see that for x c = y c = 0, the equilibrium condition is satisfied by an arbitrary value of λ, which we denote as λ a in the following tables.
Before concluding the section, we introduce a definition which will be useful in the analysis to simplify the notation. We define the function f (λ) as follows For each critical point, we also discuss its existence, whether it leads to acceleration or not and the corresponding value for the scalar field energy density Ω φ , its equation of state w φ and the effective equation of state of the system w eff . The condition for the scalar field to produce accelerated expansion is w eff < −1/3.

The Recliner Potential
Let us consider the Quintessence model for which the potential is chosen of the form of the Recliner potential in Equation (8). It is easy to show that for this model the λ and Γ parameters can be written as follows: Therefore, we find that Γ(λ) and the system (12) is indeed autonomous. The corresponding solutions for the equation dλ dN = 0, are Hence, 0 ≤ λ a ≤ λ c,2 . We include in Table 1 the list of the critical points of the system (12) and in Table 2 the analysis of the stability. Note that we implicitly assumed that −∞ < φ/M pl < +∞. This range for the scalar field in fact translates into λ c,1 ≤ λ ≤ λ c,2 . Table 1. Critical points for the dynamical system with the Recliner model. We label with λ a an arbitrary point satisfying the equilibrium condition and together with the coordinates of the critical point, we present its existence, whether it leads to an accelerated expansion phase, the scalar field energy density parameter Ω φ , the equation-of-state parameter for the scalar field and the total one, w φ and w eff , respectively.
From the definition of the variables x and y (x 2 + y 2 = Ω φ ) we easily find that 0 ≤ x 2 + y 2 ≤ 1 (as we will shortly see, several critical points lie on the boundary) and recalling that y ≥ 0, the dynamical system describes the half-cylinder with boundaries The usual phase-space for x and y is, thus, augmented by the boundaries λ = λ c,1 and λ = λ c,2 . Although the orbits in these boundaries will not be physical, we can determine, by continuity, the behaviour of the orbits lying inside the half-cylinder. Note that, by construction, therefore, the critical points will lie in the boundaries and trajectories connecting the start and end on these boundaries. Table 2. Stability analysis for the dynamical system with the Recliner model. We show the values of the eigenvalues µ i , with i ∈ (1, 2, 3) and the condition(s) for stability.

Point
(µ 1 , µ 2 , µ 3 ) Stability We now discuss the properties of the critical points in detail: • P 1 : matter point. This non-isolated critical point exists for any allowed value of w m and λ, namely, 0 ≤ w m < 1 and 0 ≤ λ ≤ λ c,2 , and it represents a matter dominated (Ω m = 1) solution. More precisely, this is a set of critical points forming a, normally hyperbolic, critical line.
Therefore, the stability conditions for this point can be analysed within linear theory. P 1 has one null eigenvalue. Indeed, we have: with µ 2 < 0 and µ 3 > 0. The eigenvectors are, respectively: From this it follows that on the eigenspace of the zero eigenvalue µ 1 , the solutions are time independent, i.e., this eigenspace is a line of equilibrium. Being the two eigenvalues µ 2 and µ 3 of opposite sign we can conclude that it will be a saddle line of equilibria. This point can represent the early time phase of matter domination. • P 2 : scalar field point. This point represents a scalar field dominated solution with Ω φ = 1 and w φ = w eff = −1.
This critical point has one null eigenvalue as in the previous case, indeed we have: with both µ 1 and µ 3 negative. The eigenvectors are, respectively: This point is, therefore, a non-hyperbolic critical point and linear theory cannot be used to establish its stability. Using the centre manifold theorem, we find that it is a stable point. Details are shown in Appendix A. • P 3,± : Stiff matter point. These are two points, labelled P 3,± = (±1, 0, 0), having the same background properties but satisfying different stability conditions. They represent scalar field dominated solutions (Ω φ = 1) with w eff = 1. In more detail, P 3,+ is a saddle point and P 3,− is an unstable one. The latter corresponds to a stiff matter solution, which is expected to be relevant at early time during inflation. P 3,− is the only point of cosmological interest between the two; • P 4,± : Stiff matter point. The phenomenology of these critical points is very similar to that of points P 3± : Ω φ = 1 and w eff = 1. They correspond to: P 4,± = (±1, 0, λ c,2 ), with different stability configurations. According to the eigenvalues, we have that P 4,− is a saddle point, P 4,+ can be either unstable for α ≥ 1/36 or a saddle point for 0 < α < 1/36. The unstable branch of P 4,+ is a stiff matter solution and as such it is relevant at early time as for P 3,− ; • P 5 : Scalar field point. This is another scalar-field-dominated critical point with x 2 + y 2 = Ω φ = 1. In this case, the equation of state of the scalar field and the total one are w φ = w eff = λ 2 c,2 /3 − 1, therefore, acceleration is possible for λ 2 c,2 < 2. However, this point is not an attractor for the system, being a saddle point. It cannot be unstable because it would require α < 1/36, thus violating the existence condition, i.e. α ≥ 1/36; • P 6 : Matter scaling point. This point is characterised by w eff = w φ = w m , which exists for λ 2 c,2 ≥ 3(1 + w m ). This point has the density of the scalar field, which scales as the matter one and as such it is said to be a scaling critical point. When λ 2 c,2 = 3(1 + w m ), we have Ω φ = 1. These points are of interest to alleviate the Coincidence problem as the initial conditions for the scalar field are fixed and set to follow the matter field. Stability analysis identifies it to be a saddle point. It is of interest to investigate in more detail the stability conditions of P 6 for particular choices of α and w m . One of these choices is: w m = 0 and α = 1/36. This combination of parameters is appropriate, as it stays in the existence domain, i.e, 0 < α ≤ 1/18. In this case, λ c,2 = √ 6 and P 6 has coordinates {1/2, 1/2, √ 6} and Ω φ = 1/2. The corresponding three eigenvalues now read {3, − 3 4 (1 ± √ 3 ı)} and the point is, obviously, a saddle, but it becomes interesting because both µ 2 and µ 3 are complex and the orbits on the plane λ = λ c,2 spiral towards the critical point. For the choice λ 2 c,2 = 3(1 + w m ) it follows f (λ c,2 ) = 1 and P 6 is still a saddle with eigenvalues {3(1 + w m ), − 3 2 (1 − w m ), 0}. Let us summarise the results of this section. We found one matter point (P 1 ) and two stiff matter points (P 3 , P 4 ). We also find one matter scaling point (P 6 ). These critical points are of interest in cosmology for early-time evolution of the system; in particular the scaling point, as it is seen as a possibility to alleviate the coincidence problem. Indeed, for this case the density of the scalar field is no longer negligibly small compared to the other matter densities at early-time and the model can be compatible with the energy scale associated with particle physics. Moreover, for scaling solutions the initial conditions (ICs) are completely fixed by the model's parameters. Two other points (P 2 , P 5 ) are scalar-field dominated points and under the conditions to have accelerated behaviour and stable configurations, only P 2 is an attractor solution and as such of cosmological interest at late-time. Viable models undergoing the transitions from P 1 to P 2 can be constructed. On the other hand, if we want to consider trajectories evolving from the P 6 scaling point towards the attractor point P 2 , one has to employ a second potential. That is because once the trajectories enter the scaling regime, the scalar field can no longer exit because the potential is too steep, therefore, a second potential is necessary to realise the cosmic acceleration. In Section 4 we will show a concrete example with this configuration.
Finally, it is of interest to compare the results of this section with the critical points and their stability in the case of the popular exponential potential (see, e.g., Ref. [34]). Except for point P 2 , the other points exist also for the exponential potential. However, the stability of the critical points for the two potentials differs. In fact, while for the exponential potential the critical points P 5 and P 6 can be stable under certain conditions, this is not the case for the Recliner potential. For example, for point P 6 , the presence of the eigenvalue 3(1 + w m ) forces the point to be a saddle, while its absence for the exponential potential makes it stable.

Broken Exponential-Law Potential
Finally, we discuss the evolution of the dynamical system defined by the Quintessence model with the potential given in Equation (9). For this model, in order to close the system of equations in (12), we use: Then, we end up solving the equation λ 2 (Γ − 1) = 0, which has two solutions: Let us notice that because the three parameters {n, p, q} are considered to be positive, it follows λ c,1 > λ c,2 . However, there is no guarantee that λ c,1 is positive. Since the background behaviour of the critical points does not depend on λ, we will have pairs of points with the same background properties but which satisfy different stability conditions according to the exact value of λ c . Additionally, the two-dimensional phase-space of the exponential potential is now augmented by the critical points for λ defining the boundaries λ = λ c,2 and λ = λ c,1 and the phase-space of the dynamical system is, therefore, the half-cylinder with boundaries {x 2 + y 2 ≤ 1, y ≥ 0, −p ≤ λ ≤ 2nq − p}, being 2nq the height of this half-cylinder. Also in this case, the scalar field can assume both positive and negative values. In fact, φ/M pl → ∞ corresponds to λ c,1 and φ/M pl → −∞ corresponds to λ c,2 . Therefore, −∞ < φ/M pl < ∞ translates to λ c,2 ≤ λ ≤ λ c,1 .
Even though we have the freedom to choose the values of the free parameters, we cannot set n = 0, otherwise Γ − 1 diverges. This case, therefore, needs a specific treatment. If we set n = 0, the broken-exponential-law potential and the relevant quantities λ and Γ − 1 read Let us also note that for λ = −p constant, we fall back to the case of the exponential potential. The system is, therefore, bidimensional and we do not study it here, as it has been studied in detail in Ref. [10]. We recover the same situation if we set q = 0. Indeed, if we redefine V 0 = V 0 /(1 + A) 2n , we obtain exactly the same expressions for n = 0. Hence, in the following analysis we will not consider the cases n = 0 nor q = 0.
In principle, the free parameters could assume any positive or even null value. In this work, we will not impose any restriction on the value of these parameters. This implies that a critical point can have different stability regions depending on the combination of these parameters. We present the study of the background properties of the critical points in Table 3 and their stability in Table 4. Table 3. Critical points for the dynamical system with the broken exponential-law. λ a represents an arbitrary value of the parameter λ. For each critical point we present its coordinates, the conditions for existence and acceleration and the corresponding values of the scalar field energy density Ω φ and equation of state w φ and the total equation of state of the system w eff . Table 4. Stability analysis for the dynamical system characterised by a broken power-law potential. We present the eigenvalues and the stability condition for each critical point.

])
Stable point for 2nq − p < 0 Saddle point 2nq − p > 0 Non-hyperbolic saddle point for p = 0 Non-hyperbolic stable point for p = 2nq Let us now discuss in more detail each critical point: • P 1 : Matter point. The point P 1 represents a matter-dominated solution with Ω m = 1 and w eff = w m , and λ a can assume any allowed value, i.e., λ c,2 ≤ λ a ≤ λ c,1 . It describes a line of hyperbolic critical points (critical line). The stability requirements identify it as a saddle point. Therefore this point is of interest in the early stage of the universe; • P 2 : Scalar field point. This point represents a scalar-field dominated solution with Ω φ = 1 and equations of state w φ = w eff = −1 which leads to an accelerated expansion. Depending on the sign of the quantity 2nq − p it is a stable or a saddle point. The configuration that is of interest for cosmological purpose is the one for which 2nq − p < 0, which leads to a late-time attractor solution. When p(2nq − p) = 0, P 2 becomes a non-hyperbolic critical point and to investigate its stability we use centre manifold method, as shown in Appendix A. We find that it is either a stable (p = 2nq) or an unstable (p = 0) point; • P 3± and P 4± : Stiff matter points. These points represent matter-dominated solutions with Ω m = 1 and equation of state w φ = w eff = 1. This is the typical behaviour of stiff matter. They are never stable points, but they can be either saddle or unstable points, depending on the values of the free parameters as shown in Table 4. In more detail, P 3,+ and P 4,− are always saddle points and P 4,+ is always unstable; P 3,− , instead, can be either unstable or a saddle, according to the value of 2nq − p, see Table 4. The saddle configurations are those of cosmological interest at early time; • P 5 : Scalar field point. This point exists for λ 2 c,1 ≤ 6. It represents another scalar-field dominated solution with Ω φ = 1 and w φ = w eff = λ 2 c,1 /3 − 1 and when λ c,1 = 0, this point coincides with P 2 . When λ 2 c,1 < 2 it defines a solution with an accelerated expansion. The stability conditions for P 5 are strongly dependant on the values of λ c,1 . It is a stable point if 0 ≤ (2nq − p) ≤ 3(1 + w m ); otherwise it is a saddle point. If λ c,1 = 2nq − p = 0 or λ c,1 = 2nq − p = 3(1 + w m ), P 5 is a stable non-hyperbolic critical point, as shown using the centre manifold theorem. The stable condition is of interest for cosmological purpose as it can define a late-time attractor; • P 6 : Scalar field point. This point exists for λ 2 c,2 ≤ 6. It is another scalar-field dominated point which satisfies the accelerated condition for p 2 < 2. When λ c,2 = 0, this point coincides with P 2 . The point is stable for 0 < p < 3(1 + w m ). The stable configuration has the requirements to be an attractor solution at late time. When p = 0 and p = 3(1 + w m ) it is a non-hyperbolic critical point. It is unstable for p = 0 and stable for p = 3(1 + w m ). • P 7 and P 8 : Matter scaling points. Both critical points represent matter scaling solutions being the energy density of the scalar field Ω φ = 3(1 + w m )/λ 2 c and w φ = w eff = w m with λ c being either λ c,1 = 2nq − p (P 7 ) or λ c,2 = −p (P 8 ). The point exists if λ 2 c ≥ 3(1 + w m ). Whilst P 8 can only be a saddle point, the phenomenology of P 7 is much richer and the point can undergo differ stability configurations. In particular it can be a saddle point for 2nq − p < − 3(1 + w m ) and a stable one when 2nq − p > 3(1 + w m ). P 7 can also be a non-hyperbolic point when 2nq − p = 3(1 + w m ). In this case P 7 is stable. For both points, it is possible to realise the situation in which f (λ c ) is complex (i.e., when the argument of Equation (13) is negative), in this case the orbits spiral towards the critical point.
In summary, the model characterised by the broken exponential-law allows for one matter point (P 1 ) and two stiff matter points (P 3 , P 4 ) that are of interest for the early time evolution of the system. Then we also find two matter scaling points (P 7 , P 8 ) which can be used to set the evolution of the scalar field to follow that of the matter component at early time. This is of cosmological interest to alleviate the coincidence problem. Three other points (P 2 , P 5 , P 6 ) are scalar field dominated solutions and have the proper features to be late-time attractors. While the critical point P 2 has the same configuration of the Quintessence model with a Recliner potential, the other two points are not in common. Under the conditions discussed in this section it is possible to realise a viable evolution for the system with trajectories starting from P 1 to any of the scalar field dominated critical points. Instead what is difficult to realise is the condition to have trajectories from one of the matter scaling points to a scalar field dominated critical point. This is the same situation present for the Recliner potential. In Section 4 we will show that it is possible to overcome this issue if one considers a second potential. In this case, while the first potential controls the scaling behaviour, the second one allows to exit such a regime and to reach the scalar field dominated critical point.

Working Examples
In this section we will show some concrete examples for both potentials we have previously analysed. In particular, we will show how it is possible to realise an evolution for the scalar field which first follows a scaling behaviour with the matter component and then gives rise to the late-time accelerated expansion. For this purpose, we have to consider a double potential, namely V(φ) = V 1 (φ) + V 2 (φ). That is because for the cases considered in this work it is hard for a single potential to have the scaling behaviour at an early time and later to source the accelerated expansion. In this regard each potential will generate one of the two evolutions. This is a typical procedure already adopted in the presence of exponential potentials [25,35]. Therefore, we introduce the following parameters: for i = 1, 2. Equation (11) then modifies as follow: The same holds for the system (12), which becomes Hereafter, we will consider w m = 0 and we will neglect the radiation contribution. This assumption does not change the results of this section. We will set ICs for the system (28) in such a way that {x i , y 1,i , λ 1,i } are chosen on the scaling critical point, while {y 2,i , λ 2,i } are selected in order to have Ω φ,0 = 0.68. The initial time is set at a i = 3.3 × 10 −4 .
For the Recliner potential case we will then consider a transition between the matter scaling point P 6 and the scalar field point P 2 . The Recliner potential V 1 with free parameter α 1 will realise the scaling regime and the second Recliner potential V 2 , with free parameter α 2 , will then boost the evolution towards the attractor solution. In Figure 3 we show the evolution of the matter density ρ m (green dotted line) and the scalar field density ρ φ (red solid line). At early time, the scalar field shows a scaling behaviour ρ φ ∝ ρ m . For this case, the field density parameter corresponding to point P 6 is Ω φ (P 6 ) = 18α 1 = 0.0018. The latter is consistent with the bound found by the Planck collaboration in 2015 using Cosmic Microwave Background (CMB) radiation measurements, Ω φ < 0.02 at 95% C.L. [64]. Around a ∼ 0.08 the scalar field exits the scaling regime and enters in the attractor one. In the same Figure 3 we show the evolution for the scalar field density for the Broken exponential-law potential case (blue dot-dashed line). In this case we are considering a transition between the matter scaling point P 7 and the scalar field point P 2 . The Broken exponential-law potential V 1 with free parameters {p 1 , n 1 , q 1 } will be responsible for the scaling regime and V 2 with free parameters {p 2 , n 2 , q 2 } for the late-time accelerated expansion. Additionally, in this case the scaling behaviour is evident from the Figure and the field density parameter Ω φ (P 7 ) = 3/(2n 1 q 1 − p 1 ) 2 = 2.52 × 10 −4 . Then, quite early, around a ∼ 2 × 10 −3 , it starts to slightly deviate from a strict scaling regime (although very close to it).  Other working examples can be realised by considering the transition P 1 → P 2 . In this case we can use for the ICs the values of x i , y 1,i , λ 1,i , corresponding to the critical point P 1 . In the case of the Recliner potential we can use x i = y 1,i = 0 and an arbitrary value for λ 1,i , which we chose to be λ 1,i = 1. The other two ICs are selected such that Ω φ,0 = 0.68 as for the previous cases: y i,2 = 9 × 10 −6 and λ i,2 = 10 −3 . The latter is the parameter which mostly affects the present day value of w φ , w 0 φ . The chosen value allows to have w 0 φ = −1, large values of λ i,2 move it to w 0 φ > −1. These ICs allow the trajectories to start on the critical point P 1 and to closely approach P 2 . We found that the values of αs do not have any significant impact on the evolution of the matter and scalar field densities for this case. We also ensure that under this condition the matter era is long enough, with a transition to the scalar field dominated era to be at a = 0.77. Let us now consider the case of the broken exponential-law potential. The ICs for the transition between P 1 and P 2 can be chosen as in the previous case, because they share the same coordinates. For this case we fix p 1 = 1, n 1 = 5, q 1 = 11. Then we consider for P 2 the configuration for which p 2 = 0 and q 2 = 0.1. For n 2 we checked that changing its value of several orders of magnitude (10 −7 < n 2 < 1) will barely change the general behaviour of the system. However, trajectories closely approach P 2 for the smaller value. For completeness, we select the ICs for y i,2 ≈ 0.89 × 10 −5 (to be slightly adjusted depending on n 2 in order to have Ω 0 φ = 0.68) and λ i,2 = −0.1. The transition between the matter domination and the scalar field domination happens around a ≈ 0.78.
We have shown how viable evolution can be obtained with a double potential of the same type. We have tested that scaling regimes realised with any of the two potentials considered in this work can then be followed by any different potential characterised by a late-time stable point, e.g., the exponential potential.

Conclusions
We have presented a thorough analysis of the phase-space dynamics of two Quintessence models characterised by either the Recliner potential, usually used to model both the late-time accelerated expansion and the inflation one, or by the broken exponential-law potential, which is a novel proposal in this work. We used the formalism of the dynamical systems, which means that the dynamics of the system under consideration needs to be described by an autonomous system of first-order differential equations. In particular, for Quintessence models, the system becomes autonomous if some assumptions on the roll (λ, proportional to the first derivative of the potential with respect to the scalar field) and tracker (Γ, proportional to the second derivative of the potential with respect to the scalar field) parameters are made. The standard procedure is to assume Γ = 1 and λ = const, which corresponds to an exponential potential and the resulting system is bi-dimensional. Otherwise, a different approach considers Γ(λ) and the system is then three-dimensional. We have considered the latter approach for our analysis, which is more general, and we proceeded in deriving the analytical relations, Γ(λ), for the two potentials considered. Let us stress that this approach can be used only if it is possible to find an analytical expression φ(λ), which inserted into Γ leads to Γ(λ) and closes the system. Nevertheless, this expression might also not be easily managed, and this explains why only a few potentials have been investigated so far.
For each potential we have then identified the critical points and by studying their stability we have determined the evolution of the system around them, and as such the evolution of the Universe governed by a Quintessence field during different cosmological phases. The important critical points for both the systems are the matter scaling solutions (P 6 for the Recliner potential and P 7 , P 8 for the broken exponential-law potential) and the scalar-field dominated solutions (P 2 for the Recliner potential and P 2 , P 5 , P 6 for the broken exponential-law potential). The former can be used to set initial conditions at early-time and as such it justifies why the scalar-field density is not negligibly small at early-time and can mitigate the coincidence problem, and the latter provides solutions for the late-time accelerated expansion. In order to have solutions for the system that go from a matter scaling point to a scalar-field dominated point one has to adopt a double potential, for which each potential drives the evolution towards the two epochs. We showed some practical examples in Section 4, where the double potential is the result of the combination of either two Recliner potentials or two broken exponential-law potentials. For these two concrete examples we also found that the early-time scalar field density is consistent with CMB constraints.
We conclude that both proposed models show very interesting features for realising a viable expansion history and alleviating the coincidence problem. Therefore, it will be of great interest to analyse the phenomenology of perturbations and the effects on cosmological observables and to set constraints on the model parameters using the most recent datasets. We leave this analysis for future work. and the matrices A ∈ R c × R c and B ∈ R s × R s having, respectively, zero real part and negative real part eigenvalues.
We define the centre of the manifold of the system in Equations (A1) and (A2) as for δ sufficiently small. For sufficiently small u ∈ R c , the reduced system readṡ From the definition y = h(x), we haveẏ = ∇h(x) ·ẋ and using (A1) and (A2) we arrive to a quasi-linear partial differential equation which h(x) has to satisfy: and once solved, we insert h(x) in Equation (A5) and we can study its properties. In general, there is no analytical solution to this equation, but its knowledge is not necessary as we only need it around the critical point. Often it suffices to assume h(x) of the form h i ( i ) and to determine the first non-trivial terms of the Taylor expansion. If they are all null, we need to consider higher order terms.
Operationally, let us consider a system in three variables. If the critical point has coordinates (x c , y c , λ c ), we define new variables such that X = x − x c , Y = Y − y c and Z = λ − λ c , then we rewrite the dynamical system in the new variables and we determine its Jacobian matrix J, evaluating it in the centre of the reference frame. We then construct the matrix T whose columns are the eigenvectors of J, evaluate its inverse T −1 and introduce new variables (u , v , w) defined as (u , v , w) T = T −1 (X , Y , Z) T , where T indicates transposition. In these new variables, the system reads The matrix D = T −1 JT is now in the desired form and we can determine the reduced system, the vector functions f and g, and the unknown function h.
In the following we will apply this technique to evaluate the stability of the non-hyperbolic critical points of the Recliner and broken-exponential law potentials.
We first construct a new set of variables, X = x, Y = y − 1 and Z = λ and construct the new dynamical system in these new variables. After this, we define a new coordinate set (u, v, w) as To apply the centre manifold theorem, we need the three functions f , g 1 , and g 2 , which read which is satisfied for The dynamics of the system restricted to the centre manifold is thus and the point is stable, being the coefficient of u 2 negative [62]. We can then conclude that the P 2 point is an attractor for the system for any value of α.
Appendix A.2. The Broken Exponential-Law Potential The phenomenology of the broken exponential-law potential is much richer and, under certain combinations of the free parameters, we can have many more non-hyperbolic critical points than the Recliner potential, as we already discussed in the main text. Indeed, we can have up to four of such points. We will examine them independently. Given their length we will not report all the expressions, but exclusively the function N and the matrix relating the coordinates (X, Y, Z) to (u, v, w) and the final function f (u).
Firstly we consider the critical point P 2 ≡ (0, 1, 0). When p(p − 2nq) = 0, the eigenvalues of the critical point become (0, −3, −3(1 + w m )). This is exactly the same configuration of the point P 2 for the Recliner model investigated in the previous section and the same variables (X, Y, Z) and (u, v, w) can be constructed. The condition p(p − 2nq) = 0 is satisfied either for p = 0 or for 2nq − p = 0. For the former case, Equation (A6) reads which is satisfied for The dynamics of the system restricted to the centre manifold is We then conclude that for p = 0, the point P 2 is a saddle as q > 0 by assumption. When p = 2nq, we find The dynamical equation is now We conclude that P 2 for p = 2nq is stable. The second point we investigate is P 5 . The phenomenology of this point is very rich and it can be a non-hyperbolic point for different values of the free parameters. It is a non-hyperbolic point either for 2nq − p = 0 or 2nq − p = 3(1 + w m ). In the first case, it coincides with P 2 and it is stable as q > 0. In fact, the restricted dynamical system reads In the second case, the eigenvalues are {0, −q 3(1 + w m ), −3(1 − w m )/2}. Since the expression for the vector N and the relevant coefficients are rather complicated, here we present the reduced dynamical system up to second order, which reads