Hybrid Loop Quantum Cosmology: An Overview

Loop Quantum Gravity is a nonperturbative and background independent program for the quantization of General Relativity. Its underlying formalism has been applied successfully to the study of cosmological spacetimes, both to test the principles and techniques of the theory and to discuss its physical consequences. These applications have opened a new area of research known as Loop Quantum Cosmology. The hybrid approach addresses the quantization of cosmological systems that include fields. This proposal combines the description of a finite number of degrees of freedom using Loop Quantum Cosmology, typically corresponding to a homogeneous background, and a Fock quantization of the field content of the model. In this review we first present a summary of the foundations of homogeneous Loop Quantum Cosmology and we then revisit the hybrid quantization approach, applying it to the study of Gowdy spacetimes with linearly polarized gravitational waves on toroidal spatial sections, and to the analysis of cosmological perturbations in preinflationary and inflationary stages of the Universe. The main challenge is to extract predictions about quantum geometry effects that eventually might be confronted with cosmological observations. This is the first extensive review of the hybrid approach in the literature on Loop Quantum Cosmology.

of cosmological spacetimes, both to test the principles and techniques of the theory and to discuss its physical consequences. These applications have opened a new area of research known as Loop Quantum Cosmology. The hybrid approach addresses the quantization of cosmological systems that include fields. This proposal combines the description of a finite number of degrees of freedom using Loop Quantum Cosmology, typically corresponding to a homogeneous background, and a Fock quantization of the field content of the model. In this review we first present a summary of the foundations of homogeneous Loop Quantum Cosmology and we then revisit the hybrid quantization approach, applying it to the study of Gowdy spacetimes with linearly polarized gravitational waves on toroidal spatial sections, and to the analysis of cosmological perturbations in preinflationary and inflationary stages of the Universe. The main challenge is to extract predictions about quantum geometry effects that eventually might be confronted with cosmological observations. This is the first extensive review of the hybrid approach in the literature on Loop Quantum Cosmology.

I. INTRODUCTION
Modern Physics has two basic pillars in Quantum Mechanics and Einstein's theory of General Relativity (GR). However, the latter is a geometric description of the gravitational field that does not incorporate the principles of Quantum Mechanics. Numerous attempts have been made to construct a quantum theory of the spacetime geometry but, at present, there is still no proposal that the scientific community accepts unanimously as fully satisfactory. One of the proposals for the quantum description of gravity that has reached more impact with a robust mathematical development is Loop Quantum Gravity (LQG) [1][2][3]. It is a quantization formalism for globally hyperbolic spacetimes, based on a canonical and non-perturbative formulation of the geometric degrees of freedom. The fundamental novelty with respect to other pre-existing canonical proposals [such as the Wheeler-DeWitt (WdW) quantization, also called quantum geometrodynamics [4,5]] lies in the use of techniques imported from Yang-Mills gauge theories [6], known by their success in explaining non-perturbative regimes of the strong and electroweak interactions. In addition, the formulation of LQG is independent of any spacetime background structure and is aimed to respect the general covariance of Einstein's theory (in its canonical formulation). To achieve this goal, LQG adopts the quantization scheme proposed by Dirac for systems with constraints [7]. In particular, in GR the Hamiltonian is a linear combination of constraints which, via Poisson brackets, generate diffeomorphism transformations, that are the fundamental symmetries of the theory. Dirac's proposal consists in requiring that those constraints are satisfied at the quantum level on the physical states of the system. In more detail, the geometric degrees of freedom in vacuo are described in LQG by pairs of canonical variables that consist of the components of a densitized triad and a gauge connection [2,3]. Their respective fluxes through surfaces and holonomies form an algebra under Poisson brackets, which is the algebra that one wants to represent quantum mechanically over a Hilbert space, where the constraints of the theory should finally be imposed.
A major obstacle that the different candidates for a theory of quantum gravity have to face, regardless of their nature, is the extreme difficulty that is found to confront them with experimental data. Most of the possible effects of a quantum spacetime are expected to occur in regimes of very high curvatures or energies. In this sense, the Universe that we observe appears to be very classical, and GR explains it almost perfectly. However, there are observational windows to regimes of the Universe in which traces of a phenomenology that exceeded Einstein's theory might be found. A relevant example is the so-called Cosmic Microwave Background (CMB). This approximately black-body radiation reaches us from regions so far away that provides information about how the Universe was like at the early times when it became transparent. Under certain circumstances, this information might as well contain some details about very previous stages of the Universe when the spacetime geometry could have experienced quantum effects [8], especially if the observable Universe had in those epoques a size of the order of the Planck scale. A second obstacle for most of the quantum gravity proposals is the complication to extract concrete predictions about those regimes where quantum effects may have been important. Therefore, from a physical point of view, it is greatly convenient to consider the specialization of those formalisms to more restricted scenarios that, even without contemplating all the phenomena that may be accounted for in the full quantum theory, are able to describe regions of the Universe of particular interest. With this motivation, approximately at the beginning of this century, it was suggested to apply LQG methods to cosmological systems that possess a finite number of degrees of freedom, owing to the presence of certain symmetries. This effort crystallized in the appearance of Loop Quantum Cosmology (LQC) [9][10][11], a branch of LQG aimed to deal with the quantum analysis of cosmological systems.
The first cosmologies that were studied in LQC were homogeneous and isotropic universes of the Friedmann-Lemaître-Robertson-Walker (FLRW) type, that classically provide a good approximation to the behavior of the observed Universe in large scales. The quantization of this type of cosmological systems was consistently completed, providing satisfactory results both from a formal and from a physical point of view. Probably the most remarkable of these results is the resolution of the Big Bang cosmological singularity, that is replaced in this formalism with a quantum bounce, usually called the Big Bounce. In addition to these investigations on homogeneous and isotropic spacetimes, other homogeneous cosmologies with a lower degree of symmetry have also been considered in LQC to discuss the role of anisotropies. In particular, a special attention has been devoted to the quantization of Bianchi I models [12][13][14][15][16].
Although homogeneity and isotropy are very successful hypotheses to describe our universe at large scales, it is necessary to give an explanation to the existence and evolution of the observed inhomogeneities. In fact, the temperature of the CMB itself presents anisotropies that contain information about the small inhomogeneities in the geometry and matter content of the primeval Universe. Such inhomogeneities should be ultimately responsible of having given rise to the structures that we observe today [17]. As we have pointed out, it has been recently proposed that the power spectrum of the CMB anisotropies might even encode information about quantum effects that were relevant in the very early stages of the Universe, if the scale of the region that we observe nowadays was of the Planck order at that time [8,18]. Other information about those epochs of the Universe with extremely high curvature might be present in non-gaussianities of the CMB, or in the spectrum of tensor cosmological perturbations. obtains for those models. We then pass to discuss the hybrid approach in LQC, studying in Sec. 3 the cosmological system that was first analyzed in this quantization scheme. In the rest of sections, we focus our attention on the more interesting case (from a physical point of view) of a perturbed homogeneous and isotropic spacetime, in order to explore how quantum gravity effects may have affected the cosmological perturbations in the primeval universe. With this aim, we first review in Sec. 4 the procedure to construct a canonical formulation for this cosmological system in terms of gauge invariants and gauge constraints for the perturbations, together with their momenta, and of suitable zero modes for the background. In Sec. 5 we consider the possible introduction of a Dirac field in the formalism. We then explain in Sec. 6 the implementation of the hybrid LQC approach in this canonical system. Next, in Sec. 7 we discuss how we can derive modified propagation equations for the perturbations starting with the quantum Hamiltonian constraint and introducing a convenient ansatz for the quantum states, as well as some plausible approximations. These modified equations for the gauge invariants can be studied to deduce predictions that ultimately might be confronted with observations. In doing this, a key piece of information are the initial conditions that one must choose, both for the background cosmology and for the perturbations, in order to fix their vacuum state. These issues are discussed in Sec. 8. Finally, in Sec. 9 we explore the possible determination, or at least restriction, of the viable choices of a vacuum for the gauge invariant perturbations that result from demanding a good behavior in the quantum Hamiltonian operators and the evolution of those perturbations, putting an emphasis on a procedure of asymptotic diagonalization of their Hamiltonians. Section 10 contains the conclusions. In the rest of this work, we adopt units such that c = = 1, where c is the speed of light in vacuo and is Planck reduced constant. Nonetheless, we maintain Newton constant explicitly in all our formulas. Owing to these conventions and to some convenient redefinitions of quantities with respect to the notation employed in previous works, special care must be taken when comparing numerical factors in our equations with those appearing in the published literature.

II. LOOP QUANTUM COSMOLOGY
Let us first introduce the formalism of LQC, applied in this work to spacetime systems that in Einstein's theory correspond to homogeneous cosmologies. We will focus our attention on a flat FLRW model, minimally coupled to a homogeneous scalar field. In this section, we will review in detail the case of a massless field, because then the quantum constraints can be solved exactly. Later in our discussion, when we consider inflationary cosmologies, we will introduce a potential in the action of the scalar field, that can be viewed as the inflaton field of the system. An FLRW spacetime is the model typically used to describe the evolution of an expanding homogeneous and isotropic universe in GR. Here, we will study only the case of flat spatial curvature. In the following section we will also generalize our analysis to globally hyperbolic spacetimes (that admit a foliation on compact Cauchy hypersurfaces) of Bianchi I type.
LQC starts from a Hamiltonian formulation of the system under consideration, selecting as canonical variables for the geometry those used in LQG. In the Arnowitt-Deser-Misner (ADM) formulation of GR, given an arbitrary Cauchy hypersurface Σ, the dynamical degrees of freedom of the spacetime metric can be captured by the spatial metric induced on Σ, h ab (we use lower case letters from the beginning of the alphabet to denote spatial indices), and its variation along the normal surface vector. This variation is called the extrinsic curvature and is given by the tensor K ab = L n h ab /2, where L n is the Lie derivative along the normal vector n. Taking the spatial metric h ab as the configuration variable, a linear function of the extrinsic curvature determines its canonically conjugate momentum. Starting from this canonical pair for the geometry and canonical pairs corresponding to the matter content, if a Legendre transformation is carried out in the Hilbert-Einstein Lagrangian (with suitable boundary terms), one obtains a Hamiltonian that is a linear combination of first-class constraints. Their coefficients are nondynamical Lagrange multipliers, provided by the lapse function N and the three components of the shift vector N a . Each of these constraints vanishes on the solutions of GR. They are the generators of the fundamental symmetries of GR: the spacetime difeomorphisms. More specifically, the constraint that is multiplied by the lapse function is called Hamiltonian or scalar constraint, and generates time reparametrizations, modulo a spatial difeomorphism. The three constraints that come multiplied by the components of the shift vector are called momentum constraints, and generate spatial difeomorphisms.

A. Ashtekar-Barbero variables: Holonomies and fluxes
The Hamiltonian description of GR can be reformulated in terms of geometric canonical variables that, involving a gauge connection, simplify the functional form of the constraints [1]. Under the quantization scheme proposed by Dirac, these variables may seem more convenient for developing the quantum theory. In addition, the introduction of a gauge connection allows the controlled use of structures that are well known in group theory, and that can facilitate the construction of a well-defined Hilbert space. These variables can be introduced as follows.
First, in the spatial sections we can make use of triads, which are defined as a local basis of vectors e a i of the considered Cauchy hypersurface, and in terms of which the spatial metric can be expressed locally as where the co-triads e i a are the inverse of e a i . Since the Kronecker delta is the Euclidean metric in three dimensions, the relationship (2.1) is invariant under the transformation of the triadic basis under three-dimensional rotations, at each point of the Cauchy sections. Therefore, the use of co-triads for the description of the spatial metric automatically introduces additional local symmetry into the theory, provided by the group SO(3). Any Cauchy hypersurface is thus supplied with a principal fiber structure of three-dimensional reference systems, with SO(3) as the gauge group [57]. We employ lower case Latin letters from the middle of the alphabet for indices corresponding to components in a local triadic basis of orthonormal frames, or equivalently in a basis of the three-dimensional Lie algebra so (3). A section of the bundle is a locally smooth assignment of an element of the group to each point of the manifold. Different choices of triads, related to each other point to point by gauge transformations, can then be understood as different sections of the bundle.
In order to define a notion of horizontality between the different fibers, as well as the associated parallel transport, one introduces a gauge connection, characterized by a one-form on the spatial hypersurfaces with components that take values in the three-dimensional Lie algebra so(3). We will call this connection Γ i a . Of all the possible connections, there is one that is uniquely determined by the densitized triad through the metricity condition where ijk is the totally antisymmetric Levi-Cività symbol (its indices are raised and lowered using the Kronecker delta), (3) ∇ b denotes the covariant derivative associated with the Levi-Cività connection compatible with h ab , and E a i = √ he a i is the densitized triad, with h equal to the determinant of the spatial metric. In adddition, taking into account that the extrinsic curvature in triadic form can be understood as a vector of so(3) with respect to the gauge transformations, as well as a one-form in spacetime, it is possible to consider, instead of Γ i a , the so-called Ashtekar-Barbero connection [58]: In this definition, γ is a real non-vanishing number, of arbitrary value in principle, which is known as the Immirzi parameter [59]. The pair formed by this connection and the densitized triad, the so-called Ashtekar-Barbero variables, turns out to be canonical for GR: where δ 3 ( x − y) is the three-dimensional Dirac delta. LQG starts with these variables in the attempt to construct a quantum theory of gravity. Actually, in order to allow the coupling to the gravitational field of matter with a half-integer spin, A i a is considered as a connection that takes values in the three-dimensional Lie algebra su (2). That is, in practice the gauge group SO(3) of the principal bundle is replaced by its double cover, SU (2).
In terms of the Ashtekar-Barbero variables, the Hamiltonian constraint H and the momentum constraints H a of GR (in vacuo) take the form [3]: where F i ab is the curvature of the Ashtekar-Barbero connection: Finally, the introduction of an additional gauge symmetry in the theory translates into the appearance of three new constraints, that generate spin rotations in SU (2) (once this is considered as the cover of the group of three-dimensional rotations): Owing to their form as a divergence, given by the covariant derivative of the triadic field with respect to the connection A i a [57] and contracted in spacetime indices, these three constraints resemble the Gauss law of electromagnetism, and, accordingly, they are usually called the Gauss constraints. For the type of spacetimes with homogeneous spatial surfaces that we want to consider, and with a suitable choice of reference system, the Gauss and the momentum constraints are automatically satisfied, therefore involving no restriction on the system. From a systematic point of view, the first step in the construction of a quantum theory of gravity, based on the introduced Hamiltonian formalism with a gauge connection, would be to find a representation, as operators on a Hilbert space, of an algebra that captures all the relevant information about the canonical pair (2.5). Now, the formulation of the quantum theory must reasonably be such that physical results turn out to be described by quantities that do not depend on the choice of SU (2) gauge. With this purpose, it is convenient that the elements of the algebra of classical variables to be quantized do not vary, or vary as little as possible, under the SU (2) transformations that change the sections of the bundle. A well-known construction in Yang-Mills theories that captures the gauge invariant information about the connection is the holonomy. Given a curveγ on a spatial hypersurface Σ, the holonomy along it of the connection A i a is defined as follows: where P denotes path ordering and τ i provides a basis of the Lie algebra su (2). Holonomies determine the parallel transport defined by the Ashtekar-Barbero connection between the SU (2) fibers that are assigned to each point of the manifold. Given any section of the principal bundle and the curveγ, the holonomy dictates how this curve should be lifted to the fiber so that its tangent vector is parallelly transported [57]. Under a change of section, the holonomy is only affected by the gauge transformation evaluated at the end points of the curve. On the other hand, it is clear that its construction does not depend at all on any fixed spacetime structure, nor on the choice of coordinate system. All these properties of the holonomies make them good candidates to be the variables represented quantum mechanically in order to capture the relevant information about the configuration space of Ashtekar-Barbero connections. In LQG, one considers holonomies along edges e, typically piecewise analytic, defined as an embedding of the interval [0, 1] in our Cauchy hypersurface [2]. The variables that represent the rest of the phase space must contain the densitized triad E a i . Since this triad is a vector density in Σ, its Hodge dual can be directly integrated over two-dimensional surfaces S. The result is a flux through them, which again does not depend on any additional spacetime structure nor on the choice of coordinates, where f i is a function that takes values on the algebra su(2) and can be treated as a vector with respect to gauge transformations. The space of holonomies h e and fluxes E(S, f ) forms an algebra under Poisson brackets that no longer possesses the distributional divergences of the canonical relations (2.5). This is the algebra chosen in LQG to represent quantum mechanically the canonical commutation relations of GR.

