Classical and quantum cosmology of minimal massive bigravity

In a Friedmann-Robertson-Walker (FRW) space-time background we study the classical cosmological models in the context of recently proposed theory of nonlinear minimal massive bigravity. We show that in the presence of perfect fluid the classical field equations acquire contribution from the massive graviton as a cosmological term which is positive or negative depending on the dynamical competition between two scale factors of bigravity metrics. We obtain the classical field equations for flat and open universes in the ordinary and Schutz representation of perfect fluid. Focusing on the Schutz representation for flat universe, we find classical solutions exhibiting singularities at early universe with vacuum equation of state. Then, in the Schutz representation, we study the quantum cosmology for flat universe and derive the Schrodinger-Wheeler-DeWitt equation. We find its exact and wave packet solutions and discuss on their properties to show that the initial singularity in the classical solutions can be avoided by quantum cosmology. Similar to the study of Hartle-Hawking no-boundary proposal in the quantum cosmology of de Rham, Gabadadze and Tolley (dRGT) massive gravity, it turns out that the mass of graviton predicted by quantum cosmology of the minimal massive bigravity is large at early universe. This is in agreement with the fact that at early universe the cosmological constant should be large.


I. INTRODUCTION
Since 1916, Einstein's general relativity theory (GR) [1] has explained the majority of the phenomena related to gravity. Newtonian dynamics is reproduced from GR in the weak field limit that makes observable predictions such as the solar system tests beside the bending of light by massive objects up to a very high precision. The scalar curvature of a metric tensor beside the Einstein-Hilbert action in GR theory, proposes a geometrical interpretation of the gravitation. In other words, Einstein's field equations represents the interplay between the space-time geometry and the matter. However, beyond the solar system scales, different astrophysical observations have raised questions that remained unanswered in the framework of GR. These open questions include the cosmological constant problem [2], the advent of an invisible dark matter component in the universe [3], and the recently observed acceleration of the universe [4]. Finding explanations for these unanswered questions convinced people to study the possibility of modifying GR theory. From a more theoretical point of view, searching for the alternatives of GR is motivated by string theory, the well-known candidate for a quantum theory of gravity. In the context of field theory, GR is presenting non-linear self interactions of a massless spin-2 particle. From this point of view, modifying GR can make the spin-2 particle massive. Constructing a consistent theory that describes a massive spin-2 particle is an old challenge. Fierz and Pauli in 1939 presented linearized massive spin-2 field fluctuation [5]. After 30 years, van Dam, Veltman and Zakharov discovered that the helicity-zero mode of the massive spin-2 field does not decouple in the zero-mass limit [6,7]. This is called vDVZ discontinuity and implies that even a very small graviton mass has serious influence on the gravitational interactions between sources. After a while, Vainshtein presented the Vainshtein mechanism showing that around massive sources the linear approximation (Fierz-Pauli mass term) loses its validity below the Vainshtein radius which leads to the requirement for a nonlinear extension of the Fierz-Pauli mass term [8]. In 1979, Boulware and Deser showed that any nonlinear extension of Fierz-Pauli theory would exhibit a ghost instability [9,10], and this finding left the massive gravity theory without any significant progress for about 40 years. Eventually, in 2010, de Rham, Gabadadze and Tolley (dRGT) found that it is possible to construct a theory of ghost-free nonlinear massive gravity [11,12]. They showed that this theory becomes free of ghost just in a certain decoupling limit which was not enough to guarantee the consistency of the theory. Later, Hassan and Rosen proved the absence of the Boulware-Deser ghost in a Hamiltonian constraint analysis [13,14]. dRGT model is describing a nonlinearly interacting massive spin-2 field in flat space in which two metric tensor components are playing roles, one as dynamical g µν and one as non-dynamical f µν (reference metric). As a result, in the cosmological context of massive gravity theory, it has been shown that the flat Friedmann-Robertson-Walker (FRW) universe dose not exist [15], but this is not supported by the recent observational results. However, open FRW solutions are allowed [16] but they are involved with the problems of strong coupling [17] and ghostlike instabilities [18]. Attempting to obviate the problems with flat FRW solutions, guided Hassan and Rosen to extend the massive gravity theory, beyond the dRGT setup, to a theory with two dynamical symmetric tensors g µν and f µν as foreground and background metrics, respectively, having a completely symmetric role. They called this theory as massive bigravity theory and showed that the corresponding Hamiltonian description is a ghost-free bimetric theory containing nonlinear interactions of a massless and massive spin-2 field in a dynamical background [19]. It is worth to mention that bimetric theory was first introduced by Isham, Salam and Strathdee [20] to describe some features of strong interactions and it was later followed and renovated by Damour and Kogan in order to address new physics scenarios [21]. Apparently, this modified model covers the massive gravity and consider two metrics thoroughly in a symmetric way that annihilate the aether-like concept of reference metric in massive gravity [22]. Massive gravity cosmology have been studied in Refs [23][24][25][26][27][28][29][30][31], meanwhile in bigravity model some regular cosmological solutions have been derived [32][33][34][35][36][37][38].
At present there is an unsolved question that whether our universe contains closed, flat or open spatial threegeometries. A flat universe looks like the ordinary three-dimensional space we experience around us. In contrast, the spatial sections of a closed universe looks like three-dimensional spheres with a very large but finite radius. An open universe looks like an infinite hyperboloid. Curvature of space is important if the universe is assumed to be created in an inflationary state. Because the closed models have finite size they have generally been thought to be most relevant for quantum cosmology as they present finite action leading to a non-vanishing nucleation probability. In fact, Atkatz and Pagels have already shown that only a closed universe can arise via quantum tunneling [39]. However, strong evidences from present cosmological observations favour a flat or open universe. Therefore, if we believe that quantum mechanics is the fundamental theory of whole nature, the study of quantum cosmology for flat or open universe is inevitable and deserves more investigation. In 1983, Hawking and Hartle developed a theory of quantum cosmology known as the "No Boundary Proposal" [40]. We know that the application of path integral to cosmology involves a sum over four dimensional geometries that have boundaries matching onto the initial and final three geometries. The Hartle-Hawking proposal simply avoids the initial three geometry and only includes four dimensional geometries that match onto the final three geometry. Therefore, path integral is interpreted as giving the probability of a universe with certain properties being created from nothing. In practice, the calculation of probabilities in quantum cosmology using the path integral is rather difficult and semiclassical approximation has to be used, where one argues that most of the four dimensional geometries considering in the path integral give rise to very small contributions to the path integral. Indeed, the path integral can be calculated by just considering a few specific geometries known as "Instantons" having considerably large contributions. An instanton describes the spontaneous appearance of a universe from nothing. People have found different types of instantons that can provide the initial conditions for the realistic universes. The first attempt to find an instanton within the context of the 'no boundary' proposal was made by Hawking and Moss [41]. The Hawking-Moss instanton describes the creation of an eternally inflating universe with closed spatial three-geometries. Later in 1987, the Coleman-De Luccia instanton was discovered to overcome the limitation of only having closed spatial three-geometries [42]. They showed that the false vacuum decay proceeds via the nucleation of bubbles whose interior is an infinite open universe in which inflation may occur. Hawking and Turok in 1998 proposed a new class of instantons that give rise to open universes, in a similar way to the instantons of Coleman and De Luccia, without requiring the existence of a false vacuum [43]. The Hawking-Turok instantons essentially make use of the fact that in de Sitter space all curvatures are equivalent. Different slicing of the 5-dimensional de Sitter hyperboloid correspond to different curvatures when considered as a 4-dimensional model [44]. But if we are to require an open or flat universe based on observations, instead of closed one, the production of such universes using Coleman-De Luccia instantons or Hawking-Turok instantons would seem a rather convoluted procedure. Indeed, since the non-closed universes can also be compact by topological identifications [45] it seems possible and plausible to study the quantum cosmology of a universe of arbitrary curvature, instead of starting with a closed universe that can later by quantum tunneling create locally open regions. Specifically, it is discussed in [46] that although one is generally interested in quantum cosmology of closed universes, but one may retain all three values of the curvature in this study. Such point of view has also been followed in detail by Coule and Martin [47]. They have shown that in the flat or open universe, the superpotential of Wheeler-DeWitt equation is significantly modified, namely the forbidden region (Euclidean region) goes away and the qualitative behavior of a typical wavefunction differs from that of closed universe. Because of the absence of Euclidean nature of the model, the smooth geometric picture of the Hartle-Hawking "no boundary proposal" seems to be lost. Restricting to the Tunneling boundary condition [48], and applying it to each of closed, flat and open universes, it is shown that the quantum cosmology actually favors the open universe. Considering the above discussions, it turns out that unlike closed models which have a forbidden Euclidean region at small scale factors, one can work with quantum cosmological models with arbitrary curvature, such as flat and open universes, whose superpotential of Wheeler-DeWitt equation does not contain a forbidden Euclidean region. Actually, there are plenty of interesting quantum cosmological models for flat and open universes in which rather than focusing on the closed universes and calculating the "creation from nothing" probability, the Wheeler-DeWitt equation is obtained for flat or open universe and those solutions of Wheeler-DeWitt equation are taken with the criteria of just having good asymptotic behaviour in minisuperspace giving rise to normalizable states or wave-packets [49]. A relevant work to the present paper has also been recently reported [50] in which the quantum cosmology for the open FRW universe was studied based on the Hamiltonian formalism for massive gravity theory and the corresponding wave packet solutions were obtained. Motivated by this recent study, in this paper, we shall study the classical cosmology for the flat and open FRW universes as well as quantum cosmology for the flat FRW universe based on the Hamiltonian formalism of the FRW cosmology for massive bigravity theory [51][52][53][54][55]. In this regard, we use the canonical formulation and quantization of the phase space variables of the minimal massive bigravity model [56]. We derive the Schrodinger-Wheeler-DeWitt equation and find its exact and wave packet solutions and then discuss on their properties.
The organization of this paper is as follows. In section 2, first we review the minimal bigravity action with two metrics g µν and f µν in the presence of perfect fluid and find the point-like Lagrangian where the graviton's mass plays the role of a cosmological constant. Then, we extract the Hamiltonian density and also the Hamiltonian constraint equations for both metrics g µν and f µν . In section 3, we obtain the classical solutions, specially in the Schutz representation of perfect fluid [59,60], and investigate their different aspects such as the effect of graviton's mass as a cosmological constant, in addition to the appearance of the singularities and also the late time accelerated expansion. In section 4, we study the quantum cosmology in Schutz representation by extracting the Schrodinger-Wheeler-DeWitt equation for g µν and f µν , and find the wave function describing the quantum behavior of the universe. Interpretation of the wavefunction leads to the result that the graviton's mass at early universe should be large and this may justify the large cosmological constant at early universe. The paper ends with a conclusion in section 5. Here, we work in units where c = = 1.

