# Qualitative Analysis of

Viscous Fluid Cosmological Models

satisfying the Israel-Stewart

theory of Irreversible Thermodynamics

###### Abstract

Isotropic and spatially homogeneous viscous fluid cosmological models are investigated using the truncated Israel-Stewart [W. Israel, Ann. Phys., 100 1976; W. Israel and J.M. Stewart, Proc. R. Soc. Lond. A., 365 1979; ibid. Ann. Phys., 118 1979] theory of irreversible thermodynamics to model the bulk viscous pressure. The governing system of differential equations is written in terms of dimensionless variables and a set of dimensionless equations of state is then utilized to complete the system. The resulting dynamical system is analyzed using geometric techniques from dynamical systems theory to find the qualitative behaviour of the Friedmann-Robertson-Walker models with bulk viscosity. In these models there exists a free parameter such that the qualitative behaviour of the models can be quite different (for certain ranges of values of this parameter) from that found in models satisfying the Eckart theory studied previously. In addition, the conditions under which the models inflate are investigated.

###### pacs:

04.20.Jb, 98.80.Hw6

[Causal Viscous Fluid FRW Models]

## 1 Introduction

Recently, spatially homogeneous imperfect fluid models have been investigated using techniques from dynamical systems theory [1, 2, 3, 4]. In these papers, dimensionless variables and a set of dimensionless equations of state were employed to analyze various spatially homogeneous imperfect fluid cosmological models. It was also assumed that the fluid is moving orthogonal to the homogeneous spatial hypersurfaces; that is, the fluid 4-velocity, , is equal to the unit normal of the spatial hypersurfaces. The energy-momentum tensor can be decomposed with respect to according to [5]:

(1) |

where and is the energy density, is the thermodynamic pressure, is the bulk viscous pressure, is the anisotropic stress, and is the heat conduction vector as measured by an observer moving with the fluid.

In papers [1, 2, 3, 4] it was also assumed that there is a linear relationship between the bulk viscous pressure and the expansion (of the model), and a linear relationship between the anisotropic stress and the shear ; that is, \numparts

(2) | |||||

(3) |

where denotes the bulk viscosity coefficient and denotes the shear viscosity coefficient. Equations (2) and (3) describe Eckart’s theory of irreversible thermodynamics. Eckart’s theory is a first order approximation of the viscous pressure and the anisotropic stress and is assumed to be valid near equilibrium [6]. Coley and van den Hoogen [1] and Coley and Dunn [3] have studied the Bianchi type V cosmological models using equations (1.2) as an approximation for the bulk viscous pressure and the anisotropic stress. They found that if the models satisfied the weak energy conditions, then the models necessarily isotropize to the future. Belinskii and Khalatnikov [7], with different assumptions on the equations of state, found similar behaviour present in the Bianchi type I models. The addition of viscosity allowed for a variety of different qualitative behaviours (different from that of the corresponding perfect fluid models). However, since the models studied in [1, 2, 3, 4] satisfy Eckart’s theory of irreversible thermodynamics [6], they suffer from the property that signals in the fluid can propagate faster than the speed of light (i.e., non-causality), and also that the equilibrium states in this theory are unstable (see Hiscock and Salmonson [8] and references therein). Therefore, a more complete theory of irreversible thermodynamics is necessary for fully analyzing cosmological models with viscosity.

Among the first to study irreversible thermodynamics were Israel [9] and Israel and Stewart [10, 11]. They included additional linear terms in the relational equations (1.2). Assuming that the universe can be modelled as a simple fluid and omitting certain divergence terms, the linear relational equations for the bulk viscous pressure, the heat conduction vector, and the shear viscous stress are [11]: \numparts

(4) | |||||

(5) | |||||

(6) |

where is the projection tensor and

Belinskii et al. [12] were the first to study cosmological models satisfying the truncated Israel-Stewart theory of irreversible thermodynamics. Using qualitative analysis, Bianchi I models were investigated with a relational equation for the bulk viscous pressure and the shear viscous stress of the form (1.3). They also assumed equations of state of the form [12]

(7) |