B. Homogeneous cosmologies: Polymer quantization
We will now summarize the methodology used in LQC for the quantization of flat FLRW cosmologies, following a strategy inspired by LQG. First, let us recall that the momentum and Gauss constraints are trivial in this hogeneous system, setting at convenience the reference system for the description of the spatial metric h ab , as well as the internal gauge of the triads that determine it [see Eq. (2.1)]. For simplicity, we will choose spatial coordinates adapted to the homogeneity of the spatial sections and homogeneous diagonal triads, proportional to the Kronecker delta δ a i . For these triads, the connection Γ i a vanishes. We will also assume that the spatial hypersurfaces of the chosen foliation are compact, with a three-torus topology (T 3 ). This compactness ensures that spatial integrations do not give rise to infrared divergences (in non-compact cases, this problem can be handled by introducing fiducial structures [10,16]). In addition, restricting the study to compact hypersurfaces is very convenient if these cosmologies provide the homogeneous sector of other more general scenarios, because the application of the hybrid strategy for their quantization would use QFT techniques that are known to be well-posed and robust in the case of compact Cauchy sections.
Taking all these considerations into account, and choosing the compactification period in T 3 of each of the orthogonal directions adapted to homogeneity equal to 2π, the geometric sector of flat FLRW cosmologies can be described using Ashtekar-Barbero variables that adopt for them the specific form For any global function of time t, the canonical variables p(t) and c(t) can be classically related with the scale factor a(t), typically used in cosmology, and with its temporal derivativeȧ(t) by the equations Note that the geometric sector of the phase space has a finite dimension (equal to two), a fact that will greatly facilitate its quantum description. Inspired by the methodology of LQG, we construct holonomies that describe the degree of freedom c characterizing the connection. Thanks to the symmetries of the spatial sections, it is sufficient to consider straight edges e a of length 2πµ, with µ ∈ R, in the three orthogonal directions adapted to the spatial homogeneity [60]. The holonomies of the connection A i a along these edges have the simple expression where I is the identity in SU (2). Similarly, spatial symmetries allow us to restrict all our considerations just to fluxes of the densitized triad through squares, formed by edges along two of the reference orthogonal directions adapted to homogeneity. These fluxes are then completely determined by the variable p, that hence describes the geometric sector of the momentum space. Holonomies, or equivalently their matrix elements, describe the rest of the geometric sector of the phase space. More specifically, the geometric configuration space consists of the algebra formed by functions that depend on the connection through finite linear combinations of the complex exponentials N µ (c) = exp(iµc/2), with µ ∈ R. On the other hand, we recall that the mater content of our FLRW cosmology is given by a homogeneous (massless) scalar field φ. This scalar field is minimally coupled to the geometry. We will call π φ its canonically conjugate momentum. The canonical algebra that we want to represent has then the following non-trivial Poisson brackets: In LQC, the quantum representation of this algebra parallels the strategy adopted in LQG. In that theory, the geometric configuration space is described by means of the so-called cylindrical functions. These are functions that depend on the Ashtekar-Barbero connection through holonomies along graphs that are formed by a finite number of edges. The algebra of cylindrical functions is completed with respect to the supreme norm, obtaining a commutative C * -algebra with identity element. Gel'fand's theory guarantees that this algebra is isomorphic to an algebra of continuous functions over a certain compact space, called the Gel'fand spectrum, that contains the smooth connections as a dense subspace. The Hilbert space for the representation of the algebra of holonomies and fluxes in LQG is then that of square integrable functions on the Gel'fand spectrum, with respect to a certain measure [2,61].
In the homogeneous and isotropic scenarios that we are considering, on the other hand, the geometric sector of the configuration space, when completed with respect to the supreme norm, is the C * -algebra of quasi-periodic functions over R. The complex exponentials that describe the holonomies, N µ : R → S 1 , where S 1 is the circumference of unit radius, form a basis in it [61]. Its Gel'fand spectrum is the Bohr compactification of the real line, R B , that contains R as a dense subspace [62]. The space R B can be characterized as the set of all homomorphisms between the additive group of real numbers and the multiplicative group of complex unit module numbers. That is, every x ∈ R B is a map x : R → S 1 such that This space R B admits a compact topological group structure [61]. All functions F µ : R B → S 1 such that F µ (x) = x(µ), for any µ ∈ R, are continuous with respect to that topology. In addition, since it is a compact group, it admits a unique Haar measure M H , which is invariant under the group action. The Hilbert space for the representation of the algebra of holonomies and fluxes in homogeneous and isotropic LQC is then L 2 (R B , M H ). It follows that the set of functions {F µ , µ ∈ R} are an orthonormal basis of this Hilbert space [61], that is therefore not separable. Using Dirac's notation, we will denote this basis as {|µ , µ ∈ R}, where µ|µ = δ µµ is the inner product on L 2 (R B , M H ). The quantum representation of the algebra (2.15) that describes the gravitational sector of the phase space is [10,61] N µ |µ = |µ + µ ,p|µ = 4πGγ 3 µ|µ . (2.17) This representation is often called polymer representation, and its Hilbert space is isomorphic to that of functions over R that are square summable with respect to the discrete measure. Making use of this isomorphism, it is clear that the states of the polymer Hilbert space must have support only on a countable number of points, and, when this number is finite, they are the direct analogue of the cylindrical functions of LQG. Besides, note that the representation of the basic operators that describe the holonomies is not continuous. As a consequence, the operator that would represent the Ashtekar-Barbero connection is not well defined, a fact that also occurs in LQG.
At this point of our discussion, it may be worth noticing that the construction of L 2 (R B , M H ) as the Hilbert space of the representation strongly depends on the choice of the Haar measure. In fact, it is possible to find another measure in R B that results in a standard Schrödinger representation for the connections and triads [61]. This alternate representation, unlike the polymer one, is continuous and is employed in the more familiar WdW quantization of this cosmological system [4]. However, owing to the discrete character of M H , the two measures, and hence their corresponding quantum theories, are not equivalent. It is therefore not surprising that LQC can provide different predictions than traditional geometrodynamics about the quantum regimes of this cosmological model.
Finally, a standard continuous Schrödinger representation is chosen for the matter sector of the phase space, that can be described by the scalar field and its conjugate momentum. The corresponding Hilbert space is L 2 (R, dφ), where the scalar field acts by multiplication and its momentum as the derivativeπ φ = −i∂ φ .

C. LQC: Hamiltonian constraint
The Hilbert space obtained by the tensor product of the polymer space and L 2 (R, dφ) does not necessarily contain the physical states of the quantum theory. They should still satisfy the Hamiltonian constraint, that is the only non-trivial constraint that exists on the system, and which should be imposedà la Dirac quantum mechanically [7]. For this reason, the elements of the considered Hilbert space are often called kinematic states. The next step in our quantization is then the representation of the Hamiltonian constraint as an operator on the kinematic Hilbert space. The gravitational part of this constraint is given by Taking into account that the lapse function is homogeneous, we have already considered the integrated version of the constraint over the three spatial directions. In terms of the Ashtekar-Barbero variables introduced before for the geometric sector, the constraint H S that we obtain for homogeneous and isotropic cosmologies, including the contribution of the homogeneous massless matter field φ, has the following form (2.19) The first evident obstruction for a polymer quantization of this constraint is the absence of an operator to represent the connection. However, this difficulty can be surpassed if the following classical identity is taken into account: where the symbol tr(·) stands for the trace and h ab is the holonomy along a certain circuit that encloses a coordinate area A . For spacetimes such as flat FLRW and Bianchi I cosmologies, one can consider a rectangular circuit in the plane formed by the directions a and b. Thus, in our specific flat FLRW case, this holonomy can be written as (2.21) and the enclosed coordinate area is A F LRW = 4π 2 µ 2 . If these holonomies are represented polymerically, the limit contained in expression (2.20) is not well defined, because neither is the connection operator. Therefore, in LQC, the enclosed coordinate area is not made to tend to zero, but instead one appeals to the existence of a minimum value, characterized by the minimum coordinate length 2πμ of the edges that enclose it: A F LRW min = 4π 2μ2 . It seems clear then that this prescription for the quantum representation of the curvature introduces a new scale. The arbitrariness in its choice can be fixed by recurring to full LQG, where the geometric area operator has a minimum non-zero eigenvalue ∆. Drawing inspiration from this fact, in LQC one postulates that this value coincides with the geometric physical area corresponding to A F LRW min . Recalling the (second) classical relation in (2.13), one concludes then that the minimum coordinate length should satisfy [63] Once we have determined the scaleμ (now turned into a dynamical variable) that sets the minimum coordinate area in the FLRW cosmological model, we have to represent the classical expression (2.20) on the polymer Hilbert space, taking in it the limit A F LRW → A F LRW min . In practice, this prescription amounts to the replacement of c with the function sin(μc)/μ in the classical expression of the Hamiltonian constraint H S before one proceeds to its quantum representation. The dependence of the constraint on the connection is thus captured by the complex exponentials N ±2μ (c), that in particular depend on p throughμ. We must then specify their representation on the polymer Hilbert space, since the classical dependence of N ±2μ (c) on c and p makes their construction ambiguous in terms of the operators that we have taken so far as basic for the FLRW geometry, namelyN µ andp. With this aim, it is useful to introduce first the following operator, constructed by means of the spectral theorem [64]: where p(µ) is the eigenvalue ofp corresponding to the eigenstate |µ , given in Eq. (2.17). The direct classical counterpart of this operator is proportional to the physical volume of the FLRW flat and compact universe. In addition, it has a Poisson bracket with b =μc equal to minus two. If we relabel the orthonormal basis {|µ , µ ∈ R} using the eigenvalues v ofv, the operatorsN ±μ are then defined so that their action is simply a constant translation, namely they simply shift the new label by a constant [63]: The square of these operators definesN ±2μ . The prescription that we have explained in order to represent the elements of holonomies that appear in the Hamiltonian of LQC is commonly called the improved dynamics scheme 1 . However, this scheme alone is not enough to complete the representation of the constraint H S . The presence of the zero eigenvalue in the spectrum of the operatorp creates problems added to those already mentioned. Indeed, the Hamiltonian constraint depends on the inverse of the geometric variable p via ratios that involve the determinant of the spatial metric. The quantum representation of this inverse cannot be defined in the kinematic Hilbert space using the spectral theorem, because zero is part of the discrete spectrum ofp. This difficulty can be circumvented again by appealing to the following classical identity, employed as well in LQG adapted to a more general context [69]: The quantum operators that correspond to inverse powers of the triadic variable p are then defined by representing, on the improved dynamics scheme, the variables appearing on the right-hand side of this equality. In particular, Poisson brackets are represented by −i times the commutator of the corresponding operators. If one follows this quantization procedure, the operator representing the homogeneous Hamiltonian constraint H S in LQC is [70] H S = 1 where we have defined the densitized operatorĤ In these definitions, we have used a prescription for the factor ordering of the involved operators proposed by Martín-Benito, Mena Marugán, and Olmedo [70], known with the initials of these authors as the MMO prescription. Its most characteristic feature is the symmetrization of the sign of the orientation of the triad with the holonomy elements in Eq. (2.28). The operatorĤ S defined in this way presents certain interesting properties thanks to its symmetric factor ordering, compared with the quantum constraint obtained with another, frequently adopted prescription, proposed by Ashtekar, Pawlowski, and Singh (APS) [60,63,71]. In particular, its action decouples the state of the polymer basis with v = 0 from its orthogonal complement. This allows for a neat densitization of the Hamiltonian constraint [70]. Besides, the action of the operator does not mix states with positive and with negative values of v [70] (namely, it does not change the orientation of the triad). We can then restrict the quantum analysis of the FLRW cosmologies to the linear subspace generated by states |v with v ∈ R + , for example. Actually, the action of the constraint leaves invariant smaller and separable subspaces that are called superselection sectors, and that are simpler with the MMO prescription than in the APS case owing to the separation between sectors of different triad orientation. In more detail, the action ofĤ S (or, equivalently, ofĤ S ) turns out to preserve every one of the linear subspaces generated by states |v with v belonging to the semilattice L + ε = {ε + 4k, k ∈ N}, entirely characterized by its smallest point ε ∈ (0, 4]. Notice that the value ofε is always strictly positive, and therefore the same is automatically true for the variable v in each of the considered superselection sectors. Finally, a comment is due about the imposition of the (densitized) Hamiltonian constraintĤ S . Although the kernel of this operator is not a proper subspace of the kinematic Hilbert space, the quantum constraint can be rigorously imposed in our representation by its adjoint action, allowing in this way generalized solutions that should provide the physical Hilbert space once they are supplied with a suitable inner product, different from the kinematic one.
In fact, it has been possible to characterize the resulting physical Hilbert space, together with complete sets of Dirac observables. The resolution of the constraint is straightforward once one completes the spectral analysis of the operatorΩ 2 0 . It has been proven that this operator has a non-degenerate absolutely continuous spectrum equal to the positive real line [14,70]. The eigenvalue equation of the operatorΩ 2 0 can be regarded as a second-order difference equation. With the MMO prescription adopted in its definition, the generalized eigenfunctions turn out to be entirely determined by their value atε, point from which they can be constructed by solving the eigenvalue equation. Besides, up to a global phase, these eigenfunctions eε δ (v) are real, because the second-order difference equation associated with the action of the operator is a real equation. With the eigenfunctions at hand, one easily obtains the solutions to the Hamiltonian constraint, which take the form (2.30) Therefore, physical states can be identified, for instance, with the positive frequency solutions ξ + (δ) exp (i √ δφ) that are square integrable over the spectral parameter δ ∈ R + [70] (negative frequency solutions provide essentially the same Hilbert space). A complete set of Dirac observables is given byπ φ and |v| φ0 , where this latter operator is defined as the action of the operatorv when the scalar field equals φ 0 . On the Hilbert space of physical states specified above, these observables are self-adjoint, property that in fact characterizes the inner product on the space of solutions described by the functions ξ + (δ).
On the other hand, a numerical analysis of the dynamics, with respect to the homogeneous scalar field, of certain families of states with a semiclassical behavior at large volumes shows that they remain sharply peaked during the quantum evolution [63]. Actually, the trajectories of their peaks coincide, in the regions of low matter density, with those dictated by Einstein's equations in the considered FLRW cosmology. However, when the matter density reaches values that are comparable to the Planck density ρ Planck , the trajectory of the peak separates from the classical solution and turns to describe a transition from a universe in contraction to an expanding one (or vice versa) [63]. In particular, the matter density reaches a critical value when it is equal to 0.41 ρ Planck (for the value of the Immirzi parameter that leads in LQG to the Bekenstein-Hawking law for the entropy of black holes [72][73][74]). As we explained in the Introduction, this phenomenon of quantum nature that eludes the cosmological singularity of the Big Bang is known by the name of Big Bounce. There is also evidence that it occurs in LQC beyond the context of homogeneous and isotropic cosmologies, with indications that it is present as well in anisotropic cosmologies such as Bianchi I models [75], or in inhomogeneous cosmologies such as the linearly polarized Gowdy model with toroidal spatial sections [76].

