5.5. Kozmologija¶

Donji kod je napravljen za Sage sustav. Možete ga izvršiti ako instalirate Sage ili online na sagemathcloud-u. Dobra je vježba konvertirati to u Jupyter/Python.

Nobelova nagrada za fiziku 2011. godine dodijeljena je Perlmutteru, Schmidtu i Riessu “za otkriće ubrzanog širenja svemira putem opažanja dalekih supernova”. U ovom odjeljku ćemo prikazati taj rezultat.

Udaljenost (luminosity distance) \(D_L\) (u megaparsecima Mpc) do supernove s prividnim sjajem \(m\) (mjeri se izravno) i apsolutnim sjajem \(M\) (poznat je za supernove tipa Ia) je dana formulom \(\mu = m-M = 5 \log_{10}D_L + 25\). Na internet adresi http://supernova.lbl.gov/Union/ nalaze se podaci o nekoliko stotina supernovih organiziranih u tablici sa stupcima:

ime supernove
\(z\) - crveni pomak
\(\mu\) - atenuacija sjaja supernove
\(\Delta\mu\) - eksperimentalna neodređenost (greška) od \(\mu\)

Učitavamo ih koristeći NumPy funkciju loadtxt() koja numeričke podatke organizirane u stupce obične datoteke pretvara u NumPy listu.

sage: import matplotlib.pyplot as plt
sage: import numpy as np
sage: import scipy
sage: import scipy.stats
sage: import urllib

sage: sns = urllib.urlopen(
....:  'http://supernova.lbl.gov/Union/figures/SCPUnion2_mu_vs_z.txt')
sage: dat = np.loadtxt(sns, usecols=(1,2,3)); dat.shape
(557, 3)

sage: dat[:2,:]
array([[  2.84880000e-02,   3.53355511e+01,   2.26143926e-01],
       [  5.00430000e-02,   3.66754416e+01,   1.67114252e-01]])

Pretvaramo sada mjerene podatke iz \(\mu\) i \(\Delta\mu\) u podatke za \(D_L\) i \(\Delta D_L\) pazeći na ispravnu propagaciju greške.

sage: var('mu DL')
(mu, DL)
sage: DL(mu) = DL.subs(solve(mu == 5*log(DL, base=10) + 25, DL)[0]); DL
mu |--> 1/100000*e^(1/5*mu*log(10))

sage: DLdat = np.array([(z, DL(valmu).n(), diff(DL)(valmu).n()*errmu)
....:                             for z,valmu,errmu in dat]); DLdat.shape
(557, 3)

Definiramo funkcije za prilagodbu, poznate iz odjeljka o prilagodbi funkcija, a dodajemo i dio za ocjenu dobrote prilagodba, vidi odjeljak o testiranju hipoteze.

sage: def dist(p, fun, data):
....:     return np.array([(y - fun(p, x))/err for (x,y,err) in data])
....:
sage: def chisq(p, fun, data):
....:     return sum( dist(p, fun, data)^2 )
....:
sage: def leastsqfit(fcn, init, data):
....:     """Fit function fcn, starting from init initial values to data."""
....:     p_final, cov_x, infodict, msg, ier = scipy.optimize.minpack.leastsq(
....:                              dist, init,  (fcn, data), full_output=True)
....:     parerrs = sqrt(scipy.diagonal(cov_x))
....:     for k in range(len(init)):
....:         print "p[%d] = %.3f +- %.3f" % (k, p_final[k], parerrs[k])
....:     chi2 = chisq(p_final, fcn, data)
....:     df = len(data)-len(init)
....:     print "chi^2/d.o.f = %.1f/%i  (p-value = %.4f)" % (chi2, df,
....:       scipy.stats.chisqprob(chi2, df))
....:     return list(p_final)

5.5.1. Određivanje \(H_0\) iz podataka o bliskim supernovama¶

Za male crvene pomake, ovisnost udaljenosti i crvenog pomaka dana je Hubbleovim zakonom

\[H_{0}D_{L} = c z\]

gdje je \(c\) brzina svjetlosti, a \(H_0\) Hubbleova konstanta. Određujemo Hubbleovu konstantu (u jedinicama km/s/Mpc) prilagodbom Hubbleovog zakona na podskup mjerenja za koja je \(z<0.12\).

sage: lowz = np.array([pt for pt in DLdat if pt[0]<0.12]); lowz.shape
(174, 3)

sage: c=3.e5
sage: def linfun(p, z):
....:     return c*z/p[0]

