Expansion-Free Dissipative Fluid Spheres: Analytical Solutions

Abstract

We search for exact analytical solutions of spherically symmetric dissipative fluid distributions satisfying the vanishing expansion condition (vanishing expansion scalar $\mathsf{\Theta }$). To accomplish this, we shall impose additional restrictions allowing integration of the field equations. The solutions are analyzed, and possible applications to astrophysical scenarios as well as alternative approaches to obtaining new solutions are discussed.

1. Introduction

Many years ago, V. Skripkin published a paper describing the evolution of a spherically symmetric distribution of incompressible non-dissipative fluid following a central explosion [1] (see also [2]). Although not explicitly stated in the Skripkin paper, his model implies the vanishing of the expansion scalar. The interest in this kind of model stems from the fact that the expansion-free condition necessarily implies the appearance of a cavity around the center, suggesting that they might be relevant for the modeling, among other phenomena, of voids observed at cosmological scales.

A general study on shearing expansion-free spherical fluid evolution (including pressure anisotropy) was carried out in [3], where the unavoidable appearance of a cavity surrounding the center in expansion-free solutions was explained as consequence of the fact that the $\mathsf{\Theta }=0$ condition requires that the innermost shell of fluid should be away from the center, initiating from the formation of a cavity.

Thus, the evolution will proceed expansion-free if the decrease (increase) in the areal radius of the outer boundary surface of the fluid distribution is compensated by a decrease (increase) in the areal radius of the boundary of the cavity (see [3] for a detailed discussion on this issue).

General results regarding expansion-free fluids may be found in [4,5], and extensions to other theories of gravitation or charged fluids as well as different kind of symmetries have been discussed in [6,7,8,9,10,11] and the references therein.

The problem of (in)stability under the expansion-free condition was addressed for the first time in [12], and afterward, this problem was addressed within the context of modified theories of gravity in [13,14,15,16,17,18] and the references therein.

Aside from the solutions presented in [3], exact solutions describing expansion-free fluids may be found in [19,20,21,22,23,24,25].

Due to the interest attracted by expansion-free fluids, it is our purpose in this work to find exact analytical solutions to the Einstein equations for dissipative expansion-free fluids. These solutions will be found by imposing additional restrictions allowing the full integration of the field equations. Some of these restrictions are endowed with a distinct physical meaning, while others are just of a heuristic nature. Among the former, the vanishing complexity factor condition and the quasi-homologous evolution stand out. A discussion on all the presented models is brought about in the last section.

2. Relevant Physical and Geometric Variables, Field Equations and Junction Conditions

We consider spherically symmetric distributions of fluid bounded on the exterior by a spherical surface ${\Sigma }^{\left(e\right)}$. The fluid was assumed to be locally anisotropic (principal stresses unequal) and undergoing dissipation in the form of heat flow (to model dissipation in the diffusion approximation).

The motivation to include dissipation was provided by the well-known fact that dissipation due to the emission of massless particles (photons or neutrinos) seems to be the only plausible mechanism to carry away the bulk of the binding energy of the collapsing star, leading to a neutron star or a black hole. Furthermore, the diffusion approximation is, in general, very sensible, since the mean free path of the particles responsible for the propagation of energy in stellar interiors is, in general, very small compared with the typical length of the object. Thus, even though in many other circumstances the mean free path of particles transporting energy may be large enough to justify the free streaming approximation, there are many physically meaningful scenarios where diffusion approximation is justified.

The fluid under consideration was anisotropic in terms of pressure. This was quite justified since local isotropy, as has been shown in recent years, is a too stringent condition which may excessively constrain the modeling of self-gravitating objects. Furthermore, local anisotropy pressure may be caused by a large variety of physical phenomena, of the kind we expect being in compact objects [26]. Moreover, as has been recently shown [27], dissipation produces anisotropic pressure, and we do not know of any physical process able to erase the acquired anisotropy during the dynamic process.

As we mentioned in the introduction, the expansion-free models present an internal vacuum cavity. We shall denote with ${\Sigma }^{\left(i\right)}$ the boundary surface between the cavity (inside which we have a Minkowski spacetime) and the fluid.

By choosing co-moving coordinates inside ${\Sigma }^{\left(e\right)}$, the general interior metric can be written as follows:

$$d{s}^2=-A^{2}d{t}^2+B^{2}d{r}^2+R^2(d{\theta }^2+{sin}^2\theta d{\varphi }^2),$$

(1)

where A, B and R are functions of t and r and are assumed to be positive. We set the coordinates to $x^0=t$, $x^1=r$, $x^2=\theta$ and $x^3=\varphi$. Observe that A and B are dimensionless, whereas R has the same dimension as r.

The energy momentum tensor in the canonical form reads as follows:

$$T_{\alpha \beta }=\mu V_{\alpha }{V}_{\beta }+P{h}_{\alpha \beta }+{\Pi }_{\alpha \beta }+q_{\beta }{V}_{\alpha }+q_{\alpha }{V}_{\beta },$$

(2)

with

$$P=\frac{P_r+2P_{\perp }}{3},h_{\alpha \beta }=g_{\alpha \beta }+V_{\alpha }{V}_{\beta },$$

$${\Pi }_{\alpha \beta }=\Pi \left(K_{\alpha }{K}_{\beta }-\frac{1}{3}{h}_{\alpha \beta }\right),\Pi =P_r-P_{\perp },$$

where $\mu$ is the energy density, $P_r$ the radial pressure, $P_{\perp }$ is the tangential pressure, $q^{\alpha }$ is the heat flux, $V^{\alpha }$ is the four-velocity of the fluid, and $K^{\alpha }$ is a four-vector unit along the radial direction. Since we are considering co-moving observers, we have

$$\begin{array}{ccc}\hfill V^{\alpha }& =& A^{-1}{\delta }_0^{\alpha },q^{\alpha }=q{K}^{\alpha },K^{\alpha }=B^{-1}{\delta }_1^{\alpha }.\hfill \end{array}$$

(3)

These quantities satisfy

$$\begin{array}{c}\hfill V^{\alpha }{V}_{\alpha }=-1,V^{\alpha }{q}_{\alpha }=0,K^{\alpha }{K}_{\alpha }=1,K^{\alpha }{V}_{\alpha }=0.\end{array}$$

(4)

It is worth noticing that we do not explicitly add the bulk, shear viscosity or dissipation in the free streaming approximation because they can be trivially introduced by redefining the radial and tangential pressures, $\mu$ and q, respectively.

The acceleration $a_{\alpha }$, the expansion $\mathsf{\Theta }$ and the shear ${\sigma }_{\alpha \beta }$ of the fluid are given by

$$a_{\alpha }=V_{\alpha ;\beta }{V}^{\beta },\mathsf{\Theta }={V^{\alpha }}_{;\alpha },$$

(5)