II. POINT-LIKE LAGRANGIAN AND HAMILTONIAN CONSTRAINT IN MINIMAL BIGRAVITY THEORY
The action of Hassan-Rosen theory named bigravity theory [19] has the following structure (1) Here, g µν and f µν are two dynamical metric tensors with R andR as the scalar curvatures for tensors g µν and f µν , respectively, and L m (g, Φ) is the matter source containing an scalar field Φ. Meanwhile, m is the mass of graviton and M eff is defined as The tensor g −1 f means The trace of this tensor as X µ µ or [X] helps us to write the following expressions for e n (X)'s More simplification motivates us to study the minimal but a non-trivial case as follows [56] S bi =M 2 where use has been made of (3) as Note that considering the non-minimal models leads to quite complicate calculations, whereas in the minimal model the interaction term of two metrics g µν and f µν is just obtained by the trace of g −1 f µ ν . Let us now obtain the equations of motion by varying the action (4) with respect to g µν and f µν , respectively as and As a consequence of the Bianchi identity and also the covariant conservation of T µν , equation (5) leads to the Bianchi constraint for g µν Similarly, equation (6) gives us the Bianchi constraint for f µν In order to extract the point-like Lagrangian of this theory, we assume two dynamical line elements as where N (t) and M (t) are the lapse functions, a(t) and b(t) are the scale factors of metrics g µν and f µν , respectively, and K is the space curvature which is assumed to be the same for both metrics [57,58].
Upon substitution of these metric coefficients and the definitions of Ricci scalars into the equations (3) and (4), besides the simplification assumption M 2 g = M 2 f = M 2 eff /2, we can obtain a point-like form for the gravitational Lagrangian in the minisuperspace {N, a, M, b} as For the metrics (9) and (10), the Bianchi constraints (7) or (8) equivalently gives which is an important result for expediting our next calculations. In order to check the correctness of the relation (11), we calculate the equations of motion of the variables M and N , respectivelẏ On the other hand, by inserting the line elements (9) and (10) in the field equations (5) and (6), defining the Hubble parameters H =ȧ N a and L =˙b b , and making use of the Bianchi constraint (12), we can find the following Friedmann equations The above equations are exactly the same equations (13) and (14) obtained by using the point-like Lagrangian (11), so this shows the correctness of the point-like Lagrangian. Comparing (15) (and (16)) with the standard Friedmann equation reveals that the squared mass of graviton plays the role of a cosmological term (provided that b a takes constant value) which is positive (negative) for b > a and negative (positive) for b < a. This is interesting because one can interpret the change from a negative to a positive cosmological term, due to the dynamical competition between two scale factors, as a phase transition from deceleration to acceleration era. In other words, it is possible to account for the recent acceleration of the universe by this model which accommodates a desirable sign and value change for the cosmological term so that the universe can experience an acceleration era after a dynamical sign and value change of the cosmological term from negative to positive. The dynamical character of the cosmological term, due to the dynamical characters of two scale factors, may also alleviate the coincidence problem in the sense that it can make it possible to have a variable dark energy density with the same order of magnitude of the variable matter density. We will not pursue these interesting issues here and leave them to be reported elsewhere.
Clearly, for the time being we have not considered the matter term in the above Lagrangian. The momentum conjugate to a and b are obtained as and respectively. The Hamiltonian can be extracted by Legendre transformation H =ȧP a +ḃP b − L, where the terms P N and P M do not appear in the Hamiltonian because they are the conjugate momenta of two non-dynamical variables N and M respectively, so we have Therefore, the Hamiltonian constraint reads as H = 0. At quantum level, the operator form of this constraint is called Wheeler-DeWitt equation which annihilates the wave function of the universe asĤΨ = 0. Now, after the study of geometric sector of the model, we turn to the sector of matter field. It should be noted that the matter part of the action is thoroughly independent of the modification of the bigravity model by mass term. Thus, we can add the matter part Hamiltonian to the geometric part Hamiltonian (19) to obtain the total Hamiltonian expression. In order to get the matter part Hamiltonian we use a perfect fluid matter with the equation of state where p is the pressure, and ρ is the energy density. According to Schutz representation for the perfect fluid [59,61], the matter part Hamiltonian can be written as where P T is the conjugate momentum corresponding to the dynamical variable T as the collective thermodynamical parameters of the perfect fluid. In fact, this representation will show its advantages in the study of quantum cosmological part of our model in that a time parameter is constructed naturally in terms of thermodynamical variables of the perfect fluid. Regarding (17) and (18), we can writė Therefore, the Hamiltonian constraint casts in the following form where it is seen that the matter is just coupled with the metric g µν as had been specified in the Lagrangian. In the next section, we are going to study the classical and quantum cosmologies of this model by considering its Lagrangian and Hamiltonian which have been extracted explicitly in this section.

