Application Center - Maplesoft

App Preview:

Alexander Friedmann's Cosmic Scenarios

You can switch back to the summary page by clicking here.

Learn about Maple
Download Application

Alexander Friedmann's Cosmic Scenarios

Frank Wang,


The Russian mathematician and physicist Alexander Friedmann (1888-1925) is well known among relativists, but his contributions to cosmology are largely misunderstood.  Even the Royal Swedish Academy of Sciences misrepresented Friedmann's work in the 2011 Nobel Prize scientific background essay.  Friedmann was the first physicist who demonstrated that Albert Einstein's general relativity admits non-static solutions, and the universe can expand, oscillate, and be born in a singularity.  Friedmann's conclusion was based on his analysis of an elliptic integral; this worksheet employs Maple's utility of handling elliptic integrals to present Friedmann's results graphically.  Friedmann's differential equation governing the evolution of the universe based on Einstein's general theory of relativity is also derived using Maple's tensor package.          


The 2011 Nobel Prize in Physics was awarded to Saul Perlmutter, Brian P. Schmidt and Adam G. Riess for "the discovery of the accelerating expansion of the Universe through observations of distant supernovae."  In the Scientific Background on the Nobel Prize in Physics 2011 compiled by the Royal Swedish Academy of Sciences, the Russian mathematician and physicist Alexander Friedmann was said to have "studied the problem of the dynamics of the Universe using essentially the same assumptions as Einstein, and found in 1922 that Einstein's steady state solution was really unstable."  As pointed out by Ari Belenkiy in a Physics Today article, this account significantly distorts Friedmann's contributions.  The Swedish Academy's oversight is just among the many examples of a wide spread misinformation of Friedmann's paper.    


Before discussing Friedmann's work, we need to review Einstein's general theory of relativity.  In 1915, Einstein obtained the celebrated field equations


R[mu*nu]-(1/2)*g[mu*nu]*R = -kappa*T[mu*nu]