III. HYBRID LQC: THE GOWDY MODEL
Historically, hybrid LQC was first developed for the Gowdy linearly polarized cosmological model, with spatial sections with the topology of a three-torus, T 3 , and after carrying out a partial gauge fixing that removes all constraints on the system except for the zero mode of the Hamiltonian constraint and of one of the momentum constraints (these zero modes are the average of those constraints on the spatial sections, modulo a constant numerical factor). Gowdy models are inhomogeneous cosmological spacetimes with compact spatial sections and two axial Killing vector fields [24,25]. The case with three-torus spatial topology is the simplest one. The linear polarization restriction on the gravitational wave content of the model amounts to require that the two Killing vectors are hypersurface orthogonal. The Killing symmetries then imply that the physical degree of freedom still present in those gravitational waves can be thought of as varying in only one spatial direction. After the mentioned partial gauge fixing, the phase space of this Gowdy model can be identified with that of a Bianchi I spacetime containing a linearly polarized wave. Furthermore, given the spatial behavior of this wave, we can describe it as a scalar field χ defined on the circle, S 1 , that corresponds to the cyclic spatial direction in which the wave propagates. In addition, we couple a massless scalar field Φ with the same symmetries as those of the geometry [26]. The hybrid quantization of this model will therefore be based on the quantization of Bianchi I cosmologies according to the LQC formalism (and the MMO prescription, see Refs. [20][21][22][23]), as well as on a suitable Fock representation of the matter field Φ and of the scalar field χ assigned to the linearly polarized gravitational wave. For convenience, we extract the zero mode φ of the matter field, that behaves as a homogeneous scalar field giving a non-trivial matter content to the Bianchi I cosmology, and assume that χ has vanishing zero mode (this assumption is only meant to simplify the discussion and involves no relevant conceptual consequences). The two sectors of the model, namely the homogeneous Bianchi I sector and the inhomogeneous scalar field sector, get mixed in a non-trivial way by the zero mode of the Hamiltonian constraint, that must be imposed on the considered system.
Let us describe the model in more detail. It is most convenient to choose coordinates {t, θ, σ, δ} adapted to the Killing symmetries, such that ∂ σ and ∂ δ are the Killing vectors. Then, the degrees of freedom of the model only depend on the time t and on the cyclic spatial coordinate θ ∈ S 1 . Starting with the canonical formulation of GR, we can then peform a symmetry reduction to take into account the Killing symmetries, as well as a partial gauge fixing that removes the momentum and Hamiltonian constraints except for the zero modes of the latter and of the momentum constraint in the θ-direction [21]. In this way, we get a reduced phase space that is formed by four canonical pairs of degrees of freedom (corresponding to the zero modes of the model), a gravitational field (that describes the linearly polarized gravitational wave of the model, and that we consider devoid of zero mode), and the inhomogeneous part of the massless scalar field. The homogeneous sector, composed of the four pairs of zero modes, coincides with the phase space of a Bianchi I cosmology coupled to a homogeneous massless scalar field φ. Besides, using a Fourier transform, we decompose the gravitational and matter scalar fields of the model, χ and Φ respectively, passing to describe them in terms of their Fourier (non-zero) modes. These modes and the corresponding modes of the canonical momenta of the two fields form the inhomogeneous sector of the system. The obtained reduced phase space is subject to two constraints, that were not eliminated in the process of partial gauge fixing. One of them is the zero mode of the momentum constraint in the θ-direction, H θ . This constraint generates rigid rotations in the circle coordinatized by θ, and imposes a restriction that affects exclusively the inhomogeneous sector. The other constraint that must be imposed is the zero mode H S of the Hamiltonian constraint of the system. This constraint H S is the sum of a homogeneous term H hom , that is the Hamiltonian constraint of the Bianchi I model, and an additional term H inh , that couples the homogeneous and inhomogeneous sectors and vanishes when the inhomogeneous sector is not present. The next step in our analysis consists in describing the Bianchi I cosmologies with three-torus spatial topology in terms of Ashtekar-Barbero variables, in order to quantize them by using LQC techniques. We can adopt a suitable internal gauge and adopt a diagonal form for the triads and connections. In this way, each of the Ashtekar-Barbero variables can be totally characterized by three homogeneous functions, that determine the diagonal components. We will call these functions p a and c a , corresponding to the densitized triad and su(2)-connection, respectively, and with a = θ, σ, δ. The only non-trivial Poisson brackets for them are {c a , p b } = 8πγGδ ab , so they form canonical pairs. We call H BI kin ⊗ L 2 (R, dφ) the kinematic Hilbert space for the Bianchi I model in LQC, where H BI kin denotes the polymer representation space of the Bianchi I geometries in LQC and L 2 (R, dφ) is the space of square integrable functions for the homogeneous scalar field, defined on the real line [16]. With this choice of Hilbert space, we adopt again a standard Schrödinger representation for the zero mode of the matter field, φ, so that its canonical conjugate momentum acts as a derivative,π φ = −i∂ φ . The construction of H BI kin , on the other hand, is similar to that explained for the FLRW geometry in LQC, except for the fact that we now have three pairs of canonically conjugated Ashtekar-Barbero degrees of freedom instead of only one. The inner product on this Hilbert space is discrete, so that the triad operatorsp a have a point spectrum equal to the real line. Defining the tensor product ⊗ a |p a (with a = θ, σ, δ) of the eigenvectors of each of the triad operators we obtain eigenstates |p θ , p σ , p δ that form an orthonormal basis of the Hilbert space H BI kin .
On the other hand, we can extend the improved dynamics proposal from the FLRW geometries to the Bianchi I model as proposed in Ref. [16], introducing in this way minimum coordinate length scalesμ a for each of the spatial directions. Explicitly, these length scales areμ with a = b = c and a, b, c ∈ {θ, σ, δ}. We then introduce the operatorsN ±μa to represent the holonomy elements N ±μa = exp(±iμ a c a /2) along an edge in the a-direction of coordinate length 2πμ a . These operators appear in the regularization of the curvature operator in the Hamiltonian constraint and are defined in a similar way as we did for the isotropic case in FLRW. The action of these holonomy operators on the states |p θ , p σ , p δ is rather complicated, since each of the length scalesμ a depends not only on p a , but also on the triad variables in the other directions. To simplify the expressions, it is convenient to relabel these states in the form |v, λ σ , λ δ , where Apart from an orientation sign, v is equal to 1/(2πGγ √ ∆) multiplied by the physical volume of the Bianchi I universe, volume that we will often call also the homogeneous volume. The action of the holonomy operatorsN ±μ θ just scale λ θ in such a way that the label v is shifted by the unit [16]. In full detail, we havê On the other hand, the holonomy operatorsN ±μσ andN ±μ δ additionally produce state-dependent scalings of λ σ and λ δ , respectively. For example, we havê To complete the ingredients that are needed for the hybrid quantization of the Gowdy model, we have to select a Fock quantization of the inhomogeneous sector. Actually, it has been proven that it is possible to single out a unique Fock quantization (given by a Fock representation and a Heisenberg dynamics for the background independent part of the fields), up to unitary equivalence, by imposing certain natural conditions, that require that the symmetry generated by H θ and the quantum evolution of the creation and annihilation operators be unitarily implementable [77,78]. In particular, this result removes the freedom of choice among the infinite number of inequivalent Fock representations, that may lead to different physics. Besides, the unitarity of the Heisenberg dynamics also imposes a concrete parametrization for the non-zero modes of both the gravitational field χ and the matter field Φ in terms of the background variables (namely the zero modes). Following these criteria, we represent the inhomogeneous sector of our hybrid model in Fock spaces F α (with α = χ, Φ) chosen in Refs. [77,78]. An orthonormal basis of each of these Fock spaces is provided by the n-particle states |n α = | · · · , n α −2 , n α −1 , n α 1 , n α 2 , · · · , where n α l denotes the occupation number of the field α in the mode l ∈ Z − {0}. Letâ (α) l andâ (α) † l denote the annihilation and creation operators of that mode, respectively. We can then reach a kinematic Hilbert space for the hybrid quantization of the Gowdy model by taking the tensor product Excluding the zero mode of the matter field temporally from our considerations, an orthonormal basis for the Hilbert space of the rest of the system is formed by the states |v, λ σ , λ δ , n χ , n Φ for all real values of the three first labels and all sets n χ and n Φ of integers with a finite number of non-vanishing elements.
Finally, we are in an adequate position to represent the constraints of the system as densely defined operators on this Hilbert space. Choosing normal ordering, the generator of the translations in the circle reads [26] (3.6) This constraint leads to the condition ∞ l=1 α=χ,Φ l(n α l − n α −l ) = 0 (3.7) on n-particle states of the inhomogeneities. Those states that satisfy the condition span a proper Fock subspace F p of F χ ⊗ F Φ . Let us now consider the quantum Hamiltonian constraint. We choose again normal ordering for the creation and annihilation operators of the inhomogeneous sector, while for the part of the constraint that acts on the homogeneous sector we choose a convenient symmetrization inspired by the MMO prescription [22]. Rational powers of the norm of the triad variables are represented adopting an algebraic symmetrization, which decouples the states of zero homogeneous volume v from their orthogonal complement. Like in the FLRW cosmology, this fact allows us to remove the states with vanishing homogeneous volume from our kinematic space, therefore eliminating the quantum kinematic analogues of the classical singularities with v = 0. Moreover, again like in the FLRW model, one finds that the action of the Hamiltonian constraint does not mix states with different orientations of any of the components of the triad or, equivalently, with different signs of the variables v, λ σ , and λ δ . Hence, as far as the constraints of the system are concerned, one can restrict all considerations, e.g., to the sector of strictly positive labels for the homogeneous geometry. Taking this into account, we redefine Λ a = ln(λ a ) for a = σ, δ so that the anisotropy labels continue to take values over the real line. The Hamiltonian constraintĤ S =Ĥ hom +Ĥ inh that one obtains with this hybrid quantization, after performing a densitization similar to that in the FLRW case, has the following form [22,23]: Here a, b ∈ {θ, σ, δ}. As we have already commented,Ĥ hom is the constraint operator for the Bianchi I model with a homogeneous massless scalar field in LQC, according to the MMO prescription [14]. On the other hand, the operator πGγ Θ a is the representation of c a p a , which is a constant of motion in the classical theory. We have defined The contribution of the inhomogeneities is captured byĤ F andĤ I . The operatorĤ F can be understood as a free-field Hamiltonian, that leaves invariant the n-particle states. It is given bŷ (3.12) The operatorĤ I may be interpreted as an interaction Hamiltonian that creates and annihilates an infinite collection of pairs of particles, while preserving the momentum constraintĤ θ . Explicitly, It is worth remarking that the inhomogeneities of both fields contribute to the constraints in exactly the same way. The action of the Hamiltonian constraint operatorĤ S does not relate all of the states with different values of v ∈ R + and Λ a ∈ R, with a = σ, δ. There are invariant subspaces in the Hilbert space spanned by those states. Each of these subspaces provides a superselection sector for the quantum theory. The superselection sectors in the homogeneous volume v are semilattices of step four, L + ε = {ε + 4n, n ∈ N}, determined by the initial pointε ∈ (0, 4], exactly as in the FLRW model. Note that, again, the homogeneous volume is bounded from below by a strictly positive quantity in each of these sectors. The superselection sectors in Λ a are more complicated. If we fix some initial data Λ * a andε, the values of Λ a in the corresponding sector (constructed by the repeated action of the Hamiltonian constraint) take the form Λ a = Λ * a + Λε, where Λε is any of the elements of a certain set Wε that is countable and dense in the real line [22]: (3.14) Here, k n m ∈ N and z ∈ Z ifε > 2, while z = 0 whenε ≤ 2. Given that the action of Θ a is considerably complicated, it has not been possible to elucidate yet whether this operator is self-adjoint. In spite of this, it is common to assume thatĤ hom is essentially self-adjoint [16] and that the same applies toĤ S [22,23]. Regardless of this, one can try to formally solve the constraints of the Gowdy model. The solutions turn out to be completely determined by the data on the section of the v-space defined by v =ε. Thanks to this fact, one can characterize the physical Hilbert space as the Hilbert space of such initial data, with an inner product that can be determined by imposing reality conditions on a complete set of observables [79,80]. In this way, one arrives to the space where H BI phys is the physical Hilbert space for Bianchi I cosmologies derived in Ref. [23].

