I A glitch in Jorrie’s Cosmo-Calculator?

[pbuk]

Hi @Jorrie, @pbuk,

When trying LightCone8, I noticed that if z_lower = 0, a duplicated line is printed at z = 0.

Thanks for that, it is an unintended side effect of fixing the range so it includes z = 0 instead of z = 1e-6 for instance. I'll fix the fix.

 

[JimJCW]

Hi @Jorrie, @pbuk,

Let’s review the current situation of Jorrie’s calculator. First, LightCone8 updated by @pbuk has eliminated the glitch discussed in this thread (see Post #1):

ICRAR is an online calculator that includes Ω_R,0 as an input parameter, just like LightCone8 (z_eq = 3370). A comparison of the calculated D_now’s using the two models is shown below:

The two results are almost indistinguishable. However, when their difference is plotted against z, the result shows,

The fact that ‘ICRAR - LightCone8’ is not a horizontal line near y = 0 suggests that there exists discrepancy between the two calculators. To look into this further, we note that ICRAR (Ω_R=0), LightCone8 (Ω_R=0), and cosmic (no Ω_R) should give the same results, but they don’t:

The figure suggests that when Ω_R=0, ICRAR and cosmic are consistent with each other, but not with LightCone8. The source of this discrepancy is unknown. It would be nice, though, to be able to identify it.

To summarize:

LightCone8 has removed the glitch discussed in this thread. However, it is noticed that there exist discrepancies between LightCone8 and some other online calculators, such as ICRAR and cosmic.​

 

[pbuk]

The figure suggests that when Ω_R=0, ICRAR and cosmic are consistent with each other, but not with LightCone8. The source of this discrepancy is unknown. It would be nice, though, to be able to identify it.

My first suspicion is the values of certain constants, particularly the conversion factor from $km\\s^{-1}\\Mpsc^{-1}$ to $Gyr^{-1}$ which is a simple $\frac 1 {978}$ in LightCone8 instead of something like $\frac{365.2421897 \times 86400 \times 1000000 \pi}{96939420213600000}$.

I've raised an issue and will look at it: https://github.com/light-cone-calc/light-cone-calc/issues/11

 

[Jorrie]

I have commented on the Github issues page.

 

[JimJCW]

Hi @Jorrie and @pbuk,

I noticed some discrepancies in results from LightCone8. Let’s use the following input values:

Equations in Tutorial, part III,

Ω_m = (Ω - Ω_Λ) S_eq / (1 + S_eq), where S_eq = z_eq + 1​

Ω_R = Ω_m / S_eq​

give,

Ω_m = 0.308908363
Ω_R = 0.000091637

These numbers are compared with output of LightCone8 below:

​	z	OmegaM​	OmegaL​	OmegaR​	OmegaT​
LightCone8​	0​	0.309000000​	6.910000000	0.00009169139​	1.000091691​
From Equations​	0​	0.308908363​	6.91​	0.000091637​	1​

Note that under Conversions in the output, LightCone8 gives,

Conclusion: The values of Ω_m and Ω_T calculated with LightCone8 are questionable.

 

[pbuk]

I noticed some discrepancies in results from LightCone8. Let’s use the following input values:
...
Conclusion: The values of Ω_m and Ω_T calculated with LightCone8 are questionable.

Thanks, I think this is related to https://github.com/light-cone-calc/light-cone-calc/issues/6
but that should have been fixed by https://github.com/light-cone-calc/light-cone-calc/commit/88c1c7324e53476bcc81fd49b23a4280fccb8755
More investigation required, probably over the weekend.

 

[JimJCW]

Hi @Jorrie and @pbuk,

Would you please let me know where in Jorrie 2021 calculator the circled number shown below is calculated and printed, for example, lines numbers in Calculation.js. I am trying to understand the discrepancy discussed in Post #105.

 

[Jorrie]

It is calculated in the main html program, not in Calculate.js

line 356 var OmegaM = (Omega-Ynow * Ynow / Yinf / Yinf) * s_eq / (1+s_eq);

I'm not sure if this should have been "z_eq / (1+z_eq) "
Also take note that it is expressly rounded to 4 decimals for that calculation.
When calculated later in the tables, it is rounded to the requested number pf decimals.

Ps: The densities for the tables are calculated in
Calculation.js line 27: var Om = (Omega - Ol) * (s_eq + 1) / (s_eq + 2);
The (s_eq + 1) / (s_eq + 2) is a mistake that has been prviously pinted out and corrected in Lightcone8 by @pbuk. There is some further processing in lines 243 to 250.
These line numbers are only relevant to the LightCone7-2021-03-12 Version

 

[JimJCW]

It is calculated in the main html program, not in Calculate.js . . .

Thanks! This is helpful. The situation is that using the same input, two different values of Ω_m,0 are obtained from the calculator, 0.3089 from the main html program and 0.3090 from Calculation.js:

This suggests that there is an error in the calculator program. @pbuk

 

[pbuk]

This suggests that there is an error in the calculator program. @pbuk

I think the error is in the UI: we should have $\Omega_{m,0} = \Omega - \Omega_{\Lambda,0} = 1 - 0.691 = 0.309$.

 

[JimJCW]

I think the error is in the UI: we should have $\Omega_{m,0} = \Omega - \Omega_{\Lambda,0} = 1 - 0.691 = 0.309$.

In Calculation.js:

var Ok = 1 - Om - Or - Ol; // curvature density par

 

[Jorrie]

I think the error is in the UI: we should have $\Omega_{m,0} = \Omega - \Omega_{\Lambda,0} = 1 - 0.691 = 0.309$.

I have checked the UI calculation and it seems to be right if radiation density is included, as it should.
$\Omega_{m,0} = 1-\Omega_{\Lambda,0} - \Omega_{r,0}= 1 - 0,691 - 0,000091637 = 0,308908363$
I obtained the value of $\Omega_{r,0}$ by substituting $\Omega_{m,0} = S_{eq} \Omega_{r,0}$ into the balancing equation.
On the other hand, if I set the columns for enough decimal places in Ligtcone8, $\Omega_{T,0}$ does not equal 1.

So the error seems to lie in calculation.js for the columns.

 

[pbuk]

I have checked the UI calculation and it seems to be right if radiation density is included, as it should.
$\Omega_{m,0} = 1-\Omega_{\Lambda,0} - \Omega_{r,0}= 1 - 0,691 - 0,000091637 = 0,308908363$
I obtained the value of $\Omega_{r,0}$ by substituting $\Omega_{m,0} = S_{eq} \Omega_{r,0}$ into the balancing equation.
On the other hand, if I set the columns for enough decimal places in Ligtcone8, $\Omega_{T,0}$ does not equal 1.

View attachment 303936
So the error seems to lie in calculation.js for the columns.

It's the same in LightCone7 (bar imprecision)

this isn't surprising as I've used essentially the same calculation as LightCone7

          rholNow = rhocritNow * Ol;
          rhol = rholNow;
          rhomNow = rhoNow - rhol;
          rhom = rhomNow *s*s*s;
          rhorNow = rhomNow/s_eq;
          rhor = rhorNow *s*s*s*s;

i.e. $\rho_{m,0} = \rho_{0} - \rho_{\Lambda,0}$. I guess you are saying this should be $\rho_{m,0} = \rho_{0} - \rho_{\Lambda,0} - \rho_{R,0}$?

There's a few other oddities I'm coming across, I think I'm going to end up rewriting the physics from scratch (which means I'm going to have to learn it first!)

 