where and are constants and and are parameters. The isotropizing effect found in the Eckart models no longer necessarily occurred in the truncated Israel-Stewart models. It was also found that the cosmological singularity still exists but is of a new type, namely one with an accumulated “visco-elastic” energy [12]. Similar to the work done by Belinskii et al. [12], Pavón et al. [13] and Chimento and Jakubi [14] studied the flat Friedmann-Robertson-Walker (FRW) models. They assumed the same equations of state as Belinskii et al. [12], namely \erefBel eqs state, but studied the models using slightly different techniques. Chimento and Jakubi [14] also found exact solutions in the exceptional case (which will be of interest later). They found that the future qualitative behaviour of the model was independent of the value of ; however, to the past, “bouncing solutions” and deflationary evolutions are possible [14].

In addition, Hiscock and Salmonson [8] investigated further generalizations of the truncated Israel-Stewart theory of irreversible thermodynamics. Namely, they included the non-linear divergence terms in the relational equations (1.3). We shall refer to such a theory as the “Israel-Stewart-Hiscock” theory. Hiscock and Salmonson used equations of state arising from the assumption that the fluid could be modelled as a Boltzmann gas. They concluded that when the Eckart equations, (1.2), or the truncated Israel-Stewart equations, (1.3), were used, inflation could occur, but if the non-linear terms of the Israel-Stewart-Hiscock theory were included, inflation was no longer present. This result led them to conclude that “inflation is a spurious effect produced from using a truncated theory” [8]. However, Zakari and Jou [15] also employed equations arising from the full Israel-Stewart-Hiscock theory of irreversible thermodynamics, but assumed equations of state of the form \erefBel eqs state and found that inflation was present in all three theories (Eckart, truncated Israel-Stewart, Israel-Stewart-Hiscock). Therefore, it appears that the equations of state chosen determine if the model will experience bulk-viscous inflation. Romano and Pavón [16] also analyzed Bianchi III models using both the truncated Israel-Stewart theory and the Israel-Stewart-Hiscock theory. They only analyzed the isotropic singular points, but concluded that the qualitative behaviour of the models in the two different theories was similar in that the anisotropy of the models dies away and the de Sitter model is a stable attractor.

In this work we use equations of state that are dimensionless, and hence we are generalizing the work in [1, 2, 3, 4] in which viscous fluid cosmological models satisfying equations (1.2) were studied. One reason for using dimensionless equations of state is that the equilibrium points of the system of differential equations describing spatially homogeneous models will represent self-similar cosmological models [17]. In addition, it could be argued that the use of dimensionless equations of state is natural in the sense that the corresponding physics is scale invariant (see also arguments in Coley [18]).

The intent of this work is to build upon the foundation laid by Belinskii et al. [12], Pavon et al. [13] and Chimento and Jakubi [14], and investigate viscous fluid cosmological models satisfying the linear relational equations (1.3). We will use dimensionless variables and dimensionless equations of state to study the qualitative properties of isotropic and spatially homogeneous cosmological models. In particular, we shall study this new “visco-elastic” singularity and we shall determine whether bulk-viscous inflation is possible. We will also determine if there is a qualitative difference between these models and the models studied by Burd and Coley [4] where the Eckart equation \erefeckart was assumed.

In section 2 we define the models and establish the resulting dynamical system. In section 3 we investigate the qualitative behaviour of the system for different values of the physical parameters. In section 4 we discuss and interpret our results and in section 5 we end with our conclusions. For simplicity we have chosen units in which .

## 2 Friedmann-Robertson-Walker Models

In this paper we assume that the spacetime is spatially homogeneous and isotropic and that the fluid is moving orthogonal to the spatial hypersurfaces. The energy-momentum tensor considered in this work is an imperfect fluid with a non-zero bulk viscosity (that is there is no heat conduction, , and no anisotropic stress, ). The Einstein field equations () and the energy conservation equations () can be written as (see Burd and Coley [4]): \numparts

(8) | |||||

(9) | |||||

(10) |

where is the expansion and is the curvature of the spatial hypersurfaces. If the curvature is negative, i.e., , then the FRW model is open, if then the model is flat, and if then the FRW model is closed. Assuming that the energy density, , is non-negative, it is easily seen from \erefgen.equality that in the open and flat FRW models the expansion is always non-negative, i.e. , but for the closed FRW models the expansion may become negative. (Great care must be taken in this case because the dimensionless quantities that we will be using become ill-defined at .)

We can obtain an evolution equation for by solving \erefisrael1 for ,

(11) |

Now the system of equations defined by equations \erefraychaudhuri, \erefrho, and \erefnewtau constitute a dynamical system of the form , where . This system of equations is invariant under the mapping (see Coley and van den Hoogen [1, 17])

(12) |