$${\sigma }_{\alpha \beta }=V_{(\alpha ;\beta )}+a_{(\alpha }{V}_{\beta )}-\frac{1}{3}\mathsf{\Theta }{h}_{\alpha \beta }.$$

(6)

From Equation (5), we have for the four-acceleration and its scalar a

$$a_{\alpha }=a{K}_{\alpha },a=\frac{A^{\prime }}{AB},$$

(7)

and for the expansion

$$\mathsf{\Theta }=\frac{1}{A}\left(\frac{\dot{B}}{B}+2\frac{\dot{R}}{R}\right),$$

(8)

where the prime stands for r differentiation and the dot stands for differentiation with respect to t.

We obtain for the shear Equation (6) its nonzero components

$${\sigma }_{11}=\frac{2}{3}{B}^2\sigma ,{\sigma }_{22}=\frac{{\sigma }_{33}}{{sin}^2\theta }=-\frac{1}{3}{R}^2\sigma ,$$

(9)

and its scalar

$${\sigma }^{\alpha \beta }{\sigma }_{\alpha \beta }=\frac{2}{3}{\sigma }^2,$$

(10)

where

$$\sigma =\frac{1}{A}\left(\frac{\dot{B}}{B}-\frac{\dot{R}}{R}\right).$$

(11)

Einstein’s field equations for the interior spacetime (Equation (1))

$$G_{\alpha \beta }=8\pi T_{\alpha \beta },$$

(12)

are given in Appendix A.

Thus, in the most general case (locally anisotropic and dissipative), we have four field equations available (Equations (A2)–(A5)) for seven variables, namely A, B, R, $\mu$, $P_r$, $P_{{}_{\perp }}$ and q. Since we are going to consider expansion-free systems, we have the additional condition $\mathsf{\Theta }=0$. Therefore, in order to find specific solutions (to close the system of equations), we need to provide additional information, which could be given in the form of constitutive equations for q, equations of state for both pressures or any other type of constraint on the physical or metric variables. In this work, we shall impose (among others) the vanishing complexity factor condition and the quasi-homologous evolution condition.

Next, the mass function $m(t,r)$ introduced by Misner and Sharp [28] (see also [29]) reads as follows:

$$m=\frac{R^3}{2}{R_{23}}^{23}=\frac{R}{2}\left[{\left(\frac{\dot{R}}{A}\right)}^2-{\left(\frac{R^{\prime }}{B}\right)}^2+1\right].$$

(13)

To study the dynamical properties of the system, it is convenient to introduce the proper time derivative $D_T$, given by

$$D_T=\frac{1}{A}\frac{\partial }{\partial t},$$

(14)

and the proper radial derivative $D_R$, given by

$$D_R=\frac{1}{R^{\prime }}\frac{\partial }{\partial r},$$

(15)

where R defines the areal radius of a spherical surface inside ${\Sigma }^{\left(e\right)}$ (as measured with its area).

Using Equation (14), we can define the velocity U of the collapsing fluid as the variation of the areal radius with respect to the proper time. In other words, we have

$$U=D_{T}R.$$

(16)

Then, Equation (13) can be rewritten as

$$E\equiv \frac{R^{\prime }}{B}={\left(1+U^2-\frac{2m}{R}\right)}^{1/2}.$$

(17)

With Equation (15), we can express Equation (A6) as

$$4\pi q=E\left[\frac{1}{3}{D}_R(\mathsf{\Theta }-\sigma )-\frac{\sigma }{R}\right].$$

(18)

Using Equations (A2)–(A4) with Equations (14) and (15), we obtain from Equation (13)

$$\begin{array}{c}\hfill D_{T}m=-4\pi \left(P_{r}U+qE\right)R^2,\end{array}$$

(19)

and

$$\begin{array}{c}\hfill D_{R}m=4\pi \left(\mu +q\frac{U}{E}\right)R^2,\end{array}$$

(20)

which implies

$$m=4\pi {\int }_0^r\left(\mu +q\frac{U}{E}\right)R^2{R}^{\prime }d\tilde{r}.$$

(21)

This assumes a regular center to the distribution such that $m\left(0\right)=0$.

2.1. The Exterior Spacetime and Junction Conditions

Outside of ${\Sigma }^{\left(e\right)}$, we assume that we have the Vaidya spacetime (i.e., we assume all outgoing radiation is massless), described by

$$d{s}^2=-\left[1-\frac{2M\left(v\right)}{\mathbf{r}}\right]d{v}^2-2d\mathbf{r}dv+{\mathbf{r}}^2(d{\theta }^2+{sin}^2\theta d{\varphi }^2),$$

(22)

where $M\left(v\right)$ denotes the total mass and v is the retarded time.

The matching of the full nonadiabatic sphere to the Vaidya spacetime on a surface where $r=r_{{\Sigma }^{\left(e\right)}}=$ constant was discussed in [26]. From the continuity of the first and second differential forms (see [26] for details), it follows that

$$m(t,r)\stackrel{{\Sigma }^{\left(e\right)}}{=}M\left(v\right),$$

(23)

and

$$\begin{array}{c}\hfill 2\left(\frac{{\dot{R}}^{\prime }}{R}-\frac{\dot{B}}{B}\frac{R^{\prime }}{R}-\frac{\dot{R}}{R}\frac{A^{\prime }}{A}\right)\stackrel{{\Sigma }^{\left(e\right)}}{=}-\frac{B}{A}\left[2\frac{\ddot{R}}{R}-\left(2\frac{\dot{A}}{A}-\frac{\dot{R}}{R}\right)\frac{\dot{R}}{R}\right]+\frac{A}{B}\left[\left(2\frac{A^{\prime }}{A}+\frac{R^{\prime }}{R}\right)\frac{R^{\prime }}{R}-{\left(\frac{B}{R}\right)}^2\right],\end{array}$$

(24)

where $\stackrel{{\Sigma }^{\left(e\right)}}{=}$ means that both sides of the equation are evaluated on ${\Sigma }^{\left(e\right)}$.

By comparing, Equation (24) with Equations (A3) and (A4), one obtains

$$q\stackrel{{\Sigma }^{\left(e\right)}}{=}{P}_r.$$

(25)

Thus, the matching of Equations (1) and (22) on ${\Sigma }^{\left(e\right)}$ implies Equations (23) and (25).

As we mentioned in the introduction, the expansion-free models present an internal vacuum cavity whose boundary surface is denoted by ${\Sigma }^{\left(i\right)}$, and then the matching of the Minkowski spacetime within the cavity to the fluid distribution implies

$$m(t,r)\stackrel{{\Sigma }^{\left(i\right)}}{=}0,$$

(26)

$$q\stackrel{{\Sigma }^{\left(i\right)}}{=}{P}_r=0.$$

(27)

If any of the above matching conditions are not satisfied, then we have to assume there is a thin shell on the corresponding boundary surface.

2.2. The Weyl Tensor and the Complexity Factor