[George Jones]

I guess you are saying this should be $\rho_{m,0} = \rho_{0} - \rho_{\Lambda,0} - \rho_{R,0}$?

Yes, at any time, $\rho_{0} = \rho_{r,0} + \rho_{m,0} + \rho_{\Lambda,0}$. In the case $\Omega = 1$, this is subject to the constraint$$\rho_{0} = \rho_\mathrm{critical} = \frac{3 H_0^2 }{8 \pi G},$$
since $\Omega := \rho/\rho_\mathrm{critical}$.

 

[Jorrie]

There's a few other oddities I'm coming across, I think I'm going to end up rewriting the physics from scratch (which means I'm going to have to learn it first!)

Yes, the same differences occur in Lightcone7, but it is calculated differently than the value in the UI. I will check that part and make it compatible, because I still feel that the value in the UI is correct.

 

[Jorrie]

Yes, the same differences occur in Lightcone7, but it is calculated differently than the value in the UI. I will check that part and make it compatible, because I still feel that the value in the UI is correct.

Minor corrections (red) in Calculation.js required to improve the issue:

var Ol = (Ynow / Yinf) * (Ynow / Yinf); // Lambda density parameter
var Om = (Omega - Ol) * s_eq / (s_eq + 1); // matter density parameter (corrected for radiation)
var Or = Om / s_eq; // Radiation density parameter (corrected for error; was /(s_eq+1)
var Ok = 1 - Om - Or - Ol; // curvature density par