where the spacetime indices mu, nu run from 1 to 4 and the constant kappa = 8*Pi*G/c^4.  (Here G is gravitational constant and c the speed of light.)  The spacetime distribution of the energy-momentum tensor T[mu*nu] determines the local geometry encoded in the metric tensor g[mu*nu], and the Ricci tensor R[mu*nu] is determined by g[mu*nu] and its spacetime derivatives.  The Ricci scalar R is the local curvature of spacetime.  (In the appendix, we will demonstrate how to use Maple's tensor package to evaluate Ricci tensor and scalar for a given metric tensor.)  To simplify calculations, the universe can be approximated as a perfect fluid of uniform density rho and pressure P.  The energy-momentum has the following elements: T[11] = -P, T[22] = -P, T[33] = -P, T[44] = c^2*rho, and all off-diagonal elements are zero.  Furthermore, we consider a matter-dominated universe so that the pressure is zero: P = 0.  


In 1917, Einstein introduced the cosmological constant lambda in the hope of finding a static cosmological solution of the field equations    


R[mu*nu]-(1/2)*g[mu*nu]*R-lambda*g[mu*nu] = -kappa*T[mu*nu] .


The geometry of Einstein's static universe can be regarded as a three-dimensional spherical surface


z[1]^2+z[2]^2+z[3]^2+z[4]^2 = a^2


embedded in the four-dimensional Euclidean space (z[1], z[2], z[3], z[4]).  

Friedmann considered a space as a 3D hypersphere similar to that of Einstein's universe, but the radius changes with respect to time.     


z[1]^2+z[2]^2+z[3]^2+z[4]^2 = a(t)^2


Swedish Academy's statement that Friedmann "[used] essentially the same assumptions as Einstein" tremendously undervalued Friedmann's originality.  With this non-static model, the time-time component of the Einstein field equations yields a first-order differential equation, known as the Friedmann equation


(diff(a(t), t))^2/c^2 = (A-a(t)+lambda*a(t)^3/(3*c^2))/a(t)


where A is proportional to the total mass of space M.  (To be precise, A = kappa*M/(6*Pi^2).  See appendix for its derivation.)  The cosmic radius a(t) in Friedmann's universe is no longer constant, and several scenarios for cosmic evolution can occur.   


From the Friedmann equation, we can obtain a(t) by inversion of an elliptic integral


c*t = Int(sqrt(x/(A-x+lambda*x^3/(3*c^2))), x = a0 .. a)+b


in which a0 and b are constants.  Let us enter this integral to Maple and see what happens.

restart: with(plots):

int(sqrt(x/(A-x+lambda*x^3/(3*c^2))), x):

Here we end the Maple input with colon to suppress the output, because it is lengthy and not too illuminating.  The reader can change the colon to semicolon and see that the result is expressed in terms of EllipticF and EllipticPi.  

The integral has physical meaning only when the cubic denominator




under the square-root sign is positive.  Friedmann denote lambda/(3*c^2) by y and consider in the (x, y) plane the family of curves of third degree:


x^3*y+A-x = 0 .  With A = 1/2, 1, 3/2, he produced this figure:

plot({(x-1/2)/x^3, (x-1)/x^3, 4/27, (x-3/2)/x^3}, x=0..7, -0.5..0.7, labels=[x,y]);


For the equation x^3*y+A-x = 0 or y = (x-A)/x^3, we can prove that the curve has a maximum at


x = 3*A*(1/2) , y = 4/(27*A^2)  .  See below.

y := (x-A)/x^3;



eq1 := diff(y,x) = 0;

1/x^3-3*(x-A)/x^4 = 0


solve(eq1, x);



eval(y, x=%);



Keep in mind that in the cubic denominator x^3*y+A-x, A is proportional to the total mass of space (A = kappa*M/(6*Pi^2)), and y is proportional to cosmological constant (y = lambda/(3*c^2)).  If y < 4/(27*A^2), or lambda < 4*c^2/(9*A^2), the cubic equation x^3*y+A-x = 0 has two positive roots.  If 4/(27*A^2) < y, or 4*c^2/(9*A^2) < lambda, the equation has no positive roots.  Let us set A = 1 and use y = 2/27, 4/27, 8/27 as representative numerical examples to plot Friedmann's cubic expression.  

cd1 := 8/27*x^3 - x + 1;



cd2 := 4/27*x^3 - x + 1;



cd3 := 2/27*x^3 - x + 1;



plot([cd1, cd2, cd3], x=0..3.5, thickness=2, legend=['y'>4/(27*A), 'y'=4/(27*A), 'y'<4/(27*A)], labels=[x, "cubic denominator"], labeldirections=[horizontal, vertical]);


We consider the first scenario that x^3*y+A-x has no positive roots by using y = 8/27 as an example.  The integration of the cosmic radius becomes

Int(sqrt(x/(8/27*x^3 - x + 1)), x);

Int((x/((8/27)*x^3-x+1))^(1/2), x)


i1 := value(%);

-6*(EllipticF(2*((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)*x*(4+2*2^(1/2))^(1/3)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)))^(1/2), ((I*(4+2*2^(1/2))^(2/3)*3^(1/2)-3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)-6)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)-(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)-2)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)))^(1/2))-EllipticPi(2*((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)*x*(4+2*2^(1/2))^(1/3)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)))^(1/2), (I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)/(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6), ((I*(4+2*2^(1/2))^(2/3)*3^(1/2)-3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)-6)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)-(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)-2)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)))^(1/2)))*2^(1/2)*(((4+2*2^(1/2))^(2/3)+2)*((3*I)*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-8*x*(4+2*2^(1/2))^(1/3)-(6*I)*3^(1/2)+6)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)))^(1/2)*(((4+2*2^(1/2))^(2/3)+2)*((3*I)*(4+2*2^(1/2))^(2/3)*3^(1/2)-3*(4+2*2^(1/2))^(2/3)+8*x*(4+2*2^(1/2))^(1/3)-(6*I)*3^(1/2)-6)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)-(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)-2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)))^(1/2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)^2*((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)*x*(4+2*2^(1/2))^(1/3)/((I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)))^(1/2)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+2)*(8*x^3-27*x+27)*(x/(8*x^3-27*x+27))^(1/2)*3^(1/2)/((4+2*2^(1/2))^(2/3)*(-2*x*(3*(4+2*2^(1/2))^(2/3)+4*x*(4+2*2^(1/2))^(1/3)+6)*((3*I)*(4+2*2^(1/2))^(2/3)*3^(1/2)-3*(4+2*2^(1/2))^(2/3)+8*x*(4+2*2^(1/2))^(1/3)-(6*I)*3^(1/2)-6)*((3*I)*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-8*x*(4+2*2^(1/2))^(1/3)-(6*I)*3^(1/2)+6)/(2+2^(1/2)))^(1/2)*(I*(4+2*2^(1/2))^(2/3)*3^(1/2)+3*(4+2*2^(1/2))^(2/3)-(2*I)*3^(1/2)+6)*((8*x^3-27*x+27)*x)^(1/2))