IV. HYBRID LQC: COSMOLOGICAL PERTURBATIONS
After testing the viability of the hybrid quantization strategy in the Gowdy model, the approach was also applied to the discussion of a much more relevant scenario in cosmology, namely the study of primordial cosmological perturbations in the very early stages of the Universe. Using that the inflationary Universe is usually described as an FLRW cosmology that plays the role of a background where the perturbations develop and propagate, the idea was to quantize this background in the framework of LQC and treat the perturbations with the techniques of QFT in a curved spacetime. The hybrid approach then transforms the curved, FLRW classical background into a quantum spacetime with which the quantum field excitations corresponding to the perturbations coexist and interact by means of the gravitational constraints. For simplicity and for a better control of the mathematical techniques of QFT, we will again assume that the spatial sections are compact, with a three-torus topology. On the other hand, in order to isolate the perturbative degrees of freedom that do not depend on a possible perturbative diffeomorphism of the FLRW spacetime, that would result in a new identification of the background geometry, we will adopt a description in terms of perturbative gauge invariants. For cosmological scalar perturbations, one can employ the invariants introduced by Mukhanov and Sasaki (MS) [40][41][42] (considered as a pair of canonical fields). Gauge invariants are also the tensor perturbations, as well as the degrees of freedom of a Dirac field if it is present [81] (treating this field entirely as a perturbation). The description of the phase space of the perturbations can be completed with suitable redefinitions of the generators of the perturbative diffeomorphisms and canonical momenta of them. For the hybrid quantization, a piece of information that is most relevant as far as the inhomogeneities are concerned is the choice of a Fock representation for the gauge invariant fields. To restrict this choice and adopt a representation with especially appealing physical properties, we will still adhere to the criterion that the Fock quantization must allow a unitary implementation of the spatial symmetries of the model and of the Heisenberg dynamics associated with the creation and annihilation operators. With these ingredients, we will proceed to construct a hybrid quantum theory for the perturbed system. On this system, we will see that the only non-trivial constraint turns out to be the zero mode of the Hamiltonian constraint. We will then discuss its quantum imposition. Moreover, we will show how to extract from it (with a convenient ansatz and plausible approximations) Schrödinger equations for the perturbations, as well as effective equations to describe the propagation of the perturbations on the FLRW geometry subject to quantum effects. These equations can be used to study modifications to the power spectra of the cosmological perturbations, originated from quantum gravitational effects. The program that we have outlined will be implemented in this and the following five sections.
We start by constructing a convenient canonical description of the system formed by the FLRW cosmology and its perturbations that contains a complete set of gauge invariants. As in the case studied in our introduction to LQC, the FLRW spacetimes that we will consider possess compact sections with the topology of a three-torus. Their geometry can be described by a scale factor a and its canonical momentum π a (or equivalently by the pair of variables c and p that determine the Ashtekar-Barbero variables in LQC). With the same choice of reference system for this cosmological background that we employed in the previous expositions about LQC, the coordinate volume of the spatial sections equals 8π 3 . As before, these spacetimes will contain a homogeneous scalar field, φ, responsible of the expansion and that consequently will play the role of an inflaton. This inflaton can be interpreted as the zero mode of a generally inhomogeneous scalar field Φ, interpretation that will be especially useful at the moment of introducing perturbations in the system. On the other hand, the main difference with respect to our previous studies is that we will now allow the possible existence of a potential V(φ) for this inflaton.
The FLRW system is subject only to a non-trivial homogeneous Hamiltonian constraint, as we have discussed above. It can be written as H 0 = 0 where 2 Let us next introduce perturbations in this system, both for the geometry and for the matter scalar field. It is also possible to introduce a Dirac field to describe fermions, regarded as perturbations of the FLRW cosmology [81]. We postpone the consideration of these fermions to the next section. We can separate the metric and inflaton perturbations into scalar, vector, and tensor, depending on their behavior under the symmetries of the spatial sections (notice that these symmetries provide the Euclidean group in the limit in which the sections become non-compact). In addition, using that the spatial Laplacian (and the Dirac operator) corresponding to the auxiliary Euclidean metric 0 h ab (with unit determinant) defined on our toroidal sections respect these symmetries, we can expand the different perturbations in eigenmodes of this differential operator. Moreover, since the spatial sections are compact, these modes are discrete. In this way, we can deal with the spatial dependence of the perturbations by considering infinite sequences of modes. For instance, choosing again (orthogonal) spatial coordinates of period equal to 2π, we expand the scalar perturbations of the metric and the matter field in a Fourier basis of sines and cosines, Here, the vector notation θ stands for the spatial coordinates (θ, σ, δ), and the Euclidean scalar product has been denoted with a dot symbol. Each mode is characterized by a wavevector k ∈ Z 3 − {0}, with strictly positive first nonvanishing component. Note that, in this way, we are not including the zero mode, that is part of the degrees of freedom considered in the FLRW background. The eigenvalue of the spatial Laplacian corresponding to k is −ω 2 k = − k · k. Scalar perturbations are described then by the corresponding Fourier coefficients of the scalar field Φ (without the zero mode, namely the inflaton), the trace and traceless scalar parts of the spatial metric h ab (without the FLRW contribution), the lapse N (without its homogeneous part), and the scalar part of the shift N a .
Similarly, tensor perturbations are described by the Fourier-like coefficients of the tensor part of the spatial metric, with two possible polarizations. These coefficients arise from the expansion in terms of the real tensor harmonics G k,ε,± ab , eigentensors of the spatial Laplacian [82]. As above, the tuple k can take here any value in Z 3 − {0}, with positive first non-vanishing component, while ε is the dichotomic label that specifies the polarization, and the superscripts ± indicate whether the harmonic is even or odd under a periodic translation of θ, as in the scalar case. Vector perturbations, on the other hand, are described by the remaining parts of the shift and the spatial metric, that can be expanded in Fourier-like coefficients in terms of eigenvectors S k a of the Laplace operator and of tensors obtained from those by spatial covariant derivatives [46]. Here, k is again any non-vanishing tuple of integers. It is convenient to parametrize all of these mode coefficients as follows: (4.7) A vertical bar stands for the spatial covariant derivative with respect to the Euclidean metric 0 h ab , and a parenthesis enclosing two spatial indices indicates symmetrization. Thus, the scalar perturbations are determined by a k,± , b k,± , g k,± , l k,± , and f k,± , whereas the tensor perturbations are described by the coefficients d k,ε,± , and the vector perturbations by c k and v k . We have normalized some of these coefficients in a convenient way to absorb several factors in the formulas that we will use in our discussion. Inserting these expressions in the Hilbert-Einstein action minimally coupled to the scalar field Φ (with suitable boundary terms) and truncating the result at quadratic order in the coefficents of the perturbations, it is possible to reach a Hamiltonian formulation for our system [46,82]. In this formulation, the above coefficients for the perturbations either play the role of Lagrange multipliers of some of the constraints, or form a canonical set together with the FLRW scale factor, the inflaton, and suitable momenta for all of them. In other words, at the order of our truncation in the action, the system formed by the FLRW cosmology and its perturbations is a totally constrained system that admits a canonical symplectic structure [34]. It is worth emphasizing that, at the considered truncation order, we are treating exactly the zero modes that determine the FLRW background, so that the perturbations that we have expressed explicitly do not contain zero modes. On the other hand, with the kind of matter content considered in our discussion, it is possible to show that the vector perturbations do not include any physical degree of freedom, but are pure gauge. To simplify our exposition, we will therefore eliminate them from our analysis in the following.
The perturbed system that we have constructed is subject to two types of constraints. On the one hand, the perturbations of the momentum and Hamiltonian constraints lead to a collection of constraints that are linear in the perturbations, and that appear accompanied by Lagrange multipliers that are also linear perturbative factors. Explicitly, the mode H k,± ↑1 of the linear perturbative momentum constraint has Lagrange multiplier given by the coefficient l k,± of the perturbations of the shift vector, while the mode H k,± 1 of the linear perturbative scalar constraint adopts as Lagrange multiplier the coefficient g k,± of the perturbation of the lapse. These constraints depend exclusively on the scalar perturbations of the metric, once the vector perturbations have been gauged away. A different type of constraint is the zero mode of the Hamiltonian constraint, which can be considered a global restriction on the system formed by the FLRW cosmology and the perturbations. This constraint has a Lagrange multiplier given by the homogeneous lapse function N 0 , and is the sum of two contributions: a term that reproduces what would have been the constraint H 0 of the FLRW cosmology in the absence of perturbations, and an additional term H 2 that contains the perturbative contribution and that is quadratic in the perturbations. This latter term is composed in turn of two parts, (s) H 2 and (T ) H 2 , respectively formed by the contributions of the scalar and the tensor perturbations. In this way, the total Hamiltonian takes the expression Moreover, the quadratic perturbative contributions to the zero mode of the Hamiltonian constraint can be decomposed as the sum of the contributions of each of the modes of the perturbations as follows: (4.9) The variables that we have chosen to describe the perturbative degrees of freedom have the drawback that they do not commute with the linear perturbative constraints, even when the FLRW cosmology is taken as a fixed entity with vanishing Poissson brackets. As a consequence, those variables would change if one performs a perturbative diffeomorphism, that would alter the form of the FLRW background without affecting the physics. To avoid this problem with the physical identification of the background, it is most convenient to use a set of variables that indeed commutes with the linear perturbative constraints when the zero modes are frozen in the computation of Poisson brackets. This leads us to consider gauge invariants for the perturbations. In the case of flat spatial sections, the gauge invariant degrees of freedom of the scalar perturbations are usually described in cosmology employing MS invariants, because they are straightforwardly related to the co-moving curvature perturbations. The variables that we have introduced for the tensor perturbations, on the other hand, are directly gauge invariant, and we will only redefine them linearly to re-express their dynamical contribution to the Hamiltonian in a convenient way.
With this motivation, we are going to introduce a change of variables for the perturbations, from the canonical set that we have been using to a new set formed by the following variables [34,43]: • The mode coefficients of the MS gauge invariant field, ν k,± . Explicitly, they are given by the formula [31,33] We notice that these coefficients mix the scalar perturbations of the metric and the perturbations of the matter scalar field.
• The mode coefficients of the tensor perturbations conveniently rescaled [37]: This rescaling simplifies the dependence of the Hamiltonian on the tensor perturbations.
• The mode coefficients π ν k,± and πd± k,ε of the canonical momenta of the above fields, defined also as gauge invariants. There exists a certain ambiguity in the specification of these momenta, since once can always add a contribution that is linear in the configuration fields, multiplied by any function of the FLRW background. A convenient criterion to fix this contribution is to require that the time derivative of each of these momenta, as dictated by Hamilton's equations, is proportional to the corresponding configuration variable. This condition amounts to demand that the Hamiltonian that generates the dynamics of the scalar and tensor perturbations should not contain cross terms between the configuration fields and their momenta, and turns out to determine the latter of these variables completely.
• An Abelianization of the linear perturbative constraints. Actually, at the order of our perturbative truncation in the action, it is possible to modify these constraints on the scalar perturbations with terms that are linear in those perturbations and such that the new constraints that one obtains commute under Poisson brackets among them, as well as with the MS field and its momentum, after freezing the zero modes. To achieve this Abelianization, it suffices to introduce the replacement This new linear perturbative scalar constraint is used together with H k,± ↑1 as additional variables in our canonical set.
• Suitable momenta of the Abelianized linear perturbative constraints. As far as those constraints generate gauge transformations consisting of perturbative diffeomorphisms, their momenta can be interpreted as variables that parametrize possible gauge fixations for the perturbations. A especially simple choice is Remarkably, the introduced change of variables for the scalar and tensor perturbations can be completed into a canonical transformation for the entire system (that is, without freezing the background), at the considered truncation order, by modifying the zero modes with terms that are quadratic in the perturbations [34,44,45]. For this, we can proceed as follows. We substitute the old perturbative variables in the Legendre term of the action (or, equivalently, in the symplectic potential) as functions of the new ones and, after several integrations by parts and convenient identifications of factors, we find new zero modes that keep the canonical form of the Legendre term up to perturbative contributions that are negligible in our truncation scheme. The new configuration variables obtained in this way adopt the generic expressionw where we have called w ι q = {a, φ} the configuration variables of the zero mode sector, {w p ι } are their momenta (ι = 1, 2), and {X k,± qm , X k,± pm } are the old variables for the scalar and tensor perturbations, each of which is given by a different value of the label m. A tilde on top of any of these canonical quantities indicates its new counterpart, defined according to the above procedure.
In the case of the momentum variables for the zero modes, the change is given by a formula of the same kind, but with a flip of sign in the term that provides the corrections quadratic in the perturbations, The availability of a canonical set for the entire perturbed system, formed by the FLRW cosmology and the perturbations, is of the greatest importance. In particular, it makes possible an easy implementation of the hybrid strategy following canonical quantization rules. But, in order to proceed to this quantization, we still have to determine the form of the zero mode of the Hamiltonian constraint (the only non-linear perturbative constraint of the system) in terms of the new canonical set, keeping the quadratic truncation order. In order to do this, we notice that, since the change of zero modes is quadratic in the perturbations, an expansion of the FLRW contribution H 0 around the new zero modes leads immediately to the desired constraint if we only include the first derivative contribution. Let us introduce the compact notation X k,± m = X k,± qm ,X k,± pm . (4.17) Then, according to our comments, the expression of the new global scalar constraint at our truncation order is [34,37] with the phase space dependence of H 0 , its derivatives, (s) H k,± 2 , and (T ) H k,ε,± 2 evaluated directly at (w ι ,X k,± m ). Namely, in formula (4.18), the evaluation of the Hamiltonian functions must be made as if one identified the old and the new set of variables. As a consequence, the contribution of each of the modes of the perturbations to the new global scalar constraint is where the symbol ∼ indicates equality modulo the Abelianized linear constraints and up to the relevant perturbative order in our truncation. Besides, we have called where the superscripts (s) and (T ) stand for the quadratic contributions of scalar and tensor nature, respectively. It is possible to prove that the sum of contributions in the left-hand side of Eq. (4.19) gives precisely the MS Hamiltonian [34], i.e. the Hamiltonian that generates the dynamical evolution (in proper time) of the MS field on the FLRW background. Likewise, the sum (T )H k,ε,± 2 of the tensor contributions to the constraint provides a dynamical Hamiltonian of harmonic oscillator type for the tensor perturbations on the FLRW cosmological background.
It is worth pointing out that, in the definition of these perturbative Hamiltonians, we can replace the squared momentum of the inflaton with π 2 φ − 16π 3 a 3 H 0 . This is so because all the new terms proportional to H 0 that are produced in this way are quadratic in the perturbations and can thus be absorbed in a redefinition of the zero mode of the lapse function up to a modification of the total scalar constraint that is at least quartic in those perturbations. Hence, such a modification is negligible at our truncation order. The freedom available in using this replacement can be fixed by requiring that the perturbative contribution to the Hamiltonian constraint be at most linear in the inflaton momentum, because one can always use the commented replacement to decrease the polynomial order in π φ by two units until one reaches either a linear contribution of the inflaton momentum or a term that is independent of it. The total Hamiltonian of the system then becomes [34] whereN 0 is the suitably redefined homogeneous lapse function that differs from the original one, N 0 , in perturbative terms that are quadratic. Similarly,g k,± andl k,± are Lagrange multipliers that differ from the original ones, g k,± and l k,± respectively, by linear perturbative contributions. Their explicit expressions can be found in Ref. [34], but they are not especially relevant for the rest of our discussion. Notice that this total Hamiltonian, imposed as a collection of constraints on the system, would include backreaction at the considered perturbative order. As we have commented, the contribution of the scalar perturbations to the global Hamiltonian constraint is nothing but the MS Hamiltonian, which is a sum of terms that are quadratic in the MS configuration variables and of quadratic terms in their momenta, but without terms that mix these two types of variables. This is a consequence of our choice of MS momentum field, as we explained when we introduced the new perturbative variables. A similar behavior is found in the contribution of the tensor perturbations to the Hamiltonian constraint. In more detail, Here, s (s) + r (s) πφ and s (T ) play the role of effective background dependent masses for the scalar and the tensor perturbative modes, respectively. The expressions of these background functions are The prime symbol denotes de derivative of the potential V with respect to the inflatonφ, and We notice that In particular, substituting Eq. (4.28), we see that s (s) = s (T ) when the inflaton potential V vanishes.
In total, after a convenient change of densitization similar to that explained in homogeneous and isotropic LQC (via multiplication by the homogeneous physical volume V = 8π 3ã3 ), we obtain a Hamiltonian constraint that can be written in the form where we have introduced the notation that explicitly separates the linear term in the inflaton momentum. Clearly, we have the identities with ϑ e = 8π 3ã2 . We note that there is no tensor contribution to Θ o . It is also worth remarking that all the ϑ-functions are independent of the particular mode that one considers. Besides, the part of the Hamiltonian constraint that contains the perturbative contribution is the same for the scalar and for the tensor perturbations except for the difference in their background dependent mass. This shows up in the appearance of the terms (s) Θ k,± o and in the difference between (s) ϑ q e and (T ) ϑ q e .

