This isn't yet a complete question; rather, I'm looking for a qual-level question and answer describing a gravitationally bound system in an expanding universe. Since it's qual level, this needs a vastly simplified model, preferably one in a Newtonian framework (if this is even possible) that at least shows the spirit of what happens. If it isn't possible, and there is some important principle that any Newtonian model will miss, that is probably the first thing to point out. Otherwise, ...
I'm thinking of a setup that starts something like this:
A body of mass $m$ orbits another of mass $M \gg m$ (taken to be the origin) in a circle of radius $R$ (and consequently with a period $T = 2 \pi \sqrt{ R^3 / GM }$). At time $t = 0$, space begins to slowly expand at a constant Hubble rate $H$ (where "slowly" means $ HT \ll 1 $, so that the expansion in one period is negligible compared to the orbital radius).
Now here is where I'm not sure how to phrase the problem. I'm thinking instead of something like the "bead on a pole" problems. Then the bead has generalized Lagrangian coordinates relative to a fixed position on the pole, but the pole is expanding according to some external function $a(t)$, a la FLRW metric (For a constant Hubble rate, $H = \dot a / a$ and so $ a(t) = e^{Ht} $). So, we can call the bead's distance from the origin in "co-moving units" the Lagrangian coordinate $q$. Then the distance of the bead from the origin is $d(t) = q(t) a(t)$, giving us a time dependent Lagrangian.
Basically, I need a 3D version of this that will work for the gravitationally bound system described. I also hope that the problem can be solved as an adiabatic process, which might mean changing coordinates to be around the original position relative to the more massive body at the origin, rather than co-moving coordinates.
EDIT 1: I won't claim that I fully understood the details of Schirmer's work, but I think one of the big takeaway points is that the cosmological expansion damages the stability of circular orbits and causes them to decay. I finished my very hand-waivey "Newtonian" model of a solar system in an expanding universe, and I don't think it captures the essence of this part of the GR solution. The model does have several sensical properties:
- There is a distance at which bound solutions are impossible and the two bodies are guaranteed to expand away from each other
- There is still a circular orbit (I haven't checked that it is stable) whose radius is modified by a term involving the Hubble constant.
- Using solar parameters, this shift is entirely negligible; it is larger (though still small) on a galactic scale.
Here is the model:
Without this scaling factor, the Lagrangian would be $$ \mathcal{L} = \frac{1}{2} m \left( \dot r^2 + r^2 \dot \theta^2 + r^2 \sin^2 \theta \, \dot \phi^2 \right) + \frac{GMm}{r} $$ Adding in the scale factor, all the kinetic terms will pick up a factor $e^{2Ht}$ and the potential a factor $e^{-Ht}$: $$ \mathcal{L} = \frac{1}{2} m \left(\dot r^2 + r^2 \dot \theta^2 + r^2 \sin^2 \theta \, \dot \phi^2 \right) e^{2Ht}+ \frac{GMm}{r} e^{-Ht} $$ (I have not added in any derivatives of the scaling factor; this would mean I was just changing coordinates and not adding in an externally controlled expansion of space). Since the problem is spherically symmetric, the total angular momentum is conserved and the motion lies in a plane (the expansion of space does not change this) at $\theta = \pi/2$. This reduces the effective Lagrangian to $$ \mathcal{L} = \frac{1}{2} m \left( \dot r^2 + r^2 \dot \phi^2 \right) e^{2Ht} + \frac{GMm}{r} e^{-Ht} $$ Now let's change coordinates to $r = \eta e^{-Ht}$. This represents a point that is stationary with respect to the origin (i.e. the massive body that the smaller body orbits). Then $\dot r e^{Ht}= \dot \eta - H \eta $, so the Lagrangian becomes $$ \mathcal{L} = \frac{1}{2} m \left((\dot \eta - H \eta)^2 + \eta^2 \dot \phi^2 \right) + \frac{GMm}{\eta} $$ This looks like a Lagrangian that has been modified slightly by the parameter $H$. After these manipulations, $\phi$ is still a cyclic coordinate, and its conjugate momentum $p_\phi = m\eta^2 \dot \phi = L$ is still conserved. The momentum conjugate to $\eta$ is $$ p_{\eta} = m (\dot \eta - H \eta)$$ So, the equation of motion for $\eta$ is $$ m \frac{d}{dt} (\dot \eta - H \eta) = m (\ddot \eta - H \dot \eta) = - mH(\dot \eta - H \eta) + m \eta \dot \phi^2 - \frac{GMm}{\eta^2} $$ Canceling terms and substituting in the equation of motion for $\phi$, this becomes $$ \ddot \eta = H^2 \eta + \frac{L^2}{m^2 \eta^3} - \frac{GM}{\eta^2} $$ So, the expansion of space shows up as a force term proportional to $\eta$. For large $\eta_0$ at time $t = 0$, there is a solution $\eta(t) = \eta_0 e^{Ht}$, where $\dot \eta = H \eta$ agrees with Hubble's law (the corresponding comoving coordinate is $r(t) = \eta_0$). There is also an a highly unstable orbital solution for large $\eta$, where gravitational force just balances cosmological expansion. This new term also shifts the location of the "stable orbit" slightly (I haven't checked that it actually is stable in this situation). The location of the stable circular orbit when $H = 0$ is $$ \eta_0 = \frac{L^2}{GMm^2} $$ Then let $$ \eta = \eta_0 + \delta \eta = \eta_0 \left( 1 + \frac{\delta \eta}{\eta_0} \right)$$ To lowest order, the shift in location of the stable orbit is $$ \frac{\delta \eta}{\eta_0} = \frac{H^2}{GM/\eta_0^3 - H^2} $$ For the orbit of the earth around the sun, this shift would be entirely negligible- about 15 pm. Taking the mass to be the mass of the Milky way and the distance to the distance of the Sun from the center, the shift would be a bit larger - a fractional amount of about 2e-7, which is several hundred AU - but still largely small.