rholNow = rhocritNow * Ol;
rhol = rholNow; // vacuum energy density remains constant
rhomNow = (rhoNow - rhol)*s_eq/(s_eq+1); // Matter energy density corrected for radiation
rhom = rhomNow *s*s*s; // Matter density at time T
rhorNow = rhomNow/s_eq; //radiation at present
rhor = rhorNow *s*s*s*s; // radiation density at time T
rhoT = rhol+rhom+rhor; // total energy density at time T
rhocrit = rhoConst*(H_t/secInGy)*(H_t/secInGy); // critical density at time T
OmegaMatterT = rhom/rhocrit; // Omegas at time T
OmegaLambdaT = rhol/rhocrit;
OmegaRadiationT = rhor/rhocrit;
OmegaTotalT = rhoT/rhocrit;
In Lightconne 7, it gives the correct matter density of 0.308908 at z=0 and a more accurate OmegaT. The drift in OmegaT may be rounding off in Lightcone 7. It will be interesting to see what it does in Lightcone8 with its improved calculation module.

 

[Jorrie]

Yes, at any time, $\rho_{0} = \rho_{r,0} + \rho_{m,0} + \rho_{\Lambda,0}$. In the case $\Omega = 1$, this is subject to the constraint$$\rho_{0} = \rho_\mathrm{critical} = \frac{3 H_0^2 }{8 \pi G},$$
since $\Omega := \rho/\rho_\mathrm{critical}$.

George, how do the equations change in the case where $\Omega \ne 1$?
I'm asking because I don't think Lightcone7 gives the correct results in such a case.

 

[George Jones]

George, how do the equations change in the case where $\Omega \ne 1$?
I'm asking because I don't think Lightcone7 gives the correct results in such a case.

I am not sure which equations you mean. $\Omega = \Omega_r + \Omega_m + \Omega_\Lambda$ is always true, as is $\Omega + \Omega_k = 1$. If $\Omega = 1$, then $\Omega_k =0$, but if f $\Omega \neq 1$, then $\Omega_k \neq 0$.

Did you have other equations in mind?

 

[pbuk]

I am not sure which equations you mean. $\Omega = \Omega_r + \Omega_m + \Omega_\Lambda$ is always true, as is $\Omega + \Omega_k = 1$. If $\Omega = 1$, then $\Omega_k =0$, but if f $\Omega \neq 1$, then $\Omega_k \neq 0$.

Did you have other equations in mind?

I suspect @Jorrie is suffering from the same problem as I am with the confusing terminology. I am working on increasing consistency for my rewrite:

Parameters at any particular time: we calculate these.

• $\Omega_{m}$ omegaM - density of mass as a fraction of the critical density.
• $\Omega_{\Lambda}$ omegaLambda - density of the cosmological constant as a fraction of the critical density.
• $\Omega_{rad}$ omegaRad - density of radiation as a fraction of the critical density.
• $\Omega$ omega - total density of [stuff] as a fraction of the critical density $\Omega = \Omega_{m} + \Omega_{\Lambda} + \Omega_{rad}$.
• $\Omega_{K}$ omegaK - Do we bother with this?

Parameters now: we measure these.

• $\Omega_{[m,0]}$ omegaM0 - density of mass now as a fraction of the critical density now.
• $\Omega_{[\Lambda,0]}$ omegaLambda0 - density of the cosmological constant now as a fraction of the critical density now.
• $\Omega_{[rad,0]}$ omegaRad0 - density of radiation now as a fraction of the critical density now.
• $\Omega_{0}$ omega0 - total density of [stuff] now as a fraction of the critical density now $\Omega_{0} = \Omega_{[m,0]} + \Omega_{[\Lambda,0]} + \Omega_{[rad,0]}$.
• $\Omega_{[K,0]}$ omegaK0 - curvature now as a fraction of the critical density now (is this actually constant in all models?) $\Omega_{[K,0]} = 1 - \Omega_{0}$.