Some of the solutions exhibited in the next section were obtained from the condition of a vanishing complexity factor. This is a scalar function intended to measure the degree of complexity of a given fluid distribution [30,31], and it is related to the so-called structure scalars [32].

In the spherically symmetric case, the magnetic part of the Weyl tensor ($C_{\alpha \beta \mu }^{\rho }$) vanishes, and accordingly, the Weyl tensor is defined by its “electric” part $E_{\gamma \nu }$ alone:

$$E_{\alpha \beta }=C_{\alpha \mu \beta \nu }{V}^{\mu }{V}^{\nu },$$

(28)

whose nontrivial components are

$$\begin{array}{ccc}\hfill E_{11}& =& \frac{2}{3}{B}^2\mathcal{E},\hfill \\ \hfill E_{22}& =& -\frac{1}{3}{R}^2\mathcal{E},\hfill \\ \hfill E_{33}& =& E_{22}{sin}^2\theta ,\hfill \end{array}$$

(29)

where

$$\begin{array}{ccc}\hfill \mathcal{E}& =& \frac{1}{2A^2}\left[\frac{\ddot{R}}{R}-\frac{\ddot{B}}{B}-\left(\frac{\dot{R}}{R}-\frac{\dot{B}}{B}\right)\left(\frac{\dot{A}}{A}+\frac{\dot{R}}{R}\right)\right]\hfill \\ & +& \frac{1}{2B^2}\left[\frac{A^{″}}{A}-\frac{R^{″}}{R}+\left(\frac{B^{\prime }}{B}+\frac{R^{\prime }}{R}\right)\left(\frac{R^{\prime }}{R}-\frac{A^{\prime }}{A}\right)\right]-\frac{1}{2R^2}.\hfill \end{array}$$

(30)

Observe that the electric part of the Weyl tensor may be written as

$$E_{\alpha \beta }=\mathcal{E}\left(K_{\alpha }{K}_{\beta }-\frac{1}{3}{h}_{\alpha \beta }\right).$$

(31)

As it is shown in [30,31], the complexity factor is identified with the scalar function $Y_{TF}$, which defines the trace-free part of the electric Riemann tensor (see [32] for details).

Thus, let us define tensor $Y_{\alpha \beta }$ as follows:

$$Y_{\alpha \beta }=R_{\alpha \gamma \beta \delta }{V}^{\gamma }{V}^{\delta },$$

(32)

which may be expressed in terms of two scalar functions $Y_T,Y_{TF}$ as

$$\begin{array}{c}\hfill Y_{\alpha \beta }=\frac{1}{3}{Y}_T{h}_{\alpha \beta }+Y_{TF}\left(K_{\alpha }{K}_{\beta }-\frac{1}{3}{h}_{\alpha \beta }\right).\end{array}$$

(33)

Then, after lengthy but simple calculations, using field equations, we obtain

$$\begin{array}{c}\hfill Y_T=4\pi (\mu +3P_r-2\Pi ),Y_{TF}=\mathcal{E}-4\pi \Pi .\end{array}$$

(34)

Next, by using Equations (A2), (A4) and (A5) with Equations (13) and (30), we obtain

$$\frac{3m}{R^3}=4\pi \left(\mu -\Pi \right)-\mathcal{E},$$

(35)

and

$$Y_{TF}=-8\pi \Pi +\frac{4\pi }{R^3}{\int }_0^r{R}^3\left(D_R\mu -3q\frac{U}{RE}\right)R^{\prime }d\tilde{r}.$$

(36)

It is worth noticing that due to a different signature, the sign of $Y_{TF}$ in the above equation differs from the sign of $Y_{TF}$ used in [30] for the static case.

For reasons explained in detail in [30], the scalar $Y_{TF}$ is the variable identified with the complexity factor. As follows from the equations above, it may be expressed through the Weyl tensor and the anisotropy of the pressure or in terms of the density inhomogeneity, the dissipative variables and the anisotropy of pressure.

In terms of the metric functions, the scalar $Y_{TF}$ reads as follows:

$$\begin{array}{c}\hfill Y_{TF}=\frac{1}{A^2}\left[\frac{\ddot{R}}{R}-\frac{\ddot{B}}{B}+\frac{\dot{A}}{A}\left(\frac{\dot{B}}{B}-\frac{\dot{R}}{R}\right)\right]+\frac{1}{B^2}\left[\frac{A^{″}}{A}-\frac{A^{\prime }}{A}\left(\frac{B^{\prime }}{B}+\frac{R^{\prime }}{R}\right)\right].\end{array}$$

(37)

3. The Transport Equation

In the dissipative case, we shall need a transport equation in order to find the temperature distribution and its evolution. Assuming a causal dissipative theory (e.g., the Israel–Stewart theory [33,34,35]), the transport equation for the heat flux reads as follows:

$$\begin{array}{ccc}\hfill \tau h^{\alpha \beta }{V}^{\gamma }{q}_{\beta ;\gamma }& +& q^{\alpha }=-\kappa h^{\alpha \beta }\left(T_{,\beta }+T{a}_{\beta }\right)\hfill \\ & -& \frac{1}{2}\kappa T^2{\left(\frac{\tau V^{\beta }}{\kappa T^2}\right)}_{;\beta }{q}^{\alpha },\hfill \end{array}$$

(38)

where $\kappa$ denotes the thermal conductivity and T and $\tau$ denote the temperature and relaxation time, respectively.

In the spherically symmetric case under consideration, the transport equation has only one independent component which may be obtained from Equation (38) by contracting with the unit’s space-like vector $K^{\alpha }$, which reads as follows:

$$\tau V^{\alpha }{q}_{,\alpha }+q=-\kappa \left(K^{\alpha }{T}_{,\alpha }+Ta\right)-\frac{1}{2}\kappa T^2{\left(\frac{\tau V^{\alpha }}{\kappa T^2}\right)}_{;\alpha }q.$$

(39)

Sometimes, the last term in Equation (38) may be neglected [36], producing the so-called truncated transport equation, which reads as follows:

$$\tau V^{\alpha }{q}_{,\alpha }+q=-\kappa \left(K^{\alpha }{T}_{,\alpha }+Ta\right).$$

(40)

4. The Homologous and Quasi-Homologous Conditions

In order to specify some of our models, we shall impose the condition of a vanishing complexity factor. However, for time-dependent systems, we also need to elucidate what the simplest pattern of evolution is for the system.

In [31], the concept of homologous evolution was introduced, in analogy with the same concept in classical astrophysics, to represent the simplest mode of evolution of the fluid distribution.

Thus, the field equation in Equation (A3), written as

$$D_R\left(\frac{U}{R}\right)=\frac{4\pi }{E}q+\frac{\sigma }{R},$$

(41)

can be easily integrated to obtain

$$U=\tilde{a}\left(t\right)R+R{\int }_0^r\left(\frac{4\pi }{E}q+\frac{\sigma }{R}\right)R^{\prime }d\tilde{r},$$

(42)