This is quite a lengthy output.  Note that we obtain time as a function of radius a (or x).  We are more interested in observing the radius as a function of time.  We use the parametric form of the plot command to do this trick.

plot([i1, x, x=0..7], thickness=2, color=black); pm1 := %:


We see that when 4*c^2/(9*A^2) < lambda, that is, when the cosmological constant exceeds some critical value that depends on density rho or total mass, the cosmos starts at t = 0 from the singularity a(t) = 0, and it expands forever.  Friedmann called this scenario "a monotonic world of the first kind."  A careful reader might notice that this curve has an inflection point.  To locate the inflection point, we need to find the second derivative diff(t(r), r, r) and solve the equation


diff(int(sqrt(x/(A-x+lambda*x^3/(3*c^2))), x), x, x) = 0 , or diff(sqrt(x/(A-x+lambda*x^3/(3*c^2))), x) = 0 .  

solve(diff(sqrt(x/(A-x+lambda/(3*c^2)*x^3)), x)=0, x);

(1/2)*12^(1/3)*(A*c^2*lambda^2)^(1/3)/lambda, -(1/4)*12^(1/3)*(A*c^2*lambda^2)^(1/3)/lambda+((1/4)*I)*3^(1/2)*12^(1/3)*(A*c^2*lambda^2)^(1/3)/lambda, -(1/4)*12^(1/3)*(A*c^2*lambda^2)^(1/3)/lambda-((1/4)*I)*3^(1/2)*12^(1/3)*(A*c^2*lambda^2)^(1/3)/lambda


Taking the real solution, we conclude that the expansion rate changes from deceleration to acceleration at a = (3*c^2*A/(2*lambda))^(1/3).  For large a, we can neglect the constant and linear terms of the cubic denominator:


t = (int(1/(sqrt(lambda/(3*c^2))*x), x))/c .


This is an elementary integration, see the Maple output below.

eqasym := t=int(1/sqrt(lambda/3)/x, x);

t = 3^(1/2)*ln(x)/lambda^(1/2)


solve(eqasym, x);



We conclude a(t) grows asymptotically like exp(sqrt((1/3)*lambda)*t).  This scenario introduced by Friedmann appears to correspond to reality, as observed by the three Nobel laureates.  


The second scenario occurs when the cosmological constant is less than some critical value (0 < lambda and 4*c^2/(9*A^2) < lambda).  Here we use y = 2/27 as a representative example to plot the cubic denominator.

plot(cd3, x=0..3.5, labels=[x, "cubic denominator"], labeldirections=[horizontal, vertical]);


In this situation, a(t) can lie in the intervals (0, x[1]) and (x[2], infinity).  For this particular numerical example, we can use Maple to find the two positive roots.

sort([solve(cd3=0, x)]);

[3, -3/2-(3/2)*3^(1/2), -3/2+(3/2)*3^(1/2)]


x1 := %[3]; x2 := %%[1];





Consider the scenario that x[2] < a(t).   