Now we have a very important equation from a computational point of view:
$$ [something] = \Omega_{[m,0]} s^3 + \Omega_{[rad,0]} s^4 + \Omega_{[\Lambda,0]} s^{3(1+w)} + \Omega_{[K,0]} s^3 $$
Yes I know that $[something] = \frac{H^2}{H_0^2}$, but we need a name for this: density parameter?

 

[George Jones]

Now we have a very important equation from a computational point of view:
$$ [something] = \Omega_{[m,0]} s^3 + \Omega_{[rad,0]} s^4 + \Omega_{[\Lambda,0]} s^{3(1+w)} + \Omega_{[K,0]} s^3 $$
Yes I know that $[something] = \frac{H^2}{H_0^2}$, but we need a name for this: density parameter?

Common practice in cosmology is to write "something" in terms of redshift $z$, i.e., to substitute $s = 1 + z$ for every $s$. Common cosmological notation is $[something] = \frac{H\left(z\right)^2}{H_0^2} =: E\left(z\right)^2$. I am not sure why.

 

[pbuk]

Common cosmological notation is $[something] = \frac{H\left(z\right)^2}{H_0^2} =: E\left(z\right)^2$.

Excellent, thank you, so adapting that slightly I can write:
$$ E(s)^2 = \Omega_{[m,0]} s^3 + \Omega_{[rad,0]} s^4 + \Omega_{[\Lambda,0]} s^{3(1+w)} + \Omega_{[K,0]} s^2 $$
[The last term in the equation above was corrected from $\Omega_{[K,0]} s^3$.]

[code lang=javascript title="Which I can implement for w = -1 as"]
const getESquaredAtS = (s) => {
const s2 = s * s;
return omegaM0 * s2 * s + omegaLambda0 + omegaRad0 * s2 * s2 + OmegaK * s2;
};
[/code]
.. and that saves a LOT of calculating square roots and then squaring, converting to different units and then back again.

 

[Jorrie]

Excellent, thank you, so adapting that slightly I can write:
$$ E(s)^2 = \Omega_{[m,0]} s^3 + \Omega_{[rad,0]} s^4 + \Omega_{[\Lambda,0]} s^{3(1+w)} + \Omega_{[K,0]} s^3 $$

I think the last term should be $\Omega_{[K,0]} s^2$, because curvature dilates by an inverse square law of the expansion factor.
I see in the code piece you have it right...

 

[Jorrie]

I am not sure which equations you mean. $\Omega = \Omega_r + \Omega_m + \Omega_\Lambda$ is always true, as is $\Omega + \Omega_k = 1$. If $\Omega = 1$, then $\Omega_k =0$, but if f $\Omega \neq 1$, then $\Omega_k \neq 0$.

Did you have other equations in mind?

No, I just wanted confirmation of my intuition, but I still need to think how to implement it in code.
Perhaps @pbuk has an implementation in mind.

 

[pbuk]

I see in the code piece you have it right...

Yes, $\LaTeX$ copy and paste error, thanks.

 

[George Jones]

No, I just wanted confirmation of my intuition, but I still need to think how to implement it in code.
Perhaps @pbuk has an implementation in mind.

I might be able to help, but I do not understand the notation. I have figured out "s" for "stretch". To start, I would like to ask about stuff in the following post.

It is calculated in the main html program, not in Calculate.js

line 356 var OmegaM = (Omega-Ynow * Ynow / Yinf / Yinf) * s_eq / (1+s_eq);

I'm not sure if this should have been "z_eq / (1+z_eq) "
Also take note that it is expressly rounded to 4 decimals for that calculation.
When calculated later in the tables, it is rounded to the requested number pf decimals.

Ps: The densities for the tables are calculated in
Calculation.js line 27: var Om = (Omega - Ol) * (s_eq + 1) / (s_eq + 2);
The (s_eq + 1) / (s_eq + 2) is a mistake that has been prviously pinted out and corrected in Lightcone8 by @pbuk. There is some further processing in lines 243 to 250.
These line numbers are only relevant to the LightCone7-2021-03-12 Version

Does OmegaM represent $\Omega_{m,0}$ or $\Omega_m \left(t\right)$?

What are Ynow and Yinf?

 

[pbuk]

Perhaps @pbuk has an implementation in mind.