where $\tilde{a}$ is an integration function or

$$U=\frac{U_{{\Sigma }^{\left(e\right)}}}{R_{{\Sigma }^{\left(e\right)}}}R-R{\int }_r^{r_{{\Sigma }^{\left(e\right)}}}\left(\frac{4\pi }{E}q+\frac{\sigma }{R}\right)R^{\prime }d\tilde{r}.$$

(43)

If the integral in the above equations vanishes, we have from Equation (42) or (43) that

$$U=\tilde{a}\left(t\right)R.$$

(44)

This relationship is characteristic of the homologous evolution in Newtonian hydrodynamics [37,38,39]. In our case, this may occur if the fluid is shear-free and non-dissipative or if the two terms in the integral cancel each other out.

In [31], the term “homologous evolution” was used to characterize relativistic systems satisfying, aside from Equation (44), the condition

$$\frac{R_I}{R_{II}}=\mathrm{constant},$$

(45)

where $R_I$ and $R_{II}$ denote the areal radii of two concentric shells ($I,II$) described by $r=r_I=\mathrm{constant}$ and $r=r_{II}=\mathrm{constant}$, respectively.

The important point we want to stress here is that Equation (44) does not imply Equation (45). Indeed, Equation (44) implies that for the two shells of fluids I and $II$, we have

$$\frac{U_I}{U_{II}}=\frac{A_{II}{\dot{R}}_I}{A_1{\dot{R}}_{II}}=\frac{R_I}{R_{II}},$$

(46)

which implies Equation (45) only if $A=A\left(t\right)$, which through a simple coordinate transformation becomes $A=\mathrm{constant}$. Thus, while in the non-relativistic regime, Equation (45) always follows from the condition that the radial velocity is proportional to the radial distance, in the relativistic regime, the condition (44) implies Equation (45) only if the fluid is geodesic.

In [40], the homologous condition was relaxed, leading to what was defined as quasi-homologous evolution restricted only by Equation (44), implying

$$\frac{4\pi }{R^{\prime }}Bq+\frac{\sigma }{R}=0.$$

(47)

5. Shearing Expansion-Free Motion

If the fluid evolves with the vanishing expansion scalar ($\mathsf{\Theta }=0$), then from Equation (8), we have

$$\frac{\dot{B}}{B}=-2\frac{\dot{R}}{R},$$

(48)

or, integrating

$$B=\frac{g\left(r\right)}{R^2},$$

(49)

where $g\left(r\right)$ is an arbitrary function of r.

By substituting Equation (48) into (A3), we obtain

$$\frac{{\dot{R}}^{\prime }}{R}+2\frac{\dot{R}}{R}\frac{R^{\prime }}{R}-\frac{\dot{R}}{R}\frac{A^{\prime }}{A}=4\pi qAB,$$

(50)

which can be integrated for $\dot{R}\ne 0$, producing

$$A=\frac{R^2\dot{R}}{{\tau }_1\left(t\right)}exp\left[-4\pi \int qAB\frac{R}{\dot{R}}dr\right],$$

(51)

where ${\tau }_1\left(t\right)$ is an arbitrary function of t. With Equations (49) and (51), the line element in Equation (1) becomes

$$\begin{array}{c}\hfill d{s}^2=-{\left\{\frac{R^2\dot{R}}{\alpha }exp\left[-4\pi \int qAB\frac{R}{\dot{R}}dr\right]\right\}}^{2}d{t}^2\\ \hfill +\frac{{\alpha }^2}{R^4}d{r}^2+R^2(d{\theta }^2+{sin}^2\theta d{\varphi }^2),\end{array}$$

(52)

which is the general metric for a shearing expansion-free anisotropic dissipative fluid, where without loss of generality (by reparametrizing r and t), we have $g={\tau }_1=\alpha$, where $\alpha$ is a unit constant with dimensions $\left[r^2\right]$.

We shall now proceed to find exact analytical solutions satisfying the expansion-free condition.

6. Solutions

As mentioned before, the expansion-free condition alone is not enough to integrate the field equations. We need to resort to additional restrictions in order to close the full system of equations. Here, two different families of solutions will be obtained. On the one hand, we shall consider non-geodesic fluids satisfying the vanishing complexity factor condition, complemented by the quasi-homologous condition or with some simple assumptions on the metric variables. On the other hand, we shall consider geodesic fluids satisfying the vanishing complexity factor condition or the quasi-homologous condition.

6.1. Non-Geodesic, $Y_{TF}=0$, Quasi-Homologous Evolution and $\mathsf{\Theta }=0$ Solutions

We shall start by assuming the quasi-homologous condition and the vanishing complexity factor condition for a non-geodesic fluid. As mentioned before, by imposing $\mathsf{\Theta }=0$, we obtain

$$\begin{array}{ccc}\hfill \frac{\dot{B}}{B}+\frac{2\dot{R}}{R}& =& 0\Rightarrow B=\frac{\alpha }{R^2}.\hfill \end{array}$$

(53)

Let us now impose the quasi-homologous condition from Equation (44):

$$\dot{R}=\tilde{a}\left(t\right)AR.$$

(54)

Then, from $\mathsf{\Theta }=0$ and $U=\tilde{a}R$, we have

$$\begin{array}{ccc}\hfill B& =& \frac{\alpha }{R^2},\hfill \end{array}$$

(55)

$$\begin{array}{ccc}\hfill A& =& \frac{\dot{R}}{\tilde{a}R}.\hfill \end{array}$$

(56)

With the above conditions, the physical variables read as follows:

$$\begin{array}{ccc}\hfill 8\pi \mu & =& -3{\tilde{a}}^2-\frac{R^4}{{\alpha }^2}\left[\frac{2R^{″}}{R}+5{\left(\frac{R^{\prime }}{R}\right)}^2-\frac{{\alpha }^2}{R^6}\right],\hfill \end{array}$$

(57)

$$\begin{array}{ccc}\hfill 4\pi q& =& \frac{3\tilde{a}}{\alpha }R{R}^{\prime },\hfill \end{array}$$

(58)

$$\begin{array}{ccc}\hfill 8\pi P_r& =& -3{\tilde{a}}^2+\frac{R^4}{{\alpha }^2}\left[\frac{2{\dot{R}}^{\prime }}{\dot{R}}\frac{R^{\prime }}{R}-{\left(\frac{R^{\prime }}{R}\right)}^2-\frac{{\alpha }^2}{R^6}\right],\hfill \end{array}$$

(59)

$$\begin{array}{ccc}\hfill 8\pi P_{\perp }& =& 3{\tilde{a}}^2+\frac{R^4}{{\alpha }^2}\left[\frac{{\dot{R}}^{″}}{\dot{R}}+\frac{{\dot{R}}^{\prime }}{\dot{R}}\frac{R^{\prime }}{R}+{\left(\frac{R^{\prime }}{R}\right)}^2\right],\hfill \end{array}$$

(60)