and this invariance implies that there exists a symmetry in the dynamical system [19]. Therefore, we introduce new dimensionless variables , and a new time variable as follows:

(13) |

and consequently the Raychaudhuri equation, \erefraychaudhuri, effectively decouples from the system.

In order to complete the system of equations we need to specify equations of state for the quantities , , and . In principle equations of state can be derived from kinetic theory, but in practice one must specify phenomenological equations of state which may or may not have any physical foundations. Following Coley [18, 20], we introduce dimensionless equations of state of the form \numparts

(14) | |||||

(15) | |||||

(16) |

where , , and are positive constants, and , , and are constant parameters ( is the dimensionless density parameter defined earlier). In the models under consideration, is positive in the open and flat FRW models, thus equations (2.6) are well defined. In the closed FRW model the expansion could become zero, in which case these equations of state break down. However, we can utilize these equations to model the asymptotic behaviour at early times, i.e., when . The most commonly used equation of state for the pressure is the barotropic equation of state , whence and (where is necessary for local mechanical stability and for the speed of sound in the fluid to be no greater than the speed of light).

We define a new constant . Using equations (2.5) and (2.6), we find that equations \erefrho and \erefnewtau reduce to \numparts

(17) | |||||

(18) |

Also, from the Friedmann equation, \erefgen.equality, we obtain

(19) |

Thus, the line divides the phase space into three invariant sets, , , and . If , then the model is necessarily a flat FRW model, if then the model is necessarily an open FRW model, and if the model is necessarily a closed FRW model.

The equilibrium points of the above system all represent self-similar cosmological models, except in the case . If , the behaviour of the equations of state, equation (2.6), at the equilibrium points, is independent of the parameters and ; namely the behaviour is

(20) |

Therefore, natural choices for and are respectively , . We note that in the exceptional case , there is a singular point which represents a de Sitter solution and is not self-similar. (This is also the case in the Eckart theory as was analyzed by Coley and van den Hoogen [1].)

To further motivate the choice of the parameter , we consider the velocity of a viscous pulse in the fluid [15],

(21) |

where corresponds to the speed of light. Using \erefnew vars and equations (2.6), we obtain

(22) |

Now, if then not only do we obtain the correct asymptotic behaviour of the equation of state for the quantity but we are also allowed to choose since then the velocity of a viscous pulse is less than the velocity of light for any value of the density parameter . Thus in the remainder of this analysis we shall choose . In order for the system of differential equations (2.7) to remain continuous everywhere, we also assume .

## 3 Qualitative Analysis

### 3.1

We now study the specific case when . In this case there are three singular points,

(23) |

where

(24) |

The point has eigenvalues

(25) |

This point is either a saddle or a source depending on the value of the parameter, ;

(26) |

If , then the point is a saddle point and if , then the point is a source. If (the bifurcation value), then the point is degenerate (discussed later).

The point has eigenvalues

(27) |

If , then the point is a saddle point and if , then the point is a source. If (the bifurcation value), then the point is degenerate (discussed later).

The singular point has eigenvalues

(28) |

This singular point is a sink for (See also \Treftable I for details).

In addition to the invariant set , there exist two other invariant sets. These are straight lines, , where

(29) |

The invariant line passes through the singular points and while the line passes through the singular points , and . These invariant sets represent the eigendirections at each of the singular points [see also A].

In order to sketch a complete phase portrait, we also need to calculate the vertical isoclines which occur whenever . From (2.7) we can see that this occurs either when or when . This straight line passes through the origin, and through the singular point if . If then the vertical isocline has a negative slope which is greater than the slope of the slope of the invariant line , [i.e., ], and when the vertical isocline has a negative slope which is less than the slope of the invariant line [i.e., ].

To complete the analysis of this model we need to analyze the points at infinity. We do this by first converting to polar coordinates and then compactifying the radial coordinate. We change to polar coordinates via

(30) |

and we derive evolution equations for and . We essentially compactify the phase space by changing our radial coordinate and our time as follows,

(31) |

that is, the plane is mapped to the interior of the unit circle, with the boundary of this circle representing points at infinity of . We have (for and general )

(32) | |||

(33) |

We easily conclude that if (or any ), then the entire circle, , is singular. Therefore, we have a non-isolated set of singular points at infinity. To determine their stability we look at the sign of as . In this case we see

(34) |

which implies that points above the line are repellors, while those points which lie below the line are attractors.