I do: by eliminating unnecessary conversions I have got the calculations down to this:
[code lang=javascript]
const getParamsAtStretch = (s: number): LcdmModelVariables => {
const eSquared = getESquaredAtStretch(s);
const s2 = s * s;
const h = h0 * Math.sqrt(eSquared);
const omegaM = (omegaM0 * s2 * s) / eSquared;
const omegaLambda = omegaLambda0 / eSquared;
const omegaRad = (omegaRad0 * s2 * s2) / eSquared;
return {
h,
omegaM,
omegaLambda,
omegaRad,
temperature: cmbTemperature * s,
rhoCrit: rhoCrit0 * eSquared,
};
};
[/code]
... and this:
[code lang=javascript]
for (let i = integrationResults.length - 1; i >= 0; --i) {
const { s, t, dNow: dUnsafe, dPar } = integrationResults;

const params = model.getParamsAtStretch(s);
const hGy = params.h * model.kmsmpscToGyr;

// Force dNow to zero at zero redshift.
const dNow = s === 1 ? 0 : dUnsafe;

results.push({
z: s - 1,
a: 1 / s,
s,
t,
dNow,
d: dNow / s,
r: 1 / hGy,
dPar,
vGen: params.h / (s * model.h0),
vNow: dNow * model.h0Gy,
v: (dNow * hGy) / s,
...params,
});
}
[/code]

I only have one variable I don't quite understand - I've called it vGen, it's called XDpar in the LightCone7 code but displayed in the table as $V_{gen} / c$ and has the description

Is there anything snappier/more descriptive than vGen?

 

[pbuk]

Does OmegaM represent $\Omega_{m,0}$ or $\Omega_m \left(t\right)$?

What are Ynow and Yinf?

Yes, good questions, but the work is done moving away from that code now as in my last post. You can see the work in progress at https://github.com/light-cone-calc/light-cone-calc/blob/refactor/src/model.ts.

My next focus is on the calculation of the current density parameters, here is what we have now ($\Omega, \Omega_{[\Lambda,0]}$ and $s_{eq} = z_{eq} + 1$ are provided as parameters e.g. from Planck 2015).

  const rhoCrit0 = rhoConst * h0Seconds * h0Seconds;
  const omegaM0 = ((omega - omegaLambda0) * seq) / (seq + 1);
  const omegaRad0 = omegaM0 / seq;
  const OmegaK0 = 1 - omegaM0 - omegaRad0 - omegaLambda0;

 

[George Jones]

My next focus is on the calculation of the current density parameters, here is what we have now ($\Omega, \Omega_{[\Lambda,0]}$ and $s_{eq} = z_{eq} + 1$ are provided as parameters e.g. from Planck 2015).

$\Omega$ or $\Omega_0$? If $\Omega_0$, then this is equivalent to taking $\Omega_{r,0}$, $\Omega_{m,0}$, and $\Omega_{\Lambda,0}$ as input parameters.

  const rhoCrit0 = rhoConst * h0Seconds * h0Seconds;
  const omegaM0 = ((omega - omegaLambda0) * seq) / (seq + 1);
  const omegaRad0 = omegaM0 / seq;
  const OmegaK0 = 1 - omegaM0 - omegaRad0 - omegaLambda0;

At the risk of outlining what folks already know ... radiation-matter equilibrium happens when $\rho_r = \rho_m$, or, equivalently, when
$$\begin{align}
\Omega_{r,eq} &= \Omega_{m,eq} \\
\Omega_{r,0} s_{eq}^4 &= \Omega_{m,0} s_{eq}^3 \\
\Omega_{r,0} &= \Omega_{m,0} / s_{eq},
\end{align}$$
which is the second-last line in the above code.

To get the second line in the above code, use $\Omega = \Omega_r + \Omega_m + \Omega_\Lambda$ to substitute for $\Omega_{r,0}$.
$$\begin{align}
\Omega_{r,0} &= \Omega_{m,0} / s_{eq} \\
\Omega_0 - \Omega_{m,0} - \Omega_{\Lambda,0} &= \Omega_{m,0} / s_{eq}. \\
\end{align}$$

Solving for $\Omega_{m,0}$ gives the second line of the code.

 

[Jorrie]

Does OmegaM represent $\Omega_{m,0}$ or $\Omega_m \left(t\right)$?

What are Ynow and Yinf?

OmegaM represents $\Omega_{m,0}$. Later on OmegaMatterT represents $\Omega_m \left(t\right)$

Ynow and Yinf are inverse Hubble values, i.e. Hubble times now and at t -> infinity
These still come from the late Marcus's inputs.

 

[Jorrie]