where we have chosen $\tilde{a}=constant$.

Next, the condition $Y_{TF}=0$ produces

$$-3{\tilde{a}}^2+\frac{R^4}{{\alpha }^2}\left[\frac{{\dot{R}}^{″}}{\dot{R}}-\frac{{\dot{R}}^{\prime }}{\dot{R}}\frac{R^{\prime }}{R}-\frac{R^{″}}{R}+{\left(\frac{R^{\prime }}{R}\right)}^2\right]=0.$$

(61)

In order to find a solution to the above equation, we shall proceed as follows. We shall write R as

$$R=F({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)\equiv F\left(z\right),$$

(62)

where F is an arbitrary function of its argument with dimensions $\left[r\right]$, ${\delta }_1$ and ${\delta }_2$ are two arbitrary constants with dimensions $[1/r]$ and ${\delta }_3$ is a dimensionless constant.

By feeding Equation (62) back into Equation (61), we obtain

$$-3{\tilde{a}}^2+\frac{{\delta }_1^2{F}^4}{{\alpha }^2}\left[\frac{\frac{{\partial }^{3}F}{\partial z^3}}{\frac{\partial F}{\partial z}}-2\frac{\frac{{\partial }^{2}F}{\partial z^2}}{F}+{\left(\frac{\frac{\partial F}{\partial z}}{F}\right)}^2\right]=0.$$

(63)

By introducing the variable $\omega \left(F\right)\equiv \frac{\partial F}{\partial z}$, we have

$$\begin{array}{ccc}\hfill \frac{\partial F}{\partial z}& =& \omega ,\hfill \end{array}$$

(64)

$$\begin{array}{ccc}\hfill \frac{{\partial }^{2}F}{\partial z^2}& =& {\omega }_F\omega ,\hfill \end{array}$$

(65)

$$\begin{array}{ccc}\hfill \frac{{\partial }^{3}F}{\partial z^3}& =& {\omega }^2{\omega }_{FF}+\omega {\omega }_F^2,\hfill \end{array}$$

(66)

with the help of which Equation (63) becomes

$$-3{\tilde{a}}^2+\frac{{\delta }_1^2{F}^4}{{\alpha }^2}\left[\omega {\omega }_{FF}+{\omega }_F^2-\frac{2\omega {\omega }_F}{F}+\frac{{\omega }^2}{F^2}\right]=0,$$

(67)

where the subscript F denotes the derivative with respect to F and whose solution reads

$$\omega =\frac{k}{F},k=\frac{\tilde{a}\alpha }{\sqrt{2}{\delta }_1}.$$

(68)

Using the above results, we may write the following for R:

$$R=\sqrt{2k({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)},$$

(69)

and the physical variables read as follows:

$$\begin{array}{ccc}\hfill 8\pi \mu & =& -\frac{9}{2}{\tilde{a}}^2+\frac{1}{2k({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)},\hfill \end{array}$$

(70)

$$\begin{array}{ccc}\hfill 4\pi q& =& \frac{3{\tilde{a}}^2}{\sqrt{2}},\hfill \end{array}$$

(71)

$$\begin{array}{ccc}\hfill 8\pi P_r& =& -\frac{9}{2}{\tilde{a}}^2-\frac{1}{2k({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)},\hfill \end{array}$$

(72)

$$\begin{array}{ccc}\hfill 8\pi P_{\perp }& =& \frac{9}{2}{\tilde{a}}^2.\hfill \end{array}$$

(73)

The corresponding expressions for the mass function and the shear are

$$m=\frac{{\tilde{a}}^2}{4}{\left[\sqrt{2k({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)}\right]}^3+\frac{1}{2}\sqrt{2k({\delta }_{1}r+{\delta }_{2}t+{\delta }_3)},$$

(74)

$$\sigma =-3\tilde{a},$$

(75)

whereas for the temperature, using Equation (40), we find

$$T=\frac{2\tilde{a}{T}_0\left(t\right)}{{\delta }_2}\left({\delta }_{1}r+{\delta }_{2}t+{\delta }_3\right)+\frac{3\tilde{a}}{8\pi k},$$

(76)

where $T_0\left(t\right)$ is a function of integration which, in principle, may be obtained from the boundary conditions on either boundary surface.

6.2. Non-Geodesic, $Y_{TF}=0$, $\mathsf{\Theta }=0$, $A=\gamma B$ and $\gamma =Constant$ Solutions

The next model will also be obtained by imposing $Y_{TF}=0$ and $\mathsf{\Theta }=0$, but instead of assuming the quasi-homologous evolution as in the preceding model, we shall assume that A and B are proportional to each other ($A=\gamma B$, with $\gamma =constant$).

Thus, from the three conditions above, we obtain

$$\frac{3}{{\gamma }^2}\frac{\ddot{R}}{R}-\frac{2R^{″}}{R}+{\left(\frac{2R^{\prime }}{R}\right)}^2=0.$$

(77)

In order to find a solution to Equation (77), we assume for R the form

$$R=F(s_{1}r+s_{2}t+s_3)\equiv F\left(z\right),$$

(78)

where $s_1$ and $s_2$ are constants with dimensions $[1/r]$ and $s_3$ is a dimensionless constant.

By replacing Equation (78) in (77), we obtain

$$\left(\frac{3s_2^2}{{\gamma }^2}-2s_1^2\right)\frac{F_{zz}}{F}+\frac{4s_1^2{F}_z^2}{F^2}=0.$$

(79)

Next, by introducing the variable y, defined by

$$y=\frac{F_z}{F},$$

(80)

we may write Equation (79) as

$$y_z+{\beta }_0{y}^2=0,with{\beta }_0=\frac{3s_2^2+2s_1^2{\gamma }^2}{3s_2^2-2s_1^2{\gamma }^2}.$$

(81)

The above equation may be easily integrated, producing

$$y=\frac{1}{{\beta }_{0}z+{\beta }_1},$$

(82)

where ${\beta }_1$ is a constant of integration.

By feeding Equation (82) back into Equation (80) and integrating once again, we obtain

$$F=R={\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{1}{{\beta }_0}},$$

(83)

where ${\beta }_2$ is a new constant of integration with dimensions $[1/r^{{\beta }_0}]$.

The physical variables for this model read as follows:

$$\begin{array}{ccc}\hfill 8\pi \mu & =& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{4}{{\beta }_0}}\left[s_1^2(2{\beta }_0-7)-\frac{3s_2^2}{{\gamma }^2}\right]}{{\alpha }^2{({\beta }_{0}z+{\beta }_1)}^2}\hfill \\ & +& {\left(\frac{{\beta }_2}{{\beta }_{0}z+{\beta }_1}\right)}^{\frac{2}{{\beta }_0}},\hfill \end{array}$$

(84)

$$\begin{array}{ccc}\hfill 4\pi q& =& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{4}{{\beta }_0}}{s}_1{s}_2(5-{\beta }_0)}{\gamma {\alpha }^2{({\beta }_{0}z+{\beta }_1)}^2},\hfill \end{array}$$

(85)

$$\begin{array}{ccc}\hfill 8\pi P_r& =& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{4}{{\beta }_0}}\left[\frac{s_2^2(2{\beta }_0-7)}{{\gamma }^2}-3s_1^2\right]}{{\alpha }^2{({\beta }_{0}z+{\beta }_1)}^2}-{\left(\frac{{\beta }_2}{{\beta }_{0}z+{\beta }_1}\right)}^{\frac{2}{{\beta }_0}},\hfill \\ \hfill \end{array}$$