For completeness, we would also like to determine the qualitative behaviour of the system at the bifurcation value where the singular points are and the line of singular points . (Note that since , .) Fortunately we are able to completely integrate the equations in this case to find

(35) |

where is an integration constant. We see that all trajectories are straight lines that pass through the point . It is straightforward to see that the line of singular points are repellors while the point is an attractor. We are now able to sketch complete phase portraits (See Figures References, References and References).

### 3.2 and

This is a case of particular interest since it represents the asymptotic behaviour of the FRW models for any and (since at the singular points the viscosity coefficient behaves like and the relaxation time like ). Note that the physical phase space is defined for , but the system is not differentiable at . In this case there are four singular points,

(36) |

where

(37) |

and and are given be equation \erefsingular points 2.

The dynamical system, (2.7) is not differentiable at the singular point . We can circumvent this problem by changing variables to and a new time variable defined by . The system then becomes \numparts

(38) | |||||

(39) |

In terms of the new variables the system is differentiable at the point , but one of the eigenvalues is zero and hence the point is not hyperbolic. Therefore in order to determine the stability of the point we change to polar coordinates, and find that the point has some saddle-like properties; however, the true determination of the stability is difficult. [We investigate the nature of this singular point numerically — see B.]

The singular point has eigenvalues

(40) |

This singular point varies both its position in phase space with its stability depending upon whether is less than, equal to or greater than one. If , then and the point is a saddle point. If , then and the point is a source. Finally, if (the bifurcation value), then and the point is degenerate (discussed later).

The stability of the points and is the same as in the privious case, see equations \erefeig1 and \erefeig2 for their eigenvalues and the corresponding text. (See also \Treftable I for details).

The vertical isoclines occur at and . This straight line is easily seen to pass through the origin and the point . If , then the vertical isocline lies below the point and if , the vertical isocline lies above the point . Finally if (the bifurcation value), the vertical isocline passes through the point .

From an analysis similar to that in the previous subsection, we conclude that there is a non-isolated set of singular points at infinity. Their qualitative behaviour is the same in this case as in the previous case; namely, points which lie above the line are repellors, while those points which lie below the line are attractors.

At the bifurcation value , the points and come together; consequently these points undergo a saddle-node bifurcation as passes through the value . The singular point is no longer hyperbolic, but the qualitative behaviour near the singular point can be determined from the fact that we know the nature of the bifurcation. Hence the singular point is a repelling node in one sector and a saddle in the others. A complete phase portrait is sketched in Figures References, References and References.

\br | ||||
---|---|---|---|---|

\mr and | source | saddle | sink | |

and | source | source | sink | |

and | saddle | source | sink | |

and | saddle | saddle | sink | source |

saddle | saddle-node | sink | saddle-node | |

and | saddle | source | sink | saddle |

\br |

.

.

.

,

These points are part of the non-isolated line singularity .

This is the situation when the points and coalesce.

## 4 Discussion

### 4.1 Exact Solutions

The exact solution of the Einstein field equations at each of the singular points represent the asymptotic solutions (both past and future) of FRW models with a causal viscous fluid source. The solution at each of the singular points represents a self similar cosmological model except in one isolated case [see the singular point ].

At the singular point we have (after a re-coordinatization) , , , , which represents the standard vaccuum Milne model.

The singular point represents a flat FRW model with a solution (after a re-coordinatization) , , , , where .

The singular point represents a flat FRW model. If then the solution is , , , , where . (Note that in this case we cannot simply change coordinates to remove the constants of integration.) The sign of depends on the sign of . If then , and if then . Thus if then is positive only in the interval and hence we can see that after a finite time , , and all approach infinity. (We will see later that the WEC is violated in this case). If then we can re-coordinatize the time so as to remove the constant of integration, , and the absolute value signs in the solution for . If then and the solution is the de Sitter model with (after a re-coordinatization) , , , . This exceptional solution is the only one that is not self-similar. It can be noted here that this is precisely the same situation that occurred in the Eckart models studied in [4] and [1].

The singular point represents either an open, flat, or closed model depending on the value of the parameter . The solution in all cases is (after a re-coordinatization) , , , .

### 4.2 Energy Conditions

The weak energy condition (WEC) states that for any timelike vector [21]. In the model under investigation this inequality reduces to and . Assuming that from here on, the WEC in dimensionless variables becomes

(41) |

The dominant energy condition (DEC) states that for every timelike , and is non-spacelike [21]. Here this inequality reduces to and which when transformed to dimensionless variables becomes

