Recent Advances on Inflation

We review recent trends on inflationary dynamics in the context of viable modified gravity theories. After providing a general overview of the inflationary paradigm emphasizing on what problems of hot Big Bang theory inflation solves, and a somewhat introductory presentation of single field inflationary theories with minimal and non-minimal couplings, we review how inflation can be realized in terms of several string motivated models of inflation, which involve Gauss-Bonnet couplings of the scalar field, higher order derivatives of the scalar field, and some subclasses of viable Horndeski theories. We also present and analyze inflation in the context of Chern-Simons theories of gravity, including various subcases and generalizations of string corrected modified gravities which also contain Chern-Simons correction terms, with the scalar field being identified with the invisible axion, which is the most viable to date dark matter candidate. We also provide a detailed account of vacuum $f(R)$ gravity inflation, and also inflation in $f(R,\phi)$ and kinetic-corrected $f(R,\phi)$ theories of gravity. In the end of the review we discuss the technique for calculating the overall effect of modified gravity on the waveform of the standard general relativistic gravitational wave form.


I. INTRODUCTION
Physicists of the current era are lucky to live in the era of precision cosmology, in which a plethora of observational data are available. Without exaggeration, nearly every 3 years a great discovery takes place. It started in 2012 with the Higgs discovery at the LHC [1], followed by the first direct observation of gravitational waves in 2015 [2], followed by the one in a million kilonova event in 2017, nowadays known as GW170817 [3]. After that, in 2020 the NANOGrav collaboration reported on the discovery of something that could be either a pulsar red noise, or a gravitational wave [4]. The latter was verified in 2023 in which year the NANOGrav reported the first detection of a stochastic gravitational wave background, verified by Hellings-Downs correlations, so the signal is clearly a gravitational wave of cosmological or astrophysical origin [5]. The chorus of observations will be further augmented by future experiments, like the stage 4 Cosmic Microwave Background (CMB) experiments [6,7], expected to commence operations in 2027, and the future gravitational wave experiments [8][9][10][11][12][13][14][15][16] like LISA and the Einstein Telescope which will commence their operation in 2035. All these new experiments will shed light to fundamental problems in cosmology and astrophysics. The most important of all, they will probe the primordial tensor modes in our Universe and give a definitive answer on the question whether inflation ever occurred. The inflationary scenario [17][18][19] is one of the most viable theoretical proposals for the early time era of our Universe. This is because inflation as a proposal solves all the basic problems of the standard hot Big Bang model and in addition it serves as a mechanism for generating matter structures at large scales in our Universe. The primordial quantum fluctuations generated during the inflationary era act as attractors on which baryons and cold dark matter are accumulated, so the large scale structure at large redshift up to z = 6 may be explained by a nearly scale invariant power spectrum of primordial scalar perturbations. Thus from a theoretical point of view, the inflationary scenario is the most appealing and vital for the viable description of the primordial evolution of the Universe, and for explaining the large scale structure of the Universe. In principle many theoretical frameworks may generate an accelerating era primordially, for a main stream of articles and reviews see . Traditionally, the inflationary scenario was firstly considered in terms of some false vacuum decay of a scalar field [17], but it was then realized that a slow-rolling canonical scalar field may appropriately describe the inflationary era [18]. Usually inflation is formalized by the use of a canonical minimally coupled scalar field or by a non-minimally coupled scalar field, however there are many shortcomings of the scalar field description of inflation that make an alternative description rather compelling. The most important shortcoming is that the single scalar field description relies on a scalar field, which must couple to all the Standard Model particles in order to reheat the Universe. Thus there are too many unknown arbitrary couplings that must be explained and also the inflaton itself must be identified or experimentally verified. The only fundamental scalar field that has ever been observed is the Higgs particle, so inflationary scenarios that use the Higgs particle as the inflaton are well motivated [95]. There exist however models that assume that the axion is the inflaton, with the axion being also a dark matter candidate. However, this could be deemed problematic, since if the axion is the inflaton, it must couple to Standard Model particles in order to reheat the Universe. This is a rather unwanted situation, since in most contexts the axion is a non-thermal relic. Regardless the theoretical shortcomings, the scalar field inflationary models are the most common and frequently used descriptions of inflation. An alternative to scalar field inflationary models comes from modifications of Einstein-Hilbert gravity which contain higher order curvature corrections. There are many modified gravity models and for some important reviews in the field, see Refs. [96][97][98][99][100]. The most important and more common models of modified gravity use f (R) gravity, with R being the Ricci scalar, see Refs. [96][97][98][99][100]. These models are simple and mainstream models and can also have an Einstein frame counterpart theory which is a minimally coupled scalar field. One important feature of this f (R) gravity framework is that it is possible to achieve a unified description of inflation with the dark energy era, see the pioneer work on this [101] and several works thereafter [43,[102][103][104][105][106][107][108][109][110]. Apart from f (R) gravity inflation, modified Gauss-Bonnet theories of gravity are also well studied in the literature of the form f (G) [111][112][113][114][115] or f (R, G) [116] where G is the Gauss-Bonnet scalar. The f (R, G) theories are plagued with ghost instabilities and degrees of freedom, thus they are less frequently used and are less appealing. In addition, quite popular theories are the Einstein-Gauss-Bonnet theories of gravity [112,113,., but these theories are plagued with a non-trivial gravitational wave speed which is distinct from the speed of light in vacuum. These theories were severely constrained after the GW170817 event since the latter indicated that the electromagnetic signal arrived almost simultaneously with the gravitational waves. However a theoretical solution for this problem was given recently in Refs. [157,158] in which case the theory resulted to a constraint between the scalar field potential and the Gauss-Bonnet coupling function. In addition, models of f (R, ϕ) gravity are also used to produce inflation and also Chern-Simons models of inflation are also used. In addition, several f (R, T ) models can also describe the inflationary era, where T is the trace of the energy momentum tensor, see the reviews [96,97]. Also two scalar field models are used for inflation see [83].
The recent NANOGrav observation of the stochastic signal of gravitational waves if interpreted as a cosmological signal, indicates that the inflationary era must have had a blue-tilt in its tensor spectral index, a fact that severely constrains the inflationary models that can reproduce such a signal. This can be achieved by Einstein-Gauss-Bonnet models, at the cost of having a low-reheating temperature, see for example [55,160], while single scalar field models cannot explain the 2023 NANOGrav signal. Thus the current epoch puts severe constraints on theoretical frameworks, which must explain the observations that come out at an incredible speedy way nowadays. Thus to stand to the challenge, the modern theoretical physicist must master the techniques of inflationary cosmology and know how to judge whether a theoretical model is viable or not. In this review we aim to provide a timely text which contains all the recent trends and techniques on inflationary dynamics, but also discussing the standard problems of Big Bang cosmology and why inflation itself describes successfully in a theoretical way the primordial era of our Universe. We shall analyze inflationary dynamics in a quantitative way for single scalar field inflation, both minimally and non-minimally coupled, and for most mainstream modified gravity theories, including those for which the tensor spectral index can be blue tilted, a necessary ingredient in order for the models to be compatible with the NANOGrav stochastic gravitational wave observation, if the cosmological description is responsible for the signal of course. Our analysis will be limited to providing the necessary tools in order for someone to be able to produce viable inflationary cosmologies, compatible with the most recent (2018) Planck constraints on inflation [161]. The scalar and the tensor power spectrum can be written in terms of a certain set of dimensionless parameters which are known as "slow-roll" parameters. When these are smaller than unity, a perturbation expansion can be performed on the scalar and tensor power spectrum and the observable quantities which quantify the inflationary era can be expressed in terms of the slowroll parameters. Regarding the observable quantities of inflation, we shall focus on the most important ones, which are the spectral index of the scalar primordial curvature perturbations, the tensor spectral index of the primordial tensor perturbations and the tensor-to-scalar ratio. Depending on the theoretical framework, the number and the complexity of the slow-roll parameters varies, so we shall emphasize on the calculation of the slow-roll parameters for various mainstream theoretical frameworks and we shall express the observational quantities of inflation in terms of the slow-roll indices needed for each theory. We shall also provide relations in closed form for all the observational quantities and the corresponding slow-roll indices so that the reader is able to reproduce the results quoted in each case. In the end of the review we provide a concrete self-contained section on the evolution of tensor perturbations in the context of modified gravity and we quantify the effect of modified gravity in terms of a single parameter. Then we analyze how the general relativistic waveform of the tensor perturbation may acquire a non-trivial multiplicative factor which contains the overall modified gravity effect from the present day back to the redshift corresponding to the mode that reentered the Hubble horizon back in our Universe's past.
We need to note that the inflationary scenario is a theoretical necessity in order to alleviate the problems of standard hot Big Bang cosmology. Even if we have no direct sign that inflation ever occurred, the observations "cry for inflation" since it is the only consistent scenario that can be compatible with the observational necessity of having a nearly scale invariant power spectrum of primordial perturbations, and it is the only consistent answer to the question how large scale matter structure was generated in the first place. Inflation and dark matter are theoretical predictions that still wait to be revealed observationally and experimentally. We may still be far away from discovering those theoretical predictions and even if we did not find them yet, it is almost certain in the minds of theorists that both play an important role in the evolution of the Universe. The situation here is the same as in the discovery of the Higgs, every theorist believed that the elusive spinless particle gives mass to the Standard Model particles, but it was never observed until 2012. We all knew it was there and only the proof of its existence remained. And we knew because the Higgs particle and the electroweak symmetry breaking was the only theoretical mechanism that could yield a mass to the Standard Model particles. The same applies with the inflationary paradigm. It must be the underlying mechanism responsible for the CMB anisotropies and the reason that large scale matter structure exists, and we need to find a proof for its existence. The road might be long and thorny till we unveil inflation, or the proof might be right at the corner of our "local time frame of inertia". Nobody knows and that is the magic of nature.
This review is organized as follows: In section II we provide an overview of the inflationary paradigm. We point out the shortcomings of the standard hot Big Bang scenario, and we explain how the inflationary paradigm theoretically solved these problems. We emphasize how important is the inflationary paradigm theoretically and why it eventually should be the correct description of nature, since it is the only scenario which provides a nearly scale invariant power spectrum of primordial scalar curvature fluctuations, which are necessary in order to explain large scale structure in our Universe as a whole. In section III we analyze in detail how the inflationary era may be generated by a single scalar field theory with minimal and non-minimal couplings. We calculate the necessary slow-roll indices, we provide and prove in detail several well-known formulas regarding single scalar field inflation. In section IV we provide a brief account of the swampland criteria, while in section V we discuss in brief the constant-roll evolution as an alternative to the standard slow-roll evolution. In section VI, we study and analyze several string motivated models of inflation, which also involve Gauss-Bonnet couplings of the scalar field, higher order derivatives of the scalar field, and some subclasses of viable Horndeski theories. In section VII, we present and analyze inflation in the context of Chern-Simons theories of gravity, presenting various subcases and generalizations of string corrected modified gravities which also contain Chern-Simons correction terms, with the scalar field being identified with the invisible axion, which is the most viable to date dark matter candidate. Section VIII is devoted to vacuum f (R) gravity inflation in its most general form, while section IX is devoted to generalized f (R, ϕ) theories of gravity, in the form of scalar field assisted f (R) gravity inflation. Section X focuses on k-essence f (R) gravity theories, while section XI focuses on kinetic-corrected f (R, ϕ) theories of gravity. Finally in section XII we provide a concrete overview of the evolution equations for the tensor perturbations in the context of modified gravity and we review how to calculate the overall effect of modified gravity on the general relativistic gravitational wave waveform. Finally the conclusions follow at the end of the review.

II. BRIEF OVERVIEW OF INFLATION
In the early period of the 1970's, physicists started to realize some problematic elements of the conventional Big Bang model and thus the first aspects of inflationary cosmology started to formulate. Different characteristics of the mechanics for inflation were discovered and the first somewhat realistic model, in which the early universe went through an inflationary de Sitter era, was proposed by Starobinsky [19]. Also an equally important point in the historical development of inflation was the model proposed by Guth [17], in which the inflationary era was a period when the Universe exponentially expanded in a super-cooled false vacuum state (a metastable state with no particles or fields, but large vacuum energy density). However due to the possible outcomes of this model, it was recognized that it is not realistic and viable even with improvements in order to explain several shortcomings of Big Bang cosmology. The solution was given by Linde [18] , who proposed the "new inflationary theory" [18], in which inflation can start either in a false vacuum or at the top of an effective potential in an unstable state and then the inflationary field slowly rolls down to the minimum of that effective potential. Various models for inflation have been constructed ever since. The main fundamental idea for all of them is that: -Inflation is an era of an abrupt accelerated expansion that took place in a very early period after the beginning of the Universe.-Some are following a canonical approach using scalar fields and others are using a description based on modified gravity for inflation, with the most interesting and promising ones being the f (R)-gravity models. Descriptions for a lot of these cases are going to be presented in this text, however we start by presenting the key problematic elements of the Big Bang model that were the reason for introducing inflation as the optimal theoretical description of our Universe's primordial era.
A. The Shortcomings of the Hot Big Bang Cosmology The Standard Big Bang (SBB) cosmology paradigm was a very successful model, since there are various successful observations about the properties of cosmic objects based on the SBB theory. Also, along with the study of the cosmic microwave background (CMB), it guided us to our first understanding of various cosmological phenomena and subsequently how the Universe evolved at the very early epochs of its existence, and how it evolved to its subsequent large-scale form. However, even though the Standard SBB scenario found success and was highly embraced, one cannot ignore that some of its fundamental features can be problematic.
In the SBB, the early Universe is adiabatically expanded and radiation-dominated. The model depends on the assumption of homogeneity and isotropy on large scales (cosmological principle), which lead to the ability to use the Friedmann-Robertson-Walker(FRW) metric for the space-time of the Universe, given by the following form, where a(t) is the scale factor that is associated with the spatial part of the metric and its evolution in time. Depending on the value of constant K, we have three different case. For K = 0 the space described by this metric is flat (flat Universe) , for K = +1 is a closed Universe that could be described with a sphere and for K = −1 the Universe is open, with the spatial part of the spacetime being some hyperbolic hypersurface, like a saddle. In Eq. (1) the coordinates used are called the comoving coordinates, meaning that as the Universe and thus space itself expands, the coordinates r, θ and ϕ are not affected by this expansion being expanded in the same way. To obtain a form of the metric with the physical distance, the scalar factor is multiplied by r, R = a(t)r. Therefore the metric is now written in the following time-dependent form [34]: where, Assuming the FRW metric and c = 1 (in natural units for example), the evolution of the Universe is mainly described by the form of the scalar factor a(t) and by transforming the form of the metric with respect to the conformal time τ we have, where, So by using this metric and these coordinates, the structure of spacetime is more easily visualized. Specifically in an isotropic Universe, the geodesics of propagating light are null ds 2 = 0. Thus the radial null geodesics for propagating light in an isotropic, homogeneous Universe in the conformal time frame are, which correspond to straight lines at angles of ±45 • in the τ − χ (conformal time-space) plane (see 1). At this point, whereȧ = da dt . The Hubble parameter is positive for an expanding Universe and negative for a collapsing one. The other important quantity is the e-foldings (or e-folds) number defined as, which represents the number of Hubble times between two epochs with scale factors a 1 and a 2 . The Hubble time is H −1 , which represents the time period that takes for the Universe to expand substantially. The e-folds number can also be given as N = t2 t1 Hdt. The range of the e-folds number for a successful inflationary model that can effectively solve the flatness and horizon problems, which are going to be explained later, is from N ≈ 45 to N tot min ≃ 60, in the context of most inflationary cosmologies, with the estimations being based on single scalar field inflationary models. However this value is model dependent and for some specific inflationary models , it is possible that N has a higher value than 60 e-folds, a feature that mainly depends on the total equation of state parameter of the Universe at the end of the inflationary era [162], and in others N tot ≫ 60, with N tot having no upper bound that is specifically correlated with the idea of "eternal" inflation [163] . [2pt] In general, any physical system and its behavior can be described, in the context of the fundamental laws of physics, with the action, where L is the Lagrangian of the system, with an allowed explicit dependence on t, [8]. According to Hamilton's principle by taking the variation of the action δS = 0, the path such as to maximize or minimize the action integral can be derived as, which leads to the Euler-Lagrange equations of motions for the system, Now accordingly in the context of field theory, a Universe as a physical system can be described by the scalar field action with generic coordinates, where g is the determinant of a metric g µν and L = L(ϕ, ∂ µ ϕ) is the Lagrangian density that has dimensions of [m] 4 in natural units, it is Lorentz invariant and related to the Lagrangian of the fields by L = d 3 x L. The equations of motion and basically the evolution of the field can be determined by the same variation principle for the action with respect to the field ϕ and they are called the field equations . Relating to the Lagrangian density L, the energy-momentum tensor T µν is defined as, and it is conserved with, Eq. (13) can be thought as expressing two continuity equations, where in Eq. (14) corresponds to the energy continuity equation and Eq. (15) corresponds to the momentum continuity equation in curved spacetime, also known as the Euler equation as well.
On another note, related to the Hubble rate H, by inserting Eq. (2) into the Einstein equations, we can derive the following two equations also known as the Friedmann Equations, where , G is the gravitational constant, M P l is the reduced Planck mass and c the speed of light. When Eq. (16) and (17) are combined they give rise to the continuity equation, or the energy conservation equation, By defining the equation-of-state parameter w, and integrating Eq.(18) we get, which indicates how the scale factor evolves with time [34]. The value of the equation-of-state parameter depends on what perfect fluid energy density dominates the scale factor of the flat Universe. Possible values of w could be: w = 0 for non-relativistic non-baryonic matter, w = −1 for vacuum energy (for a cosmological constant), w = 1/3 for radiation or relativistic matter and w = +1 for stiff matter. Additionally, causality demands that w has to be strictly less or at the very least equal to 1 , while no lower-bound appears. For negative values w < − 1 3 , the equation-of-state parameter relates to the presence of a dark fluid. For a quintessence and quintessential evolution of the Universe, w ̸ = −1 and −1 < w < −1/3, while also data from Planck+WP, SNLS and BAO suggest the limits of −1.3 < w < −0.81 (95%C.L.) [165]. For w < −1, relates to the existence of hypothetical phantom dark energy, and it can lead to a Big Rip singularity. While the equation-of-state parameter can be taken as a constant, it can also be a function of redshift or time. It is also possible to have significant contributions from various kinds of matter to the total energy density ρ and pressure p, thus they can be given by the sum of all the individual components as, For each of the components there is also the ratio of the observed energy density to the critical density for the present time t 0 , where ρ crit = 3H 2 0 /κ 2 is the density required to make the expansion of the Universe slow down, stop and reverse, meaning the existence of a closed Universe. Additionally for the curvature contribution at t 0 , and using Eq. (16) and setting the scalar factor as a(t 0 ) = a 0 = 1, Therefore, Eq.(25) for the present time t 0 leads to, Generally as a function of the cosmic time t the curvature parameter is, and by assuming the simplest case of the expansion being dominated by a form of matter with an equation of state parameter w and with the scale factor given by a(t) = t 2/3(1+w) and also by taking the first derivative with respect to time, we obtain the following equation, From Eq. (28), it is deduced that the value w = − 1 3 is an equilibrium point and an unstable one, more specifically. This indicates that if the strong energy condition is satisfied, when w deviates a bit from that unstable equilibrium point, for w > −1/3 with Ω K > 0, the parameter Ω K keeps growing and with Ω K < 0, Ω K keeps decreasing away from zero. From experimental data and observations of large-scale structures and the CMB it is concluded that −0.0179 < Ω K < 0.0081 (95%CL). Meaning that the presently observed Universe is very flat and therefore it used to be even more flat in the past with the value of parameter Ω K reaching exceptionally small values, Ω K ∼ 0. This could be the case, if we assumed that K = 0 precisely, for the Universe in its very early initial state. However, it seems quite peculiar to have such a precise value of K and no explanation why this would be the case.
And this is where the first signs of how some characteristics of the SBB could be problematic and a new approach is needed. The conventional SBB model describes the Universe as very homogeneous, isotropic and flat, which are features that do not emerge from a fundamental mechanism that explain how the Universe came to be or evolved in that fine-tuned way, but only from the initial assumptions of the construction of the model. Apart from the flatness problem, there exist several other shortcomings of the SBB model, and inflation, provides concrete theoretical explanations for these problems [? ? ] and there are presented one by one in the following subsections.

The Horizon problem
Before continuing, we will provide definitions of some quantities which are important for better understanding of the concepts described here. In Eq. (7), the Hubble parameter was defined. There are also the quantities of the Hubble distance, or horizon, given as H −1 (t) and the particle horizon R, which is the distance that light could have travelled since the beginning when a = 0 and regions separated by distance more that the particle horizon are "causally disconnected", meaning they can never communicate with each other. The particle horizon is also called the comoving horizon and is given by the following relation, where (aH) −1 is the comoving Hubble radius and if some regions are separated by a distance greater that the coming Hubble radius, they cannot communicate with each other at that time. By this definition, for a Universe dominated with a fluid with an equation of state parameter w, the Hubble radius is, So by the conventional SBB, where w ≥ 0, the Hubble radius grows monotonically and the comoving horizon R increases with time. This means that comoving regions that are entering the horizon in the present time, and thus becoming causally connected now, used to be outside of the limits of the horizon and therefore causally disconnected in the past and especially in the primordial era of the Universe. Consequently, the high homogeneity and thermal equilibrium of far distant regions of the Universe observed in the CMB, that were supposed to be causally disconnected in the early times, is quite problematic. This is also sometimes mentioned as the Horizon problem. So if homogeneity is not assumed as an initial feature, there must be some mechanism that made the Universe to evolve that way. However the SBB model cannot provide an explanation on how the Universe could be that homogeneous without fine-tuning the initial conditions. Inflation solves this problem by introducing a period with a decreasing comoving Hubble radius (aH) −1 . The exponential growth of the scale factor a during inflation and the relatively constant H allows the Hubble radius to be decreasing while Inflation is happening. So in the very early period, the comoving Hubble radius was much larger than the comoving horizon R and all the scales-regions of the Universe (relevant to cosmological observation) used to be small enough to be inside the Hubble radius and thus they were in fact causally connected. Then during inflation, the radius started to decrease, became smaller than the horizon while the volume of the horizon itself expanded. After Inflation ended, the Hubble radius started increasing its size again and in consequence, at present time, the horizon is larger than the Hubble radius and regions that were causally connected, now seem that they are causally disconnected, see Fig. (2). According to this mechanism, different regions of the Universe, which used to be closer together and inside the horizon in the early epoch, were able to communicate, become homogeneous with respect to each other and established thermal equilibrium. When the Hubble radius decreased and the Universe rapidly expanded they became causally disconnected, as we observe them to be today.

The Flatness Problem
We briefly discussed the Flatness Problem, now we shall further elaborate on this. The second issue that arises is the Flatness problem, where the conventional Big Bang model points to the fact that the Universe is very flat, with, which is verified by observational data. However, in the conventional model, this value is an unstable fixed point of the equations and there is no significant reason why the Universe should be at that unstable point, except in the case in which the initial conditions were fine-tuned that way. Combining Eqs. (26) and (27), we obtain, Therefore in the duration of the inflationary period, since there is the decrease of the Hubble radius (aH) −1 , the parameter Ω K decreases towards zero as well and thus the Universe evolves towards its flatness naturally, which is in At the beginning, all the Universe, whose ends are represented by the edge of the comoving horizon, is inside the Hubble radius. Thus all the regions of the Universe are causally connected and during this period there was enough time for them to become homogeneous and isotropic. As inflation begins the Hubble radius (red line) decreases in size and it keeps decreasing as inflation is in process. On the contrary, the comoving horizon (blue patch) is increasing in size, as the Universe is expanding. At some point the comoving horizon becomes greater in size than the Hubble radius and thus some regions of the Universe become causally disconnected from the regions of the Universe that are still inside the Hubble radius, therefore impossible for those different parts to communicate with each other. After inflation the Hubble radius starts to increase again and therefore those regions re-enter the Hubble radius and become visible again.
accordance with the experimental observations. It also justifies the disregard of the cases of a closed or open Universe since if the Universe was not flat, and for example we had a closed Universe with an intrinsic spatial curvature K = +1, no phase transition can occur and even if inflation washes out the effect of curvature such that the density parameter Ω K ≃ 0, K can never transition from +1 to 0, or in the case of an open Universe from −1 to 0.

The Primordial Relics Problem
In a very early period after the "Big Bang", the Universe can be described by Grand Unified Theories (GUT) or string theory or some Standard Model extensions. In a GUT scenario the Universe used to have a temperature of the order of the GUT temperature, which is about T GU T ∼ 10 28 K and the electromagnetic, weak and strong forces were unified. According to GUT, the Universe went through a phase transition when the temperature of the Universe dropped below the T GU T , and during that transition there was a production of primordial relics (e.g. domain walls, magnetic monopoles or other topological defects etc.), which are described as point-like topological defects in the scheme of GUT. So at the time of their creation, the number and energy density of primordial relics like magnetic monopoles should had been large, but still smaller than the ones for radiation during the GUT period and thus the Universe in that early stage is still radiation-dominated. Later on, as those relics should have been quite large, they would have quickly become non-relativistic, since ρ M ∝ a −3 for magnetic monopoles (massive particles) and ρ r ∝ a −4 for radiation, thus they would have dominated over radiation and ordinary matter until the present time. However, from observations there is no evidence for the existence of such primordial relics and definitely no signs that they dominate the Universe, with research setting an upper limit to their number density today of n M (t 0 ) ∼ 10 −19 cm −3 . This very small number density and the lack of this kind of observations today in the Universe compared to the predictions of the Standard Big Bang scenario along with particle physics, is the Primordial Relics Problem.
In the case of inflation, once the monopoles were created before or during the inflationary period, thereafter their number density would have decreased significantly during the rapid exponential expansion of the Universe. After the Universe expanded, the magnetic monopoles were basically so outspread in space that their number density n M (t 0 ) today would have reached such a small value, which renders them nearly impossible to detect.

B. Conditions for Inflation
Having presented the issues that motivated the construction of the inflationary paradigm, it is time to introduce the conditions under which inflation takes place. During the inflationary period the Hubble radius decreases with time, as mentioned in the previous paragraphs, therefore, d dt Since the Hubble rate is given by Eq. (7), two more equivalent conditions can be deduced for inflation, From condition Eq. (35), if there is a strong inequality, |Ḣ| ≪ H 2 , then the Hubble rate H is almost constant for many Hubble times and there is approximately an exponential expansion with a(t) ∝ e Ht . For H exactly constant over many Hubble times, inflation is described by a de Sitter expansion. If however the Hubble rate contains linear or higher order time dependencies, like for example H(t) ∼ H 0 − H i t, then this evolution is called quasi-de Sitter evolution. Thus inflation is described as an early era of rapid nearly exponential expansion of the Universe. From the conditions of Eq. (34), it is also derived that during inflation the scale factor a increases fast as a function of the cosmic time t. This condition is also achieved by an expansion with a power-law scale factor of a(t) ∼ t n with n > 1, however the tensor-to-scalar ratio for these power-law inflation models is larger than the limits set by Plank data, thus they are not appropriate to describe the inflationary era. Furthermore, from Eq. (21) for an exponential form of the scale factor, the equation-of-state parameter is w = −1, which corresponds to a cosmological constant and from Eq. (20), the density is ρ ∼ constant, which can also be approximately taken for a dominant dark energy era (w < −1) with w ∼ −1. There is also an inflationary condition for the pressure that arises if Eq. (17) and Eq. (35) are combined, where Λ is the cosmological constant. For Λ = 0 or absorbed into ρ and p, ρ + 3p Again here, it can be seen that for a dark energy dominated era w < −1 or for w ∼ −1, from the definition of the equation-of-state parameter Eq. (19), this condition if fulfilled.
In the next sections, we shall present the standard models that are used in the literature to describe inflation, namely scalar field inflation models and modified gravity models. While scalar field inflationary models are more customary, nowadays it seems that single scalar field inflation might be insufficient for describing the primordial era of our Universe. The reason is two-fold, firstly it is theoretically unappealing to have to explain the couplings of the inflaton scalar field to all the Standard Model particles in order to reheat the Universe and secondly, after the recent NANOGrav detection of the stochastic gravitational wave background [5], the cosmological perspective of such a background complexifies the single scalar field description of inflation, while leaving room for modified gravity descriptions [160].

A. Minimally-coupled Scalar field Inflation
For the scalar field models, the Lagrangian density for N real fields is assumed to have the form, where V is a function of all the fields and the fields with this Lagrangian are said to be canonically normalized. For the simplest case of a single scalar field in flat spacetime, the Lagrangian is given as, where the first is a kinematic term, V (ϕ) is the scalar field potential and both terms separately have mass dimensions of [m] 4 in natural units. In this section we will focus on the single scalar field inflation models with minimal coupling to gravity. In this case the inflationary period is associated with a single scalar field, dubbed inflaton and the minimal coupling to gravity conveys that the action for the inflationary field ϕ is not coupled with the scalar curvature in any way but through a term of the Lorentz invariant √ −gd 4 x from the metric. The energy density of the inflaton is dominant compared to the rest of the matter fields for this period and thus, no additional field emerges. We also consider a flat Universe, K = 0, described by a FRW metric. In the literature there are many models that may consider the effect of both a scalar and gauge fields simultaneously or model with intrinsic curvature, however this is not relevant to the context of this text.
So the action for minimally-coupled Single scalar field inflation considering Eq.(1) is, and g is the determinant of the FRW metric g µν for K = 0 [96]. The first term is the gravitational Einstein-Hilbert action including the Ricci scalar R, which in terms of the metric connection it has a form of R = g µν R µν with R µν being the Ricci tensor. The last two terms in the action of the scalar field are the canonical kinetic term and the scalar field potential V (ϕ) which takes into account the self-interactions of the scalar field. It is worth mentioning that quantum fluctuations can have as a result perturbations in the inflaton field and in the metric tensor with, ϕ → ϕ 0 + δϕ , g µν →ḡ µν + δg µν (41) whereḡ µν is the FRW metric and ϕ 0 the classical solution for homogeneous, isotropic evolution in the inflationary era. For the scalar field action of Eq. (40), through its variation, the scalar field equation is, where ∇ 2 is determined with comoving coordinates and for a homogeneous field, it takes the form, where V ′ (ϕ) is the derivative of the potential with respect to the field. Accordingly, the energy-momentum tensor here is defined as, and therefore the pressure and the energy density of the Universe can be obtained from T i j = −P δ i j and T 0 0 respectively, and In order for a scalar field to be able to produce a viable model for inflation, it needs to slow-roll under the influence of the inflationary potential V (ϕ) as it evolves towards the minimum of the potential, and also in a sufficiently slow manner, in order for the scale factor to increase enough to be able to resolve the standard Big Bang cosmology problems. The condition that is imposed for this effect to be secured is called the slow-roll assumption, which indicates that the kinetic term is sub-leading and becomes really small compared to the potential. Also under this condition the Friedmann equation (16) takes the form, and from the derivative of the condition (47), it is indicated thatφ ≪ V ′ , and in effect the scalar field equation of (43) results in, and from the Raychaudhuri equation (17), this condition, implies thatḢ ≪ H 2 , which is in agreement with the condition for the Hubble radius for homogeneity. In order to quantify the consequences of this assumption, the slowroll indices are introduced and for the single scalar field inflation, they can be expressed with respect to the Hubble parameter as follows, and therefore the results of the slow-roll assumption for the canonical field of (47) are translated as the following conditions for the slow-roll indices, for the time when inflation is taking place [96], with ϵ 1 ≪ 1 ensuring the occurrence of the inflationary era and ϵ 2 ≪ 1 ensuring that inflation lasts a sufficient amount of time so that the scalar field slowly evolves with respect to cosmic time t for a large number of e-folds, and the density parameter for curvature vanishes from the background equations. It is worth mentioning that the condition on the ϵ slow-roll parameter is not just a mathematical boundary of a specific mechanism in an inflationary model, but it is linked to the natural process of inflation. By elaborating on the mathematical form of (50) we get that ϵ = −Ḣ H 2 = 1 −ä ȧ a 2 . Independently of any specific mathematical formulation of a particular model, generally a key feature of inflation is that during the inflationary era the Hubble radius decreases and so d dt 1 aH < 0 → −ä a 2 < 0 and thus −ä ȧ a 2 < 0 → 1 −ä ȧ a 2 < 1 → ϵ 1 < 1. The same applies to the conditions for the end of inflation, since the Hubble radius rate with respect to time is equal to zero the moment that inflation ends. Another more familiar expression for the slow-roll parameters is the one with respect to the canonical scalar field potential V (ϕ), where V ′ , V ′′ are the first and second derivatives of the potential with respect to the field ϕ and the form of the potential V (ϕ) can take various forms depending on the model. Also in (52), the sign of η is insignificant and only its order of magnitude matters. The connection between the two representations of the slow-roll parameters is [96], where ϵ 1 = V ′ (ϕ) 2 /6H 2 V (ϕ) and ϵ 2 =φ/Hφ, with the derivative of (49). The end of the inflationary era occurs when basically the inflationary condition is violated, and the slow-roll perturbative expansion for the power spectrum breaks down, which is equivalent to the condition, After the end of inflation, the reheating era follows and the field oscillates about the minimum value of its potential. The quantum fluctuations δϕ of the scalar field basically generate the CMB fluctuations nearly 60 e-folds before the end of inflation. Lastly, two of the most important observational quantities, the spectral index of the primordial scalar curvature perturbations n S and the tensor-to-scalar ratio r, can be also expressed with respect to these slow-roll parameters as, or as n S = 1 − 4ϵ 1 − 2ϵ 2 , r = 16ϵ 1 for the slow-roll indices of (50). The are various models that consider different forms for the inflationary potential V (ϕ). Some examples are presented in the following subsections.

Massive Scalar Field
This is a simple case of single scalar field inflation, called the Chaotic Inflation, driven by a mass term. The field starts at a large value and rolls down towards the origin of the scalar potential and the form of the potential is, Since the potential is axial symmetric it is expected that either positive or negative values for the scalar potential could work and therefore the results are independent of the sign of ϕ k , which is the value of the field at the horizon crossing, something which is hinted by the slow-roll indices as well. The choice of sign can however in principle affect the evolution of the scalar field with respect to time. From (52) the slow-roll parameters are, and the end of inflation occurs when, The integral definition of the e-folds number N with respect to the filed ϕ is, By taking the slow-roll approximation, with (48) and (49), the integral for this specific potential form is, By replacing the ϕ end , in (60), the value of the field at horizon crossing is determined as, and inserting this value into (57) and then in (55), the spectral index and the tensor-to-scalar ration can be determined in the horizon crossing as, According to the relations of Eq.(69), for a number of e-folds of N ≃ 60, the spectral index and the tensor-to-scalar ratio have values of n s ≃ 0.9889 and r ≃ 0.044 respectively. From observations of the CMB by the Plank collaboration, the constrains on the values of the spectral index and of the tensor-to-scalar ratio are determined to be, n S = 0.9649 ± 0.0042 (68%CL) and r < 0.064 (95%CL).
Therefore we see that the values of n S and r deviate enough from the one determined by observation, thus this potential of the power-law form of Eq. (56) cannot be viable to generate an inflationary period. It is noted, as previously mentioned, that "power-law" models like the examples 1, 2 and 3 are not viable due to the range of values for the quantities n S and r, simultaneously compared to the limits from observational data from Planck. However this is specifically for the case of minimally-coupled scalar field inflation model and for the case of non-minimally coupled theories those models could possibly be proven viable under certain circumstances.

Self Interacting Scalar Field
In this case, the potential has a quadratic form with respect to the field, which has mass dimension of [m] 4 , since λ is a dimensionless coupling constant and ϕ 4 has a dimension of [m] 4 . The slow-roll parameters here are, therefore for the end of inflation, Similarly by taking the integral form of N and the slow-roll approximation, with (48) and (49), the integral for this specific potential form is, where . So by replacing the ϕ end in (67), the value of the scalar field ϕ k can be determined as, where M P l is the reduced Planck mass and N the number of e-folds. By inserting this value in (57) and then in (55), the scalar spectral index and the tensor-to-scalar ratio can be determined at the horizon crossing as, According to the relations of Eq. (69), for a number of e-folds of N ≃ 60, the spectral index and the tensor-to-scalar ratio have values of n S ≃ 0.95 and r ≃ 0.26 respectively. So based on these values in comparison again with the Planck values of (63), a self interacting scalar field with the potential of (64) of a power-law form, cannot be viable and is not fit to generate an inflationary era according to observational data.

Natural Inflation
Let us now consider the natural inflation model, usually considered in axion field contexts. In this case the inflationary field, inflaton, is represented by a pseudo-Nambu-Goldstone boson(PNGB), which could be e.g. an axion . The potential of the inflaton is, where L and f are two mass scales related to the height and the width of the potential respectively and they are of . The inflationary era corresponds in the region of 0 < ϕ < πf and it occurs as the inflaton evolves towards the potential minimum at ϕ = πf . The slow-roll parameters for this potential are, as seen, they depend solely on f ∼ M P l and not on L. Following the same process as before, it can be determined that, In order for this model to be viable for inflation and to obtain a number of e-folds N ≳ 70, it is required that the initial value of the inflaton is ϕ k ≲ 0.1M P l . It is also worth mentioning that the Natural inflation model, or axion model, can safely produce the power-law models examined before by simply assuming that ϕ f ≪ 1 and performing a Taylor expansion, which is connected to the kinetic axion model.

B. Observable quantities in the inflationary paradigm
In the previous sections, two important observable quantities, that of the spectral index of the scalar perturbations n S and the tensor-to-scalar ratio r, have already been introduced in the context of the single scalar field inflation. This section, concentrates more on the origin of these quantities and on the introduction of the tensor spectral index n T as well.
Even though the primordial Universe is considered to be homogeneous, CMB observations have proved that this is not entirely the case and it has anisotropies of lower order ∼ 10 −5 than the homogeneous background. Inflation can explain sufficiently these anisotropies, with the existence of quantum fluctuations in sub-horizon scales during the early periods of the inflationary epoch. So during inflation, perturbations are defined around the homogeneous background solutionsφ(t) of the inflaton and the metricḡ µν , as also similarly seen in (41), Specifically, during the inflationary era, the comoving Hubble radius decreases as the Universe expands and it becomes smaller than the comoving wavelength (horizon), 2. So when these fluctuations exit the horizon, they become causally disconnected and they remain frozen until the end of inflation, when the physical horizon expands again and they gradually reenter as classical density perturbations. During the time of inflation, the stress-energy tensor contributions are heavily dominated by the energy of the inflaton and therefore the perturbations of the inflationary field have some effect on the geometry of the spacetime through the field equations. Also, since the background spacetime is considered fairly symmetric, justified by being spatially flat, homogeneous and isotropic, the decomposition of the metric and stress-energy perturbations into independent scalar, vector and tensor components is possible. This approach is called the SVT decomposition, it can be described in the Fourier space and each type is able to evolve independently and treated separately. For Vector perturbations, it can be seen from the decomposition of metric perturbations that they are not created by inflation and nevertheless they are resolved while the Universe expands. So the focus is going to be on the scalar and tensor perturbations, which are observable as density fluctuations and gravitational waves. Depending on the comoving wavelength k of a mode it can be characterized as super-horizon when k < αH and sub-horizon for k > αH, while the sub-horizon modes also satisfy k ≫ αH when inflation is considered to be in its vacuum state and thus the fluctuations are produced at very scales inside the horizon. So after a mode has exited the horizon during its contraction, they can be described by a classical probability distribution, whose invariance is determined by the power spectrum at horizon crossing. Typically the condition reads c S k = aH however the model at hand predicts a sound wave velocity equal to the speed of light. k = αH. This is true for a single scalar field theory in a homogeneous flat background like the FRW spacetime. For scalar perturbations the power spectrum is expressed as, which relates to the primordial scalar curvature perturbations and for tensor perturbations, which corresponds to the power spectrum of primordial gravitational waves. The dependence of the power spectra on the scale is described through the scalar and tensor spectral indices respectively with, and Additionally, the tensor-to-scalar ratio is defined as, in order to correlate the amplitudes of scalar and tensor fluctuations and be able to compare them. It can be noted that these quantities can be considered scale invariant, since their values remain essentially unaffected under the change of scale k. Under the assumption of the slow-roll approximation, for canonical-scalar field models the power spectra can be expressed solely with respect to the potential of the field V (ϕ) and considering the definitions of the slow-roll parameters in (52) and the relations between the two representations (53), we arrive in the relations of (55) [96], For more complicated models than the ones presented previously, there is the introduction of extra slow-roll parameters, which are going to be included in the presentation of each case in later sections. Also for modified gravity models, the sound speed and the propagation speed of primordial gravitational waves is non-trivial too.
Observations that can contribute to the evaluation of these quantities play a detrimental role to obtaining valuable insight for physics in the primordial Universe. The tensor-to-scalar ratio r is an auxiliary parameter and as previously mentioned, it quantifies the ratio of the amplitude of tensor over scalar perturbations and it is evaluated at the CMB pivot scale k = 0.002 Mpc −1 . Some of the parameters may be scale dependent, for example in some models the scalar spectral index may have a non-trivial scale dependence, called "running", but we shall not consider such issues here. For the scalar spectral index n S , in principle it is predicted that for a completely homogeneous Universe n S = 1. However perturbations, which are, quantified by the power spectrum of (74), result in an apparent deviation observed in the CMB as mentioned in (63). Lastly, in contrast to the scalar spectral index n S , the tensor spectral index n T has not been computed yet, due to the lack of B-modes (curl modes) in the CMB. B-modes, which are a specific mode of polarization, can arise from the conversion of the E-mode polarization modes to B-modes that occurs at late times or on small angular scales, or from primordial tensor perturbations which are the inflationary tensor modes. So a detection of such B-modes directly give a verification for the existence of the inflationary era.

C. Non-minimally coupled Scalar field Inflation
There are a lot of models that can be used to describe the inflationary era, which have a more complicated theoretical background than the canonical minimally-coupled single scalar field that was mentioned previously. Some models may include further curvature correction terms with respect to the Ricci scalar, for the coupling to gravity described as f (R) gravity corrected canonical scalar field models or include multiple fields. A more general class of inflationary models can be described by the following action [96], where f (R, ϕ) is a smooth function of R and ϕ indicating the non-minimal coupling to gravity and the kinetic term ω(ϕ), for which if ω(ϕ) ̸ = 1 it refers to a non-canonical scalar field. In this section, we focus on the canonical non-minimal coupled model for scalar field inflation. In this case, the scalar curvature is no longer coupled with gravity only through the Lorentz invariant term √ −gd 4 x, but there is also another term that couples the field with the scalar curvature of the form f (R, ϕ), f (R, ϕ) = f (ϕ)R. This is a sub-case of the more general class from (80), described by the following action [96], with f (ϕ) being a dimensionless scalar coupling function. In principle when ϕ reaches its vacuum expectation value, it can become equal to unity and generate Einstein's gravity at late times. It has to be specified that the action of (81) is in the Jordan frame and one may choose to work in the Einstein frame, which means rewriting the action so that a linear term of the Ricci scalar appears as the sole contribution of curvature while higher order curvature corrections are described by means of a scalar field by performing a conformal transformation. Performing a conformal transformation in order to change the scale is not forbidden since general relativity has not an exclusive scale. Nevertheless, the description should in essence be the same between the two frames, when conformal invariant quantities are considered, while also taking into consideration the differences when imposing the slow-roll conditions. By varying the action (81) with respect to the metric and the scalar field ϕ, assuming the flat FRW metric, we get the equations of motion [96],φ The slow roll indices in this case are defined as [96], where E is a function defined as The parameters ϵ 1 and ϵ 2 are the slow-roll parameters also used previously in the minimally coupled scalar field and the two new parameters ϵ 3 and ϵ 4 were added in light of the additional functional degree of freedom f (ϕ) introduced in (81) for this case. By taking the slow-roll assumption holds true and the slow-roll condition that ϵ i ≪ 1, i = 1, 2, 3, 4, then the observational quantities can be expressed with respect to these parameters as [96], where Q s is defined as the function of, Also from the imposed slow-roll condition to the parameters, the equations of motion from (82), (83) and (84) take the following form, Taking the slow-roll condition and (88), (89) and (90) into consideration, the parameter Q s can also be approximated as, therefore from (86), the tensor-to-scalar ratio r and the scalar spectral index n s during the slow-roll era are, Now since the observational quantities of r and n S are expressed with the forms of (92) and (93) respectively, specifically under the slow-roll condition, we are able to analytically compute them for the duration of the slow-roll era for any given function f (ϕ) [96]. For the impact of the non-minimal coupling to the quantity of tensor spectral index n T , included also in the context of string corrections.

Example for Specific Form of the Function f (ϕ)
To implement the above formalism, a specific example, which can also be found in [166], for the form of the function f (ϕ) is considered as follows [166], where β is a constant. This form of f (ϕ) has a special symmetry and β → 1 β for β > 1 and for further simplification ξ = 1, f (ϕ) can be approximated to the form of, Additionally, the most simple form for the potential V (ϕ) is assumed with V (ϕ) = Λ, where Λ is a positive constant parameter. Thus by considering the slow-roll approximation and κ 2 = 1 for simplicity, from (89) we can derive the expression forφ as,φ and also forḟ as,ḟ Thus by taking the formulas for the slow-roll parameters from (85) and considering the relations of (96) and (97), it is determined that, Additionally, by substituting these relations into (93) and (92), the scalar spectral index n S and the tensor-to-scalar ratio r can be computed as, Now by using the integral (59), the number of e-folds N can be determined as, where ϕ end is the value of the field when inflation ends and ϕ k at horizon crossing. Also for the result of (102) with respect to ϕ, the approximation of ϕ k ≫ ϕ end was used, which is justified for the duration of the slow-roll era. So lastly, the observable quantities of n S and r take the following forms with respect to the e-folds number N , An interesting result of this case, is when we select the value n = 2 √ 3 , in which case the observational quantities become n S = 1 − 2 N and r ≃ 12β 2 N 2 , which is the same relations derived by α-attractor models [167]. The difference is that in the α-attractor models β ≪ 1 and so the same attractor behavior can be exhibited by a non-minimal theory with correctly chosen parameters. It is also worth noting that here we have a different formalism than that of a strongly coupled non-minimal theory and in that case the parameter ξ is large and f (ϕ) is chosen under a different criteria.

IV. A BRIEF ACCOUNT OF THE SWAMPLAND CRITERIA
As a sidenote, we shall briefly discuss, but shall not cover explicitly in every example, the completeness of models through the prism of the Swampland criteria. The interested reader may check for example [168,169]. for a wide range of applications. In short, the gravitational action that one may chose to work with can in principle be regarded as a low-energy effective model. In order to distinguish theories based on their UV completeness and therefore ascertain whether they serve indeed as effective modes or not, a set of criteria can be investigated. Let us showcase them explicitly.
The first criterion is the Swampland distance conjecture, This condition states that the field range of the scalar field must not be arbitrary during inflation but in principle is smaller than or equal to the Planck mass. As shown, the condition is the same irrespective of the sign therefore the scalar field could increase in value as time flows by. The second criterion is the de Sitter conjecture, It is applied at the start of inflation and suggests that the slope of the scalar potential has a lower bound. The same criterion can be written in a different form as − V ′′ κ 2 V ≥ O(1). This in turn implies that for a positive scalar potential, its form is specific with −V ′′ < 0 and it also has a lower bound. These conditions can be applied to several models. The reader should also keep in mind that the aforementioned criteria can in principle be satisfied as separate conditions and not simultaneously. In fact, if one can ensure that a single criterion is satisfied then in consequence the model belongs to the Swampland and serves an effective model which is UV incomplete in the high energy regime. A characteristic example is the power-law model where it was shown that while the de Sitter conjecture is indeed satisfied, the equivalent condition is not. A similar example will be covered subsequently in the following sections however the Swampland criteria shall not be covered in this review but are nonetheless mentioned here for the sake of completeness.

V. EVADING THE SLOW-ROLL EVOLUTION: THE CONSTANT-ROLL EVOLUTION
Here we shall briefly discuss a different approach on the scalar field evolution that has interesting phenomenological implications for the inflationary era. Previously, it was shown that under the slow-roll assumption, the scalar field evolves slowly and therefore issues like the apparent flatness and the horizon problem can be explained properly. In this approach, for potential driven inflation, it is shown that the necessary condition is the dominance of the scalar potential over the kinetic term, i.e. 1 2φ 2 ≪ V . This condition, along with the continuity equation, can be used in order to derive the additional slow-roll conditionḢ ≪ H 2 andφ ≪ Hφ, where it should be stated that these inequalities are indicative of the order of magnitude of the respective object and not its sign. These two conditions are not postulated but are derived from the assumption, or from a different perspective, the necessity of the dominance of the scalar potential. In this approach, the scalar field is said to slowly evolve with respect to time, or the e-foldings number.
Another assumption that can be made about the dynamics of the scalar field is known as the constant roll condition. In this case, the scalar field evolves approximately under the condition, where β is an auxiliary dimensionless parameter which is not necessarily constant for a ϕ dependent constant roll condition. This evolution rate can be used as an approximation during the inflationary era and in principle can be used along with the slow-roll condition for β ≪ 1, however it is not necessary. One of the advantages of the constant-roll condition is that the contribution of the second order derivative of the scalar field is now considered in the continuity equation of the scalar field however now the degrees of freedom are increased. In addition, the value of the second slow-roll index ϵ 2 =φ Hφ is specified completely by the aforementioned parameter and it is constant in the constant-roll case. The evolution of the scalar field under the constant-roll condition could be a leading factor in the production of primordial scalar non-Gaussianities in the CMB, however it is not the only factor that leaves a non-Gaussian imprints, as subsequent cosmological eras could leave an imprint as well. In short, anisotropies in the CMB are not only expected, but they are also quantified by the scalar spectral index. Such anisotropies must obey a specific distribution pattern. Information about such pattern can be derived by examining curvature perturbations of uniform density hypersurfaces. Theoretically, the origin of curvature perturbations can be quantum fluctuations in cold inflationary models or thermal fluctuations in warm inflationary models, if not both, on superhorizon scales. Regardless of their origin, the distribution pattern of CMB anisotropies can be studied through the bispectrum, meaning the Fourier transform of the three-point correlation function. As it is known, a Gaussian distribution implies that even correlations can be written as combinations of lower but nonetheless even correlations, i.e. the four point correlation can be written through combinations of the two point correlation function and so on. It is therefore expected that the bispectrum is zero in usual Gaussian distributions. Since secondary anisotropies are not always linear in physical systems, one can introduce a nonlinear parameter f N L in order to quantify the deviation of CMB anisotropies from a Gaussian distribution. In principle, the distribution pattern, or equivalently the numerical value of such nonlinear parameter, differs if the observer used different wavelengths in order to perform the measurement, for instance the equilateral nonlinear term f eq N L , in which the wavelengths in momentum space are equal. In consequence, the three point correlation function is dominated by the scalar field dynamics primordially and is connected to the numerical value of auxiliary parameters such as the slow-roll indices during the first horizon crossing. In the literature, there exist several studies that analyze such patterns for several scalar-tensor models of gravity. The main result is that the amount of non-Gaussianities in the CMB, or in other words the deviation of the distribution of the CMB anisotropies from a Gaussian distribution, is negligible as f N L ∼ O(10 −2 ) and only more involved models that predict a propagation velocity of scalar perturbations that clearly deviates from the speed of light can produce a larger value, closer to the upper bounds currently available. Obviously, the same applies to tensor perturbations as well where information can be extracted by examining the three-point correlation function for gravitons.

VI. STRING INSPIRED MODELS OF GRAVITY
In this section we expand on the previously presented canonical scalar field theory by including additional terms related to the scalar field, originating from string corrections of the scalar field Lagrangian. In general, the four dimensional scalar field Lagrangian which contains at most two derivatives has the following form, Note that when the scalar fields are considered in their vacuum configuration, the scalar field has to be either conformally or minimally coupled. When quantum corrections of the local effective action are considered, with the quantum corrections being consistent with diffeomorphism invariance of the action and also contain up to fourth order derivatives, the scalar field action is generalized to [170], with the parameters Λ i , i = 1, 2, ..., 6 being appropriate dimensionful constants. For the purposes of this section, the gravitational action of the model is defined as [117], , similar to the non-minimal case, is an arbitrary dimensionless function depending on the scalar field ϕ while ξ(ϕ) is an arbitrary function of the scalar field with parameters c i 's being auxiliary parameters with mass dimensions of eV −i+1 for the sake of consistency. For generality, an additional dimensionless parameter ω is introduced so that one can distinguish between the canonical (ω = 1) and phantom (ω = −1) case however for the time being it shall be considered as a dynamical variable depending solely on the scalar field. This is the most general string inspired model that can be introduced where for simplicity the same coupling function ξ(ϕ) is considered however this is not mandatory. Indeed, the user may feel free to change the coupling function that accompanies each c i factor in subsequent computations. Referring to the contributions themselves, the first term describes the Gauss-Bonnet density G which serves as a non minimal coupling between curvature and the scalar field. The introduction of the coupling ξ(ϕ) is in fact quite important because, due to the nature of the Gauss-Bonnet density, it does not participate in the background equations as a total derivative if it is introduced linearly in the gravitational action, therefore the coupling function ξ(ϕ) keeps in place the Gauss-Bonnet density. Of course in the literature there exist several extensions that do not require an arbitrary coupling such as f (G) gravity, and even a linear model αG in D dimensions which, upon rescaling the auxiliary parameter α as α → α D−4 and taking the limit D → 4, the Gauss-Bonnet density indeed participates in the equations of motion, however we shall not consider these examples in this brief review. The Gauss-Bonnet model is commonly known as a low-energy effective string model, and in an essence, when introduced, it affects accordingly not only the background equations but also the behavior of scalar and tensor perturbations respectively for as long as the scalar coupling function evolves dynamically. Now the second term that is introduced in the gravitational action (109) is the kinetic coupling. As the name stands, it serves as a coupling between curvature and the kinetic term of the scalar field and serves as an effective corrective term. It should be stated that the inclusion of only these two terms manages to affect the propagation velocity of tensor perturbations, therefore the model may be at variance with recent observations such as the GW170817 event, however as we shall showcase subsequently, there exists a way that the model can be rectified. The kinetic coupling manages to effective shift the contribution of the kinetic term of the scalar field and in principle it does not require a dynamical coupling function in front, in order to participate but is nonetheless introduced for the sake of generality. In fact, as long as the scalar field evolves dynamically then the kinetic coupling has an active role on the field equations. The third contribution in (109) is the Galilean model [171,172] which is a type of higher order coupling between the scalar field and its kinetic term and finally the last term can be interpreted as a coupling between the scalar field with the square of its kinetic term. These corrections refer to the kinetic term as well, however they are more involved as they are more intricate. The Galilean term serves as a non-minimal coupling between first and second order variations of the scalar field and are introduced in a nonlinear way while the final term, reminiscing of the k-essence models, serves exactly as a higher power of the kinetic term and can be treated either as an important contribution or higher order corrections based on the occasion. The reason why this action is considered to be a general case is due to the fact that all these additions in the action, even though they result in several inclusions in the background equations which are not only nontrivial but also nonlinear, the continuity equation of the scalar field still remains a second order differential equation. Let us show this explicitly by working on the background equations. By performing a similar work to the previous sections, one can easily see that the field equations for gravity read [117], where as usual the energy-stress tensor due to the presence of string corrections is defined as T where L string is the Lagrangian density of the ξ dependent part in action (109). Also, the continuity equation of the scalar field reads, which as mentioned before is a second order differential equation with respect to the scalar field. In principle, the inclusion of additional scalar terms, either minimally or non-minimally coupled to curvature, is not forbidden however the model at hand is the most general case that yields second order nonlinear differential equations. Now in this context, the contribution of the string correction terms included in (109) is quite lengthy and is showcased below, where for simplicity, ϕ ,µ = ∇ µ ϕ, while the contribution in the continuity equation is, Here, it becomes abundantly clear that their contribution is important only if the scalar field evolves dynamically, therefore de Sitter solutions are not affected by the inclusion of string corrections. At this point it should also be stated that the above expressions are valid even for the case of a nonzero spatial curvature K however hereafter we shall limit our work to only the flat case for the sake of simplicity. The generalization to nonzero spatial curvature is relatively straightforward at the level of background equations however it becomes tedious when linear perturbations are considered. Let us now focus on the equations of motion for the case of vanishing curvature. In this case, the temporal and spatial components of the field equations are written as [117], where according to the results of (112), whereas (113) is identically equal to, In this approach, a homogeneous scalar field was once again considered. Hence, it becomes clear that the inclusion of additional string corrective terms, apart from introducing a new degree of freedom if the non-minimal coupling function ξ(ϕ) is indeed dynamical, something which is mandatory for the Gauss-Bonnet model at the very least, it results in the appearance of several terms that evolve dynamically with respect to the scalar field. In particular, each c i factor generates a term proportional toφ i therefore string corrective terms seem to introduce corrections to the kinetic term of the canonical scalar field which in the context of higher order gravity are well motivated. Their contribution as showcased appears in a highly nonlinear manner as nowφ is also coupled to not only ξ and its derivatives but also appropriate powers ofφ. Obviously the same thing applies to the case of higher powers of the kinetic term, not just quadratic, that may be inserted in the gravitational action. The same term also results in the appearance of the second time derivative of the scalar field in continuity equation, as stated by (118) hence the reason why they were selected. In the literature, their contribution has been thoroughly investigated, mainly in separated models however it is not obligatory to consider only one correction at a time. As an example, one could consider the k-essence model that was presented previously, with or without the inclusion of the scalar potential and combine it with a kinetic coupling, meaning that c 2 is the only nonzero parameter in (109). In this approach, the kinetic coupling introduces an additionalφ 2 in the background equations therefore phenomenologically speaking, it acts as a shift in the kinetic term. String corrections, although they have a similar behavior in the continuity equation of the scalar field, meaning that they actively affect the evolution of the scalar field, their contribution is quite different at the level of perturbations as we shall showcase explicitly.
Let us now see how the inflationary era is described in this context. In order to do so, we shall follow similar steps as in the previous sections given that the previous results are obviously subcases. Firstly, the slow-roll indices are defined as [117], where in this case a total of five indices are introduced. Here, it should be stated that each index describes a specific part. In particular, the first slow-roll index is specified by all components of the gravitational action (109) and is completely specified by the background equations (114)- (115). It is mainly used in order to quantify the duration of inflation and in consequence extract information that is crucial for the computation of the observables. The second index can be extracted from the continuity equation of the scalar field (111) and, while it is affected by a possible nonminimal coupling f (ϕ), it mostly captures the contribution of the potential and the string corrective terms however every scalar coupling function participates in principle. In an essence, it carries information about the dynamical evolution of the scalar field which is important for the study of scalar perturbations. The third index addresses the contribution of the non-minimal coupling, if existent, while the last two indices contain both contributions. Index ϵ 3 , similar to the previous non-minimal case, is indicative of the dynamical evolution of the respective non-minimal scalar coupling function and is used in both scalar and tensor perturbations. In addition, ϵ 4 shows the apparent dominance of the scalar functions of the model in a compact manner and is mainly used in order to address scalar perturbations but in principle it can be written with respect to the rest indices. Finally the last index in indicative of the running Planck mass, where Q t is the modified Planck mass squared due to the presence of string corrections. It should be stated that the fifth index is identical to the third in the limit of vanishing string corrections therefore ϵ 3 can carry information about the running of the Planck mass but only for the non-minimal coupling. The auxiliary parameters in this context read, which can be used in order to study scalar and tensor perturbations. As shown, dynamical evolution of the scalar field is crucial otherwise a finite value, at best nonzero, is obtained and it thus results in vanishing indices. Note that the previously defined indices, similar to the quite simpler canonical scalar field model, can be specified and in the end, both the scalar and tensor spectral indices as well as the tensor to scalar ratio can be extracted. Before we proceed however with the specification of such indices, let us study briefly the impact that string corrections have on scalar and tensor perturbations. Firstly, consider that the FRW metric is perturbed as, where α(η) is the cosmic scale factor depending on conformal time η specified by the condition c dt = αdη. Here, the conformal time is used for the sake of simplicity. Furthermore, φ,α, β and γ describe scalar perturbations whereas h ij tensor perturbations and in particular traceless g ij h ij = 0 and transverse, ∂ i h ij = 0. The latter is connected to gravitational waves and in the presence of certain string corrections, the propagation velocity of tensor perturbations is affected as we shall showcase subsequently. Typically, α is the lapse function which specifies the connection between proper time τ and conformal time η while β ,j is the shift function, indicative of the apparent velocity between the threading and the worldlines which are orthogonal to the chosen slicing. By studying scalar perturbations using the above perturbed metric, one can see that scalar modes propagate with a nontrivial velocity [117], which obviously differs from the sound wave velocity of a perfect fluid c 2 s =Ṗρ . This is the most general expression that the sound wave velocity has and as shown, the inclusion of string corrections in the gravitational action has a major impact on the propagation velocity of scalar perturbations. In the end, the above expression should be well behaved for a viable inflationary model, meaning that its numerical value during the first horizon crossing where modes start becoming superhorizon should be in the range 0 < c S ≤ c. The fact that it is strictly real suggests that no ghost instabilities are initially present while the upper bound ensures that the model is in agreement with causality. Physically speaking, while the sound wave velocity can be equal to the speed of light, consider for instance the canonical scalar field that was presented initially or a non-minimal coupling between the scalar field and the Ricci scalar, the lower-bound is purely mathematical since in a sense it implies no propagation. Now for the tensor mode, once can show that it is affected by string corrections as, where a M =Q t HQt is an auxiliary dimensionless parameter which specifies the running Planck mass [173] and is introduced for simplicity since in subsequent sections a brief comment on the energy spectrum of primordial gravitational waves shall be made. In this case it becomes apparent that the non-minimal coupling f (ϕ) and string corrections have a major impact on the behavior of tensor perturbations. In particular, the Gauss-Bonnet density ξ(ϕ)G and the kinetic coupling ξ(ϕ)G µν ∇ µ ϕ∇ ν ϕ are the only string corrective terms in the gravitational action (109) that influence tensor perturbations. This can be seen not only from the definition of a M but also from the propagation velocity of tensor perturbations which is defined as, where similar to the case of the sound wave velocity c S previously discussed, it must also satisfy the relation 0 < c T ≤ 1. This is a really interesting result and is worth discussing further. Recently, the GW170817 event which was a signal from the merging of two Neutron Stars (NS) has made it abundantly clear that gravitational waves propagate through spacetime with the velocity of light. As a result, the above subclasses are in peril since they predict a propagation velocity for tensor perturbations which deviates from that of light's. However, imposing the condition Q f = 0 restores compatibility and the models can in fact be salvaged. This condition imposes a rather strong constraint on the Gauss-Bonnet scalar coupling function and has a significant implication on the inflationary phenomenology. For instance, the degrees of freedom are in fact decreased and one scalar coupling function can be extracted from the continuity equation (111) once the rest have been specified. In turn, the continuity equation can be treated either as a first order differential equation with respect to the scalar potential or a second order differential equation with respect to the coupling function ξ(ϕ), this will become more clear in the following. Note also that this constraint is discussed as a plausible scenario even in the early era due to the fact that a primordial propagation velocity c T which deviates from the speed of light predicts massive gravitons which should not be the case for this subclass of models since they can be considered as low-energy effective string theory models. Nevertheless, for the sake of generality, several results shall be presented where the aforementioned constraint is either implemented or not. The distinction will be clear. At this stage, the observed indices can be computed by making use of (119). It should be stated that these indices are not necessarily slow-roll indices, meaning that it is not mandatory for their values to be approximately of order O(10 −3 ) and lesser even if the slow-roll assumption is imposed but their combination should be quite small in order to obtain results compatible with the latest (2018) Planck data. The only conditions that need to be respected are, for the extraction of the scalar and tensor spectral index respectively. These conditions do not require that (119) are actually slow-roll conditions in order to obtain a viable inflationary era, where the aforementioned indices are of order O(1) and greater. In addition, in order to extract the expressions for the spectral indices, it is common to use the conditionε i = 0 however in Ref. [174] it was shown that the same results are indeed extracted for slow-varying variables, that isε i Hϵi ≤ O(10 −3 ) approximately. In the end, the spectral indices and the tensor-to-scalar ratio obtain the following forms [117], where everything is written as a function of the auxiliary variables. It should be stated that this is the most inclusive expression in terms of string corrections as it contains all the known expressions in the literature for the observed quantities. As an example, for the case of a minimally coupled canonical scalar field, only indices ϵ 1 and ϵ 2 are needed, which are showed are connected to the usual indices ϵ and η and therefore n S = 1 − 6ϵ + 2η, n T = −2ϵ and r = 16|ϵ| as shown before. As another example, consider the k-essence models obtained for c i = 0 for i = 1, 2, 3 and ξ(ϕ)c 4 = c4 2 , for which after a few calculations, it becomes clear that the tensor to scalar ratio as usual reads r = 16c S ϵ 1 . These forms also showcase the previous statements made on the influence of indices ϵ i on perturbations as now indices ϵ 2 -ϵ 4 affect scalar perturbations while ϵ 5 tensor perturbations. On the other hand, all auxiliary parameters needed participate in the tensor-to-scalar ratio. Its form is general however under the slow-roll assumption one can obtain a more simplified expression. At this point, it should also be stated that the analysis of scalar and tensor perturbations suggests that the respective power-spectra obtain the following forms, where ε = ϵ 3 f κ 2 Qt + Qa 2HQt , ν = 4−n S 2 and ν T = 3−n T 2 . These are the expressions from which the spectral indices are essentially derived however they are presented here not only for the sake of generality but also in order to impose constraints on the free parameters of given models given that the amplitude of scalar curvature perturbations A S is known. Since string corrections affect tensor perturbations as well, it stands to reason that the amplitude of tensor perturbations A T could also impose further constraints provided that it is observed in the future.
Here, it is worth making a comment on the impact of string corrections and in particular on couplings between gravity and the scalar field. Previously, it was shown that the Gauss-Bonnet term ξ(ϕ)G and the kinetic coupling ξ(ϕ)G µν ∇ µ ϕ∇ ν ϕ manage to affect the propagation velocity of tensor perturbations. Recall that the latter does not require a dynamical coupling ξ(ϕ), in contrast to the first. These string corrections are the only additions in Eq. (109) that result in a nontrivial running Planck mass a M =Q t HQt , something which can easily be seen from the fact that the slow-roll index ϵ 5 does not coincide with ϵ 3 in those cases. In turn, the tensor spectral index is affected by string corrections not only indirectly through the first slow-roll index ϵ 1 , but also directly thought the additional index. This realization has interesting phenomenological implication on string inspired models since a blue spectrum can easily be manifested provided that the condition ϵ 5 < −ϵ 1 is now satisfied. This is an rather interesting statement since the canonical scalar field model is incapable of producing a blue-tilted tensor spectral index however the inclusion of well-motivated string corrections can in fact affect tensor perturbations to quite an extend. Indeed, in Ref. [158] it was shown that constrained Gauss-Bonnet models that satisfy the conditionξ = Hξ can manifest a blue spectrum if certain criteria are met. As a result, the blue spectrum leaves a major impact on the energy spectrum of gravitational waves and in particular, in the high frequency regime. This is because high frequency modes are the first to reenter the horizon, therefore they become sub-horizon in the early Universe and can thus carry information about such mysterious eras. The fact that a blue-tilted spectrum is a possibility implies that the energy spectrum of high frequency modes is amplified, compared to the general relativistic description, and therefore it becomes easier to spot them. Indeed, this will be shown explicitly in subsequent sections but it is worth mentioning here since string corrections can leave a major imprint on the energy spectrum of primordial gravitational waves.
Let us now see a few examples of string inspired models that have been presented in the literature. For simplicity, in order to avoid presenting a plethora of models, we shall restrict the analysis on three models, one unconstrained and two constrained models, with the first two revolving around the Einstein-Gauss-Bonnet model while the last containing in addition a kinetic coupling. This does not mean however that the rest string corrections are not important of course.

A. Gauss Bonnet Model
For the first example, let us consider the dilaton model presented in [137], where f is an auxiliary parameter with mass dimensions of eV whereas V 0 and ξ 0 serve as the potential amplitude and the Gauss-Bonnet scalar coupling function amplitude respectively with the former having mass dimensions of eV 4 while the latter is dimensionless for the sake of consistency. In this approach, both functions are exponential however they behave differently in the limit ϕ f ≪ 1 since in this case the potential vanishes while the Gauss-Bonnet coupling dominates. Note that the product V (ϕ)ξ(ϕ) = V 0 ξ 0 has a constant value, that is V (ϕ)ξ(ϕ) = V 0 ξ 0 and is specified completely by the amplitudes. These types of models in which the Gauss-Bonnet scales inverse with the potential amplitude are used in order to unify early and late-time era, therefore it is interesting to present the inflationary phenomenology. Note that in this approach the propagation velocity of tensor perturbations is assumed to be unconstrained. Another reason that makes the exponential potential an interesting model is the fact that the first two slow-roll index in Eq. (119) is now ϕ dependent in contrast to the canonical scalar field case where ϵ 1 and ϵ 2 are constant and cannot describe a graceful exit from inflation. For the case at hand, ϵ 4 can be used as an index that quantifies inflation if the slow-roll assumption is imposed therefore the inflationary era can be properly studied. In consequence, by using index ϵ 4 , one can speculate that the inflationary era ceases when the condition ϵ 4 = 1 is satisfied. In this case, apart from the slow-roll conditions, the following assumptions are also made, which are motivated by the slow evolution of the scalar fieldφ ≪ Hφ. In principle, the slow-roll condition is not mandatory in order to study the inflationary era however a slow-varying scalar field, as mentioned before, manages to solve both the flatness and horizon issues as a large duration can be obtained. As a result, the scalar spectral index and the tensor-to-scalar ratio are given by the following expressions,

B. Constrained Einstein-Gauss-Bonnet Model
For the second example we shall consider that the propagation velocity of tensor perturbations is in fact identical to that of light's. This condition can easily be satisfied if the Gauss-Bonnet scalar coupling function satisfies the differential equation,ξ = Hξ , regardless of the cosmological era that is studied. Therefore, even in primordial eras, a concrete description that predicts primordial massless gravitons can be realized. Now for simplicity we shall assume that the scalar coupling functions of the model are given by the following expressions, where f 0 is an auxiliary dimensionless parameter, n is the power-law exponent which is not necessarily integer and φ is a free parameter of the model with mass dimensions of [φ] =eV. These models were initially considered for two reasons. Firstly, due to their respective forms, one can easily see that f ′ f = n ϕ , f ′′ f = n(n−1) ϕ 2 and ξ ′′ ξ ′ = − 1 φ therefore the inflationary phenomenology should be quite straightforward since the aforementioned rations participate in the slowroll indices. The second reason that this model is considered is because the exponential Gauss-Bonnet scalar coupling function on its own was proven to be inconsistent with observations under the assumption (132) since, depending on the magnitude of φ, it may lead to eternal or no inflation [157]. As a result, it is intriguing to examine whether it can actually become viable if it is paired with the simplest choice of a Ricci scalar coupling. Furthermore, the constraint in the propagation velocity of primordial tensor perturbations manages to decrease the degrees of freedom. In particular, assuming that the second slow-roll index ϵ 2 is actually quite lesser than unity, Eq. (132) yields, hence the reason why the exponential Gauss-Bonnet scalar coupling function is regarded as a convenient replacement, sinceφ and in consequence every time derivative is simplified to a great extend. Note that Eq. (134) is valid only for ϵ 2 ≪ 1. In turn, the ϕ dependence ofφ is provided solely from the Hubble rate expansion and since we assume a potential driven inflation, the extraction of the potential is provided by the continuity equation of the scalar field (111). As a result, the scalar potential reads, with [V 0 ] =eV 4 being the potential amplitude. As shown, the potential is now a combination of a power-law and an exponential function. For a specific exponent n and small exponent values, it stands to reason that the potential is in fact the linear combination of power-laws however their exponents is not necessarily integer as well. Non integer exponents are not a new inclusion but they have been considered in cosmology. Now choosing the values φ = 0.01M P , n = 1 2 and f 0 = 100 for N = 60, implies that the exponential model accompanied by a power-law non minimal coupling is in fact a viable inflationary model as now n S = 0.967045, r = 0.000551 and n T = 0.000069 which are in agreement with the latest observations. Interestingly enough, the inclusion of string corrections can in fact result in the manifestation of a blue tilted tensor spectral index. Obviously, this is not limited to the non-minimal case, see for instance Ref. [157] where it becomes clear that having a positive parameter λ < 1 where λ = 4κ 2 ξ ′′ V 3 suffices. The results also suggest that during the first horizon crossing, ϕ k = 0.6025M P and therefore the exponential factor in the potential can be expanded. This leads in the inclusion of additional power-law factors such as a linear and ϕ 3 2 . This is not limited to the above designation of free parameters and in fact several other choices can lead to a viable inflationary era.

C. Constrained Viable Horndeski Theories
For the final model we shall consider that the action for the model contains both terms that affect the propagation velocity of tensor perturbations as shown below, where c is an auxiliary parameter with mass dimensions of eV −2 . Here, we consider both terms from Eq. (112) that affect the propagation velocity of primordial tensor perturbations, see equations (120) and (124). In principle one can include only the kinetic coupling however the theory cannot become compatible with the GW170817 event in that case. For the sake of consistency we shall assume that Q f = 0 in this case. This in turn implies that the differential equation thatφ must satisfy is given by the following expression, where upon implementing the constant-roll condition ϵ 2 = β, it generates the following solution, Obviously, in the limit c, β → 0 it becomes abundantly clear that the previous form in Eq. (134) is recovered, as it should. Here we shall consider the constant-roll condition in order to avoid presenting the same phenomenology over and over again. Let us also assume that the Gauss-Bonnet scalar coupling function has a linear form, that is, where [f ] =eV. This choice is extremely convenient as now ξ ′′ = 0 and thus the time derivative of the scalar field scales asφ ∼ 1 ϕ . In turn, the continuity equation of the scalar field is simplified and the resulting expression of the scalar potential reads, where similar to the previous case, [V 0 ] =eV 4 . Here, the potential seems to have a power-law form with a fixed exponent, which once again is not necessarily an integer. In particular, having c = − 4 M 2 P , β = 0.009, N = 60 and f = 10 8 M P , which is a set of parameters that produces results that are compatible with observations, suggests that the exponent of the potential is approximately 1. Regarding the results, we report that the spectral indices are n S = 0.968604 and n T = −0.00666 while the tensor to scalar ratio has the value of r = 0.0531572. Therefore the model predicts a red spectrum and is in agreement with observations. In this approach, it becomes clear that the potential has an almost linear form, that is V ∼ ϕ. In other words, this model belongs to the category ξ(ϕ) = λκ 4 V .

VII. CHERN-SIMONS AXIONIC GRAVITY
In this section of the review we shall consider a really interesting model that makes quite unique predictions for the inflationary dynamics. Let us introduce a gravitational Chern-Simons model of the form, where RR = η µνρσ R αβ µν R ρσαβ is the Chern-Pontryagin density and for the sake of generality, an arbitrary coupling between the scalar field and the Ricci scalar is introduced. In an essence, it is reminiscing of the product F * µν F µν which appears in fibre bundles and is frequently named a Chern-Simons term because it is connected to a 3D Chern-Simons term cohomologically, that is ν(ϕ)RR = d(Chern − Simons), therefore the above characterization shall be used. Both scalar coupling functions f (ϕ) and ν(ϕ) in this scenario as treated as dimensionless functions. It should be stated that while the inclusion of f (ϕ) is not mandatory and was considered only for illustrative purposes, the Chern-Simons scalar coupling function ν(ϕ) is of paramount importance, similar to the case of the Gauss-Bonnet coupling. Here however, the Chern-Simons term does not vanish identically in D = 4 only but due to the fact that it is a parity odd term, its contribution on the energy momentum tensor vanishes as a whole. Truthfully, even under the assumption that an arbitrary coupling between the scalar field and the gravitational Chern-Simons term exists, the background equations of motion remain unaffected and thus the equations of motion remain exactly the same as in the case of the non-minimally coupled canonical scalar field. This can easily be ascertained from the definition of the energy stress tensor T which in this case the variation of the aforementioned contribution yields [117], where it becomes abundantly clear that indeed if the scalar coupling function ν(ϕ) is not dynamical then it does not contribute. Note also that in this approach the field equations obtain the following form, For the sake of simplicity, let us assume that the spatial curvature is equal to zero similar to the previous cases. In consequence, by taking into consideration that the metric corresponds to a flat, isotropic and homogeneous background as recalling that the scalar field is assumed to be time dependent only, then the energy-stress tensor obtains the following form, where recall that h ij (t, x) describes traceless and transverse tensor modes (121) in the 3D surface and −□h ij = h ij + 3Hḣ ij − ∇ 2 a 2 h ij due to the FRW metric. Now owning to the fact that the totally antisymmetric Levi-Civita tensor in 3D ϵ ijk is present, it can easily be inferred that the aforementioned tensor does not contribute in the background equations of motion since the diagonal parts are needed for the Friedmann and Raychaudhuri equations and of course the same applies to the continuity equation of the scalar field. This is a really interesting result since the inclusion of a new scalar coupling function, namely ν(ϕ), does not influence the evolution of the scalar field and in a sense serves as a new degree of freedom. Of course, it could be possible to implement certain constraints on the free parameters of the scalar coupling function ν(ϕ) if the amplitude of tensor perturbations A T (k) were to be computed at the pivot scale, or in other words if the tensor spectral index was computed, similar to how the potential amplitude in potential driven inflation is constrained from the amplitude of scalar curvature perturbations A S . In the end, the equations of motion coincide with the non-minimal case [117], where similar to the previous models, the prime denotes differentiation with respect to the scalar field. The inclusion of the gravitational Chern-Simons term however has a major impact on tensor perturbations. In this case, by assuming that the background is perturbed as, one can shown that the tensor mode must satisfy the following equation in configuration space [117], which as shown, it is greatly affected by the Chern-Simons scalar coupling function. This can be shown explicitly by performing a Fourier expansion, in order to move to momentum space and thus one can easily see that Eq. (149) obtains the following form, where Q t,l = f κ 2 + 2λ lν k a is an auxiliary parameter denoting the shifted squared Planck mass but obviously differing from the previous string corrections definition. Here, e (l) ij is the circular polarization of tensor perturbations which in this approach is either left or right handed. Furthermore, h lk (t) describes tensor perturbations in momentum space while c T which stands for the propagation velocity of tensor perturbations is identically equal to the speed of light therefore the model is in agreement with the GW170817 event. Finally, parameter λ l is used in order to distinguish between the two different polarization states with λ L = −1 and λ R = +1. Hence, in the gravitational Chern-Simons model, tensor perturbations satisfy different differential equations in momentum space depending on their polarization state and they are essentially distinguished based on their chirality [175]. This parity-odd term can also be used in order to study gravitational leptogenesis. As a result, the different behavior of modes based on their chirality has a major impact on the tensor spectral index and the tensor-to-scalar ratio. In particular, it turns out that [117], where, with E = f κ 2 1 + 3ḟ 2 2κ 2 fφ 2 . It should be stated that the factor of 1 2 in front of ϵ 5 is for averaging the contribution of the two polarization states. It is worth making a comment on the expression of the tensor spectral index n T . Due to the fact that index ϵ 5 participates, a blue tilted tensor spectral index can be obtained only if ϵ 5 < −ϵ 1 which can result in an amplification of the energy spectrum of primordial gravitational waves.
Before we proceed with any examples, it is worth making two comments on the Chern-Simons model. Firstly, as demonstrated before, the parity odd terms breaks chirality and therefore modes with different circular polarization have a different behavior. This difference is of course model dependent and is connected to the Chern-Simons scalar coupling functions and the apparent dominance of the factor 2ν k a over the Planck mass squared, recall the new expression for Q t . This realization, combined with the fact that a blue spectrum can in principle be manifested seems to have interesting phenomenological implications. Consider for instance high frequency modes that become sub-horizon after inflation. Since the scalar field evolves dynamically with respect to cosmic time, then chirality is broken and therefore the difference between circular polarization states should leave an imprint on the energy spectrum of primordial gravitational waves. The fact that a blue-tilted tensor spectral index is obtained simply implies that the amplitude is enhanced compared to the GR predictions therefore it is easier, in principle, to observe it however the same distinction is present even if the spectrum is red, as most scalar-tensor theories predict. Obviously, low frequency modes should not be affected since chirality is restored provided that the scalar field has reached its vacuum expectation value. Therefore, for the Chern-Simons model and in the high frequency regime, one should expect for a given frequency, two signals due to the fact that chirality is broken and this is quite different compared to the previous string inspired models. It can be shown that even though the tensor spectral index is one, since it is averaged over the circular polarizations exactly as the tensor-to-scalar ratio, chirality is still broken and since different evolution laws must be satisfied based on the circular polarization, the running Planck mass a M,l =Q t,l HQ t,l differs between modes meaning that the enhancement factor, which we shall showcase explicitly in the following sections, also differs.
Another issue that one should keep in mind when working on Chern-Simons models is the possibility of producing ghost instabilities. In principle, the effective Chern-Simons mass scale is a parameter which can be used in order to extract information about ghosts. The aforementioned mass is defined as, and to no ones surprise is connected to the arbitrary scalar coupling function ν(ϕ), or more specifically its dynamical evolution. An observant reader might notice that the denominator appears in the definition of the auxiliary parameter Q t,l which is indeed the case and serves as an effective shift to the Planck mass, depending on the polarization of the mode. Primordially, the scalar field evolves with respect to cosmic time, assuming that it is only homogeneous, therefore the Chern-Simons mass scale has a finite value. As the Universe expands and subsequently cools down, reaching equilibrium, the scalar field itself decreases in magnitude, compared to the near Planck scale values that it originally had, and tries to reach its vacuum expectation value. When the scalar field reaches its minimum value, it no longer evolves dynamically with respect to time thereforeν vanishes identically. Note that quantum fluctuations around the vacuum expectation value are still present however this is different from dynamical evolution. Hence, since the vacuum expectation value has been reached, the denominator in Eq. (154) tends to zero and in consequence the effective Chern-Simons mass scale explodes. Phenomenologically speaking, this means that the impact of the parity odd term on the model is lifted, therefore one obtains the GR description, or any other modified expression that was accompanied by the gravitational Chern-Simons term. This is in agreement with the previous statement about low frequency modes on the energy spectrum of gravitational waves. In order to effectively avoid the appearance of ghost instabilities, once should make sure that the numerical value of the Chern-Simons mass scale is in fact greater than a given threshold during the inflationary era. Subsequent eras are also affected but provided thatφ decreases in magnitude with respect to time until it reaches zero, the above mass scale increases in turn until it reaches infinity. Any viable models should respect these thresholds otherwise it is in peril as ghost instabilities are manifested. Let us now focus on a toy models to see what is the impact of the Chern-Simons term in the inflationary era that we are interested in. The simplest example and most interesting from a phenomenological point of view is the case of chaotic inflation. Let us assume that in the gravitational action (141) is described by f (ϕ) = α and a power-law potential V (ϕ) = V 0 (κϕ) n . The power-law model is known for being incompatible with observations because there exists no pair of exponent n and viable e-folding number N that result in acceptable values for the scalar spectral index and the tensor-to-scalar ratio simultaneously, as it was also presented in the canonical case as well. The gravitational Chern-Simons model however has the advantage that only tensor perturbations are affected by its inclusion. Therefore, for a viable e-folding number N = 60 and the Chaotic model n = 2, irrespective of the potential amplitude V 0 , the new degree of freedom, meaning the scalar coupling function ν(ϕ), along with the auxiliary parameter α, can be chosen such that the a viable tensor-to-scalar ratio is manifested. For instance, having a linear coupling function ν(ϕ) = ϕ/f with f = M P , V 0 = M 4 P and α = 1 suggests that the chaotic model becomes viable as the values n S = 0.96667 and r = 0.00639 seem to be in agreement with the latest Planck data. For the sake of generality, it should be mentioned that for the exact same set of parameters but a quadratic Chern-Simons coupling, the tensor spectral index is slightly positively defined as n T ≃ 4 × 10 −8 and increases with the increase of the exponent of the Chern-Simons coupling. On the other hand the tensor-to-scalar ratio decreases. Overall, it can easily be inferred that the gravitational Chern-Simons model is quite useful since due to the fact that tensor perturbations behave differently depending on their polarization state, recall Eq. (151), previously discarded models can now be rectified.

VIII. VACUUM f (R) GRAVITY INFLATION
Now let us consider the description of inflationary dynamics in the context of vacuum f (R) gravity, for details we refer to Ref. [176]. In the most general vacuum scenario with a modified gravity of the f (R) type, the gravitational action reads, and by varying with respect to the metric one obtains the gravitational field equations of motion, where f R = ∂f ∂R and □ = ∇ a ∇ a . The Ricci tensor is expressed in terms of the Christoffel symbols as, We are using the flat FRW metric, for which the non-zero Ricci tensor components are, The Ricci scalar is obtained by contracting the Ricci tensor with the metric tensor, where, is the Hubble rate and its first derivative with respect to time is found to be, Hence, the field equations (156) take the form, Our ultimate goal is to find an expression for the spectral index n S and the tensor-to-scalar ratio r in terms of the slow-roll indices. For a theory without matter components the latter read, The observational parameters n S and r are given by, and r = 48 Now, the relations (165)- (168) hold true for any model. If we further assume that the we have slow-roll inflation, the conditions,Ḣ correspond to ϵ i ≪ 1, i = 1, 3, 4 and we can use them to simplify the expressions of the slow-roll indices and observational parameters. To this end, ϵ 1 can be approximated as, Hence, the scalar spectral index can be written as, and the tensor-to-scalar ratio, Now let us try to simplify the expression of ϵ 4 , which is, where the dot denotes derivative with respect to time. The time derivative of the Ricci scalar is, where we have made use of the slow-roll approximation. Hence, ϵ 4 becomes, Butε 1 as seen by (165) is,ε So we can further approximate ϵ 4 as, where we named the dimensionless parameter in front of ϵ 1 , x/2, Hence, by combining (167), (168), (177) we find, Let us take a moment to discuss this result. The dimensionless parameter x = 48f RRR H 2 f RR does not have to be constant, as it depends on the Hubble rate and its derivatives both directly and through the Ricci scalar R appearing in f (R) and the derivatives, so its evolution can be calculated if one can solve the Friedmann equations for the evolution of the Hubble rate in the case of vacuum f (R) gravity. Using the definition of e-foldings number N = t end t k H(t)dt and inverting it, having found the evolution of H from the Friedmann equations we can express the initial horizon crossing time instance in terms of the e-foldings N number and the time instance t end that corresponds to the end of inflation. However, the condition that marks the end of the inflationary era ϵ 1 (t end ) = 1 can be used to extract the final time t end in terms of the model's free parameters, which we call σ. So, following this path we end up with an expression of t k as a function of t end (σ) and the e-foldings number N , which we can plug in to the solution of the Hubble rate evolution and following to the expression of x, so we will end up with x = x (N, σ), an expression dependent of the model's free parameters and the e-foldings number which we know for inflation has to be around 60. The "recipe" described above is not so easy to follow whatsoever, since solving the Friedmann equations is usually a difficult task for the majority of cases even assuming slow-rolling inflation. So instead of calculating the analytical form of x, let us discuss about its estimated behavior in some limiting cases, specifically, The first case we will examine is x ≪ 1, which also includes the case x = 0. In this scenario the expression for the tensor-to-scalar ratio (180) can by a series expansion as, which to leading order is given by, This is also the expression obtained directly from setting x = 0, so at leading order the expansion of the x ≪ 1 case gives the same r − n S relation as the R 2 model (f RRR = 0). Now, in the case where x is in the order of magnitude of 1, then we cannot approximate the r − n S relation further than (180), so profoundly this case is more complex and also blows up for f (R) theories which yield x ∼ 4, which are obviously non-viable. However, there may be found cases for which x is in the vicinity of 4 and yield a viable phenomenology. For x ≫ 1 we can approximate (180) as, There is one problem in this case, recall that ϵ 4 = − x 2 ϵ 1 − ϵ 1 and from the slow-roll conditions we have ϵ i ≪ 1, so if x ≫ 1 then the slow roll condition for ϵ 3 does not hold true anymore and the formulation we have described so far is not valid in this case, although for a spectral index respecting the Planck 2018 constraints x ≫ 1 would probably yield a very small tensor-to-scalar ratio and satisfy the Planck 2018 constraint. If x is constant, hence it does not depend on the e-foldings number which quantifies the evolution, then the r − n S relation takes the form, where α = 16 (4−x) 2 , which is the exact relation for α-attractors, f (R) = αR and some other string theory motivated f (R) gravities. However, this similarity does not render the x =constant case viable, since if x ≫ 1 yields a breakdown of the slow-roll conditions for the reason explained in the previous paragraph. However, it can produce a viable phenomenology for x ≪ 1 and for some x ∼ 4. Nonetheless, the results are model dependent in any case. Let us explore the x =constant case a little further. The Planck 2018 data constraint the scalar spectral index and the tensor-to-scalar ratio as 0.962514 ± 0.00406408 , r < 0.064, and the slow-roll index ϵ 1 ∼ O(10 −3 ). Additionally, for the slow-roll condition ϵ 4 ≪ 1 to hold true, the maximum value of ϵ 4 can be ∼ O(10 −2 ), so (177) constraints x to be x ∼ O(10) in order of magnitude. Combing the restrictions for n S and r and using (184) there can be found the ranges of x =constant that can yield viable slow-roll inflation era with N ∼ 50 − 60.

IX. SCALAR FIELD ASSISTED F(R) GRAVITY
In this section we present the case of an inflationary era generated by a canonical scalar field in the presence of a f (R) gravity. The action contains a gravitational term, the quadratic kinetic term of the scalar field and its potential, From the variation principle of the action δS = 0, with respect to the metric we can obtain the Einstein field equations and with respect to the field ϕ we obtain the equation of motion for the field. The Einstein field equations for the action (185) read, where, Now, taking into account (160), R = 12H 2 + 6Ḣ, and using it in the Einstein field equations, the latter reduce to the Friedmann equations which along the field equation of motion read, where the dot denotes differentiation with respect to the cosmic time and the prime with respect to the field ϕ. We assume slow-roll inflation, hence,Ḣ but also the slow-rolling of the field demands that its potential term is significantly larger than its kinetic term, but also the kinetic part of the field does not change rapidly so that inflation can last a sufficient amount of time, thus, which immediately reduces the field equation of motion to, Now that we have a complete set of equations, we can use the Friedmann equations (188), (189) to express the Hubble rate and its derivatives in terms of the potential V (ϕ) and the fields derivativesφ,φ, which we can also express in terms of the potential through the equation of motion (193), so we ultimately have the Hubble rate and its derivatives expressed only as a function of ϕ, the potential's free parameters and the function's f (R) free parameters. Having these, we can calculate the expression for the slow-roll indices, where, which quantify the slow-roll inflationary evolution and since we assume the slow-roll approximation, then the slow-roll indices have to satisfyε i ≪ 1, i = 1, 2, 3, 4. So now one can find the scalar spectral index n S and the tensor-to-scalar ratio r for such models, which are given by, The end of the inflationary era is characterized by ϵ 1 (ϕ end ) = 1, because at that instance the derivative of the Hubble rate is zero, so the Hubble horizon stops decreasing and will begin increasing again, marking the end of this cosmic era. From ϵ 1 (ϕ end ) = 1 we can then express ϕ end in terms of the model's free parameters, but also using the definition of the e-foldings number, we can write it as a function of the potential and the field and change the integration variable from cosmic time to the scalar field ϕ, and use it to solve for the initial value of the field at the beginning of inflation in terms of the e-foldings number N and ϕ end , so ultimately have the initial ϕ k only as a function of the model's free parameters, which can then use to find the slow-roll indices and observational parameters for this ϕ only depending on these free parameters.
Then, demanding confrontation with the Planck 2018 constraints we can find for which values of the parameter space these models are viable. We mention two examples in this section, one is a string-theory motivated corrected gravity, where M is a free parameter, and the latter is an effective modified Einstein-Hilbert action, motivated by the Lagrangian contained in the action, with α = 1 − λ, but in the large curvature regime R → ∞, at leading order is, Let us begin with the first example, For the field equations we need the derivatives of f (R) with respect to the Ricci scalar and the derivatives of the latter with respect to the cosmic time. They read, Combining (160) We now make use of the slow-roll approximations (191),(192) and we further assume that, which is an assumption we are making to simplify the formulation , based on the fact that for usual quasi-de Sitter expansion this holds true, but has to constantly be tested throughout calculations. So the field equations become, We can solve the second equation (210) forḢ, and we get, and by plugging it in (209) we get, At this point we will make another assumption, that the following approximation also holds true, the validity of which should also be tested as we move ahead to calculations. Finally, we Taylor expand the term 1 − 2κ 2φ2 M 2 and we arrive at the final forms of the Friedmann equations and field equation of motion, which is our system of dynamical equations that describe this model, Making use of all the ingredients, (204), (205),(214), (215) and after some algebra, one can find the expression of the slow-roll indices for the model at hand, R 2 minimally coupled scalar field inflation.
We quote the expressions for E andĖ, but omit the full expression of ϵ 4 since it is very long. The e-foldings number (199) for this model can easily be found to be, We pick a potential to perform a realistic test of this model and apply our framework. The potential we choose is a simple power-law potential, sometimes referred to as chaotic inflation, which does not yield a viable phenomenology when considered in the context of usual Einstein-Hilbert gravity, where V 0 is a model's dimensionless free parameter. Using (217)-(221) we find the expressions of the slow-roll indices, and from (197), (198) the expressions for the observational quantities can also be extracted. We omit the latter as well as the expression for ϵ 4 due to their length. From ϵ 1 (ϕ end ) = 1 we find ϕ end = and employing (222), we calculate the integral and solve with respect to ϕ k and find ϕ k = . It is obvious that we have expressed ϕ k only in terms of the model's free parameters and the e-foldings number, which for inflation is N ∼ 60. This model is a viable model, since there can be found ranges of values of the free parameters V 0 and β, where β = κM , and we set κ = 1, that satisfy the constraints on n s and r imposed by the Planck 2018 mission, but also comply with the observation that the amplitude of the primordial scalar perturbations is, The amplitude is given by, where z = (φk) E f R /κ 2 H 2 (ϵ3+1) . Furthermore, the free parameters should also respect the validity of our assumptionḢ 2 M 2 ≪ H 2 ,Ḣ 2 M 2 ≪ V (ϕ) and 2κ 2φ2 M 2 ≪ 1. Let us demonstrate an example, for V 0 = 9.37 × 10 −13 , β = 6.8 × 10 −6 and N = 60 we find n S = 0.96611, r = 0.063968 and P ζ (k) = 2.19216 × 10 −9 . One noteworthy feature of this model is that it sets an upper bound on the e-foldings number at N = 67, since for greater values of N there cannot be found values for the free parameters such that the constraints on the spectral index n S , on the amplitude P ζ and the approximation P ζ are satisfied simultaneously.
Let us now apply our framework to yet another case, that of where α is a dimensionless parameter 0 < α ≤ 1. Using the slow-roll approximation again we obtain the field equation of motion (193), while the Einstein field equations at leading order are, Substituting our f (R) in (200) and (201) we find the Friedmann and Raychaudhuri equations for this case, 3α Using the slow-roll approximation (191) the first Friedmann equation is simplified to, Now we have the full set of dynamical equations (233),(232), (193) that describe the inflationary era generated by a minimally coupled scalar field in the presence of a αR gravity, so we charge forward to calculating the slow-roll indices. This case is rather simple and the only non-zero slow-roll indices are ϵ 1 , ϵ 2 , The e-foldings number N can be found using (199) and (233) and is, Last, we need the expressions of the observational parameters, which are, Let us mention a brief and simple paradigm. We study the potential, and the slow-roll indices are, From ϵ 1 (ϕ) = 1 marking the end of inflation we find ϕ f = ln( 1 2 ( √ 2 √ αq+2)) κq and from (237) we get N = κϕ αq − e κqϕ αq 2 . We keep the exponential term only since ϕ is of order of 10 − 20M pl so this is the leading term. Solving for ϕ i we find expression for ϕ i = ln( 1 . We can then express the slow-roll indices as well as n s ,r in the initial state, which we quote since they are brief, The expressions are rather simple in this case, so we can even solve analytically the constraint expressions 0.962514 ± 0.00406408 , r < 0.064 and we find that for the free parameters (α, q) that satisfy, the latter are satisfied.

X. INFLATION IN THE CASE OF K-ESSENCE f (R) GRAVITY
The usual inflationary scenarios employ a scalar field ϕ whose potential V (ϕ) dominates over its kinetic energyφ forcing the Universe to an exponential type of accelerating expansion. In this section we will describe a class of models which yield an inflationary evolution in the absence of a potential, but in the presence of non quadratic kinetic terms along with the vacuum f (R) gravity. So, the action looks like, where κ 2 = 8πG c 4 and X = 1 2 ∂ µ ϕ∂ µ ϕ, and if the scalar field is homogeneous, thus depends only on the cosmic time t, then X = − 1 2φ 2 . The field equations and equation of motion for the metric (158) read, where f R = ∂f (R) ∂R , G X = ∂G ∂X . We set κ 2 = 1 for simplicity and present the case of slow-roll inflation, hence, we work in the regime where the slow-roll condition holds true,φ ≪ Hφ .
Let us note that in the case of f (R) = R it reduces to the case of k-essence inflation in the familiar Einstein-Hilbert gravity framework. Now, the slow-roll indices which quantify the evolution during the slow-rolling inflation era are given by, where, We also need to express the observational quantities, i.e. the spectral index of the scalar primordial curvature perturbations, n S , and the ratio of tensor to scalar ratio, r, in terms of the slow-roll indices. The latter are, where, Now the path to be followed in order to determine the viability of a model of this class is one. After choosing a model, thus a f (R) gravity and a form for the function G(X), one can use the Friedmann equations and equation of motion (246)-(248) to calculate the evolution of the Hubble rate, employing this and also using (160) can then find the evolution of the slow-roll indices and observational parameters. Demanding that the latter satisfy the Planck 2018 constraints and also the slow-roll conditions ϵ i ≪ 1 hold true during Inflation, the values of the model's parameter space for which it is viable can be determined. As examples we will present two models, In the first case, with f 1 > 0, so that the theory is free of ghost degrees of freedom. For this action (248) becomes, We can simplify this expression using the slow-roll condition (249), which allows us to discard terms containing higher powers of the first derivative and the second derivative of the scalar field. Finally we have, and the solution is simply found to be, The fact that we can obtain explicitly the evolution of the inflaton as a function of cosmic time, using the slow-roll condition, will assist in the calculation of the slow-roll indices significantly. To move ahead we need to solve the Friedmann equations to specify the evolution of Hubble rate, to this end we need to choose the f (R) model we want to study. For the purpose of this review we will present the case of a f (R) gravity which does not produce the desired results when considered alone, to investigate the difference made by the presence of the k-essence term G(X) in the Lagrangian. The model we study is, where α > 0 and n ∈ [ 1+ √ 3 2 , 2), so that we have an inflationary acceleration and not superacceleration. For n = 2 we have the case of the so called Starobinsky model which provides a viable phenomenological description of the inflationary era, in contrast to the case of n ∈ [ 1+ √ 3 approximation,Ḧ ≪ HḢ, in essence dismissing terms containing the second time derivative of the Hubble rate and using the result for the evolution of the scalar field ϕ(t) (260) the first Friedmann equation (246) reads, The last term can be dismissed since it is subleading, considering the fact that n ∈ [ 1+ √ 3 2 , 2). Also, recall the relation between the Ricci scalar and the Hubble rate is R = 12H 2 + 6Ḣ and to leading orderṘ = 24HḢ, substituting these in (262) and solving the differential equation one obtains the evolution of the Hubble rate with cosmic time, where, and t i is an initial time instance. Using these results (260), (263), we can now express the slow-roll parameters ϵ i (250) and subsequently the observational quantities n S , r in term of the model's parameters t i , f 1 , m, n, α. In the following we quote the expression for ϵ 1 , ϵ 2 , ϵ 3 , but omit the ones for ϵ 4 as well as n S , r since they are very lengthy, but they can be found in closed form.
Testing this model in terms of satisfaction of the Planck 2018 constraints 0.962514 ± 0.00406408 , r < 0.064, one can recover viable phenomenology for a range of parameter space values. Let us provide an example. For Using again the slow-roll approximation on the corresponding field equation of motion (248) we get, whose solution is the same as in the ghost-free case (260). The slow-roll indices obtained in this case are the same as (265) except ϵ 4 and thereafter calculate the expressions of the scalar spectral index and the tensor to scalar ratio. For this model, viable phenomenology can be yielded, but only with extreme fine tuning of the model's parameters. We mention one example: for (f 1 , α, n, m, t i )=(10 −40 , 6.751 × 10 43 , 1.36602, 1.4, 10 −10 ) then, n S = 0.966 and r = 0.0613, which are compatible with the Planck 2018 constraints. One detail about this model is that by default, for f 1 > 0, due to the sign of the + sign of the kinetic term, ghosts occur, so even the viable theories that can be described by this formulation, can only be effective theories, since the physical theory has to be free of ghosts.
We will now present another approach to examining k-essence f (R) gravity inflationary theories, that of reconstruction, which is a procedure that can be used for inflationary evolution realization in the slow-roll approximation in other cases as well. So, in this approach we know the action at hand, for us (245) and we have yielded the field equations and equations of motion (246)-(248), as well as the evolution of the inflaton with cosmic time (260), which is the same for both scenarios of our action. What we want to do now is reconstruct the form of the f (R) gravity and examine for which parts of the corresponding parameter space the whole theory confronts with the Planck 2018 constraints. The desired evolution of the Hubble rate in terms of the e-foldings number we would like to realize is, where β and γ are free parameters. We set H 2 (N ) = A(N ), so the Ricci scalar is given by, and combining (268) and (269) we find, Now, we have expressed the Hubble rate in terms of the e-foldings number and the latter in terms of the Ricci scalar, in a backwards sort of procedure, so next we can use this result in the first Friedmann equation to express it only in terms of the Ricci scalar and derivatives of the function f (R) with respect to R, solve the differential equation and obtain the desired form of f (R). The first Friedmann equation (246) by plugging in the solution for the scalar field evolution (260), which takes into account the slow-roll condition, reads, with the + sign corresponding to the phantom case (266) and the − sign to the ghost-free case (256). The Friedmann equation the yields the differential equation, whose solution is, where C 1 , C 2 are constants coming from integration and the parameters µ, ν are, Since we know the expression of the Hubble rate evolution we want (268), we can easily express the Ricci scalar in terms of the e-foldings number using (269) and we also know the explicit form of our f (R) function, which can in turn express in terms of N ., we can then use this to calculate the slow-roll indices and as well as the spectral index and the tensor to scalar ratio. For this purpose we note that since, the slow-roll indices can be expressed in terms of the e-foldings number N as, and, with where, 2) we get n S = 0.964894 and r = 0.0179065. However, for the ghost-free model there cannot be found values of the free parameters that satisfy the constraints on n S ,r simultaneously, so we cannot realize a viable phenomenology with Hubble rate evolution of (268) with this model.
XI. GENERALIZED f (R, ϕ, X) GRAVITY At this stage, it becomes abundantly clear that the overall phenomenology described previously can be written in a compact manner. This realization simplifies the analysis to quite an extend and also paves the wave to other theoretical models that are not presented here. Let us consider that the gravitational action is written as, where f (R, ϕ, X) is an arbitrary function of the Ricci scalar, the scalar field and the kinetic term X = 1 2 g µν ∇ µ ϕ∇ ν ϕ while L (c) is the Lagrangian density for corrections such as string corrections and axionic strings as presented before. This choice of action corresponds to the following set of equations, As an example, one can obtain the canonical scalar field case for f (R, ϕ, X) = R − 2κ 2 (X + V ), the k-essence case for f (R, ϕ, X) = R − 2κ 2 (X + ωX 2 ). Also, for the Chern-Simons model, one simply needs to replace T (c) = 0 = T (c)µ µ . In this unified way, the background equations read, and, where the components of the energy stress tensor for string corrections are given by the expression T and it could be either the contribution from string corrections studied in (112) or from the gravitational Chern-Simons case (142). In the latter case, the equation of motion remain unaffected obviously and only tensor perturbations are affected as showcased before. In this approach, the f (R, ϕ, X) model can be used in order to study non canonical theories such as Dirac-Born-Infeld or the inclusion of higher order kinetic terms. This description is convenient when studying primordial scalar non Gaussianities under the constant-roll assumption since the sound wave velocity in this case, which reads, could obtain a quite small value and thus the equilateral nonlinear term f eq N L can be enhanced to the point where it may be detectable in subsequent experiments. For an arbitrary model, the sound wave velocity, as mentioned before, should be well-behaved, meaning that it satisfies the condition 0 < c A ≤ 1. In addition, the slow-roll indices that are required in order to study the inflationary era are, where all the auxiliary parameters can be taken from Eq. (120) with the only change being f → f ,R . Note that for the f (R, ϕ, X) model without string corrections, the auxiliary parameter E that participates in the 4-th slow-roll index is given by the following expression, which, for the case of f (R, ϕ, X) = f (R, ϕ) − 2κ 2 (X + ωX 2 + V ), it obviously coincides with the result obtained from action (109) for c 1 = 0 = c 2 = c 3 and ξ(ϕ)c 4 → ω 4 , with ω being ϕ independent. Furthermore, if the Chern-Simons model is considered, then the slow-roll index ϵ 5 should be replaced as ϵ 5 = 1 In the following we shall present certain examples based on this formalism, describe in detail the results obtained and showcase certain differences between the models at hand.

A. Kinetic Axion f (R) Model
Let us start by presenting the results from the kinetic axion f (R) model. This model has an intriguing phenomenology due to the axion dynamics and its impact on both the inflationary era and subsequent cosmological eras. For this model, the gravitational action reads, where the functions of the model are specified as, with M being a mass scale indicative of the energy scale beyond which the R 2 term dominates, m α stands for the axion mass and finally f α serves as the axion decay constant. In this approach, it is assumed that the f (R) gravity dominates primordially which is a reasonable assumption considering that the Ricci scalar R = 12H 2 + 6Ḣ has a quite large value. In turn, the canonical scalar field can safely be considered as subleading and thus the background equations of motion (284) and (285) are equivalent to the ones obtained in the pure f (R) case, that is, Therefore, the background equations coincide with the vacuum R 2 model which is known for admitting the following solution for the Hubble rate expansion, Here, it should be stated that apart from the linear dependence on cosmic time t, which in turn implies an exponential expansion (note that the conditionä > 0 is satisfied), it also implies that the (quasi) de Sitter expansion is fully determined by the f (R) gravity and in particular the effective mass scale. In order to be in agreement with observations, the mass scale is considered to be approximately M = 150 2N 10 −5 M P where N stands for the e-folding number. In addition, parameter H I is assumed to be approximately of order O(10 13 ) GeV so that a viable inflationary era is obtained. Regarding the scalar field, for the kinetic axion model, it is assumed that ϕ, which at this stage evolves dynamically, is inferior to the axion decay constant, that is ϕ ≪ f α . As a result, the canonical potential can be well approximated by a quadratic potential, which, according to the latest estimates for the axion mass which predict an upper bound m α ≤ O(10 −12 ) eV, it can easily be inferred that the scalar potential primordially is negligible. Note that while it would be acceptable to further approximate the potential with the inclusion of a quartic term, since it would describe fundamental interactions between the axions as mentioned in the canonical scalar field case previously, the interaction coupling in this case is given by the ratio mα fα 2 which is extremely small, given that the axion decay constant is quite large, thus the mass term suffices. As a result, the Klein-Gordon equation (43) can be approximated as, which in turn implies that the axion, primordially, behaves as a stiff matter perfect fluid. This can easily be inferred from the EoS which reads, therefore the energy density scales as ρ α ∼ a −6 (t). As mentioned before, this value for the EoS is the maximally allowed value that is in agreement with causality and in short states that pressure is directly proportional to the energy density, that is P = ρc 2 , therefore the deceleration parameter in this case has the largest value possible, that is q = 2. The reader should also keep in mind that while this is indeed the behavior of the axion primordially, the canonical potential becomes important subsequently when the scalar field approaches its vacuum expectation value φ = θ α f α , with θ α being the misalignment angle ranging between 0 < θ α < 1, and this in turn implies that the axion redshifts as dark matter since ρ α = 1 2φ 2 + V ∼ a −3 (t), making it a perfect candidate for successfully describing cold dark matter but for the time being, we shall focus only on the primordial phenomenology. Since the kinetic term dominates and the effective EoS at the end of inflation is specified by a stiff matter perfect fluid, the inflationary era is in this case further prolonged depending on the reheating temperature. Before we showcase this however, an important feature should be addressed. Due to the fact that the axion mass is approximately m α ∼ O(10 −12 )eV, thus rendering the scalar potential negligible, the second index is completely specified by the continuity equation of the scalar field (297) and thus it becomes of order O(1), that is ϵ 2 ≃ −3. This is a prime example as to why the indices in Eq. (288) should not be regarded blindly as slow-roll indices since their numerical value may be important. This may seem as a problem given that a quite large value for ϵ 2 accompanied by typical values for the rest indices ϵ i produce a scalar spectral index which is at variance with reality, as it was shown in Refs. [157], however for the case at hand, it is the f (R) part that actually saves the model given that at leading order, ϵ 4 = −ϵ 1 − ϵ 2 , thus resulting in an elegant cancellation between the large indices. As a result, for N = 60 and M = 1.25 × 10 −5 M P , the results obtained coincide with the vacuum R 2 model. The only constraint that needs to be considered is the initial value of the scalar field during the first horizon crossing relative to the axion decay constant. More specifically, in order for the expansion (296) to be valid, one requires f α ≫ ϕ k . In addition, the fact that ϵ 2 ≃ −3 primordially implies that during the first horizon crossing the scalar field obeys a constant-roll condition. This is connected to the appearance of scalar non-Gaussianities in the CMB however due to the fact that the R 2 term dominates, the usual result is obtained which implies that the equilateral nonlinear term is of order O(10 −2 ) so no significant prediction is made.
Let us return to the kinetic axion dynamics and focus particularly on the EoS. As briefly mentioned before, the dominance of the kinetic term of the axion compared to its scalar potential implies that the axion behaves as stiff matter which increases the duration of inflation depending on the reheating temperature. More specifically, in Ref. [177], it was shown that the duration of inflation has the following form, with k * = 0.05Mpc −1 being the pivot scale and subscripts k, end and reh are used in order to denote the moment where inflation starts and modes become superhorizon, the ending stage of inflation and the ending of the reheating era or in other words the start of the radiation domination era respectively. Note that at the end of the reheating era, the radiation fluid is assumed to be in thermal equilibrium therefore its energy density is specified by the relation ρ reh = π 2 30 g * T 4 reh with g * being the relativistic degrees of freedom at that instance hence the reason why we stated that extension of inflation is affected by the reheating temperature. Note that in principle the degrees of freedom differ between different cosmological eras therefore in a sense they also scale with temperature however for the case at hand one can safely assume that g * ∼ O(100). In a sense, the vacuum R 2 model dominates inflation however when it reaches its end, the kinetic term of the axion dominates, therefore an intermediate era of stiff matter manifests between the inflationary era and the radiation domination era, meaning during the reheating era. One can understand this as a result of the stability of the de Sitter fixed point that emerges when an autonomous dynamical analysis is performed for the power-law scalar field assisted R 2 models. In short, the de Sitter fixed point is nonhyperbolic and after a few e-folds have elapsed, depending on the initial conditions for the dynamical parameters used, the trajectories on the phase space are driven away however instead of reaching the expected EoS of ω = 1 3 corresponding to radiation, the kinetic axion kicks in and results in an intermediate stiff matter era. Now if the EoS in Eq. (299) is replaced as ω = 1 and the reheating temperature is specified, then depending on the value of the temperature, the e-folding number is increased by O(1). The change may seem not extensive however the impact of this change has a significant effect on the observational quantities. We report that according to the findings of Ref. [162], a really high reheating temperature of approximately T R ∼ O(10 12 )GeV results in an increase of the e-folding number by approximately 5 e-folds, making it equal to N = 65.3439 and in consequence the scalar spectral index now reads n S = 0.969393 while the tensor-to-scalar ratio becomes equal to r = 0.00281042. While the tensor-to-scalar ratio can be decreased by approximately 15%, the scalar spectral index lies outside the area of viability, recall that n S = 0.9649 ± 0.0042 with a 68% C.L. This implies that extremely large reheating temperatures of order O(10 12 )GeV and beyond should not be considered in the kinetic axion f (R) model since the following predictions are not in agreement with observations. On the other hand, a typical value of T reh = 10 7 GeV affects mildly the results as now n S = 0.967483 and r = 0.00317206. Overall, large reheating temperatures increase the scalar spectral index of primordial curvature perturbations while the tensor-to-scalar ratio, and as a result the tensor spectral index, decrease. Obviously, the opposite applies to the case of small reheating temperatures but a lower bound for the reheating temperature around the MeV scale should be considered since subsequent cosmological eras cannot follow for smaller values of T reh .

B. Kinetic Axion Gauss-Bonnet Model
Having introduced a phenomenologically interesting scalar field assisted f (R) model, it stands to reason that higher contribution of the scalar field as presented before can affect the inflationary era. For the second model let us consider that the gravitational action contains also string corrective terms and in particular a non-minimal coupling between the scalar field and curvature through the inclusion of the Gauss-Bonnet topological density [162], where now Eq. (292) is still valid, along with the expansion (296) and in addition, the Gauss-Bonnet scalar coupling function is assumed to be linear, that is ξ(ϕ) = ϕ f . Of course, all the previous statements are still valid, meaning that the R 2 dominates the background equations and the scalar field affects the observational quantities though the second and fourth indices ϵ 2 and ϵ 4 , not to be considered as slow-roll indices. The model is chosen for two reason. Firstly, the linear Gauss-Bonnet scalar coupling function is obviously the simplest model that one can consider. Secondly, it can be shown that if the constraint on the propagation velocity of tensor perturbationsξ = Hξ is imposed on the case of linear coupling, the constant-roll conditionφ = Hφ that emerges is in fact too strong and spoils the viability of the model regardless of the scalar potential that is chosen. In the previous example however we showcased explicitly that while ϵ 2 is quite large, ϵ 4 turns out to be proportional to −ϵ 2 and thus they cancel when the scalar spectral index is computed, leaving us only with the rest slow-roll indices which can be considered as slow-roll given that they are approximately of order O(10 −3 ). In particular, without imposing any approximations, one can show that for the R 2 model, index ϵ 4 reads, where ϵ 5 = ϵ3+y(1−ϵ1) 1+2y and for simplicity, the auxiliary variables that were introduced are equal to x = κ 2φ2 6f ,R H 2 and y = κ 2 Qa 2f ,R H . Obviously, in the limit of y → 0 and for a negligible x, the previous result ϵ 4 = −ϵ 1 − ϵ 2 is extracted. It should also be stated that a viable inflationary model should predict quite small values for the auxiliary parameters x and y, something that can easily be inferred from the background equations (284)-(285) since the f (R) term cannot dominate the equations of motion if x and y are comparable to it. In consequence, even though string corrections are present, the condition ϵ 5 < −ϵ 1 cannot be satisfied and thus the model at hand is incapable of producing a blue tilted tensor spectral index unless string corrections dominate over the f (R) part. Tensor perturbations are most likely to be affected solely by the f (R) part, meaning that the consistency relation r = −8n T is still satisfied exactly as in the canonical scalar field case. In particular, one can show that the tensor-to-scalar ratio for the constrained Gauss-Bonnet model has the following form, however for the constrained Gauss-Bonnet model, c S ≃ c even in the presence of an f (R) gravity and provided that y ≪ ϵ 1 , the vacuum R 2 result is extracted. Let us now proceed with the derivation ofφ. According to the continuity equation of the scalar field (111), assuming that the canonical potential is inferior and that the scalar field satisfies the constant-roll conditionφ = Hφ, one can easily see that, from which one can see that in this approach,φ(t end ) = 0 when inflation stops. As shown, the overall phenomenology is quite different compared to the previous case not only due to the inclusion of the Gauss-Bonnet density but also due to the constant-roll condition imposed from the propagation velocity of tensor perturbations. Let us see the impact of this result on the inflationary era. For typical values for the free parameters as M = 1.25 × 10 −5 M P , N = 60 and f = 10 −5 M P , one can show that n S = 0.96686, r = 0.00327847 and n T = −0.000413283 which are obviously in agreement with experimental data. In addition, the tensor spectral index and the tensor-to-scalar ratio coincide with the vacuum R 2 result which in turn implies that parameters x and y are not so important which should be the case in order for the R 2 to dominate the background equations. On the other hand, when x, depending on the value of f , increases quite close to ϵ 2 3 , then the scalar spectral index can be affected though index ϵ 4 . Note also that parameter f is not necessarily the axion decay, here it was simply treated as an additional degree of freedom however if one identifies it as the axion decay rate, then this in turn implies that f α ∼ O(10 13 )GeV and thus the value of the scalar field during the first horizon crossing needs to be lesser than this in order for the expansion of the scalar potential to be valid.
Concerning the form ofφ in Eq. (303), one can easily deduce that at the ending stage of inflation where ϵ 1 = 1, the kinetic term of the axion vanishes identically. Obviously this is not exactly true but the idea that the kinetic term now becomes inferior is correct. To showcase this, let us include the scalar potential in the continuity equation and go to the limit where ϵ 1 = 1. In this case, the time derivative of the scalar field reads, and thus, the ratio between the kinetic term and the scalar potential at that time instance is, and thus no extension of the inflationary era occurs. This is a drastic change between the canonical case and the Gauss-Bonnet model. It should be stated that the constraint is really powerful since it alters the continuity equation of the scalar field from a second order to a first order differential equation and thus the solutionφ(ϕ) can be found algebraically. In fact, including additional string corrective terms does not seem to affect the outcome since either cξ(ϕ)G µν ∇ µ ϕ∇ ν ϕ or ω 1 2 g µν ∇ µ ϕ∇ ν ϕ 2 will still result in the conditionφ(t end ) = 0. Obviously the phenomenology is completely altered if the constraint on the propagation velocity of tensor modes is not initially imposed.

C. Kinetic Axion Chern-Simons Model
For the final model we shall consider the kinetic axion model studied in the previous two cases however now a parity odd term will be included. Suppose that the gravitational action reads [162], where the Chern-Simons scalar coupling function ν(ϕ) is assumed to be quadratic, that is ν(ϕ) = ϕ f 2 with f being an auxiliary parameter with mass dimensions of [f ] =eV, not to be confused with the axion decay constant necessarily. In this approach, phenomenologically speaking, there exist not many differences. In fact, regarding index ϵ 4 , the correct expression is extracted from Eq. (301) in the limit of y → 0. In addition, the scalar spectral index is expected to be the same as in the previous cases if N = 60 is selected. The only element that is affected non trivially is the description for tensor perturbations and more specifically the tensor spectral index and the tensor-to-scalar ratio. For simplicity, let us designate two auxiliary dimensionless parameters as y l = 2λ l κ 2ν H F and z = ν ′′φ Hν ′ . In consequence, the fifth index (153) now reads, and thus the tensor observable, which for the model at hand are given by the following expressions, Let us now see how the results can differ from the previous cases. For the same set of values as before and for f = 10 −9 M P , ϕ k = 10 −10 M P , one finds that while the scalar spectral index and the tensor-to-scalar ratio are not affected and in fact coincide with the vacuum R 2 model, the tensor spectral index obtains the value n T = 0.0132771 thus the Chern-Simons model can in fact yield a blue-tilted tensor spectral index, in contrast to the previous two cases. This is because the new degree of freedom does not participate in the background equations and thus a large value ofν during the first horizon crossing is possible without spoiling the dominance of the f (R) gravity at the level of the equations of motion. Furthermore, due to the fact that the gravitational Chern-Simons term does not influence the background equations, recall Eq. (144), the continuity equation of the scalar field suggests that the apparent dominance of the kinetic term of the axion over its scalar potential results in the appearance of a stiff matter era during reheating, therefore the duration of the inflationary era is prolonged depending on the reheating temperature, exactly as was the case with the kinetic axion f (R) model. In addition, since the energy density during the first horizon crossing and at the end of inflation is not affected by the gravitational Chern-Simons term, the increase in the e-folds is exactly the same. We report that for high values of the reheating temperature as T reh = 10 12 GeV and further, where it was shown that the scalar spectral index is at variance with observations, the tensor spectral index decreases in value and now reads n T = 0.0117172, roughly speaking a 12% decrease while typical values for the reheating temperature as T reh = 10 7 GeV suggest that n T = 0.0128038,, only a 3.6% decrease. In short, this model proves that the tensor spectral index can easily be manipulated by the assumption that the aforementioned parity odd term participates in the gravitational action (306) without spoiling scalar perturbations or even the background equations, in other words the phase space has exactly the same fixed points whether the Chern-Simons term is considered or not, the only change lies with the behaviour of tensor perturbations depending on their polarization state. This is also an interesting outcome since it showcases that scalar-field assisted f (R) gravity models can in fact result in the amplification of the energy spectrum of primordial gravitational waves, which is an exciting result that may explain a possible signal from a SGWB provided by LISA in the next decade or so.

XII. ENERGY SPECTRUM OF PRIMORDIAL GRAVITATIONAL WAVES
In the last section of this review, we shall briefly discuss the impact of modified theories of gravities in the energy spectrum of primordial gravitational waves. This is because certain models manage to predict quite different results compared to the general relativistic description, at least in the high frequency regime. Hence, making comparisons with a future signal detected by a third generation detector, provided that it is attributed to a stochastic gravitational wave background, can result in further constraints on parameters used to quantify inflation. Even to date, the NANOGrav 2023 detection [5] of a stochastic signal of primordial gravitational waves can be attributed to cosmological sources and even to some modified gravities with blue tilted tensor spectral index [55,160] In order to showcase the impact that modified theories of gravity have on the energy spectrum of gravitational waves, it is preferable to firstly study the general relativistic predictions.
Let us commence our study by considering only tensor perturbations in the perturbed metric Eq. (121), that is, where the conformal time τ is used for simplicity and x i describes the spatial coordinates in a comoving frame. Here, h ij denotes the gauge-invariant perturbed metric which describes symmetric and traceless transverse models, meaning that h ij satisfies the following conditions, These conditions are implemented in order to properly describe gravitational waves. Now for the sake of simplicity, let us work with Einstein-Hilbert gravity such that, with G being the gravitational constant, c stands for the speed of light and L matter denoting the Lagrangian density for matter. Using the perturbed metric in Eq. (309), the gravitational action can be written up to second order as, where the anisotropic stress tensor is defined as Π i j = T i j − pδ i j while it satisfies the traceless and transverse conditions δ ij Π ij = 0 = ∂ j Π ij . Now by varying the perturbed action (312) with respect to the gauge-invariant field h ij , one finds that, where hereafter the prime is used in order to denote differentiation with respect to conformal time τ . Now in order to proceed, it is convenient to work in momentum space by performing a Fourier transformation and also treat h ij and Π ij as operators describing conjugate variables. In particular, let the Fourier expansion of the fields be, where e r ij is the polarization operator for tensor perturbations describing states + and ×. For the sake of consistency, the polarization operator must satisfy the traceless and transverse conditions δ ij e ij = 0 = k i e ij since the same requirements appear in momentum space as well. Here,ĥ r ij (τ, x) is treated as an operator since perturbations in momentum space can be interpreted as a creation of a mode with momentum −k or an annihilation of a mode with momentum k, that isĥ r k (τ ) = h k (τ )â r k + h * k (τ )â r † −k in order forĥ r ij to be a Hermitian operator. Also, in order to quantize the theory, the following equal time commutation relations are imposed, [ĥ r k (τ ),π s k ′ (τ )] = iδ rs δ (3) (k − k ′ ) , [ĥ r k (τ ),ĥ s k ′ (τ )] = 0 , [π r k (τ ),π s k ′ (τ )] = 0 , whereπ r k (τ ) = a 2 (τ )ĥ r ′ −k (τ ), or using the creation/annihilation operators, Note that the factor of (2π) 3 appears in the commutation relation between creation and annihilation operator of tensor modes for consistency only because the Fourier transformation considered for tensor perturbations in Eq. (314) has however due to the fact that a nearly scale invariant result is expected according to the latest Planck observations, one can focus on the first order running if not the pivot scale tensor spectral index. In other words one can assume that, n T (k) = n T (k * ) + a T (k * ) 2 ln k k * , and substitute in Eq. (322). In addition, referring to the amplitude of tensor perturbations, it is known that the amplitude is connected to the tensor-to-scalar ratio as, and thus, by combing all the above expressions, one can show that the energy spectrum for gravitational waves is [178], This result applies to general relativity and it should be stated that due to the fact that for a canonical scalar field that dominates during inflation the tensor spectrum is red, it can easily be inferred that a suppression of modes is expected at high frequencies. Let us now see how modified theories of gravity can in principle affect the overall procedure. By following the same reasoning as before, one needs to extract information about the behavior of tensor perturbations. In previous sections, a detailed analysis was provided, see for instance Eq. (123) for string corrections.
Overall, the main difference between GR and modified gravity lies with the running of the Planck mass as now, for an arbitrary modified scalar-tensor theory, traceless and transverse modes satisfy the following equation [178], where a M =Q t HQt as mentioned before is the running Planck mass and parameter Q t is the shifted Planck mass square and has a unique form for every model. Note that not all modified theories of gravity affect the aforementioned mode equation, see for instance k-essence models. In this case, modifications simply alter numerically the tensor-to-scalar ratio and the tensor spectral index by means of the first slow-roll index ϵ 1 as stated before. Now in order to solve this equation, one can incorporate the WKB approximation, therefore the solution can be written with respect to the GR solution as, where h ij = he ij and the exponent is equal to D = 1 2 z 0 dz ′ a M 1+z ′ . In the end, since the tensor power spectrum (319) is proportional to the vacuum expectation value of the contraction of a transverse traceless mode with itself, it becomes clear that the gravitational wave energy spectrum is affected in a similar manner by modified theories of gravity through this exponent and in the end, the final expression that takes into consideration the effect of modified theories of gravity reads [178], In principle, it should be stated that the effects of modified theories do not lie only on the aforementioned exponent D since different considerations can in principle impose certain constraints on the reheating era for instance and therefore the reheating temperature which participates in the second transfer function affects the energy spectrum, see for instance the kinetic axion Chern-Simons model. In this approach, models with a blue spectrum are actually favored since they can result in an apparent amplification of the energy spectrum at high frequencies compared to the result of GR. In consequence, a potential signal in the near future from second and third generation detectors of gravitational waves can be explained by modified gravity models that predict a blue-tilted tensor spectral index. On the other hand, if a signal is detected near the GR prediction then this could still be explained by modified theories, for instance power-law f (R) models that a small, insignificant enhancement is extracted. Finally, models with a heavy red spectrum can suppress the signal thus making the search of primordial gravitational waves more difficult. Overall, modified theories seem to perplex the results however they stand strong in explaining a potential signal that is attributed to a stochastic gravitational wave background. It should also be stated that as long as the running of the Model aM f (R, ϕ) (F,RṘ + F ,ϕφ )/(HF ) f (R, ϕ) Gauss-Bonnet (F,RṘ + F ,ϕφ − 8κ 2 (Hξ +Ḣξ))/(H(F − 8κ 2ξ H)) f (R, ϕ) Constrained Gauss-Bonnet (F,RṘ + F ,ϕφ + 8κ 2ξ H 2 q)/(H(F − 8κ 2ξ H)) f (R, ϕ) String Corrections (F,RṘ + F ,ϕφ − 8c1κ 2ξ H 2 (ξ Hξ +Ḣ H 2 ) + c2κ 4 ξHφ 2 (ξ Hξ + 2φ Hφ ))/(H(F − 8κ 2 c1ξH + c2κ 4 ξφ 2 )) f (R, ϕ) Constrained String Corrections (F,RṘ + F ,ϕφ − 8c1κ 2ξ H 2 + 2c2κ 4 ξHφ 4 (1 +ξ 2Hξ +φ Hφ ))/(H(F − 8κ 2 c1ξH + c2κ 4 ξφ 2 )) f (R, ϕ) Chern-Simons ((F,RṘ + F ,ϕφ )/(HF ) + (2λ l κ 2ν /F )(ν Hν − 1)k/a)/(1 + 2λ l κ 2ν Planck mass a M for a given model is nontrivial, then the amplitude of tensor perturbations in the energy spectrum differs from the GR prediction, hence the reason why the high frequency regime is studied since in this regime the inflaton has yet to reach its vacuum expectation value. This however does not mean that high frequencies are only important since the case of an f (R) gravity could in principle have interesting implications in low frequencies as well.
As a final note, let us present in detailed the running of the Planck mass for various models that were previously studied. For the models considered in the present article, it should be stated that only a few affect the energy spectrum of gravitational waves and are presented in Table. XII. In order to understand which models affect the spectrum, one needs to consider the mode equation and the effect that terms included in the gravitational action have on the running Planck mass. The simplest example is the f (R) gravity for which it can easily be inferred that a M =Ḟ HF where F = df dR for simplicity. The same expression applies to general f (R, ϕ) models since only the derivative with respect to the Ricci scalar appears. In consequence, the numerator that contains information about the time derivative gives rise to a term proportional toφ and affects the results for as long as the scalar field evolves dynamically. In other words, when the vacuum expectation value is reached, only the contribution from the Ricci scalar remains, as stated previously. This description is of course model dependent. The second model that we considered here is the Einstein-Gauss-Bonnet gravity. As shown, it is one of the terms that affect the propagation velocity of tensor perturbations and in general the behavior of tensor modes. By recalling the form of Q t from (120) for the case of a general f (R, ϕ) model, it becomes clear that a M is proportional to Hξ +Ḣξ. This term, similar to the previous case, appears when the scalar field evolves dynamically with respect to time. It is interesting to mention that under the assumption that tensor perturbations propagate through spacetime with the velocity of light, the Gauss-Bonnet scalar coupling function satisfies the differential equationξ = Hξ which upon solving, it becomes clear that the contribution of the Gauss-Bonnet term is proportional to the scale factor and thusξ = λ 1+z with λ being an auxiliary parameter with mass dimensions of eV. In addition, for ω ef f = − 1 3 , or in other wordsä = 0, the contribution ofξH +Ḣξ vanishes identically. Hence, regardless of the coupling function, the scaling with redshift is the same and therefore the results are in fact universal. Now in Eq. (120), it becomes clear that the running of the Planck mass is affected also by the inclusion of the kinetic coupling ξ(ϕ)G µν ∇ µ ϕ∇ ν ϕ however now the results are not universal, in contrast to the Gauss-Bonnet case, as a different coupling results in a different evolution forφ. Note that is the kinetic coupling is considered on its own without the Gauss-Bonnet density, then while it is possible to extract information about the running Planck mass and in consequence the conditions under which a blue spectrum is extracted, if they are not strictly model dependent, compatibility with the GW170817 event cannot be extracted as now the description for primordial massless gravitons requires a vanishing kinetic coupling. Finally, the Chern-Simons model takes into consideration the different circular polarizations as shown in Eq. (151). This model is quite interesting since for a given frequency, two signals are expected with a different amplitude. This can be inferred from the fact that modes satisfy a different differential equation based on their circular polarization and therefore a difference is expected, depending on the choice of the Chern-Simons scalar coupling function it may be insignificant or on the contrary, quite important. This result appears regardless of the sign of the tensor spectral index and is indicative of the chirality that gravitational waves posses in the Chern-Simons model since it violates parity. For this model it is expected that for a viable inflationary era, high frequency modes are separated in the energy spectrum based on their chirality with the difference tending to zero as the frequency decreases. This is similar to having the scalar field reach its vacuum expectation value. When the minimum is reached and no dynamical evolution is present, chirality is restored in the limit of the Chern-Simons mass scale reaching infinity and therefore no split is expected between modes any longer.

XIII. CONCLUSIONS
In this review we presented the most recent trends for inflationary dynamics in terms of modified gravity theories. The motivation for using modified gravity theories in order to describe inflation is multi-fold since the standard single scalar field description of inflation has many shortcomings, among which, there must be a large number of couplings to the standard model particles in order for the Universe to be reheated, or if inflation can explain the 2023 NANOGrav stochastic gravitational background observation the tensor spectral index has to be significantly blue tilted and so on. Apart from the inflationary era motivation, the late-time era also must be described by modified gravity since the possibility that the total background equation of state parameter is slightly phantom, cannot be described by the standard single scalar field description of general relativity without resorting to phantom scalar fields. Thus we provided a modern text that describes the most timely extensions of general relativity for describing the inflationary era. Note that the general relativistic descriptions of the inflationary era are basically provided by single scalar field and k-essence descriptions. We provided an informative overview of the inflationary paradigm in which we pointed out the shortcomings of the standard hot Big Bang scenario, and we explained how the inflationary paradigm may theoretically solve these problems. Also we emphasized how important is the inflationary paradigm theoretically and why it should eventually be the correct description of nature regarding the early time era, since it is the only scenario which provides a nearly scale invariant power spectrum of primordial scalar curvature fluctuations, which are necessary ingredients for a successful large scale matter structures existence. We provided a detailed overview of how the inflationary era may be generated by a single scalar field theory with minimal and non-minimal couplings. We calculated in some detail the necessary slow-roll indices, and we provided and proved in some detail several well-known formulas regarding the single scalar field inflationary theories. After we briefly discussed the swampland criteria, and the constant-roll evolution which serves as an alternative to the standard slow-roll evolution, we studied and analyzed several string motivated models of inflation, which involve Gauss-Bonnet couplings of the scalar field, higher order derivatives of the scalar field, and some subclasses of viable Horndeski theories. We also presented and analyzed inflation in the context of Chern-Simons theories of gravity, including various subcases and generalizations of string corrected modified gravities which also contain Chern-Simons correction terms, with the scalar field being identified with the invisible axion, which is the most viable dark matter candidate to date. We also provided a detailed account of vacuum f (R) gravity inflation, and also inflation in f (R, ϕ) and kinetic-corrected f (R, ϕ) theories of gravity. In the end of the review we discussed how the gravitational waves evolve in the context of modified gravity and we quantified the effect of modified gravity on the general relativistic waveform, quantifying the overall effect in a single parameter which must be evaluated for all the evolutionary eras of our Universe.
In the next decade, several experiments will provide concrete evidence of inflation, like the stage 4 CMB experiments in 2027 and the interferometric gravitational wave experiments like LISA and the Einstein telescope in 2035. The smoking gun for the existence of inflation is the actual observation of CMB B-modes (curl modes). In such a case, the interferometric experiments will provide hints on which model may be the correct description of the inflationary era. Indeed, if the NANOGrav signal originates from an inflationary era, then the theories that can describe such an era have specific characteristics. If the NANOGrav signal is combined with other future gravitational waves experiments, this could be useful for determining the actual theory behind inflation. For example a signal detectable in more than two distinct frequency ranges, could for example indicate a flat energy spectrum for gravitational waves, or if the signal is detectable by some experiments, this could point out specific characteristics of the theory that describes inflation and the subsequent reheating era. Thus the necessity for a multi-frequency study of gravitational waves is compelling. In addition, the observation of scalar modes of gravitational waves will provide direct evidence for a modified gravity theory describing nature, although such a task is challenging due to the sensitivity of the detectors. In all cases, the next decade belongs to the quest of the inflationary era and the early Universe probes.