I only have one variable I don't quite understand - I've called it vGen, it's called XDpar in the LightCone7 code but displayed in the table as $V_{gen} / c$ and has the description
View attachment 304033
Is there anything snappier/more descriptive than vGen?

AFAICR, the late Marcus and I were thinking of the worldlines of generic particles that presently cross our Hubble sphere. The variable XDpar appeared in the discussion as the ratio Vthen/Vnow and it was specifically included in the columns for purposes of plotting that interesting wordline mentioned. This is how Vgen[eric] came to be. It shows "dramatically" how a particle could enter the Hubble sphere early on and then be accelerated to cross it again, in the generic case, right now.
We used it to work out the crossing and inflection points.

 

[pbuk]

$\Omega$ or $\Omega_0$? If $\Omega_0$, then this is equivalent to taking $\Omega_{r,0}$, $\Omega_{m,0}$, and $\Omega_{\Lambda,0}$ as input parameters.

Good catch - yes this should be $\Omega_0$, I will change the variable name.

At the risk of outlining what folks already know
...
Solving for $\Omega_{m,0}$ gives the second line of the code.

Many thanks, it's good to have it confirmed.

 

[George Jones]

OmegaM represents $\Omega_{m,0}$. Later on OmegaMatterT represents $\Omega_m \left(t\right)$

Ynow and Yinf are inverse Hubble values, i.e. Hubble times now and at t -> infinity
These still come from the late Marcus's inputs.

Then, in

line 356 var OmegaM = (Omega-Ynow * Ynow / Yinf / Yinf) * s_eq / (1+s_eq);

$$\frac{Y_{now}^2}{Y_\infty ^2} = \frac{H_\infty ^2}{H_0^2} = \Omega_{r,0} s^4 + \Omega_{m,0} s^3 + \Omega_{\Lambda,0} + \Omega_{k,0} s^2$$
when $w = -1$ for dark energy.

If $R$ is used for the universe's scale factor, then $s := R_0 /R$. In universes close to our own, $R \rightarrow \infty$ as $t \rightarrow \infty$. Consequently, $s \rightarrow 0$ and $H\infty ^2 /H_0^2 \rightarrow \Omega_{\Lambda,0}$ as $t \rightarrow \infty$.

Thus, in universes close to our own,

line 356 var OmegaM = (Omega-Ynow * Ynow / Yinf / Yinf) * s_eq / (1+s_eq);

is equivalent to

const omegaM0 = ((omega - omegaLambda0) * seq) / (seq + 1);

The version with $\Omega_{\Lambda,0}$ should also work in universes that are not close to our own.

 

[pbuk]

I think I'm all done now, I've updated LightCone8 at https://light-cone-calc.github.io/, everything should work including unit conversions and all. One final thing: LightCone7 uses a CMB temperature of 2.75 K, I have been using 2.7548 K in LightCone8, hope that is OK.

The recession rate chart is exactly the same as LightCone7.

 

[pbuk]

The version with $\Omega_{\Lambda,0}$ should also work in universes that are not close to our own.

That's good as that's what I've used and the UI allows entering the parameters to suit different universes!

 

[Jorrie]

think I'm all done now, I've updated LightCone8 at https://light-cone-calc.github.io/, everything should work including unit conversions and all. One final thing: LightCone7 uses a CMB temperature of 2.75 K, I have been using 2.7548 K in LightCone8, hope that is OK.

Great Job! Now we have the UI and the tables agreeing on $\Omega_m$ and also $\Omega$ staying at precisely 1. 00... I have also checked that for $\Omega \ne 1$ the behavior seems to be as expected (and interesting!).

$\Omega$ had to start very much closer to 1 and diverged to the requested 1.1 at present, but then the cosmological constant will drive it back towards 1 (for $w=-1$ at least)

 

[JimJCW]

Minor corrections (red) in Calculation.js required to improve the issue:

var Ol = (Ynow / Yinf) * (Ynow / Yinf); // Lambda density parameter
var Om = (Omega - Ol) * s_eq / (s_eq + 1); // matter density parameter (corrected for radiation)
var Or = Om / s_eq; // Radiation density parameter (corrected for error; was /(s_eq+1)
var Ok = 1 - Om - Or - Ol; // curvature density par