Int(sqrt(x/(1-x+2/27*x^3)), x=x2..r);

Int((x/((2/27)*x^3-x+1))^(1/2), x = 3 .. r)


i2 := value(%);

int((x/((2/27)*x^3-x+1))^(1/2), x = 3 .. r)


plot([i2, r, r=3..5], thickness=2, color=blue); pm2 := %:


We see that the expansion of the universe starts from a nonzero value x[2] and expands forever.  Friedmann called it a monotonic world of the second kind.  

Friedmann cited Weierstrass' work and demonstrated that when lambda lies in the interval (-infinity, 4*c^2/(9*A^2)), a(t) is a periodic function of t.  The period is given through

2*Int(sqrt(x/(1-x+2/27*x^3)), x=0..x1);

2*(Int((x/((2/27)*x^3-x+1))^(1/2), x = 0 .. -3/2+(3/2)*3^(1/2)))


tp := value(%);

-24*3^(1/2)*(-EllipticK(I*(3^(1/2)-1)^(1/2)*(3+3^(1/2))^(1/2)/((3-3^(1/2))^(1/2)*(1+3^(1/2))^(1/2)))+EllipticPi((3^(1/2)-1)/(-3+3^(1/2)), I*(3^(1/2)-1)^(1/2)*(3+3^(1/2))^(1/2)/((3-3^(1/2))^(1/2)*(1+3^(1/2))^(1/2))))*(4*3^(1/2)-7)*2^(1/2)/((3-3^(1/2))^(1/2)*(1+3^(1/2))^(1/2)*(3^(1/2)-1)^2*(-2+3^(1/2)))





Int(sqrt(x/(1-x+2/27*x^3)), x);

Int((x/((2/27)*x^3-x+1))^(1/2), x)


ip := value(%);

-2*3^(1/2)*(x/(2*x^3-27*x+27))^(1/2)*(2*x^3-27*x+27)*(3^(1/2)-1)*(-(3-3^(1/2))*x/((3^(1/2)-1)*(x-3)))^(1/2)*(x-3)^2*(-(6*x+9+9*3^(1/2))/((1+3^(1/2))*(x-3)))^(1/2)*(-(-6*x-9+9*3^(1/2))/((3^(1/2)-1)*(x-3)))^(1/2)*2^(1/2)*(EllipticF((-(3-3^(1/2))*x/((3^(1/2)-1)*(x-3)))^(1/2), I*((9/2+(3/2)*3^(1/2))*(-3/2+(3/2)*3^(1/2))/((3/2+(3/2)*3^(1/2))*(9/2-(3/2)*3^(1/2))))^(1/2))-EllipticPi((-(3-3^(1/2))*x/((3^(1/2)-1)*(x-3)))^(1/2), (3^(1/2)-1)/(-3+3^(1/2)), I*((9/2+(3/2)*3^(1/2))*(-3/2+(3/2)*3^(1/2))/((3/2+(3/2)*3^(1/2))*(9/2-(3/2)*3^(1/2))))^(1/2)))/(((2*x^3-27*x+27)*x)^(1/2)*(-3+3^(1/2))*(-x*(x-3)*(2*x+3+3*3^(1/2))*(-2*x-3+3*3^(1/2)))^(1/2))


p1 := plot([ip, x, x=0..x1], thickness=2, color=red):

p2 := plot([tp-ip, x, x=0..x1], thickness=2, color=red):

display([p1, p2]); pp := %:


Friedmann called this world "the periodic world": expansion oscillates between a(t) = 0 and a(t) = x[1].  

When lambda has precisely the critical value lambda = 4*c^2/(9*A^2), then the cubic denominator is degenerate: it has a double positive root.  This scenario corresponds to Einstein's static universe.


Here is a plot of three possible scenarios that we have discussed.

display([pm1, pm2, pp], labels=["time", "cosmic radius"], labeldirections=[horizontal, vertical], labelfont=[TIMES, ROMAN, 14]);