III. CLASSICAL COSMOLOGICAL DYNAMICS
Let us take the classical point of view for the cosmological dynamics of the massive bigravity model including the matter part. We consider the Friedmann equations including the energy density ρ as Here, we should remind that although there are several related works about the classical cosmological studies in bigravity that provide an expanded spectrum of cosmological solutions, we are still motivated to seek the behavior of the scale factor in the particular case of massive bigravity theory including the line elements (9) and (10) in a more transparent way. In this regard, we define the following relations which in the minimal massive bigravity leads to some specific cosmological solutions [62]. Using (26), the Friedmann equations (24) and (25) can be written as Subtracting two above equations to eliminate H 2 and K a 2 gives the following equatioñ Therefor, we can exactly solve this equation to find, among other two additional complex solutions, the following real expression for λ It is fruitful to notice, in the late time behavior of an expanding universe, that ρ generally becomes constant so-called ρ vac which means the vacuum energy density. Obviously, the dimensionless parameter λ in (30) approaches a constant value at late times, as well. Thus, we can interpret the Friedmann equations (27) and (28) as the equations which clearly are implying that the late time universe always approaches asymptotically towards a de Sitter or anti-de Sitter universe. For the early universe limit withρ → ∞ we again use the equations (30) which shows the result λ → 0.
Having considered the Friedmann equation for g µν (27), we are ended up with an evolution for H 2 dominated by the very large value ofρ and a rather large cosmological constant, similar to the general relativity. As we mentioned, besides the solution (30) we have two more complex solutions which diverge and so are omitted from the explicit physical class of solutions in the limitρ → ∞. Thus, we can deduce that the solution (30) which vanishes in this limit, is the physical one because we may have ordinary general relativity in the very early expansion history of the universe.
In summary, the Friedmann equations in this model support the universe subject to the ordinary general relativistic Friedmann equations with a cosmological constant proportional to the squared mass of gravitons m 2 , and depending on the parameter of the theory they can describe asymptotically de Sitter or anti-de Sitter universe. Now, we study the behavior of the physical scale factor a(t) or the corresponding Hubble parameter H to show the details of the above discussion. In doing so, we may insert the λ expression (30) in (27) to extract the behavior of the scale factor with respect to time. Defining λ * = 2 (1 − λ) −ρ as an specified expression which only includesρ, we can writeȧ In order to find physical non-oscillating scale factor, we assume the gauge N = 1 at late times with constant λ * < 0. Under this assumption, we have the following cases for two values of K: We get the following form of modified Friedmann equation in bigravity model which has two following branch of solutions For a limited time interval at late time, they are two branches of solutions describing de Sitter and anti de Sitter-like solutions. Solutions for a flat space with two opposite behaviors of the scale factor, without any singularity in this limited time interval, is shown in Figure 1.
It is worth mentioning that the conformal time τ is equivalent to the cosmic time t since we have chosen the lapse function N = 1. Anyway, these solutions for an open universe are treated in a manner that with t 0 ≥ 0 , t > 0 and a(t) ≥ 0, they are divided into two branches in which a + (t) and a − (t) are valid for t ≥ −t 0 and t ≤ t 0 respectively, such that a ± (t 0 ) = 0. In general, the a + (t) evolution begins with a big-bang-like singularity at t = −t 0 and then goes on its evolution exponentially; at late time of cosmic evolution that is exactly our case, the mass term m 2 shows itself as a cosmological constant. The evolution is opposite for the case, a − (t), in a way that the universe gradually decreases its scale factor from a large initial size at t = +∞ towards a zero size at t = t 0 . Nevertheless, the solutions for K = −1 generally differ from the previous ones for K = 0, because in the present solutions the initial singularities appear, whereas the de Sitter or (Ads) solution does not carry any initial singularity; in other words, the present solutions are similar to the previous ones just at late time. It should be noted that the closed universe case K = +1 dose not have any preferable physical solution.
We can follow another approach such that the matter part is considered as perfect fluid, similar to the previous section. For simplicity in the calculations of classical and quantum cosmology (see bellow), we consider the Schutz representation [59,61] for the perfect fluid which is coupled with the gravity. In this case, the Hamiltonian (21) describes the dynamics of the system. Thus, we have the following equations of motion for T and P Ṫ Considering the above equations, we conclude that choosing N = a 3ω will help us to write showing the time role of the variable T . Before going through the Friedmann equations, we set up the following assumption where γ is a function to be obtained later. This assumption gives the opportunity to have an analytical calculation and getting some explicit results for the scale factor a(t). Taking P T = P 0 = const, with some similar calculation of the equations (27) and (28) we can writė Considering Eqs. (38), (39), we deduce that γ should be at most a function of the scale factor a which according to (12) means that the lapse coefficient M should be treated as a function of the scale factor a, too. In fact, we are interested in this behavior because it can describe a new set of cosmological solutions differing from the de-Sitter or Ads universes. Subtracting equations (38) and (39) gives For K = 0 we can classify the results in two cases: • ω = − 1 3 : cosmic string By using (40), we obtain γ = 1 + 4m 2 a 2 + P 0 2 1 3 12m 6 a 6 + 3m 4 a 4 P 0 + m 6 a 6 (m 2 a 2 − P 0 ) (4m 2 a 2 + P 0 ) 12m 6 a 6 + 3m 4 a 4 P 0 + m 6 a 6 (m 2 a 2 − P 0 ) (4m 2 a 2 + P 0 ) Inserting the above result in the equation (38) helps us to extract the time evolution of the scale factor by means of solving the following reduced equation 12m 6 a 6 + 3m 4 a 4 P 0 + m 6 a 6 (m 2 a 2 − P 0 ) (4m 2 a 2 + P 0 ) 12m 6 a 6 + 3m 4 a 4 P 0 + m 6 a 6 (m 2 a 2 − P 0 ) (4m 2 a 2 + P 0 ) The behavior of scale factor with respect to the cosmic time is plotted in Figure 2, as the solution for this equation. Clearly, this is describing a universe in which we have a big-bang singularity at t = 0 and a positive-valued scale factor increasing monotonically.