(42) |

The strong energy condition (SEC) states that [21]. Here this inequality reduces to and which when transformed to dimensionless variables becomes

(43) |

If we assume that the WEC is satisfied throughout the evolution of these models then we find that there are five distinct situations. If , then and the line intersects the line at a point . If , then and the line intersects the line at the point . If then can be of any sign or zero, but the line intersects the line at a point . If the WEC condition is satisfied throughout the evolution of these models then the possible asymptotic behaviour of the models is greatly restricted.

### 4.3 Asymptotic Behaviour

The qualitative behaviour depends on the values of and . If the parameter is different from unity then there is an additional singular point. This property is also present in the Eckart models studied by Burd and Coley [4] and Coley and van den Hoogen [1]. The value corresponds to the case when the dynamical system, (2.7), is polynomial (the only other value of that exhibits this property is ). The value is of particular interest as it represents the asymptotic behaviour of all the viscous fluid FRW models, and also, this is the case when the equation of state for is independent of (i.e., ). The parameter, , plays a role similar to the parameter found in both Burd and Coley [4] and Coley and van den Hoogen [1]. The value of the parameter determines the stability and global behaviour of the system.

One of the goals of this paper is to determine the generic behaviour of the system of equations (2.7). Using the above energy conditions, and in particular the WEC, and the phase-portraits (Figures References–References), we can determine the generic and exceptional behaviour of all the viscous fluid models satisfying the WEC. We are primarily interested in the generic asymptotic behaviour of the FRW model with viscosity: If we consider the dynamical system, (2.7), as where , are the variables, and are the free parameters then generic behaviour occurs in sets of non-zero measure with respect to the set (except for the flat models in which case the state space is a subset of ). For example, the case is a set of measure zero with respect to the set . All behaviour is summarized in \Treftable II.

\brparameters | m | models | generic behaviour | exceptional behaviour |
---|---|---|---|---|

\mr, | open | |||

flat | ||||

closed | , | |||

, | open | , | ||

flat | ||||

closed | ||||

, | ||||

, | m=1 | open | ||

flat | ||||

closed | ||||

m=1/2 | open | |||

flat | ||||

closed | ||||

, | m=1 | open | ||

flat | ||||

closed | ||||

m=1/2 | open | , | ||

flat | ||||

closed | ||||

\br |

These are exceptional trajectories and consequently do not represent typical or generic behaviour.

In the case .

Typically, if , then the open models evolve from the big-bang visco-elastic singularity at and evolve to the Milne model at [if ] or the non-vacuum open model at [if ]. If , then the closed models evolve from the big-bang visco-elastic singularity at to points at infinity. These particular points at infinity correspond to the points where (point of maximum expansion) and the various dimensionless variables breakdown.

Typical behaviour of models with depends upon the sign of . If , then all trajectories for the open models will violate the WEC and if , then the open models evolve from the big-bang visco-elastic singularity at , become open models and then evolve towards an inflationary flat FRW model at the point . A visco-elastic singularity is a singularity in which a signifigant portion of the initial total energy density is viscous elastic energy, that is . Concerning the closed models when , if then models evolve from the big-bang visco-elastic singularity at to points at infinity. However, if , then the closed models again evolve from the big-bang visco-elastic singularity at but now have two different typical behaviours. There is a class of models which approach points at infinity and do not inflate and there is a class of models which evolve towards the inflationary flat FRW model at the point .

The flat FRW models consist of a subset of measure zero of the total state space . However, the flat models are of a special interest in that the flat models represent the past asymptotic behaviour of both the open and closed models. If , then the flat models evolve from the visco-elastic singularity at to points at infinity or to the flat model located at . If , and , then the flat models evolve from the big-bang visco-elastic singularity at to points at infinity. And if and , then models evolve from the visco-elastic singularity at to points at infinity (non-inflationary) or to the inflationary model at the point .

Note, that if the WEC is dropped (i.e., let ) then a class of very interesting models occurs. There will exist models that will evolve from the visco-elastic singularity at with and start inflating at some point and at a finite time after will start expanding at increasing rates, that is , and will eventually evolve towards the point . (This is the special case mentioned in the previous subsection.) What this means in terms of the open and flat models is that they will expand with decreasing rates of expansion, start to inflate, and then continue to expand with increasing rates of expansion. For the closed models, the models will expand with decreasing rates of expansion, start to inflate, and then continue to expand with increasing rates of expansion, these models will not recollapse.