V. HYBRID LQC: INCLUSION OF FERMIONS
In the matter content of our cosmological system, we can also include fermionic fields, e.g. a Dirac field. Their pressence does not modify much our treatment if we consider them as perturbations, including the possible fermionic zero modes, so that they do not alter the dynamics of the homogeneous background cosmology in the linearized theory. Since the Dirac action is quadratic in the fermionic field, as a perturbation it couples directly only to the background FLRW geometry, but not to the perturbations of the metric, nor to the matter scalar field [81]. Moreover, for the same reason, the fermionic field does not contribute to the linearized perturbative constraints, that arise from the perturbation of the Hamiltonian and momentum constraints. As a consequence, the fermionic field can be treated as a gauge invariant perturbation at the considered order of truncation. This simplifies the formulation considerably.
If we adopt a Weyl representation [83], the Dirac field can be described by a pair of two-component spinors of definite chirality. We will call ϕ A and χ * A the respective left-handed and right-handed spinors associated with the field. Capital Latin letters from the beginning of the alphabet, both primed and unprimed, take values equal to 1 or 2, corresponding to the two components of the chiral spinors. These indices will be raised and lowered using the antisymmetric symbols AB and AB (with e.g. 12 = 1) , as well as their counterparts for right-handed chirality. It is most convenient to adopt an internal gauge such that the spatial part of the tetrad has vanishing temporal Lorentz components, namely e a 0 = 0. As a consequence of this gauge fixation on the spin structure in four dimensions, the two-component spinors of the Dirac field can be viewed as families of cross-sections of a spinor bundle defined on the compact spatial sections. On the other hand, the Hamiltonian formalism of the Dirac field is initially complicated by the existence of second-class constraints that relate the field with its momentum. Nevertheless, one can eliminate these constraints and capture the canonical anticommutation relations of the Dirac field in anticommutators of its two-component spinors. To take into full account this anticommuting character, we will treat these components as Grassmann variables [84].
In a similar way as we did with the metric and the scalar field perturbations, we can decompose the spinors of the Dirac field in modes. Since the spatial differential operator that appears naturally in the dynamical equation of our fermionic field is the Dirac operator constructed with the auxiliary Euclidean triad 0 e a i on the toroidal spatial sections of our model (with 0 e i a 0 e j b δ ij = 0 h ab being the Euclidean metric introduced above), it is logical to treat the spatial dependence of the field by an expansion in eigenmodes of this Dirac operator. The spectrum of this operator is discrete, owing to the compactness of the sections. The eigenvalues are ±ω k , where ω 2 k = k · k and k ∈ Z 3 is any tuple of integers. We are assuming a trivial spin structure on the spatial sections. Otherwise, the definition of ω k would include a constant displacement of k characteristic of the specific spin structure chosen for the fermions [85].
Using these modes, we can express the two-component spinors of the Dirac field in the form Here w k,(±) A are the left-handed Dirac eigenspinors with respective eigenvalue equal to ±ω k . With our choice of the auxiliary Euclidean triad, and recalling that we have assumed a trivial spin structure, these eigenspinors take the expression Each of the coefficients m k , s k , t k , and r k forms a Grassmann canonical pair with its respective complex conjugate [83]. Furthermore, in this sense they provide a canonical set together with the variables introduced in the previous section for the metric and scalar field perturbations and for the FLRW cosmology, once we have adopted a description of the cosmological perturbations in terms of gauge invariants [81]. For convenience, in the following we will employ the notation (x k , y k ) to refer to any of the ordered pairs of coefficients (m k , s k ) or (t k , r k ).
As we have already commented, there is no fermionic term in the linear perturbative constraints of our system, so that the only contribution of the Dirac field to the total Hamiltonian is included in the zero mode of the Hamiltonian constraint. This contribution is given by the Dirac Hamiltonian, evaluated at the variables for the cosmological zero modes defined in the previous section and at the fermionic perturbations determined by the variables (x k , y k ), as far as the respective dependence on the FLRW cosmology and the Dirac field is concerned. Motivated by previous works on this subject [83], the first approach to the treatment of fermions in hybrid LQC was to carry out a change of fermionic variables that produces a diagonalization of the Dirac Hamiltonian. To reach this diagonalization, the change of variables must depend on the FLRW geometry, a situation that is similar to that studied when we introduced gauge invariants to describe the relevant degrees of freedom of the scalar perturbations. This has two consequences, as we know. First, the zero modes have to be corrected to maintain the canonical structure in the set of variables that describe the whole of the cosmological system. Second, the fermionic contribution to the Hamiltonian constraint gets an additional term, up to quadratic order in the perturbations, owing to the background dependence of the change of variables. Even if this change was designed to diagonalize the Dirac Hamiltonian, it will generally not diagonalize the new fermionic contribution, and therefore the resulting fermionic Hamiltonian will still contain interacting terms. Actually, the new fermionic variables have Poisson brackets of the creation-annihilation type, so that we can view these fermionic interactions as the creation or annihilation of pairs of particles. Finally, we will treat the fermionic zero modes on their own, keeping their description in terms of the original variables (x 0 , y 0 ) to avoid problems with the particularization of our formulas to a vanishing Dirac eigenvalue (i.e., for ω k = 0). Explicitly, and leaving aside those zero modes, the new variables are given by with M denoting the mass of the Dirac field. Notice that the sum of the square modulus of the coefficients in each of the above linear combinations of the variables (x k , y k ) equals the unit. This ensures that the transformation is canonical in the fermionic phase space [81]. In a Fock representation with a standard interpretation, the operators representing a (x,y) k and b (x,y) k would annihilate particles and antiparticles, respectively, while their adjoints (representing the complex conjugate variables) would create them.
The fact that our change of fermionic variables depends only on the scale factor implies that we only need to modify the momentum of that background variable in order to recover a canonical set for the entire cosmological system. The modification of the momentum πã consists in adding to it the following terms that are quadratic in the fermionic perturbations, obtained in a similar way as it was explained in the previous section for the scalar and tensor perturbations [81]: For simplicity in our notation, we will denote the new momentum of the scale factor with the same symbol as before. From the context, it must be clear in our discussion whether we are referring to the original momentum or to the momentum that has been changed with the addition of fermionic contributions. On the other hand, we also notice that the variables for the scalar and tensor perturbations need not be altered at this stage, because our change of fermionic variables is independent of them. In terms of this new canonical set, the total Hamiltonian has the same expression (4.22) as before except for two things. First, its dependence on πã must be evaluated at the new momentum of the FLRW scale factor, which includes the fermionic modification. And second, the zero mode of the Hamiltonian constraint includes one additional contribution (F ) H 2 which is due to fermions, In consonance with our comments above, this fermionic contribution is given by the sum of the Dirac Hamiltonian

VI. HYBRID QUANTIZATION OF COSMOLOGICAL PERTURBATIONS: IMPLEMENTATION
Once we have at our disposal a canonical set of variables for the description of our perturbed cosmological model in which the variables that describe the perturbations are either gauge invariants, perturbative gauge generators, or associated gauge degrees of freedom, we are in an appropriate situation to face the quantization of the system. We will carry out this quantization adopting the hybrid approach within the framework of LQC. As we have already commented, this hybrid strategy is based on the hypothesis that the most relevant effects of quantum geometry for cosmology are those that affect the FLRW substrate, namely the behavior of the scale factor, while the purely quantum geometric effects on the perturbations can be approximately ignored, and handle the quantum description of those anisotropies and inhomogeneities using techniques directly related with the formalism of QFT in a curved spacetime, generalized to the case in which such a spacetime is quantum mechanical as well. In practice, we quantize the FLRW cosmology using the methods of LQC and the perturbations (essentially) with Fock quantization methods, then combine both types of quantum descriptions by adopting a tensor product representation space for the system, and finally impose on it the diffeomorphism constraints that are present in the system. We will see that, among them, the only intrincate constraint is the zero mode of the Hamiltonian constraint, that relates in a complicated way the FLRW cosmology with the scalar and tensor gauge invariant perturbations (as well as with the fermionic ones if we also consider a Dirac field). Thus, the hybrid quantization is non-trivial precisely because of the imposition of this constraint.
The linear perturbative constraints obtained from the Abelianization of the perturbations of the diffeomorphism constraints can be imposed straightforwardly by representing them as derivative operators with respect to their canonically conjugate degrees of freedom (or, if one considers an integrated version of these constraints, as operators that displace the values of such canonically conjugate degrees of freedom, resulting in transformations that should be symmetries). With this representation, the states that satisfy such constraintsà la Dirac are simply those that are independent of the gauge degrees of freedom πH k,± 1 and π H k,±