FIG. 2:
The time evolutionary behavior of the scale factor a(t) according to (42). We have considered the numerical values m = 1 and P0 = −10 6 .
• ω = −1: cosmological constant By using (40), we obtain a constant value γ 0 and so which has the following solution a ± (t) = 9 2 By means of the relation between the conformal time and the cosmic time dτ = a(t) −1 dt we make the above result more comprehensive as Considering the condition a(τ ) > 0, we definitely can ignore the case 2 (1 + γ 0 ) m 2 + P 0 > 0 and also we should emphasize that the condition 2 (1 + γ 0 ) m 2 + P 0 < 0 imposes some specific limitation on the parameters m and P 0 which we disregard here further discussion about its details. Plotting the evolution of the scale factor for this case in Figure 3, we can see a big-bang-like singularity at τ = τ 0 which starts from the singular point and grows gradually with time.

IV. QUANTUM COSMOLOGICAL DYNAMICS
Let us study the quantization of the present model by means of the canonical quantization. To this end, we obtain the Wheeler-DeWitt equationĤΨ = 0, as the operator version of the Hamiltonian constraint H = 0, and solve it to find the universe wave function Ψ. By means of (23), the Wheeler-DeWitt equation reads aŝ We sperate the variables in Ψ (a, b, T ) as follows where E is a constant quantity. This gives rise tô Let us assume, for simplicity, that (b/a) = σ = cte. Then, the Bianchi constraint (12) takes the following form Using (17), (18) and (49), we find that P b = σP b which leads to the following relation between two momentum operatorsP b = σP a .
As a result, we have the final form of the Wheeler-DeWitt equation as This is a direct consequence of the coupling between the matter source and the gravitation g µν which states that the dynamics of two metrics g µν and f µν are identified and so the two-variable system is mostly simplified in favor of a one-variable system including just one scale factor a. In other words, we are allowed to take Ψ (a, b) ≡ Ψ (a, σ) ≡ ψ (a).
Note that this Schutz representation is related directly to the well-known time problem in quantum cosmology context [63]. Clearly, in (46) the conjugate momentum associated with T appears linearly in the Hamiltonian of the model. By means of the canonical quantization, we obtain a Schrodinger-Wheeler-De Witt (SWD) equation in which the matter variable T plays the role of time. It is noticeable that this approach dose not solve completely the time problem in massive bigravity quantum cosmology in a primordial way. However, applying a perfect fluid instead of a field is somehow leads to the emergence of a Schrodinger-type equation with a useful time parameter T . The particular solution of this Schrodinger-type equation is going to be fined in the rest. Here we should consider the time ordering rule of the operatorsâ andP a by means of the relationP 2 a = −a −p ∂ ∂a a p ∂ ∂a , thus (51) becomes The generic solutions of the above equation for the case ω = −1 and K = 0 (Taking K = 0 is due to the fact that it is in agreement with the current observations) is as follows By the choice of the factor ordering value p = −2 (the factor p denotes the uncertainty in the choice of operator ordering ) [64,65] we have Imposing the boundary condition on these solutions as ψ (a = 0) = 0, we are led to the following result with c 1 = 0.
ψ Eσ (a) = ic 2 sinh It is obvious that the wave function ψ Eσ (a) does not satisfy the squared integrability condition unless we transform it into a "Sin" function. In doing so, we may impose the constraint 3E + 2m 2 3 − 4σ + σ 4 < 0 so that the wave function becomes oscillatory The eigenfunctions of the Wheeler-DeWitt equation may be written as In Figure 4, we have plotted the square of the eigenfunction ψ Eσ (a) (55) for typical values of the parameters. It seems from this figure that we have a well-defined wave function in the vicinity of a = 0 describing a universe emerging out of nothing without any tunneling. This is because there is no potential barrier in the classical Hamiltonian (52) for tunneling effect. It is seen that the squared wave function has peaks, in the vicinity of specific and countless values of growing scale factors, which are rapidly decaying as the scale factor grows. If we consider this as an instant spatial configuration of probability for observing a universe, we realize that the universe is observed with most probability, corresponding to the dominant contribution, at the smallest (nonvanishing) scale factor among other possible scale factors. This means that although there is no tunnelling from nothing scenario for this universe, however, similar to the tunnelling scenario the universe is born from nothing with most probability, close to a = 0. In other words, the quantum cosmology can mimic the tunneling scenario and avoid the classical initial singularity by supporting the probability for creation of universe from nothing at a desired small nonvanishing scale factor. By extremizing |ψ Eσ (a)| 2 , using (56), it is easily found that the size of initial scale factor corresponding to the dominant contribution of |ψ Eσ (a)| 2 is determined by the mass of gravitons, and that a small size newborn universe requires a large mass for gravitons. Such large mass for gravitons at early universe has already been proposed in the study of Hartle-Hawking no-boundary proposal in dRGT massive gravity [66]. They have shown that the contribution from the massive gravity sector can substantially enhance the probability of a large number of e-folding for a sufficiently large graviton mass comparable to the Hubble parameter during inflation, namely m 10 12 GeV , and illustrated possible models to trigger such a large graviton mass at early universe while it is negligibly small in the present universe.
We may reach to the most precise description of the wave function, provided that we construct a wave packet for this wave function. In order to obtain the general solution as a wave packet we may have superposition of the eigenfunctions Ψ Eσ (a, T ) with a suitable Gaussian weight function e −γE 2 as The Gaussian weight function e −γE 2 and other examples (quasi-Gaussian and shifted Gaussian) are extensively used in quantum mechanics as a tool to obtain localized states. These types of weight functions are concentrated about a specific value of their argument and fall off swiftly away from that central point. Multiplying the Gaussian weight function by the obtained wave function e iET ψ Eσ (a), we finally find from (58) a Gaussian-like behavior for the wave packet which is localized about some special values of its argument a. Since the above integral is too complicate to extract an analytical closed form for the wave function, we have just plotted the squared wave function |Ψ σ0 (a, T ) | 2 in Figure 5, where σ 0 = 4.7. The figure indicates that the wave packet spreads out along both the spatial axes a and temporal axes T .

