Plot Structures for Friedman-Robertson-Walker metric, k=-1

Before anything else, we know that we will need some special routines that Maple uses for plotting solutions of

ordinary differential equations which are not always in memory; therefore, we must load them in:

> with(plots):

A few numbers that will be useful later on are the following:

cMpc is the speed of light in Megaparsecs per eon (an eon is a billion years)

gMpc is the factor to multiply by to change 1/(Megaparsecs)^2 to grams/cm^3

> cMpc:=306.937; gMpc:= 1/(7.0554*10^20);

We may now begin by stating the values of the constants that will be used here;

those values are explained and/or justified in a different place. We do this one with constants so that

the lifecycle of the universe is infinite in time.

We have need for constants, A, B, k, and Lambda; we denote the value for B^2 by B2 and the value for

Lambda/3 by Lambda3. A is in Mpc, B2 in Mpc squared and Lambda in inverse Mpc squared.

> A:=10772; B2:=24240; k:=-1; Lambda3:=-1/8226^2;

We may then write down the fourth-order polynomial appropriate for this problem, where I use the

symbol x for the usual arc-time variable, eta:

> R1:=A*R(x)+B2-k*(R(x))^2-Lambda3*(R(x))^4;

and then we may use that to describe the differential equations that need to be solved. Although there are in fact two branches for the square root involved, we may simply use only the positive branch.

The reason for this is that the polynomial above has no positive roots, having all positive coefficients.

> sys1:=diff(R(x),x)=R1^(0.5), diff(t(x),x)=R(x);

To actually solve the differential equations we must append initial conditions; these are straightforward

since the universe is beginning: we begin with R(0)=0 and t(0)=0:

> rw1:=dsolve({sys1,R(0)=0,t(0)=0},[R(x),t(x)],type=numeric):

> r1:=subs(R(x)=r,R1);

We now plot this function to get a feel for its shape, knowing that we want it to begin at R=0 and go until the polynomial

has its first, real root. We change the number 7000 in the plot command below two or three times until we get it right.

> plot(r1,r=0..7000);

As we see, it has no intention of having any positive roots.

> RW1:=odeplot(rw1,[x,R(x)],0..2.2,color=red,numpoints=500):

> RWt1:=odeplot(rw1,[x,t(x)],0..2.2,color=green,numpoints=500):

> display([RW1,RWt1], title=`R(eta) and t(eta), first half`);

We have decided to choose the current radius to be the following, in Megaparsecs.

> R0:=6500;

We want to find the associated values of eta=x and of time, which we read off the datapoints for the

graph above. This is done by first typing in the RW-commands above with a semicolon, so that Maple

prints out the entire list of points, and reading off the ones desired, and then changing it back to a colon,

so that long list is not printed out unnecessarily. The value of eta will then be 1.38, and that of time is

2800 Mpc. This can be changed to eons, i.e., billions of years by using the conversion factor above.

It is worth noting that the number is TOO SMALL, being less than 10 eons; therefore, it would have been

reasonable to have changed the associated parameters already chosen. Nonetheless, this does show the

technique and I will take that as enough for now.

> eta0:=1.38; T0:=2800; t0:=T0/cMpc;

To find the density at maximum expansion, we ignore the radiation density, since we know it has been getting

smaller faster, and calculate the matter density at that time from its R-dependence: the first number is in

inverse Megaparsecs squared, the second is in grams/cubic centimeter.

> rhomax:=(3/evalf(8*Pi))*A/Rmax^3; Rhomax:=rhomax*gMpc;

The above graph gives us an idea what this first half looks like, with the red curve being R, and the green one t.

Yet a different alternative is to actually plot R versus t:

> RWrt1:=odeplot(rw1,[t(x),R(x)],0..2.2,numpoints=500,title=`R=R(t)`):

> display(RWrt1,title=`R=R(t), first half`);

We are also asked for various other quantities.

The first of these is the time when the matter and radiation densities were the same, i.e., when

A/R^3 = B2/R^4, which means that R was equal to B2/A. The answer for R is given in Megaparsecs.

Looking that value up in the tables created above give us the corresponding values of eta=x and of time.

The first answer for time is also in Megaparsecs, but it is then changed to years by dividing by the speed of

light and multiplying by 1 billion:

> Requal:=evalf(B2/A); etaequal:=.01197; Tequal:=.01312; tequal:=Tequal*10^9/cMpc;

We also want to know the current value of the Hubble constant, in km/sec/Megaparsec

the current rate of increase of the size of the universe, i.e., dR/dt, in Megaparsecs/eon,

do note that this number is greater than the speed of light measured locally,

and the current value for q0, which is the second derivative, appropriately normalized:

q0 = -(d^2R/dt^2)(R/H0^2), which is dimensionless:

> H0:=sqrt(evalf(A/R0^3+B2/R0^4 -k/R0^2-Lambda3))*(2.9979*10^5);

> Rdot:=H0*R0*cMpc/(2.9979*10^5);

> qnum:=A/(2*R)+B2/R^2+Lambda3*R^2;

> qnum0:=evalf(subs(R=R0,qnum));

> qdenom:=r1/r^2;

> qdenom0:=evalf(subs(r=R0,qdenom));

> q0:=qnum0/qdenom0;

We also want to know about the source that is observed on Earth today with a redshift of Z=5. This implies

immediately that it was emitted when R=R(t) was R0/(Z+1) = 1083.3 Mpc. Again looking at the tables, this

gives us a value for eta of 0.596 and for t of approximately 0.48 Mpc, which is 1.56 billions of years after

the initial Big Bang, or 9.12 billions of years ago.

> RZ:=evalf(R0/6); etaZ:=0.596; TZ:=0.48; tZ:=TZ/cMpc; DeltaTZ:=t0-tZ;

To determine the current metric distance of that source from us, we multiply the current value of R(t),

namely R0, by the difference in eta between the emission event and the reception event, the result being

given in Megaparsecs. We see that this source is now 89% of the distance R0 away from us.

> LZ:=R0*(eta0-etaZ); fract:=LZ/R0;

Lastly we are asked about the particle horizon. There we pretend that light could now be received that was

originally emitted at (or almost at) the time of the Big Bang, so that it has had an arc-time, eta, to travel of the

current value of eta, called eta0 above. The current metric distance to a galaxy at that horizon is then given by

the following calculation, where we see that it is 140% of the distance R0 away from us, at this time.

> LPart:=R0*(eta0-0); fractPart:=LPart/R0;

All answers to questions about a maximum radius are irrelevant here, since this particular universe has no

maximum, but continues to expand forever.