(86)

$$\begin{array}{ccc}\hfill 8\pi P_{\perp }& =& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{4}{{\beta }_0}}(1+{\beta }_0)\left(s_1^2-\frac{s_2^2}{{\gamma }^2}\right)}{{\alpha }^2{({\beta }_{0}z+{\beta }_1)}^2},\hfill \end{array}$$

(87)

$$m=\frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{7}{{\beta }_0}}\left(\frac{s_2^2}{{\gamma }^2}-s_1^2\right)}{2{\alpha }^2{({\beta }_{0}z+{\beta }_1)}^2}+\frac{1}{2}{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{1}{{\beta }_0}},$$

(88)

$$\sigma =-\frac{3s_2}{\gamma \alpha ({\beta }_{0}z+{\beta }_1)}{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{2}{{\beta }_0}},$$

(89)

$$\begin{array}{ccc}\hfill T& =& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{2}{{\beta }_0}}}{\gamma \alpha }\left[T_0\left(t\right)-\frac{\tau s_2^2(5-{\beta }_0)(2-{\beta }_0){\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{2}{{\beta }_0}}}{4\pi \kappa \gamma \alpha (1-{\beta }_0){({\beta }_{0}z+{\beta }_1)}^2}\right]\hfill \\ & +& \frac{{\left(\frac{{\beta }_{0}z+{\beta }_1}{{\beta }_2}\right)}^{\frac{2}{{\beta }_0}}{s}_2(5-{\beta }_0)}{4\pi \gamma \kappa {\beta }_0\alpha ({\beta }_{0}z+{\beta }_1)}.\hfill \end{array}$$

(90)

As in the preceding models, the temperature is obtained using the truncated transport equation.

6.3. Non-Geodesic, $Y_{TF}=0$, $\mathsf{\Theta }=0$, $A=A\left(r\right)$ and $R=R_1\left(t\right)R_2\left(r\right)$ Solutions

We shall now find a solution satisfying the conditions $\mathsf{\Theta }=0$ and $Y_{TF}=0$, as well as $A=A\left(r\right)$ and the condition that R is a separable function (i.e., $R=R_1\left(t\right)R_2\left(r\right)$).

From $\mathsf{\Theta }=0$ and $Y_{TF}=0$, we may write

$$\frac{3}{A^2}\left(\frac{\ddot{R}}{R}-\frac{2{\dot{R}}^2}{R^2}-\frac{\dot{A}}{A}\frac{\dot{R}}{R}\right)+\frac{1}{B^2}\left(\frac{A^{″}}{A}+\frac{A^{\prime }}{A}\frac{R^{\prime }}{R}\right)=0.$$

(91)

By imposing the further conditions $A=A\left(r\right)$ and $R=R_1\left(t\right)R_2\left(r\right)$, we see that a simple solution to Equation (91) reads as follows:

$$\begin{array}{ccc}\hfill R_1\left(t\right)& =& \frac{{\nu }_0}{t+{\nu }_1},\hfill \end{array}$$

(92)

$$\begin{array}{ccc}\hfill R_2\left(r\right)& =& {\nu }_2{A}^{{\nu }_3-1},\hfill \end{array}$$

(93)

$$\begin{array}{ccc}\hfill A& =& {\nu }_4{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{1}{{\nu }_3}},\hfill \end{array}$$

(94)

where ${\nu }_0$ and ${\nu }_3$ are dimensionless constants and ${\nu }_1,{\nu }_2,{\nu }_4$ and ${\nu }_5$ are constants with dimensions $\left[r\right],\left[r^2\right],[1/\left(r^{1/{\nu }_3}\right)]$ and $\left[r\right]$, respectively.

The physical variables for this model read as follows:

$$\begin{array}{ccc}\hfill 8\pi \mu & =& \frac{{(t+{\nu }_1)}^2{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{-2{\nu }_3+2}{{\nu }_3}}}{{\nu }_0^2{\nu }_2^2{\nu }_4^{2{\nu }_3-2}}\hfill \\ & -& \frac{{\nu }_0^4{\nu }_4^{4{\nu }_3-4}{\nu }_2^4(5{\nu }_3^2-12{\nu }_3+7){\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{2({\nu }_3-2)}{{\nu }_3}}}{{\alpha }^2{(t+{\nu }_1)}^4}\hfill \\ & -& \frac{3{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{-2}{{\nu }_3}}}{{\nu }_4^2{(t+{\nu }_1)}^2},\hfill \end{array}$$

(95)

$$\begin{array}{ccc}\hfill 4\pi q& =& \frac{{\nu }_2^2{\nu }_0^2(4-3{\nu }_3){\nu }_4^{2{\nu }_3-3}{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{{\nu }_3-3}{{\nu }_3}}}{\alpha {(t+{\nu }_1)}^3},\hfill \end{array}$$

(96)

$$\begin{array}{ccc}\hfill 8\pi P_r& =& -\frac{5{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{-2}{{\nu }_3}}}{{\nu }_4^2{(t+{\nu }_1)}^2}-\frac{{(t+{\nu }_1)}^2{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{-2{\nu }_3+2}{{\nu }_3}}}{{\nu }_0^2{\nu }_2^2{\nu }_4^{2{\nu }_3-2}}\hfill \\ & +& \frac{{\nu }_0^4{\nu }_4^{4{\nu }_3-4}{\nu }_2^4({\nu }_3^2-1){\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{2({\nu }_3-2)}{{\nu }_3}}}{{\alpha }^2{(t+{\nu }_1)}^4},\hfill \end{array}$$

(97)

$$\begin{array}{ccc}\hfill 8\pi P_{\perp }& =& \frac{{\nu }_0^4{\nu }_4^{4{\nu }_3-4}{\nu }_2^4(2{\nu }_3^2-3{\nu }_3+1){\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{2({\nu }_3-2)}{{\nu }_3}}}{{\alpha }^2{(t+{\nu }_1)}^4}\hfill \\ & -& \frac{2{\left({\nu }_{3}r+{\nu }_5\right)}^{\frac{-2}{{\nu }_3}}}{{\nu }_4^2{(t+{\nu }_1)}^2},\hfill \end{array}$$

(98)

$$\sigma =\frac{3}{{\nu }_4(t+{\nu }_1){({\nu }_{3}r+{\nu }_5)}^{\frac{1}{{\nu }_3}}},$$

