Loading [MathJax]/jax/output/HTML-CSS/jax.js

First time here? Check out the FAQ!

Ask Your Question
1

Can't substitute H=dot(a)/a in SageMath

asked 5 years ago

cduston gravatar image

Hello all, I'm working with the Friedman equations, and I've gotten them down to the common form presented in terms of the scale factor. MWE incoming:

M = Manifold(4, 'M', structure='Lorentzian')
fr.<t,r,th,ph> = M.chart(r't r:[0,+oo) th:[0,pi]:\theta ph:[0,2*pi):\phi')
var('G, Lambda, k', domain='real')
a = M.scalar_field(function('a')(t), name='a')
rho = M.scalar_field(function('rho')(t), name='rho')
p = M.scalar_field(function('p')(t), name='p')
g = M.metric()
g[0,0] = -1
g[1,1] = a*a/(1 - k*r^2)
g[2,2] = a*a*r^2
g[3,3] = a*a*(r*sin(th))^2
nabla = g.connection()
Ricci = nabla.ricci()
Ricci_scalar = g.ricci_scalar()
E1 = Ricci - Ricci_scalar/2*g + Lambda*g
print("First Friedmann equation:\n")
E1[0,0].expr().expand() == 0

With the result being

-Lambda + 3*diff(a(t), t)^2/a(t)^2 + 3*k/a(t)^2 == 0

Which is great, but I really want this to be in terms of the Hubble paramter H=˙a/a. I can't make that substitution happen:

var('h')
eq1.substitute(diff(a.expr(),t)/a.expr()==h)

just spits out the same thing,

-Lambda + 3*diff(a(t), t)^2/a(t)^2 + 3*k/a(t)^2 == 0

Any ideas?

Preview: (hide)

2 Answers

Sort by » oldest newest most voted
1

answered 5 years ago

eric_g gravatar image

Assuming that eq1 is E1[0,0].expr().expand() == 0, what about

sage: eq1.subs({diff(a.expr(),t): h*a.expr()})
3*h^2 - Lambda + 3*k/a(t)^2 == 0

Alternatively, you may want to keep h being a function of t:

sage: eq1.subs({diff(a.expr(),t): function('h')(t)*a.expr()})
3*h(t)^2 - Lambda + 3*k/a(t)^2 == 0
Preview: (hide)
link
1

answered 5 years ago

rburing gravatar image

I hope you appreciate that this is not such a simple operation to perform in general.

The simplest way is probably to replace ˙a by Ha.

We can use the class I gave as an answer to Substituting a differential equation into an expression:

sage: eq1 = E1[0,0].expr().expand() == 0
sage: var('h')
sage: pde = diff(a.expr(), t) == h*a.expr()
sage: s = SubstituteEvolutionaryPDE(pde, t)
sage: s(eq1)
3*h^2 - Lambda + 3*k/a(t)^2 == 0

This also works when there are second or mixed derivatives, e.g.:

sage: eq2 = diff(a.expr(), t, t) + 13 == 0
sage: s(eq2)
h^2*a(t) + 13 == 0
Preview: (hide)
link

Comments

Well, since this is my major frustration with Sage, I know that it is not a simple operation, but I can't for the life of me understand why it should be. Although your answer does the job, I'm going to select the other one as it doesn't require a new class definition. Although please, disagree and convince me!

cduston gravatar imagecduston ( 5 years ago )

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 5 years ago

Seen: 391 times

Last updated: Feb 14 '20