V. CONCLUSIONS
In this paper, we have applied the recently proposed nonlinear massive bigravity theory to a FRW cosmological model. Using the Hamiltonian formalism for bigravity theory we have studied the classical and quantum cosmological behaviors of the particular massive bigravity model so called minimal bigravity. It is shown that the classical field equations receive contribution from the massive gravitons as a cosmological term whose value can change dynamically from negative to positive value, depending on the competition between two scale factors of bigravity metrics. Such a dynamical cosmological term is capable of accounting for the recent acceleration of the universe as well as resolving the coincidence problem, to which we shall address elsewhere. We have presented the classical cosmological solutions for the cases where the universe contains a perfect fluid with the equation of state parameters ω = −1 (cosmological constant), and ω = − 1 3 (cosmic string). For the general energy density ρ, we have found some solutions of two Friedmann equations in flat (K = 0) and open (K = −1) cases. In the flat case, we have obtained two branches of solutions in which there are no singularities. Similarly, in the open universe we have obtained two branches of solutions which are analyzed in details leading to the similar results as those of the flat case except that in this case we have initial singularities. For the energy density with Schutz representation of perfect fluid, applied for ω = − 1 3 (cosmic string), we have found a universe in which there is a singularity, and then we plotted its scale factor evolution in Figure 2. We have also studied another case ω = −1 (cosmological constant) and found that the situation is similar to the previous case except for the scale factor which evolves more rapidly than that of the cosmic string, and then we plotted its evolution in Figure 3.
We have investigated the quantization of the minimal bigravity model via the method of canonical quantization in a two variable minisuperspace including two scale factors. Due to the presence of interaction term between two scale factors in the massive bigravity action and using the Bianchi constraint, we could reduce the two variable minisuperspace into one variable minisuperspace including just one scale factor, which is preferred to be "a(t)". Having used Schutz representation for the perfect fluid, beside a particular gauge choice N = a −3 we have introduced a parameter playing the role of time in a Schrodinger-Wheeler-DeWitt equation. In the presence of matter for a particular vacuum case ω = −1 and K = 0, corresponding to early universe, we have solved the Schrodinger-Wheeler-DeWitt equation exactly and constructed the wave packet corresponding to this solution, numerically. We have found that although there is no tunnelling from nothing scenario for this one variable quantum cosmology, however, this quantum cosmology can avoid the classical initial singularity by supporting the probability for creation of universe from nothing at a nonvanishing scale factor. Moreover, it is shown that the demand for a small initial scale factor for the newborn universe, requires a large mass for gravitons. Such a large mass of graviton (comparable to the Hubble parameter during inflation, namely m ∼ 10 12 GeV) at early universe is in agreement with the fact that the cosmological constant should be large at early universe. Actually, the same result has already been proposed in the study of Hartle-Hawking no-boundary proposal for the quantum cosmology of de Rham, Gabadadze and Tolley (dRGT) massive gravity [66] where the authors have given two reasons for why the graviton can have large value at early universe while it is negligible today. We may add another reason based the application of uncertainty principle on the universe as a single quantum system. The very small size of universe at early times confines the range of gravitons to be very small, at most, of the size of universe. According to the uncertainty principle ∆P ∆X ∼ , such mediating particles with a small range ∆X should have a large momentum ∆P , and this may justify the large mass for gravitons at early universe. As the universe expands, the range of gravitons becomes larger and, according to the uncertainty principle, the mass of such mediating gravitons becomes smaller. In other words, if we assume the universe with the age ∆t ∼ H −1 , the mass of gravitons, according to the uncertainty principle ∆E∆t ∼ , should be of the size m ∼ H/c 2 which is in agreement with the present prediction on the graviton's mass and also the prediction of [66] where it is shown that the mass of graviton at early universe was comparable (in the natural units = c = 1) to the Hubble parameter during inflation. The identification of graviton's mass with the cosmological term and the reduction of graviton's mass by time evolution of the universe may be considered as a solution of the cosmological constant problem.