rholNow = rhocritNow * Ol;
rhol = rholNow; // vacuum energy density remains constant
rhomNow = (rhoNow - rhol)*s_eq/(s_eq+1); // Matter energy density corrected for radiation
rhom = rhomNow *s*s*s; // Matter density at time T

I think these changes are incorrect. The following statements near Line 27 and Line 325 in Calculation.js,

correctly represent the following equations (note that s_eq = z_eq in LightCone7),

Ω_m = (Ω – Ω_Λ) (z_eq + 1) / (z_eq + 2)​

Ω_R = Ω_m / (z_eq + 1) .​

They are consistent with equations in Tutorial, part III:

Ω_m = (Ω - Ω_Λ) S_eq / (1 + S_eq)​

Ω_R = Ω_m / S_eq ,​

where S_eq = z_eq + 1.

However, the following statements near Line 261,

should be changed to,

rhomNow = (rhoNow - rhol) * (s_eq + 1) / (s_eq + 2);​

rhom = rhomNow *s*s*s;​

rhorNow = rhomNow / (s_eq + 1);​

I believe this is the source of the error discussed in Post #109. @pbuk

 

[pbuk]

I think these changes are incorrect.
...
should be changed to,

rhomNow = (rhoNow - rhol) * (s_eq + 1) / (s_eq + 2);​

rhom = rhomNow *s*s*s;​

rhorNow = rhomNow / (s_eq + 1);​

I believe this is the source of the error discussed in Post #109. @pbuk

Yes you are probably right, however LightCone8 now calculates everything differently as discussed above. By going back to square one on the Physics, aided by @George Jones cross-check to the equations used, I have eliminated the intermediate calculations involving the various densities as below:

https://github.com/cosmic-expansion...bc785935a20b65cda6beb725897/src/model.ts#L201

      const eSquared = getESquaredAtStretch(s);
      const s2 = s * s;
      const h = h0 * Math.sqrt(eSquared);
      const omegaM = (omegaM0 * s2 * s) / eSquared;
      const omegaLambda = omegaLambda0 / eSquared;
      const omegaRad = (omegaRad0 * s2 * s2) / eSquared;
      return {
        h,
        omegaM,
        omegaLambda,
        omegaRad,
        omega: omegaM + omegaLambda + omegaRad,
        temperature: temperature0 * s,
        rhoCrit: rhoCrit0 * eSquared,
      };

 

[George Jones]

aided by @George Jones cross-check to the equations used

@Jorrie and @pbuk have done tremendous work on LightCone; I haven't done that much. All the stuff I have done in this thread has been performed at coffee shops while "on vacation" at my in-laws' place (11 loud people in the house), and while 3300km from my books. Flying home today.

Sent from my Tim Hortons at Major MacK and Jane.

 

[pbuk]

@Jorrie and @pbuk have done tremendous work on LightCone; I haven't done that much.

That's exactly why your contribution has been valuable - it is always important to have eyes on the detail from outside the core team.

 

[Jorrie]

I have updated the Lightcone Tutorial sections to reflect the present status. Will the other members of the core team please comment, specifically on what we should have in Part III, the calculation section.

 

[JimJCW]

Some observations:

(1) The suggested correction to LightCone7 discussed in Post #136,

where, (s_eq = z_eq), is consistent with the following statements in LightCone8:

(2) It is noticed that whenever z_lower in LightCone8 is set to be 0, there is extra line printed. For example,

This is a little bit strange, but inconsequential.

@Jorrie, @pbuk

 

[pbuk]

Yes it is annoying: it is a consequence of ensuring that ranges spanning z = 0 always include z = 0 exactly (otherwise lines which are not smooth there look wrong). I'll add it to the bottom of the list

 

[JimJCW]

Yes it is annoying: it is a consequence of ensuring that ranges spanning z = 0 always include z = 0 exactly (otherwise lines which are not smooth there look wrong). I'll add it to the bottom of the list

I am sorry if I was a little pushy. I misunderstood your writing in Post #133:

“I think I'm all done now, . . .”​

 

[Jorrie]

Yes it is annoying: it is a consequence of ensuring that ranges spanning z = 0 always include z = 0 exactly (otherwise lines which are not smooth there look wrong). I'll add it to the bottom of the list

I have looked at all the code, tested as far as I could and I am very happy with it. Your model.ts code is a little outside of my capabilities as a coder, so I cannot really assist with the doubling of the z=0 line issue.