## 5 Conclusion

The only models that can possibly satisfy the WEC and inflate are those models with and . Therefore, we can conclude that bulk-viscous inflation is possible in the truncated Israel-Stewart models. However, in the models studied by Hiscock and Salmonson [8] inflation did not occur (note, the equations of state assumed in [8] are derived from assuming that the universe could be modelled as a Boltzmann gas), while inflation does occur in the models studied by Zakari and Jou [15] who utilized different equations of state. Furthermore, Maartens [22] has also analyzed models arising from the Israel-Stewart-Hiscock theory of irreversible thermodynamics. Maartens assumes an equation of state for the temperature of the from and finds that the inflationary attractor is unstable in the case [22]. In our truncated model we choose dimensionless equations of state and find that inflation is sometimes possible. The question of which equations of state are most appropriate remains unanswered, and clearly the possibility of inflation depends critically upon the equations of state utilized [15].

This work improves over previous work on viscous cosmology using the non-causal and unstable first order thermodynamics of Eckart [6] and differs from the work of Belinskii et al. [12] in that dimensionless equations of state are utilized. From the previous discussion we can conclude that the visco-elastic singularity at the point is a dominant feature in our truncated models. This singular point remains the typical past asymptotic attractor for various values of the parameters and . This agrees with the results of Belinskii et al. [12]. The future asymptotic behaviour depends upon both the values of and . If , then the open models tend to the Milne model at [] or to the open model at [], and if , then the open models tend to the inflationary model at the point or are unphysical. If , then the closed models tend to points at infinity, and if , then the closed models tend to the inflationary model at the point . The future asymptotic behaviour of the flat models is that they either tend to points at infinity or tend to the point , in agreement with the exact solution given in [14].

Belinskii et al. [12] utilized the physical variables , , and in their analysis and assumed non-dimensionless equations of state, and they found a singularity in which the expansion was zero but the metric coefficients were neither infinite nor zero; the authors passed over this observation stating that in a more realistic theory this undesirable asymptotic behaviour would not occur. We note that this behaviour did not occur in our analysis in which dimensionless variables and and a set of dimensionless equations of state were utilized.

The behaviour of the Eckart models in Burd and Coley [4] with is very similar to the behaviour of the truncated Israel-Stewart models studied here in the case . This result also agrees with the conclusions of Zakari and Jou [15]. However, when various new possibilities can occur; for instance, there exist open and closed models that asymptotically approach a flat FRW model both to the past and to the future. Interestingly enough this is also the case in which the future asymptotic endpoint is an inflationary attractor. This type of behaviour does not occur in the Eckart theory.

While this work (and the work of Belinskii et al. [12]) employs a causal and stable second order relativistic theory of thermodynamics, only the ‘truncated’ version of the theory has been utilized [22], rather than the full theory of Israel-Stewart-Hiscock [8, 15, 11, 23]. However, the ‘truncated’ equations can give rise to very different behaviour than the full equations; indeed, Maartens [22] argues that generally the ‘truncated’ theory will lead to pathological behaviour, particularly in the behaviour of the temperature. It is possible to reconcile the truncated and full theories, at least near to equilibrium, by using a “generalized” non-equilibrium temperature and pressure defined via a modified Gibbs equation in which dissipation enters explicitly [24]. Moreover, it is expected that the truncated theory studied here is applicable at least in the early universe. We have neglected the divergence terms here as a first approximation. If one includes the divergence terms then the system of equations describing the model requires an additional equation of state for the temperature . A reasonable assumption in this case might be to assume a dimensionless equation of state for .

Notwithstanding these comments it is clear that the next step in the study of dissipative processes in the universe is to utilize the full (non-truncated) causal theory [8, 11, 15, 23], and this will be done in subsequent work. In particular, we shall find that the work here will be useful in understanding the full theory and is thus a necessary first step in the analysis. In addition, anisotropic Bianchi type V models which include shear viscosity and heat conduction will also be investigated.

This research was funded by the Natural Sciences and Engineering Research Council of Canada and a Killiam scholarship awarded to RVDH. The authors would like to thank Roy Maartens for helpful discussions and for making available to us recent work prior to publication. The authors would also like to thank Des McManus for reading the manuscript.

## Appendix A First Integrals

We will use Darboux’s theorem [25] to find an algebraic first integral of the system in the case by first finding a number of a algebraic invariant manifolds. An algebraic invariant manifold, , is a manifold such that , where is a polynomial. The following are invariant manifolds of the system:

where

(44) |

Calculating , we find

Using Darboux’s Theorem, an algebraic first integral can be found by setting and then determining what values of satisfy the equation . Solving the resulting algebraic system we find the following algebraic first integral of the dynamical system, (2.7), in the case :

(45) |

where is a free parameter and and must satisfy

(46) |

This first integral determines the integral curves of the phase portraits in Figures References, References and References, where the value determines which integral curve(s) is being described. For example, if , the integral curves are and . Also, we can see that if , and then the integral curve describes an ellipse; however, these closed curves necessarily pass through the points and , thereby nullifying the possible existence of closed orbits.

## Appendix B Numerical Analysis

If , the dynamical system (2.7) is only defined for , but more importantly it is not differentiable at . In this section we will use numerical techniques to analyze the integral curves in the neighborhood of the singular point in the case . The integration and plotting was done using Maple V release 3. From the qualitative analysis we find that the behaviour depends on the parameter . In the first of these two plots we choose , , and , so that (see \FrefFigure 7). In the second plot we choose , , and so that (see \FrefFigure 8). From the numerical plots we can conclude that the point has a saddle-point like nature, in agreement with preliminary remarks made in the text in section 3.2.

## References

## References

- [1] Coley A.A. and van den Hoogen R.J. 1994 J. Math. Phys. 35(8) 4117–4144.
- [2] Abolghasem G. and Coley A.A. 1994 Int. J. Theor. Phys. 33(3) 695–713.
- [3] Coley A. and Dunn K. 1992 J. Math. Phys. 33(5) 1772–1779.
- [4] Burd A. and Coley A. 1994 Class. Quantum Grav. 11 83–105.
- [5] MacCallum M. A. H. Cargèse Lectures in Physics, edited by E. Schatzman. Gordon and Breach New York 1973.
- [6] Eckart C. 1940 Phys. Rev. 58 919–924.
- [7] Belinskii V.A. and Khalatnikov I.M. 1976 Sov. Phys. JETP 42(2) 205–210.
- [8] Hiscock W.A. and Salmonson J. 1991 Phys. Rev. D 43(10) 3249–3258.
- [9] Israel Werner 1976 Ann. Phys. 100 310–331.
- [10] Israel W. and Stewart J. M. 1979 Proc. R. Soc. Lond. A. 365 43–52.
- [11] Israel W. and Stewart J. M. 1979 Ann. Phys. 118 341–372.
- [12] Belinskii V.A. , Nikomarov E.S. , and Khalatnikov I.M. 1979 Sov. Phys. JETP 50(2) 213–221.
- [13] Pavón D. , Bafaluy J. , and Jou D. 1991 Class. Quantum Grav. 8 347–360.
- [14] Chimento L.P. and Jakubi A.S. 1993 Class. Quantum Grav. 10 2047–2058.
- [15] Zakari M. and Jou D. 1993 Phys. Rev. D 48(4) 1597–1601.
- [16] Romano V. and Pavón D. 1994 Phys. Rev. D 50(4) 2572–2580.
- [17] Coley A.A. and van den Hoogen R.J. Self-similar Asymptotic Solutions of Einstein’s Equations. NATO ASI 332B, Plenum New York 1994.
- [18] Coley A.A. 1990 J. Math. Phys. 31(7) 1698–1703.
- [19] Bluman G.W. and Kumei S. Symmetries and Differential Equations. Springer-Verlag New York 1989.
- [20] Coley A.A. 1990 Gen. Rel. Grav. 22(1) 3–18.
- [21] Hawking S. W. and Ellis G.F.R. The large scale structure of space-time. Cambridge Univ. Press Cambridge 1973.
- [22] Maartens Roy 1994 preprint.
- [23] Pavón D. , Jou D. , and Casas-Vázquez J. 1982 Ann. Inst. H. Poincaré A 36 79–88.
- [24] Gariel J. and Denmat G. Le 1994 Phys. Rev. D 50(4) 2560–2566.
- [25] Schomliuk Dana 1993 Trans. Amer. Soc. 338(2) 799–841.

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and with . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and with . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models with bulk viscous pressure in the case and with . The arrows in the figure denote increasing -time () or decreasing -time ().

The phase portrait describes the qualitative behavior of the FRW models in the case and .

The phase portrait describes the qualitative behavior of the FRW models in the case and .