↑1
[see Eq. (4.13)]. In other words, physical states depend only on (a complete set of compatible) zero modes and gauge invariants. Note that this result is obtained without the need to introduce any perturbative gauge fixing. Physical states must still satisfy one constraint, that is the only one remaining at this stage, namely the zero mode of the Hamiltonian constraint, as we anticipated.
The desired quantum formulation is then reached by choosing the representation space of the improved dynamics scheme of homogeneous and isotropic LQC, H grav LQC , for the perturbatively corrected volume and its momentum (namely, the zero modes of the FLRW geometry once they have been suitably modified with terms that are quadratic in the perturbations in order to maintain the canonical symplectic structure of the entire cosmological system). For the perturbatively corrected inflaton and its momentum, we use a standard Schrödinger representation L 2 (R, dφ). On the other hand, for the MS and tensor gauge invariants, we adopt Fock representations F s and F T , chosen within a unique privileged family of unitarily equivalent representations that are characterized by [86,87]: • The invariance of the vacuum under the symmetries of the spatial hypersurfaces.
• A unitarily implementable Heisenberg evolution of the creation and annihilation operators, in the context of QFT in a curved background. This Heisenberg evolution is determined by the dynamics of the gauge invariant modes that we have picked out for the description of the perturbations.
In addition, if the system contains a Dirac field, viewed as a fermionic perturbation, we employ for it a Fock representation F D in the equivalence class of the one that is naturally associated with the previously introduced choice of creation and annihilationlike variables proposed by D'Eath and Halliwell [83]. This again belongs to a uniquely distinguished class of unitarily equivalent representations characterized by the same two conditions that we have listed above, together with the requirement of recovering a standard notion of particles and antiparticles. It is worth emphasizing that the choice of Fock representation (or of a family of unitarily equivalent representations) does not determine a concrete choice of vacuum state. Any Fock state in our representation is valid for this purpose. Therefore, in order to fix a unique vacuum state, more restrictions are needed, either in the form of additional requirements about the physical properties of the subsequent quantum theory or in the form of conditions on a particular spatial section able to specify the state there. Let us continue with our hybrid approach, thus adopting as representation space the tensor product H grav LQC ⊗ L 2 (R, dφ) ⊗ F s ⊗ F T ⊗ F D . We construct our quantum representation so that the zero modes commute with the perturbations, as it happens under Poisson brackets in the classical theory, and so that all functions of the inflatoñ φ act by multiplication. As we have commented, the zero mode of the Hamiltonian constraint results in a nontrivial coupling between the various factors of our tensor product. The quantization proposed for this constraint is based on the representation adopted in homogeneous LQC 3 . For the FLRW contribution H 0 we adopt the same prescription as in LQC. In particular, we adhere to the improved dynamics proposal, so that quantities that depend on the momentum of the scale factor are represented in terms of holonomies defined employing squares with a fiducial length that guarantees that the physical area enclosed by them coincides with the area gap ∆, determined by the area spectrum of LQG. Using the resulting basic operators of homogeneous LQC, as well as the homogeneous physical volume operatorV = 2πGγ √ ∆|v| and the regularized inverse volume operator obtained from them, we get We recall that V(φ) is the inflaton potential, and thatΩ 0 was defined in Eq. (2.28).
As for the functions of zero modes of the FLRW cosmology that appear in the perturbative part of the constraint, we adopt a symmetric factor ordering that tries and respects, as far as possible, the assignation of representation from homogeneous and isotropic LQC. In more detail, we adopt the following rules for their quantum representation: • We symmetrizeà la Weyl the representation of the product Θ o π φ , to deal with the presence of functions of φ in Θ o that do not commute with the inflaton momentum.
• We adopt an algebraic symmetrization for factors of the form V r g(b), that are promoted to the operatorŝ V r/2ĝ LQCV r/2 , where r is any real number and g(b) a function of the variable b, withĝ LQC its operator counterpart in the improved dynamics scheme of LQC. This algebraic symmetric factor ordering is adopted as well for powers of the inverse volume.
• Even powers of −ãπã 4πG/3 are promoted to even powers of the operatorΩ 0 , which represents this quantity in LQC, whereas odd powers, let's say, of order 2z + 1, with z any integer, are represented as |Ω 0 | zΛ 0 |Ω 0 | z . Here |Ω 0 | is the square root of the positive operatorΩ 2 0 , andΛ 0 is defined exactly asΩ 0 , but with holonomies of double length. The result can be obtained by dividing the right-hand side of expression (2.28) by 2, and replacing b in that expression with 2b. The operatorΛ 0 defined in this way only shifts v in multiples of four units, and hence preserves the superselection sectors of the homogeneous and isotropic geometry.
With these prescriptions, we arrive at the following operator representation of the functions (4.36) that appear in the densitized Hamiltonian constraint: According to these formulas, the counterpart of relation (4.29) between the scalar and tensor background dependent masses is It is worth remarking that (s)θ o is proportional to the derivative of the inflaton potential, so that one expects its contribution to be negligible when the dependence of the potential on the inflaton is not important. The operators representing the phase space functions (4.31) can be constructed with the above operators and the Fock representation adopted for the modes of the MS and the tensor gauge invariants. In a completely similar manner, one can construct an operator representation for the fermionic contribution H F to the densitized zero mode of the Hamiltonian constraint (obtained from the original one by multiplication with the homogeneous volume), that depends only on the FLRW geometry and the Dirac field, but not on the inflaton nor on its momentum. For more details about this fermionic contribution, we refer the reader to Refs. [34,81]. In this way, we finally get