Alexander Friedmann's pioneer work is summarized in this plot.  He for the first time showed that Einstein's general relativity admits non-static solutions.  He called the black curve "a monotonic world of the first kind."  The cosmos expands at a decelerating rate from a zero-radius singularity until an inflection point, after which the expansion accelerates.  This scenario provides the basis for our current view of the Big Bang and the accelerating universe.  The blue curve is named as "a monotonic world of the second kind" in which the universe expands from a nonzero initial radius.  The periodic scenario is shown in the red curve: the universe evolves from and back to zero radius.   




John Wheeler summarized generality in nontechnical terms as "geometry tells matter how to move, and matter tells geometry how to curve."  This statement captures the essence of Einstein's field equations


R[mu*nu]-(1/2)*g[mu*nu]*R = -kappa*T[mu*nu] .


See main text for further elaboration of the meaning of the symbols.  The right-hand side of the equation is proportional to the energy-momentum tensor, which is related to matter density rho.  The left-hand side of the equation (consisting of Ricci tensor and scalar) is the Einstein tensor.  The square of the spacetime interval defines the metric g[mu*nu]:


ds^2 = g[mu*nu]*dx^mu*dx^nu


Here we use Einstein summation convention: whenever a Greek index is repeated in a term, summer over this index from 1 to 4 is implied.  


The 3D hypersphere considered by Friedmann was


z[1]^2+z[2]^2+z[3]^2+z[4]^2 = a(t)^2 .