sage: H0, = leastsqfit(linfun, (60,), lowz)
p[0] = 68.012 +- 0.385
chi^2/d.o.f = 179.1/173  (p-value = 0.3585)
sage: print "\nHubbleova constanta H0 = %.1f  km/s/Mpc" % H0
<BLANKLINE>
Hubbleova constanta H0 = 68.0  km/s/Mpc

sage: fig, ax = plt.subplots(figsize=[5,3])
sage: ax.errorbar(lowz[:,0], lowz[:,1], yerr=lowz[:,2], marker='.',
....:                                               linestyle='None')
sage: zvals = np.linspace(0.01, 0.12)
sage: ax.plot(zvals, [linfun([H0], z) for z in zvals],'r-')
sage: ax.set_xlabel('$z$')
sage: ax.set_ylabel('$D_L \; {\\rm [Mpc]}$')
sage: fig.tight_layout()
sage: fig.savefig('fig')

5.5.2. Određivanje ubrzanja ekspanzije svemira¶

Za veće crvene pomake, Hubbleov zakon se modificira i u prvoj aproksimaciji se može pisati kao

\[H_{0}D_{L} = c z \big(1 +(1 - q_0)\frac{z}{2} \big)\]

gdje je \(q_0\) parametar deceleracije. Definiran je tako da pozitivni \(q_0\) odgovara svemiru čija se brzina širenja usporava (kako se očekivalo prije ovih mjerenja). Određujemo taj parametar prilagodbom ove formule na cijeli skup mjerenja.

sage: def qfun(p, z):
....:     return c*z/p[0] * (1  +  0.5*(1 - p[1])*z)

sage: qpars = leastsqfit(qfun, (50, 0.1), DLdat)
p[0] = 69.184 +- 0.368
p[1] = -0.153 +- 0.039
chi^2/d.o.f = 553.7/555  (p-value = 0.5081)

sage: fig, ax = plt.subplots(figsize=[5,3])
sage: ax.errorbar(DLdat[:,0], DLdat[:,1]/1e3, yerr=DLdat[:,2]/1e3,
....:                                     marker='.', linestyle='None')
sage: zvals = np.linspace(0.01, 1.8)
sage: ax.plot(zvals, qfun(qpars, zvals)/1e3,'r')
sage: ax.set_xlabel('$z$')
sage: ax.set_ylabel('$D_L \; {\\rm [Gpc]}$')
sage: fig.tight_layout()
sage: fig.savefig('fig')

Iako iz mjernih podataka nije očito da postoji odstupanje od linearnog Hubbleovog zakona, podataka ima dovoljno da statistička analiza nedvosmisleno pokazuje da je parametar deceleracije negativan i različit od nule:

\[q_0 = -0.15 \pm 0.04\]

5.5.3. Određivanje gustoće materije i tamne energije¶

Gornja formula ne vrijedi za \(z\approx 1\), dakle ovakav račun nije sasvim korektan i može se rabiti samo u pedagoške svrhe. Ispravan pristup zahtijeva integraciju propagacije zrake svjetlosti od supernove do promatrača kroz svemir u širenju. Ta propagacija je određena sastavom svemira (obična i tamna materija različito djeluju od tzv. tamne energije koja uzrokuje ubrzanje širenja svemira). Hubbleov zakon postaje

\[H_{0}D_{L} = c (1+z) Z(z, \Omega_M, \Omega_\Lambda)\]

\[Z(z, \Omega_M, \Omega_\Lambda) = \int_{1/(1+z)}^{1}\frac{da}{a\sqrt{X(a, \Omega_M, \Omega_\Lambda)}}\]

\[X(a, \Omega_M, \Omega_\Lambda) = \frac{\Omega_M}{a} + \Omega_\Lambda a^2\]

\[\Omega_M + \Omega_\Lambda = 1\]

gdje je \(\Omega_M\) udio materije (obične i tamne), a \(\Omega_\Lambda\) udio tamne energije u ukupnoj energiji svemira.Još općenitije formule, koje ne rade pretpostavku \(\Omega_M + \Omega_\Lambda = 1\), mogu se naći na stranicama Neda Wrighta.

Provjeravamo egzaktne formule razvojem u red po z:

sage: var('z a Om Ov Z X q0'); assume(z>0, a>0, Om>0, Ov>0);
(z, a, Om, Ov, Z, X, q0)
sage: X(a, Om, Ov) = Om/a + Ov*a**2 + (1-Om-Ov)
sage: DLH0 = taylor(((1+z)^2 * (Z/(1+z)) * (1 +(1-Om-Ov)*Z^2/6)).subs(
....:           Z==taylor(integral(1/(a*sqrt(X(a, Om, Ov))),
....:                           a, 1/(1+z), 1), z, 0, 2)), z, 0, 2)
sage: solve(DLH0 == z*(1  +  (1 - q0)*z/2), q0)[0]
q0 == 1/2*Om - Ov