VII. HYBRID LQC: MODIFIED PERTURBATION EQUATIONS
Although we have been able to handle all the constraints of our perturbed model except the zero mode of the Hamiltonian constraint, this constraint is still so intrincate that, in the presented form, it does not seem possible to obtain its general solution analytically. In order to investigate the properties of the physical states, we will now introduce an ansatz that contemplates a situation of special interest. We will consider states in which the dependence on the FLRW geometry and on each of the gauge invariant fields can be separated. In this separation, all parts are allowed to depend on the inflaton. In more detail, from now on we analyze states of the form where we have adopted the abstract notation N s , N T , and N F to denote the dependence on the degrees of freedom of the corresponding Fock space, via a set of occupation numbers in the respective basis of n-particle states. In addition, ξ(v,φ) designates a state in the kinematic Hilbert space of homogeneous and isotropic LQC, such that it is normalized and evolves unitarily with respect toφ as whereÛ is an evolution operator with generatorˆ H 0 that is close to the unperturbed one, determined byĤ (2) 0 . This last condition can be understood as the requirement that the action ofĤ be at most of the order of the perturbative contributions when imposing the Hamiltonian constraint. Moreover, in the following, for simplicity, we will assume that this term is actually negligible in the action of the Hamiltonian constraint on the considered state, assumption that can always be checked for consistency once ξ(v,φ) is specified.
On this family of states, we still must imposeà la Dirac the Hamiltonian constraint operatorĤ S , that couples the FLRW background cosmology with the gauge invariant perturbations. To get solutions in physically relevant regimes, we can employ certain approximations that facilitate the resolution of the constraint. First, we consider regimes in which the transitions in the FLRW geometry mediated by the Hamiltonian constraint can be ignored as negligible on ξ(v,φ). In this situation, the relevant part of the Hamiltonian constraint is provided by its expectation value on ξ(v,φ) over the FLRW geometry (with the integration measure of the inner product of LQC). This expectation value provides a constraint equation on the gauge invariant perturbations of the form In what follows, we will neglect the perturbative operator Θ o ξ when compared to ˆ H 0 ξ on the left-hand side of this equation, according to our perturbative scheme.
Suppose for the moment that in Eq. (7.3) we can also neglect the first term, equal to the second derivative of the wave function of the perturbations with respect to the inflaton. This can be regarded as a kind of Born-Oppenheimer approximation, in the sense that one neglects the variation of certain degrees of freedom of the considered quantum state in comparison with the variation of others. Explicitly, we ignore the variation of the perturbations with respect to the inflaton in favor of the variation of the FLRW state, that is given in average by the expectation value ˆ H 0 ξ . Additionally, it is worth noticing that the last term in the constraint equation (7.3) affects only the scalar perturbations, becauseΘ o depends only on them. Let us assume that this term for the scalar perturbations is negligibly small. Taking into account that, in our representation,πφ acts as the derivative with respect to the explicit dependence on the inflatonφ (multiplied by −i), and thatˆ H 0 has been chosen to be close to the generator of the homogeneous and isotropic quantum dynamics with respect to the inflaton, the term under consideration can be understood as the total derivative of the operatorΘ o with respect to the inflaton, both in its explicit and in its implicit dependence. Thus, we expect that the analyzed contribution to the scalar perturbations can be ignored when the variation with respect to the inflaton is not significantly relevant. With these two approximations, the studied constraint amounts to the sum of a set of Schrödinger equations, one for each of the considered perturbations (scalar, tensor, and fermionic).
Specifically, we get the following equations for the gauge invariant perturbations: Note that the separation of variables can actually be made mode by mode in each of the gauge invariant perturbations, since these modes are not coupled by the Hamiltonian constraint.
Had we not neglected the contribution ofĤ , but considered instead that its action on the wave function of the FLRW geometry is of the same order as that of the perturbative contributions, we would have obtained an equation similar to the constraint (7.3) although with an additional term, given by the expectation value on ξ(v,φ) of the discussed difference of operators. Then, we should have added to the right-hand side of each Schrödinger equation a backreaction term, which could only depend on the inflaton. The balance between these backreaction terms C (ξ) (φ) would require that where the subscript on the backreaction tells us whether the term corresponds to the scalar (s), tensor (T ), or fermionic (F ) contribution. From this balance, we see that the sum of all the backreaction terms gives us information, in mean value and within our approximations, about how much the state ξ(v,φ) departs from an exact solution of the unperturbed homogeneous and isotropic cosmology in LQC. Moreover, let us return to Eq. (7.3) and let us assume now only that the gauge invariant perturbations admit a direct (effective) counterpart of the Heisenberg dynamics that results for their operator analogs from this Hamiltonian constraint equation, something that seems reasonable because the considered Hamiltonian is quadratic in the perturbative variables. Then, it is immediate to realize that we get a set of modified propagation equations for the MS modes, the tensor perturbations, and the fermionic perturbations. For instance, the modified MS equations are where the conformal time η ξ is defined by the equation Therefore, this time depends on the state ξ(v,φ) of the FLRW geometry. Similarly, for the modes of the tensor perturbations we obtain The conformal time is the same as for the scalar perturbations, thanks to the fact that the operatorθ e that multiplies the squared momenta in the Hamiltonian constraint (and that represents the squared scale factor, up to a constant) coincides both for tensor and scalar gauge invariants and, furthermore, for all the modes of these perturbations. In a similar way, an equation with quantum geometry corrections and in the same conformal time can be obtained as well for the femionic perturbations (see Ref. [81]). In the above propagation equations, the ratio of expectation values on the right-hand side gives the quantum corrected mass for the specific gauge invariant perturbation under consideration. We notice that this corrected mass is actually mode independent, because this is the case for the corresponding operators. Also, note that the equations contain no dissipative term. Much more important, the deduced effective equations are hyperbolic in the ultraviolet regime, regardless of the concrete behavior of the quantum state for the FLRW geometry, provided that our approximations are valid.
In order to extract predictions from the above equations about quantum geometry effects on the primordial perturbations, one needs to compute the expectation values that give the corrected masses for the MS and tensor perturbations. There are several possible strategies to reach this goal. Let us list three of these strategies, in decreasing order of accuracy but increasingly easier to implement. First, one could compute the quantum expectation values numerically. For this, one may try and ignore the backreaction (checking the validity of this approximation afterwards) and integrate numerically the quantum evolution of the FLRW state with respect to the inflaton. With the FLRW state obtained in this way, one can calculate with numerical methods the desired expectation values at each given value ofφ. The more difficult part of this program is the integration of the FLRW dynamics in the presence of non-trivial inflaton potentials. Second, taking into account the commented complication that the potential introduces, one can compute the evolution of the FLRW state not numerically, but in an interaction picture in which the potential (or part of it) is regarded as an interaction added to the homogeneous and isotropic Hamiltonian of LQC [35], and treated as a perturbation of that Hamiltonian using a Dyson series expansion [90]. And third, for suitable FLRW states, one can directly adhere to the effective dynamics description of LQC, integrating numerically only the trajectory of the peak of the state, rather that the quantum dynamics strictly speaking. Furthermore, this integration can be simplified by identifying regimes with universal behavior in the evolution from the bounce for the background solutions of interest in LQC [18,[91][92][93]. For instance, the most interesting situations to get quantum geometry corrections on the primordial spectra that can be observed nowadays are found for background solutions that are kinetically dominated around the bounce, so that the potential there has little influence. This allows us to introduce further simplifications in the integration of the FLRW trajectories that, at the end of the day, facilitate the calculation of the quantites that determine the studied masses of the perturbations.
Most of the work in the literature has indeed been done assuming an effective dynamics for the description of the FLRW cosmology in LQC. Even if, with this approximation, the problem of computing the evolution of the primordial perturbations is handleable, the results (and hence the predictions obtained from them) depend critically on the initial conditions that one chooses for the FLRW background in this effective dynamics, as well as on the initial conditions that determine the state of the perturbations subject to the propagation equations that we have derived. We will discuss these issues in the next section.
Let us point out that, adopting this effective dynamics for the description of the FLRW states, it has been proven [94] that the corrected mass that appears in the modified propagation equations for the scalar and tensor perturbations is positive around the Big Bounce, at least for the most interesting ranges of energy density contribution of the inflaton potential. Since the Big Bounce is precisely the region where the quantum effects on the geometry are more significant, one would expect that the largest departures from the classical situation described by GR cosmology happen there. This positivity of the quantum corrected mass is important to be able to define adiabatic vacua as initial states around the bounce for all the perturbative modes [91,94,95]. A negative mass involves a breakdown of the adiabatic approximation around the bounce at least for values of ω k that are not sufficiently large, invalidating the construction of adiabatic states as natural candidates for a vacuum at frequencies that can be of physical interest, for instance because they cover part of the observed spectrum in the CMB. Moreover, the positivity of the mass at the bounce is not shared by other proposals for the quantization of cosmological perturbations within the framework of LQC, like the so-called dressed metric approach that has been put forward by Agullo, Ashtekar, and Nelson [96][97][98].
Finally, in our discussion above, and owing to the compactness of the spatial sections, the modes that we have considered possessed a discrete spectrum of Laplace eigenvalues, ω k , that play the role of frequencies in the modified propagation equations for the perturbations even after the introduction of quantum geometry corrections. Nonetheless, it is possible to reach the continuum limit for this set of frequencies in the following form. One first extracts a length scale of reference from the scale factor. All observable quantities are defined with respect to this reference scale, that becomes physically irrelevant. One may choose as such scale the value of the scale factor today, or at the moment of the bounce, for instance. Then, the desired continuum limit is reached as the limit in which we make the reference scale tend to infinity. We refer the reader to Ref. [99] for more details.

VIII. INITIAL CONDITIONS
As we have commented, even if we have succeeded in deriving propagation equations for the primordial perturbations that contain modifications caused by quantum geometry effects and even if we assume FLRW states that can be described within the effective dynamics approach to LQC, in order to extract predictions about the primordial cosmological perturbations we need to specify the particular FLRW effective solution that plays the role of a background and, in addition, the vacuum state that determines the conditions on the perturbations. Both pieces of information can be supplied by giving convenient initial data on a certain spatial section. An appealing possibility is to choose this section precisely at the Big Bounce. We will concentrate our discussion on this case. Other possibilities are equally valid, for instance a section in the asymptotic past, if the effective dynamical evolution previous to the bounce connects with a manageable asymptotic region [100].
Let us consider first the initial conditions for the FLRW background, solution to the effective dynamics of LQC. The FLRW cosmology is described by two pairs of canonical zero modes, namely four variables. But we have chosen to impose initial conditions at the bounce, where the time derivative of the scale factor vanishes, reducing the liberty in one degree of freedom 4 . In addition, the Hamiltonian constraint associated with the effective dynamics reduces the degrees of freedom in one more variable. Moreover, we have also commented that we can employ the value of the scale factor at the bounce as a reference scale, depriving it of physical relevance. In practice, this allows us to set that value equal to the unit, for instance. In total, we see that only one variable must be fixed at the bounce by the initial conditions there. We choose the value of the inflaton as this piece of initial data. On the other hand, we can consider also as free data the parameters that determine the inflaton potential. Focusing our attention on the most studied case of a quadratic potential, we find only one parameter, given by the inflaton mass. From this perspective, the FLRW effective background turns out to be completely fixed if we provide the value of the inflaton at the bounce and the value of the inflaton mass.
Actually, we are only interested in effective solutions that lead to power spectra for the perturbations that are compatible with the observations, but that still retain some quantum geometry corrections. One expects that, if these corrections have survived, they should be present in the region of large angular scales or its nearby region, because it is only in this region that the agreement between GR and observations may not be completely solid [101,102]. This requirement determines a relatively narrow interval of values for the initial condition on the inflaton φ B and the inflaton mass m, around φ B = 0.97 and m = 1.2 × 10 −6 (in Planck units). For this latter choice of specific values, we show in Fig. 1 the evolution of the Hubble parameter H multiplied by the scale factor. This rescaled Hubble parameter aH vanishes at the bounce and then increases in a very short superinflationary epoch in which H grows to a value of the Planck order. This happens so fast that the scale factor remains almost constant in the process. Since we have taken the scale factor at the bounce equal to one, then aH at its maximum should be of the order of one as well in Planck units (like H). This maximum sets a scale, that we denote K LQC in terms of wavenumbers, and that should be of Planck order according to our previous arguments. From that moment on, the rescaled Hubble parameter starts to decrease until it reaches a minimum. Besides, the quantum corrections in the effective dynamical equations of the background become negligibly small, and the effective trajectory gets totally adapted to a GR cosmological solution.
For solutions that allow for quantum geometry corrections in the spectra of the perturbations at large scales, the inflaton dynamics around the bounce is dominated by its kinetic energy density, which is of Planck order, with an ignorable contribution of the potential. Since the potential is almost negligible, the effective solution behaves as if the scalar field were massless, situation in which the inflaton momentum is a constant of motion and the kinetic energy density decreases rapidly, as a −6 . The kinetic energy density continues to diminish until it becomes of the order of the potential. Given thatȧ increases when the potential drives the evolution of the scale factor (both in GR and in the effective dynamics of LQC), the coincidence between the kinetic energy density and the potential of the inflaton occurs approximately when aH reaches its minimum. On the other hand, during the evolution from the bounce to this minimum of aH, the inflaton typically increases only by a few orders of magnitude. As a result, the potential, quadratic in the inflaton, varies as well only in a few orders. Taking this into account, and since the potential at the bounce is m 2 φ 2 B /2, with values around 10 −12 in Planck units, when the kinetic and potential energies coincide we expect a density in the range [10 −12 , 10 −9 ]. This energy density determines yet another scale in the system, that we call K K−P expressed as a wavenumber.
In total, the influence of the effective background solution on the perturbations is characterized by two (wavenumber) scales, that we have already mentioned, K LQC and K K−P . In a first approximation to the problem, these scales determine the regions where the quantum geometry effects may cause departures from the standard model of inflation in GR. As we have argued, K LQC is of the order of the unit, because it is related to quantum gravity phenomena. To estimate K K−P in our solutions, note that during the epoch of kinetic dominance, the energy density decreases as a −6 from a value of the Planck order to values in the interval [10 −12 , 10 −9 ] as we have pointed out. Recalling the Hamiltonian constraint of effective LQC (or of FLRW cosmology in GR, once one is away from the immediate vecinity of the bounce), one concludes that H 2 must be proportional to the discussed energy density in the considered region, and thus decrease during kinetic dominance also as a −6 . Consequently, aH must evolve as a −2 , decreasing from the Planck order as the cubic root of the energy density, and hence reaching values in the range [10 −4 , 10 −3 ]. In Fig. 1 we see that wavenumbers larger than K LQC only intersect aH once in the evolution. Note that this intersection is the moment when the associated physical length a/k coincides with the Hubble scale 1/H, and therefore can be taken as the moment of horizon crossing. For modes between K LQC and K K−P , there exist three intersections. Essentially, the modes exit the horizon immediately after the bounce, reenter in the phase of kinetic dominance, and exit again during the inflationary expansion. We expect these modes to be severely affected by the quantum geometry effects around the bounce, and that their evolution differs considerably from that experienced in GR for a background solution with the same behavior in the classical region. Finally, modes with wavenumbers below K K−P do not exit the horizon during inflation but much before, a fact from which one may expect important departures from the predictions of standard inflation.
After the potential equals the kinetic energy density, the latter rapidly decreases while the potential becomes essentially constant and generates inflation. We are interested in effective background solutions such that the modes that experience quantum geometry effects (roughly, those with wavenumbers between K LQC and K K−P ) are entering the horizon today. If they had entered much long ago, they would correspond to scales of the power spectra where there is no discrepancy with GR, raising the problem of explaining the absence of departures from the Einsteinian predictions or implying that quantum gravity effects are too tiny to be observable in those circumstances. On the other hand, if they had not entered the horizon again (nor were about to do it), they could not be observed in the power spectra. The interesting situation is when the modes are entering the horizon at present, as we have said. But the effective FLRW backgrounds for which this happens turn out to experience a short-lived inflation. As a consequence, during the first moments of the inflationary expansion, there is still some influence of the kinetic energy density, producing departures from a genuine slow-roll behavior. This will affect the power spectra if the modes that were exiting the horizon at those moments are observed today (see e.g. the discussion of Ref. [103] in the framework of GR). Therefore, the slow-roll approximation will not be good, at least, for modes that exited the horizon during those first stages of inflation [103]. Such modes are precisely those close to the scale K K−P . Hence, those modes will experience two types of corrections from a standard inflationary scenario with slow-roll in GR: quantum geometry effects, accumulated around the bounce, and short-lived inflation effects. One of the most important challenges for LQC nowadays is to be able to separate these two kinds of effects and prove that it is possible to identify and falsify the quantum modifications in cosmological observations.
In summary, the LQC modifications in the FLRW background with respect to the standard inflationary solutions of GR may have a relevant influence on modes between the typical scale of the quantum geometry effects and the scale K K−P , close to the onset of inflation. If these include the modes that are now re-entering the horizon, so that the scale of the Universe that we observe today was at the very early stages in the range affected by the quantum effects, some traces of those quantum modifications may have survived in the CMB in spite of the later inflationary expansion, and they might be observable. The fact that the background which those modes feel effectively differs substantially from a de Sitter expanding phase should imply that the natural vacuum for them ought to differ from the standard Bunch-Davies vacuum [104]. As a result, the power spectra of the perturbations at those scales changes from the conventional predictions based on the choice of a Bunch-Davies state. Suppose that the new vacuum state is related to the standard one by a Bogoliubov transformation that does not mix modes with different values of ω k , something that is ensured if the invariance under the symmetries of the spatial sections is respected. Let us call α k and β k the coefficients of this Bogoliubov transformation, with |α k | 2 − |β k | 2 = 1. Recall that the beta-coefficient determines the antilinear part of the Bogoliubov transformation, i.e. the part that mixes creation and annihilation operators. These coefficients can be determined, e.g., at the initial time chosen in our analysis, if we know there the initial data that specify the two bases of solutions of the gauge invariant field equations, {μ k } and {µ k }, that characterize respectively the new and the old vacua. Then, if the primordial power spectrum of the standard vacuum is P R (k), the power spectrum of the new vacuum state does not need to be calculated from scratch: it is given by the formulaP Here, φ α k and φ β k are the phases of the respective Bogoliubov coefficients, treated as complex numbers, and φ µ k is the phase of the solution µ k evaluated at the time of computation of the power spectrum (typically by the end of inflation).
The second problem related with the choice of initial data is, therefore, the selection of conditions that determine the vacuum state of the perturbations. Clearly, from the above formula, a change of vacuum state may result in a radical variation of the power spectrum. The predictive power of the formalism is lost unless we have a way to select a vacuum as the preferred state for the gauge invariant perturbations. While, in situations like de Sitter inflation, the high degree of symmetry of the background can help us in picking out a unique state, invariant under the symmetries and with a local Minkowskian behavior, this does not seem possible in more general situations, like those experienced by the modes affected by quantum geometry effects in the kind of effective backgrounds that appear in LQC. In these circumstances, several proposals have been suggested in order to single out a unique Fock state that could then be viewed as privileged in the system. Among these proposals, the attempt to use adiabatic states [105,106] has received a considerable attention [18,95,97,98]. Nonetheless, their construction may find some obstructions, especially if the effective mass in the propagation equations of the perturbations becomes negative. We have seen that this does not occur in the hybrid approach in the region of important quantum geometry effects, at least for (effective) solutions with kinetic dominance in the energy balance of the inflaton [94]. Nonetheless, this is not the case for other approaches like the dressed metric quantization if one considers scales that are not sufficiently small [94]. Besides, the power spectra of adiabatic states, computed numerically, often present large oscillations, and even if these oscillations are averaged, they typically result in an increase of power that does not seem to fit properly with observations if the scales affected by the quantum geometry effects are inside the Hubble horizon today [18].
Ashtekar and Gupt have put forward a different proposal for the vacuum state [107,108]. In the region with relevant LQC effects, they have required that the quantum Weyl curvature satisfy a bound which is the lowest value compatible with the uncertainty principle and stable under evolution. This condition selects a ball of states. Among them, the vacuum of the perturbations is chosen by imposing another condition at the end of inflation, ensuring that the dispersion in the field operators be minimized [108]. In the dressed metric approach, this proposal has been shown to lead to primordial power spectra that, though still highly oscillatory, seem in very good agreement with observations after being averaged [19,108]. Nonetheless, the direct relation of this vacuum with adiabatic states is not known.
Another interesting proposal for a vacuum state is the so-called non-oscillating vacuum, suggested by Martín-de Blas and Olmedo [95]. The proposal is to select the state that minimizes the integral in a certain interval of conformal time, usually the interval from the time of the initial spatial section to a time well inside the inflationary regime. For instance, in our typical class of effective backgrounds, this can be a time when the kinetic energy density of the inflaton becomes so negligible that the inflationary expansion is completely driven by the potential. Since the primordial power spectrum for each mode is proportional to the square norm of the associated mode solution µ k , the proposal picks out a state that minimizes the power oscillations in a definite sense. In general, the determination of this vacuum state needs numerical methods, since the criterion of choice is posed as a variational problem that involves the calculation of an integral. For simple cases, the proposal can be handled analitycally and has been proven to provide a conventional choice of vacuum state. Thus, it selects the Poincaré vacuum for flat spacetime in the presence of a scalar field, either massless or with a quadratic potential. In addition, for de Sitter spacetime, the proposal selects the Bunch-Davies vacuum [95]. The primordial and angular power spectra for this vacuum state have been calculated numerically, both for scalar and tensor perturbations [36,37]. The results are compatible with observations, and they even open the possibility of explaining some of the features of the spectra that perhaps may be in tension with GR [36,91], at large angular scales or for multipoles around l = 20 [101,102].

IX. CHOICE OF VACUUM STATE FOR THE PERTURBATIONS: SPLITTING OF PHASE SPACE VARIABLES
The problem of selecting a vacuum for the perturbations appears in our formalism because the requirements that we have imposed to determine the Fock representation of the gauge invariant perturbations in the hybrid approach at most select a family of representations that are unitarily equivalent, but not a privileged state. Any of those representations, or equivalently any Fock state in the considered family, can be chosen as the vacuum. This leaves a large freedom in the selection of a vacuum for the perturbations, and so in the initial conditions that characterize it. The proposals that we have commented at the end of the previous section are some of the attempts to fix this freedom, but there is yet no general consensus about how to settle the question. Besides, although some of those proposals lead to power spectra that are compatible with observations, they happen to rely on numerical and/or minimization techniques.
In fact, we can consider families of representations related among them by unitary transformations which depend on the background. The Heisenberg dynamics of the creation and annihilation operators associated with each of these representations, even if unitarily implementable as provided by a composition of unitary transformations, would differ between them, given that part of the evolution is removed by assigning it to the background sector of phase space. What is more, by means of this type of unitary transformations with dependence on the background, we can change the splitting of the phase space degrees of freedom between the zero modes that describe the background and the modes of the gauge invariant perturbations. Actually, there exist many ways of separating the phase space into a homogeneous sector and an inhomogeneous one using canonical transformations that mix them. The specific splitting that one adopts determines the properties of the resulting quantization. In particular, the representation of the Hamiltonian constraint and its ultraviolet features strongly depend on this choice.
This fact can be employed to improve the behavior of the field operators that represent the perturbative terms, ameliorating the need for the introduction of regularization procedures typical of QFT. Indeed, as they stand, the actions of the MS, tensor, and fermionic Hamiltonians that appear in the constraint (6.7) are ill defined with a standard choice of Fock representation for the corresponding perturbations. Moreover, the backreaction is in general divergent. For instance, by constructing the unitary operator that implements the Heisenberg dynamics of the fermionic variables (in the context of QFT in a quantum mechanically corrected background), one can compute the backreaction of the fermions, C (ξ) F (φ), and show that it is not absolutely convergent [81]. These problems can be solved, or at least alleviated, by introducing new gauge invariants prior to quantization, defined by using canonical transformations that depend on the zero modes. Considering the system as a whole, we have freedom in • Changing the dynamical separation between the FLRW geometry and the gauge invariant perturbations via canonical transformations.
• Chosing the Fock vacuum for the perturbations, within the hybrid scheme, regarding this vacuum as the state from which one defines the Fock representation as a cyclic one.
All this ambiguity can be encoded in choices of the form a k,± = f k,± (ã, πã,φ, πφ)ν k,± + g k,± (ã, πã,φ, πφ)π ν k,± (9.1) for the MS annihilationlike variables. Here, so that the introduced variables satisfy, when one freezes the background, canonical commutation relations with the corresponding MS creationlike variables, defined by the complex conjugate of relation (9.1). For the tensor variables, on the other hand, one is led to consider analogous families of creation and annihilationlike variables, characterized by two functions f k,ε,± and g k,ε,± that satisfy a condition similar to Eq. (9.2). Generally, one is interested exclusively in canonical transformations of the gauge invariant variables that depend (apart that on the cosmological zero modes) only on the frequency ω k of the mode, but not on other details about the wavevector k, nor on the sine or cosine character of the Fourier mode or the polarization of the tensor mode. For those cases, we would adopt the simpler notation f k and g k for the functions that define the creation and annihilationlike variables.
In the case of fermions, the ambiguity is captured in the freedom to define annihilationlike variables for particles and creationlike variables for antiparticles as follows: b (x,y) k * = g k 1 (ã, πã,φ, πφ)x k + g k 2 (ã, πã,φ, πφ)y * − k , In the same spirit that we have commented above, one is usually interested only in cases in which the functions f k 1 , g k 1 , g k 1 , and g k 2 depend on k only via ω k . Notice that the creation and annihilationlike variables (5.9) that we used for the fermions in Sec. 5 were of this kind. We will restrict to this type of cases in the following.
As we already know, a change from the gauge invariant variables that we have adopted for our system to any of the above sets of creation and annihilationlike variables for the perturbations can be completed into a canonical set for the full cosmological model. It suffices to correct again the zero modes with contributions that are quadratic in perturbations in the way that we discussed in Sec. 4. In addition, in terms of the new canonical set, the resulting MS, tensor, and fermionic Hamiltonians are the old ones plus some known corrections. These new contributions contain, in general, both diagonal products of annihilation and creationlike variables, and terms that are responsible for the creation and destruction of pairs. The asymptotic behavior of these latter interaction terms when ω k → ∞ is what tells us if the quantization of the Hamiltonians is well defined on the vacuum, assuming normal ordering. In all cases f k , g k , f k 1 , f k 2 , g k 1 , and g k 2 can be chosen so that the dominant powers of ω k in the interaction terms that prevent a nice behavior of the Hamiltonian operators on Fock space are eliminated.
Moreover, it is possible to remove, order by order in inverse powers of ω k , all the asymptotic contribution to the interaction terms in the Hamiltonians. For example, let us consider the scalar perturbations. The MS Hamiltonian gets asymptotically diagonalized with [109] The functions γ n are determined by the recursion relation where γ 0 = s (s) + r (s) πφ is just the background dependent mass for the MS field. Creation and annihilationlike variables are then asymptotically fixed, up to a phase, since from the canonical commutation relations it generally follows that Similar asymptotic characterizations to diagonalize the field Hamiltonians can be obtained for the tensor perturbations and for the Dirac field [110]. Actually, in all of these cases, the first few terms in the asymptotic expansion are enough to construct variables with well-defined Hamiltonians (and finite backreaction contributions to the quantum constraint).
On the other hand, the phases that still remain free in the creation and annihilationlike variables can be determined univocally by means of further physical considerations. Specifically, it seems natural to demand that the background dependence extracted from the dynamics of the original perturbations by our choice of those phases is the minimum allowed, and that the resulting asymptotically diagonal Hamiltonians are positive, as functions of the background.
Given that our analysis has been carried out asymptotically for modes with large wavenumbers, the question arises of what happens for other kinds of modes and, in particular, if the asymptotic expansions provided by the Hamiltonian diagonalization for large ω k can uniquely specify a choice of creation and annihilationlike variables for all the modes. Let us consider, e.g., the scalar perturbations. In fact, the interaction terms in the Hamiltonian for each possible value of ω k are completely eliminated if and only if This is a semilinear partial differential equation for which the complex solutions satisfy where the prime stands for the operation of taking the Poisson bracketã{., H 0 }. It is worth commenting that, in the linearized context of QFT in curved spacetimes, our asymptotic characterization above can be shown to lead in a unique way to the Minkowski vacuum in the case of constant mass, and to the Bunch-Davies vacuum when the homogeneous background is taken as the de Sitter solution [109]. Thus, in these linearized contexts, the procedure of asymptotic diagonalization is able to uniquely fix a solution to Eq. (9.9) for all wavenumbers, and this solution reproduces the natural choice of vacuum state in the considered scenarios. Furthermore, the corresponding asymptotic expansions for the fermionic creation and annihilationlike variables that diagonalize the Hamiltonian have been proven to determine as well a unique choice for all scales in the linearized de Sitter context, even if it is known that those expansions have zero radius of convergence in this case [111]. In summary, the asymptotic diagonalization of the Hamiltonian of the perturbations may provide a procedure to determine a vacuum state, and therefore to fix initial conditions for the primordial perturbations in such a way that they are optimally adapted to the dynamics dictated by the Hamiltonian constraint of the total system. Moreover, recent investigations [112] support a close analytical relation between the vacuum state that would be selected in this manner and the NO vauum proposed by Martín-de Blas and Olmedo [95], at least in the context of hybrid LQC.

X. CONCLUSIONS
In this work, we have reviewed the hybrid approach to LQC. This approach to the quantum description of gravitational systems with local degrees of freedom within the framework of the loop quantization program tries to provide, in a controlled way, a formalism for the study of inhomogeneous cosmological scenarios that, yet, display some symmetries that simplify the physics, or in which the inhomogeneities can be described in a perturbative way over a highly symmetric background. In this way we have been able to analyze linearly polarized gravitational waves in Gowdy cosmologies with toroidal compact sections, and scalar, tensor, and fermionic perturbations at quadratic order in the action around an FLRW spacetime in the presence of an inflaton. In particular, for these cosmological perturbations and at the considered truncation order, we have found a canonical set for the full system composed of gauge invariant perturbations (including the MS field), linear perturbative constraints and gauge variables conjugate to them, and zero modes that contain the relevant information about the background FLRW cosmology. In a hybrid quantization of this canonical system, physical states depend only on the quantum FLRW background and on gauge invariant perturbations. Starting from the zero mode of the Hamiltonian constraint, that couples these perturbations with the FLRW background, and adopting a suitable ansatz for the quantum states of interest, we have been able to derive propagation equations for the perturbations in the primeval stages of the Universe. These equations differ slightly from those of GR by the inclusion of quantum corrections, corrections that we have succeeded to explicitly derive with our hybrid strategy taking fully into account the quantum behavior of the FLRW substrate, and therefore beyond the level of an effective description of this background within homogeneous and isotropic LQC.
In order to quantize differently the system within the framework of LQC, but still adhering to the idea of developing a QFT for the perturbations on a quantum spacetime, one can follow the so-called dressed metric approach, put forward in Refs. [97,98,113], instead of the hybrid approach. Indeed, as in the hybrid proposal, the dressed metric approach adopts also the philosophy of combining a loop representation for the homogeneous sector of the (truncated) phase space and a Fock representation for the tensor and MS perturbations (and possible fermionic perturbations, if they are present). Again, in the dressed metric approach one also introduces an ansatz for the quantum states of cosmological interest in which the dependence on the homogeneous geometry and on the perturbations factorizes. In this ansatz, all partial wavefunctions are allowed to depend on the inflaton field φ. However, in the dressed metric case there is no Hamiltonian constraint that affects the perturbations, since the whole of the truncated cosmology is not treated as a constrained symplectic system. Instead, one has the Hamiltonian constraint of the homogeneous FLRW model, and the Hamiltonian functions (4.23) and (4.24) that, classically, generate the dynamics of the perturbations. Consequently, the approach requires that the homogeneous part of the states be an exact solution of the FLRW model in LQC, and then uses this solution to define the quantum dynamics on the phase space of the gauge invariant perturbations [97,98]. In this way, the perturbations behave as test fields that see a dressed metric determined by certain expectation values of operators of the homogeneous geometry, which incorporate the most relevant quantum effects.
In spite of the similarities between the hybrid and the dressed metric approaches, the effective equations that they provide for the propagation of the gauge invariant perturbations are somewhat different even if backreaction is neglected. The discrepancy appears only in the term of the time dependent mass in the propagation equations [94]. At the end of the day, this can be traced back to the differences in the treatment of the phase space of the perturbed FLRW cosmologies in the hybrid and the dressed metric proposals. As we have emphasized, in the hybrid case the whole phase space is treated as a symplectic manifold, and accordingly it is described in terms of canonical variables. This applies, in particular, to the expression deduced for the time dependent mass. On the contrary, in the dressed metric formalism, one evaluates the time dependent mass directly on the FLRW metric dressed with quantum corrections. For states such that this metric satisfies the effective dynamics of LQC, the time derivatives involved in the corresponding expression of the time dependent mass are then computed along an effective trajectory of homogeneous and isotropic LQC. The difference then arises because of the departure of the standard classical relation between the time derivatives of the scale factor and its canonical momentum (inherent to the hybrid approach) with respect to the alternative effective relation in LQC (employed in the dressed metric case) [94]. Remarkably, this difference is specially important around the bounce, precisely the region where the quantum corrections on the propagation of the perturbations are expected to be relevant.
Several other approaches have also been suggested for the investigation of cosmological perturbations within LQC. For a comprehensive summary of such approaches, we refer the reader to the reviews listed in Refs. [48][49][50][51][52][53][54][55]. They include the deformed constraint algebra approach [114][115][116][117], the group field theory models [118][119][120], and the quantum reduced loop gravity scheme [121,122]. Additionally, different ways of addressing backreaction effects of the perturbations on the background within canonical quantum cosmology have been recently explored using techniques from space adiabatic perturbation theory [123]. Our attention here has been exclusively put on the hybrid approach in order to fill a gap in the literature, as this is the first extensive review of this proposal that includes a detailed description of the application to primordial perturbations.
A remarkable fact of the hybrid quantization is that, while inhomogeneities and background degrees of freedom are treated as parts of a single constrained system, the imposition of the quantum constraints is consistent and does not give rise to anomalies. This statement holds both in the Gowdy model and for cosmological perturbations. The precise relation of these constraints with the full set of four-dimensional spacetime diffeomorphisms is a different issue that calls for more detailed investigations. As presented in Sec. III, the Gowdy model is not only a symmetry reduction of Einstein gravity, but it is also a partially gauge-fixed system in which only two global constraints remain, namely the zero mode of the Hamiltonian constraint and the zero mode of the momentum constraint in the angular direction on which the metric fields depend 5 . It is worth emphasizing that these are only two constraints, and not two constraints per point (neither of the spatial section nor in the considered angular direction). The aforementioned momentum constraint generates rigid translations in the corresponding angle, while the Hamiltonian one generates global time reparameterizations. These two constraints of the model actually display vanishing Poisson brackets between them and, with the adopted quantization, their corresponding operators commute. For cosmological perturbations, the constraint algebra has to be consistent just up to the order of the perturbative truncation used in our treatment. We have shown that the linear perturbative diffeomorphisms and Hamiltonian constraints admit an Abelianization at this truncation order, and we have represented them directly as part as our canonical elementary variables. The only remaining constraint in the system is a global one, given by the zero mode of the Hamiltonian constraint, that includes contributions from the background and from the perturbations. Notably, its only dependence on the perturbations is via gauge invariants, and therefore commutes with the linear perturbative constraints both classically and in the quantum theory. Indeed, we recall that in the hybrid quantization the Mukhanov-Sasaki and tensor perturbations are represented as operators that commute with the linear perturbative constraints. In other words, at the level of our perturbative truncation and with our hybrid strategy, the algebra of the quantum constraints of our perturbed system does not present anomalies.
The physical relation of these constraints with the four-dimensional diffeomorphisms algebra and the extent to which recent claims about problems with general covariance in LQC [56,126,127] affect the system at the considered perturbative order deserve further study. These claims have been inspired in part by the deformed constraint algebra approach, which in particular predicts processes of effective signature change in high curvature regimes [128][129][130]. In this respect, let us point out, for instance, that some of the perturbative canonical variables used in the hybrid approach are defined with fields that are non-local functions of the spatial metric, inasmuch as they can only be obtained by taking inverse derivatives. This is a common situation even in standard cosmological perturbation theory [38][39][40][41][42][43]. Moreover, in terms of the background variables employed in our formulation, the metric functions include corrections that are quadratic in the perturbations already at the studied truncation order. A representation of these metric components as quantum operators has yet to be constructed, but it is clear that now issues such as the non-degenerate Lorentzian character of the metric become intrincate questions from a quantum perspective. Even the square scale factor, that in absence of pertubations is strictly positive in each superselection sector of homogeneous LQC with the quantization prescription adopted here 6 , might in principle turn negative by the effect of the perturbations. Nonetheless, none of these unexplored questions on the quantum geometric structure changes the hyperbolic ultraviolet behavior that we have found for the propagation equations of the perturbative modes.
Even if we have succeeded in deriving such mode equations, that rule the evolution of the primordial perturbations in the hybrid approach, we have seen that this is not yet enough to extract predictions that can be confronted with observations. For this purpose, we also need two types of initial data, namely initial values to fix the FLRW background and conditions to choose a unique vacuum state for the perturbations. With respect to the FLRW cosmology, we have seen that it suffices to provide, e.g., the value of the inflaton at the bounce, apart from the parameters that determine the inflaton potential. In the case of a quadratic potential, we have found values for the inflaton on the bounce section and for the inflaton mass such that the modes affected by quantum geometry effects are those that are re-entering the Hubble horizon nowadays, situation that is the most interesting possibility in terms of observational plausibility in the CMB. Concerning the vacuum state of the perturbations, we have commented on various proposals to select it that lead to power spectra that seem compatible with the observational data.
To go beyond those proposals and find a criterion to select the vacuum that is rooted on the hybrid strategy, that combines loop and Fock representations, we have put an additional emphasis on the choice of splitting between the homogeneous and isotropic sector of phase space and the gauge invariant perturbations. This freedom can be employed to reach a Hamiltonian constraint with nice properties, at least as far as its action on the perturbations is concerned. Requiring such good physical and mathematical properties turns out to restrict the possible quantum dynamics of the perturbative gauge invariants, as well as the Fock representation chosen for them. In turn, this can be regarded as a limitation in the admissible choices of vacuum state. In particular, we have shown that a criterion such as the diagonalization of the Hamiltonian of the gauge invariant perturbations, based on its asymptotic structure, might be able to provide a unique vacuum state, to which one may partcularize in the future the discussion of the effects of quantum geometry in cosmology to extract concrete and distinctive predictions.