(In Einstein's universe, a is independent of time.)  It is convenient to transform from Cartesian to polar coordinates:


z[1] = a(t)*sin(chi)*sin(theta)*cos(phi)

z[2] = a(t)*sin(chi)*sin(theta)*sin(phi)

z[3] = a(t)*sin(chi)*cos(theta)

z[4] = a(t)*cos(chi)


The spacetime geometry is then described by the metric


ds^2 = a(t)^2*(d*chi^2+sin(chi)^2*(d*theta^2+sin(theta)^2*d*phi^2))-c*dt^2


Providing g[mu*nu] to Maple's tensor package, one easily obtains Ricci tensor, Ricci scalar, and Einstein tensor.  



coords := [chi, theta, phi, t]:
g := array(symmetric, sparse, 1..4, 1..4):
g[4,4] := c^2: g[1,1] := -(a(t))^2: g[2,2] := -(a(t)*sin(chi))^2: g[3,3] := -(a(t)*sin(chi)*sin(theta))^2:
metric := create([-1,-1], eval(g));

metric := table([index_char = [-1, -1], compts = (Matrix(4, 4, {(1, 1) = -a(t)^2, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (2, 1) = 0, (2, 2) = -a(t)^2*sin(chi)^2, (2, 3) = 0, (2, 4) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = -a(t)^2*sin(chi)^2*sin(theta)^2, (3, 4) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 0, (4, 4) = c^2}))])


tensorsGR(coords, metric, contra_metric, det_met, C1, C2, Rm, Rc, R, G, C):

displayGR(Einstein, G);


`The Einstein Tensor`


`non-zero components :`


` G11` = ((diff(a(t), t))^2*cos(chi)^2*cos(theta)^2+2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2*cos(theta)^2+cos(chi)^2*cos(theta)^2*c^2-(diff(a(t), t))^2*cos(chi)^2-(diff(a(t), t))^2*cos(theta)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(theta)^2-cos(chi)^2*c^2-cos(theta)^2*c^2+(diff(a(t), t))^2+2*a(t)*(diff(diff(a(t), t), t))+c^2)/(sin(chi)^2*sin(theta)^2*c^2)


` G22` = ((diff(a(t), t))^2*cos(chi)^2*cos(theta)^2+2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2*cos(theta)^2+cos(chi)^2*cos(theta)^2*c^2-(diff(a(t), t))^2*cos(chi)^2-(diff(a(t), t))^2*cos(theta)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(theta)^2-cos(chi)^2*c^2-cos(theta)^2*c^2+(diff(a(t), t))^2+2*a(t)*(diff(diff(a(t), t), t))+c^2)/(sin(theta)^2*c^2)


` G33` = ((diff(a(t), t))^2*cos(chi)^2*cos(theta)^2+2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2*cos(theta)^2+cos(chi)^2*cos(theta)^2*c^2-(diff(a(t), t))^2*cos(chi)^2-(diff(a(t), t))^2*cos(theta)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(chi)^2-2*a(t)*(diff(diff(a(t), t), t))*cos(theta)^2-cos(chi)^2*c^2-cos(theta)^2*c^2+(diff(a(t), t))^2+2*a(t)*(diff(diff(a(t), t), t))+c^2)/c^2


` G44` = -3*((diff(a(t), t))^2*cos(chi)^2*cos(theta)^2+cos(chi)^2*cos(theta)^2*c^2-(diff(a(t), t))^2*cos(chi)^2-(diff(a(t), t))^2*cos(theta)^2-cos(chi)^2*c^2-cos(theta)^2*c^2+(diff(a(t), t))^2+c^2)/(a(t)^2*sin(chi)^2*sin(theta)^2)



`character : [-1, -1]`


Below we simplify the time-time component of the Einstein tensor by copying G44 from Maple's output.

Gtt := combine( -3*(-cos(theta)^2*c^2+c^2-cos(chi)^2*c^2+cos(chi)^2*c^2*cos(theta)^2+diff(a(t),t)^2-diff(a(t),t)^2*cos(theta)^2-diff(a(t),t)^2*cos(chi)^2+diff(a(t),t)^2*cos(chi)^2*cos(theta)^2)/a(t)^2/sin(chi)^2/sin(theta)^2);

(-3*(diff(a(t), t))^2-3*c^2)/a(t)^2


The Friedman equation is the time-time component of the field equations:

Eq1 := Gtt + lambda = -kappa*rho*c^2;

(-3*(diff(a(t), t))^2-3*c^2)/a(t)^2+lambda = -kappa*rho*c^2


Eq2 := isolate(Eq1, diff(a(t),t)^2);

(diff(a(t), t))^2 = -(1/3)*(-c^2*kappa*rho-lambda)*a(t)^2-c^2





The volume of a 3D hypersphere is obtained to the integral below.

Int(Int(Int(a(t)^3*sin(chi)^2*sin(theta), chi=0..Pi), theta=0..Pi), phi=0..2*Pi);

Int(Int(Int(a(t)^3*sin(chi)^2*sin(theta), chi = 0 .. Pi), theta = 0 .. Pi), phi = 0 .. 2*Pi)


V := value(%);



The total mass is the product of this volume and the density.

M := rho*V;



With A defined as A = kappa*M/(6*Pi^2), the Fridmann equation is written as

(diff(a(t), t))^2/c^2 = (A-a(t)+lambda*a(t)^3/(3*c^2))/a(t) .  

Acknowledgements:  In December 2014, I revised the worksheet after two readers, Ari Belenkiy and Frank Bensch, posted their comments online.  I thank them for their remarks.  The solution order from Maple's solve output might vary in different versions, and the sort command was added to address this issue.  Thanks to Eithne Murray, Maplesoft's Project Manager, for offering valuable suggestions to me.  My statement about the degenerate case was oversimplified.  For a more complete analysis, please consult Professor Belenkiy's paper, particularly Figure 6 and the paragraph after it.  




1. Royal Swedish Academy of Sciences, "The Accelerating Universe," (2011), available at  

2. Ari Belenkiy, "Alexander Friedmann and the origins of modern cosmology," Physics Today, October 2012, 38-43.  

3. A. Friedmann, "Ueber die Kruemmung des Raumes," Zeitschrift fuer Physik 10, 377-386 (1922).  English translation in A. Friedman, "On the Curvature of Space," General Relativity and Gravitation, 31, 2991-2000 (1999).  

4. C. W. Misner, K. S. Thorne and J. A. Wheeler, Gravitation, Freeman, 1973.

5. Frank Wang, Physics with Maple: The Computer Algebra Resource for Mathematical Methods in Physics, Wiley-VCH, 2006.