Ova korespondencija s parametrom deceleracije se slaže s (Peacock 3.34).

sage: def J(x):
....:     if x<0:
....:         return sin(sqrt(-x))/sqrt(-x)
....:     elif x==0:
....:         return 1
....:     else:
....:         return sinh(sqrt(x))/sqrt(x)

sage: def Xfun(a, Om, Ov):
....:     return Om/a + Ov*a**2 + (1-Om-Ov)

sage: def Zfun(z, Om, Ov):
....:     return numerical_integral(lambda a: 1/(a*sqrt(X(a, Om, Ov))),
....:                                                     1/(1+z), 1)[0]

sage: def DLfun(p, z):
....:     """DL(z) with pars p=[H0, Omega_m, Omega_Lambda]."""
....:     H0 = p[0]
....:     Om = p[1]
....:     Ov = p[2]
....:     Otot = Om + Ov
....:     Zval = Zfun(z, Om, Ov)
....:     return c/H0 * (1+z) * Zval * J((1-Otot)*Zval**2)

sage: def DLfunflat(p, z):
....:     """DL(z) with pars p=[H0, Omega_Lambda]. Flat universe."""
....:     H0 = p[0]
....:     Ov = p[1]
....:     Om = 1-Ov
....:     Zval = Zfun(z, Om, Ov)
....:     return c/H0 * (1+z) * Zval

Prilagodbom parametara \(H_0\), \(\Omega_M\) i \(\Omega_\Lambda\) određujemo udjele pojedinih komponenti svemira

sage: fullpars = leastsqfit(DLfun, (70, 0.3, 0.), DLdat) # long time
p[0] = 70.717 +- 0.466
p[1] = 0.368 +- 0.075
p[2] = 0.834 +- 0.122
chi^2/d.o.f = 527.7/554  (p-value = 0.7833)

Ukoliko pak fiksiramo ravnu geometriju ( \(\Omega_M + \Omega_\Lambda = 1\)) dobivamo uobičajenu jednu trećinu materije i dvije trečine tamne energije:

sage: flatpars = leastsqfit(DLfunflat, (70, 0.3), DLdat) # long time
sage: print "Om = 1-Ov = %8.3f" % (1-flatpars[1],)
p[0] = 70.412 +- 0.363
p[1] = 0.708 +- 0.021
chi^2/d.o.f = 528.7/555  (p-value = 0.7827)
Om = 1-Ov =    0.292

Dakle,

\[\begin{split}\Omega_M = 0.29 \\ \Omega_\Lambda = 0.71\end{split}\]

Crtamo ovaj rezultat (crvena linija) zajedno s nekim drugim karakterističnim scenarijima.

sage: fig, ax = plt.subplots(figsize=[8,6])
sage: ax.errorbar(DLdat[:,0], DLdat[:,1], yerr=DLdat[:,2], marker='.',
....:             linestyle='None')
sage: zvals = np.linspace(0.01, 1.8)
sage: ax.plot(zvals, [DLfunflat(flatpars, z) for z in zvals],'r-',
....:                  label='$\\Omega_m=0.29, \\Omega_\\Lambda=0.71$')
sage: ax.plot(zvals, [DLfun([H0, 0.5, 0.5], z) for z in zvals],'b:',
....:                     label='$\\Omega_m=0.5, \\Omega_\\Lambda=0.5$')
sage: ax.plot(zvals, [DLfun([H0, 1, 0], z) for z in zvals],'k-.',
....:                         label='$\\Omega_m=1, \\Omega_\\Lambda=0$')
sage: ax.plot(zvals, [DLfun([H0, 0, 1], z) for z in zvals],'g--.',
....:                         label='$\\Omega_m=0, \\Omega_\\Lambda=1$')
sage: ax.plot(zvals, [DLfun([H0, 0.3, 0.], z) for z in zvals],'k--',
....:                       label='$\\Omega_m=0.3, \\Omega_\\Lambda=0$')
sage: ax.set_ylim(0, 1.5e4)
sage: ax.set_xlabel('$z$')
sage: ax.set_ylabel('$D_L \; {\\rm [Mpc]}$')
sage: ax.legend(loc='upper left').draw_frame(0)
sage: fig.tight_layout()
sage: fig.savefig('fig')

5.5. Kozmologija¶

5.5.1. Određivanje \(H_0\) iz podataka o bliskim supernovama¶

5.5.2. Određivanje ubrzanja ekspanzije svemira¶

5.5.3. Određivanje gustoće materije i tamne energije¶

Pregled sadržaja

Prijašnja tema

Slijedeća tema