I'm happy to look after the UI and the Tutorial, so if anyone has any suggestion, we can discuss and implement if approved by the core team. The equations in Tutorial part III still represent the LCDM model accurately. Do you think we must put more in there to indicate how you have implemented the model?

 

[JimJCW]

Hi @pbuk and @Jorrie,

I am getting some puzzling results concerning output of LightCone8. Let’s use PLANCK 2015 Data and calculate OmegaM and OmegaR with three different methods:

1. LightCone8
2. Correct equations,

3. Wrong equations,

omegaM0 = (1 - OmegaL0) * zEq / (zEq + 1)​

omegaR0 = OmegaM0 / zEq​

The results are shown in the table below:

Is it possible that the LightCone8 result is calculated with the wrong equations?

 

[pbuk]

Is it possible that the LightCone8 result is calculated with the wrong equations?

No, I have looked into it and it is more subtle than that. The problem is that LightCone8 uses the same mapping to the UI as LightCone7 which does not pass the input value of $\Omega_\Lambda$ directly to the caculations, it actually passes
$$ \left ( \frac{\frac{978}{H_0}}{\frac{978}{H_0\sqrt{\Omega_\Lambda}}} \right )^2 !$$
This obviously introduces roundoff error, and because of the sqaring with a 53 bit mantissa I'd expect only about 8 decimal digits of accuracy - your calculation shows ~7.5. The fact that a similar result occurs when taking $z_{eq}$ instead of $s_{eq}$ is coincidence.

I can easily fix this by changing the inputs to the back-end to read directly from the input fields.

 

[pbuk]

I can easily fix this by changing the inputs to the back-end to read directly from the input fields.

Well I did that, which was very worthwhile but it didn't fix the problem, which was actually less subtle: the input field that is labelled $z_{eq}$ was actually

<input name="s_eq" />

so while it looked as though the conversion $s_{eq} = z_{eq} + 1$ had been done, it never was!

I have now fixed this, and tidied up the credits and a couple of other things, enjoy.

 

[JimJCW]

Well I did that, which was very worthwhile but it didn't fix the problem, which was actually less subtle: the input field that is labelled $z_{eq}$ was actually

<input name="s_eq" />

so while it looked as though the conversion $s_{eq} = z_{eq} + 1$ had been done, it never was!

I have now fixed this, and tidied up the credits and a couple of other things, enjoy.

Great! Now we have,

	z​	OmegaM​	OmegaL​	OmegaR​	OmegaT​
LightCone8	0​	0.3089083630000​	0.691​	0.0000916370107​	1​
Wrong Eqs	0​	0.3089083358054​	0.691​	0.0000916641946​	1​
Correct Eqs	0​	0.3089083629893​	0.691​	0.0000916370107​	1​

 

[pbuk]

Thanks, of course the UI only displays a maximum 9 digits, the results at full precision are:

• OmegaMatterT: 0.3089083629893239
• OmegaLambdaT: 0.691
• OmegaRadiationT: 0.0000916370106761566
• OmegaTotalT: 1

A useful feature for the UI would be display and/or download in CSV and/or JSON format of the full precision results.

 

[Jorrie]

Thanks, of course the UI only displays a maximum 9 digits, the results at full precision are:

----

A useful feature for the UI would be display and/or download in CSV and/or JSON format of the full precision results.

Could be nice to have for model comparison, but we know the cosmological inputs to no better than the 4th digit, with $H_0$ the main constraint. A more important UI change would be to have the latest (2018/2021) Planck data release as an option.

Also, historically, LightCone has concentrated on the $H_0$ and $\Omega_\Lambda$ as base inputs. This was partially because the late @marcus had a lot of threads explaining the usefulness of having $\Omega_\Lambda$ a parameter.

However, the Planck base group of sampled data does not include $\Omega_\Lambda$, but rather the two matter densities: baryonic and dark - see Table 2 on page 16 in the paper (partially shown below) - which they then statistically add together to get $\Omega_m$. So $\Omega_\Lambda$ is a derived parameter in their books.

In the past data releases, we took the LightCone inputs from the last row of the second group of derived inputs in Table 2. I suggest that we keep on doing that, except for giving priority to $\Omega_m$ rather than $\Omega_\Lambda$, the latter becoming a derived parameter.

Before we change the mapping between the UI and the calculation module, I think we should decide on this.

BTW, since the "glitch" has now been fixed, I think we should start a new thread on LightCone Improvement and close this one.