(99)

$$\begin{array}{ccc}\hfill m& =& \frac{{\nu }_0^3{\nu }_2^3{\nu }_4^{3{\nu }_3-5}{({\nu }_{3}r+{\nu }_5)}^{\frac{3{\nu }_3-5}{{\nu }_3}}}{2{(t+{\nu }_1)}^5}+\frac{{\nu }_0{\nu }_2{\nu }_4^{{\nu }_3-1}{({\nu }_{3}r+{\nu }_5)}^{\frac{{\nu }_3-1}{{\nu }_3}}}{2(t+{\nu }_1)}\hfill \\ & -& \frac{{\nu }_0^7{\nu }_2^7{\nu }_4^{7({\nu }_3-1)}{({\nu }_3-1)}^2{({\nu }_{3}r+{\nu }_5)}^{\frac{5{\nu }_3-7}{{\nu }_3}}}{2{\alpha }^2{(t+{\nu }_1)}^7},\hfill \end{array}$$

(100)

$$\begin{array}{ccc}\hfill T& =& \frac{{({\nu }_{3}r+{\nu }_5)}^{\frac{-1}{{\nu }_3}}}{{\nu }_4}\left[T_0\left(t\right)-\frac{3\tau (4-3{\nu }_3)}{4\pi \kappa {\nu }_4{(t+{\nu }_1)}^2{({\nu }_{3}r+{\nu }_5)}^{\frac{1}{{\nu }_3}}}\right.\hfill \\ & -& \left.\frac{(4-3{\nu }_3)}{4\pi \kappa {\nu }_3(t+{\nu }_1)}ln({\nu }_{3}r+{\nu }_5)\right].\hfill \end{array}$$

(101)

6.4. Geodesic Models

We shall now consider geodesic fluids, for which we have

$$A(t,r)=1.$$

(102)

That aside, from the expansion-free condition, we have

$$B(t,r)=\frac{\alpha }{R^2}.$$

(103)

From the above, it follows that the general expressions for the physical variables in this case read as follows:

$$8\pi \mu =-\frac{3{\dot{R}}^2}{R^2}-\frac{R^4}{{\alpha }^2}\left[2\frac{R^{″}}{R}+5{\left(\frac{R^{\prime }}{R}\right)}^2\right]+\frac{1}{R^2},$$

(104)

$$4\pi q=\frac{R{\dot{R}}^{\prime }}{\alpha }+\frac{2R^{\prime }\dot{R}}{\alpha },$$

(105)

$$8\pi P_r=-\frac{2\ddot{R}}{R}-\frac{{\dot{R}}^2}{R^2}+{\left(\frac{R{R}^{\prime }}{\alpha }\right)}^2-\frac{1}{R^2},$$

(106)

$$8\pi P_{\perp }=\frac{\ddot{R}}{R}-\frac{4{\dot{R}}^2}{R^2}+\frac{R^4}{{\alpha }^2}\left[\frac{R^{″}}{R}+2{\left(\frac{R^{\prime }}{R}\right)}^2\right],$$

(107)

producing

$$2\pi (\mu +P_r+2P_{\perp })=-3{\left(\frac{\dot{R}}{R}\right)}^2.$$

(108)

The first model will be obtained from the vanishing complexity factor condition. Thus, from the condition $Y_{TF}=0$, we obtain

$$\frac{\ddot{R}}{R}-\frac{2{\dot{R}}^2}{R^2}=0\Rightarrow R=\frac{1}{b_1\left(r\right)t+b_2\left(r\right)}\equiv \frac{1}{b_1\left(r\right)\left[t+\frac{b_2\left(r\right)}{b_1\left(r\right)}\right]},$$

(109)

where $b_1$ and $b_2$ are two arbitrary functions of their argument with dimensions $[1/r^2]$ and $[1/r]$, respectively.

By feeding Equation (109) back into Equation (103), we see that by reparameterizing r, we may choose without loss of generality $b_1=1/\alpha$. Thus, our metric variables become

$$R=\frac{\alpha }{\left[t+\alpha b_2\left(r\right)\right]},B=\frac{{\left[t+\alpha b_2\left(r\right)\right]}^2}{\alpha }.$$

(110)

For this metric, the physical variables, the mass function and the shear read as follows:

$$\begin{array}{ccc}\hfill 8\pi \mu & =& -\frac{3}{{\left(t+\alpha b_2\right)}^2}+\frac{2{\alpha }^3{b}_2^{″}}{{\left(t+\alpha b_2\right)}^5}-\frac{9{\alpha }^4{\left(b_2^{\prime }\right)}^2}{{\left(t+\alpha b_2\right)}^6}\hfill \\ & +& \frac{{\left(t+\alpha b_2\right)}^2}{{\alpha }^2},\hfill \end{array}$$

(111)

$$\begin{array}{c}\hfill 4\pi q=\frac{4{\alpha }^2{b}_2^{\prime }}{{\left(t+\alpha b_2\right)}^4},\end{array}$$

(112)

$$8\pi P_r=-\frac{5}{{\left(t+\alpha b_2\right)}^2}+\frac{{\alpha }^4{\left(b_2^{\prime }\right)}^2}{{\left(t+\alpha b_2\right)}^6}-\frac{{\left(t+\alpha b_2\right)}^2}{{\alpha }^2},$$

(113)

$$8\pi P_{\perp }=-\frac{2}{{\left(t+\alpha b_2\right)}^2}+\frac{4{\alpha }^4{\left(b_2^{\prime }\right)}^2}{{\left(t+\alpha b_2\right)}^6}-\frac{{\alpha }^3{b}_2^{″}}{{\left(t+\alpha b_2\right)}^5},$$

(114)

$$m=\frac{{\alpha }^3}{2{(t+\alpha b_2)}^5}-\frac{{\alpha }^7{\left(b_2^{\prime }\right)}^2}{2{(t+\alpha b_2)}^9}+\frac{\alpha }{2(t+\alpha b_2)},$$

(115)

$$\sigma =\frac{3}{t+\alpha b_2},$$

(116)

$$T=T_0\left(t\right)+\frac{1}{\pi \kappa (t+\alpha b_2)}\left[1-\frac{2\tau }{(t+\alpha b_2)}\right].$$

(117)

A second geodesic model will be obtained from the quasi-homologous condition (Equation (44)) which, as discussed before, implies Equation (45) in the geodesic case, implying in turn that R is a separable function. For this case, we obtain

$$R=\frac{g\left(r\right)}{t},B=\frac{t^2}{\alpha },$$

(118)

where $g\left(r\right)$ is an arbitrary function of r. Since R has a dimension of $\left[r\right]$ (as t), then the dimension (units) of g should be $\left[r^2\right]$.

The physical variables for this model read as follows:

$$\begin{array}{c}\hfill 8\pi \mu =-\frac{3}{t^2}-\frac{{\alpha }^2}{t^4}\left[\frac{2g^{″}}{g}+{\left(\frac{g^{\prime }}{g}\right)}^2\right]+\frac{t^2}{g^2},\end{array}$$

(119)

$$\begin{array}{c}\hfill 4\pi q=-\frac{3\alpha g^{\prime }}{t^{3}g},\end{array}$$

(120)

$$\begin{array}{c}\hfill 8\pi P_r=-\frac{5}{t^2}+\frac{{\alpha }^2}{t^4}{\left(\frac{g^{\prime }}{g}\right)}^2-\frac{t^2}{g^2},\end{array}$$

(121)

$$\begin{array}{c}\hfill 8\pi P_{\perp }=-\frac{2}{t^2}+\frac{{\alpha }^2}{t^4}\frac{g^{″}}{g},\end{array}$$

(122)

$$m=\frac{g}{2t}\left[\frac{g^2}{t^4}-{\left(\frac{\alpha g^{\prime }}{t^3}\right)}^2+1\right],$$

(123)

$$\sigma =\frac{3}{t},$$

(124)

$$T=T_0\left(t\right)+\frac{3}{4\pi \kappa t}\left(1-\frac{3\tau }{t}\right)lng.$$

(125)

7. Discusion

The main lesson we can extract from this work is that the expansion-free condition allows for a wide range of models for the evolution of spherically symmetric self-gravitating systems, including dissipative fluids with anisotropic pressure.

As mentioned in the Introduction, one of the most interesting features of expansion-free models is the appearance of a vacuum cavity within the fluid distribution. Whether or not such models may be used to describe the formation of voids observed at cosmological scales (see [41,42] and the references therein) is still an open question. We skipped over this issue in a hope of a resolution a posteriori.

Let us now analyze in some detail the obtained solutions.

The first model satisfies the quasi-homologous and vanishing complexity factor conditions and is described by Equations (69)–(76). By choosing $k,{\delta }_1,{\delta }_2$ and ${\delta }_3>0$, we ensure the positivity of the expression within the square root in Equation (69), which implies that because of Equation (68) that $\tilde{a}>0$, and therefore all fluid elements are moving outward. In the limit $t\to \infty$, the areal radii of all fluid elements tend toward infinity, and the fluid distribution becomes a shell. In addition, with the above choice, we ensure that $R^{\prime }>0$, thereby avoiding the appearance of shells crossing singularities. In this same limit, $8\pi \mu =8\pi P_r=-8\pi P_{\perp }=-\frac{9{\tilde{a}}^2}{2}$, which means that the inertial mass density $(\mu +P_r)$ is negative. It is worth noticing that the inertial mass density is always negative and not only in the limit $t\to \infty$. On the other hand, we see that the expression within the square bracket in the “gravitational term” in Equation (A9) (the first term on the right of Equation (A9)) is negative as $t\to \infty$, producing a positive $D_{T}U$ (i.e., such a term acts as a repulsive force). For sufficiently small values of $\tilde{a}$, the other parameters of the solution may be chosen such that for some finite time interval, the energy density is positive. Since the heat flux is constant, no contributions from the transient period (terms proportional to $\tau$) appear in the expression of the temperature. Thus, this solution might be used to model expansion-free evolution only for a limited time interval.

The second model also satisfies the vanishing complexity factor condition, but instead of the quasi-homologous evolution, we assumed that the metric functions A and B were proportional. Its evolution is described by Equations (83)–(90). By choosing for all the parameters of the solution to be positive (thereby avoiding shells crossing singularities), then in the limit $t\to \infty$, the areal radii of all fluid elements tend toward infinity, and the fluid distribution becomes a shell. However, in this case, depending on the specific values of the parameter ${\beta }_0$, the behavior of the model may be very different. Indeed, as follows from Equations (84)–(87) in the limit $t\to \infty$, the physical variables tend toward zero if ${\beta }_0>2$ and diverge to infinity if ${\beta }_0<2$. If ${\beta }_0=2$, then the model has a static limit described by the equation of the state $8\pi \mu =8\pi P_r=-\frac{3}{{\beta }_2^2{\alpha }^2}\left(s_1^2+\frac{s^2}{{\gamma }^2}\right)$ with a negative energy density and radial pressure. In this case (${\beta }_0=2$), both q and $P_{\perp }$ are constant at all times, implying that the transient effects in temperature vanish, as is apparent from Equation (90).

The third solution satisfies the vanishing complexity factor condition, the metric function A only depends on r, and R is a separable function. The full description of this model is provided by Equations (92)–(101). They describe a collapsing fluid for which, in the limit $t\to \infty$, the energy density and the radial pressure diverge and satisfy the equation of the state $\mu =-P_r>0$, whereas the heat flux vector and the tangential pressure vanish. In this limit, the transient effects vanish too. An appropriate choice of the parameters allows one to obtain well-behaved physical variables, at least for a finite time interval.

All three solutions described above are non-geodesic. The next two models instead have a vanishing four-acceleration.

The first one satisfies the vanishing complexity factor condition and is described by Equations (109)–(117). This model depicts a collapsing fluid for which as $t\to \infty$, the energy density and the radial pressure diverge and satisfy the equation of the state $\mu =-P_r$, whereas the heat flux vector and the tangential pressure vanish, and the temperature tends toward $T_0$. For sufficiently large (but finite) values of t, the energy density is positive, the radial pressure is negative, and the fluid evolves almost adiabatically.

Finally, the second geodesic model is described by Equations (118)–(125). In this model, the vanishing complexity factor condition is replaced by the quasi-homologous condition. As in the previous model, this one depicts a collapsing fluid for which as $t\to \infty$, the energy density and the radial pressure diverge and satisfy the equation of the state $\mu =-P_r$, the heat flux vector and the tangential pressure vanish, and the temperature tends toward $T_0$. Additionally, for sufficiently large (but finite) values of t, the energy density is positive, the radial pressure is negative, and the fluid evolves almost adiabatically.

To summarize, the five models presented here might describe some physical realistic situations for finite time intervals. We noticed that neither of them satisfy the Darmois conditions on either boundary surface, implying that these are thin shells.

We would like to conclude with the following remarks:

• The analytical models presented here have the main advantage of simplicity, which allows one to use them as test models for describing the evolution of voids. However, they were obtained under specific restrictions, some of which are of a purely heuristic nature. In order to get closer to a physically meaningful scenario, one should use some observational data as input for solving the field equations. At this point, the best candidate for that purpose appears to be the luminosity profile produced by the dissipative processes within the fluid. Afterward, it seems unavoidable to resort to a numerical approach in order to solve the field equations.
• In the first two models, the vanishing complexity factor condition leads to two differential equations (Equations (61) and (77)) which have been solved analytically, resorting to the heuristic ansatz in Equations (62) and (78), respectively. Of course, a much more satisfactory procedure would be to solve those equations using numerical methods. However, this would be out of the scope of this work.