Skip to content

Angular diameter distance between redshifts and negative curvature #4661

Closed

Description

@joergdietrich

angular_diameter_distance_z1z2 contains

        # does not work for negative curvature
        Ok0 = self._Ok0
        if Ok0 < 0:
            raise CosmologyError('Ok0 must be >= 0 to use this method.')

presumably because

            out = ((dm2 * np.sqrt(1. + Ok0 * dm1 ** 2 / dh_2) -
                    dm1 * np.sqrt(1. + Ok0 * dm2 ** 2 / dh_2)) /
                   (1. + z2))

comes from Hogg (1999, arXiv:astro-ph/9905116), which claims "Unfortunately,
the above formula is not correct for Ωk < 0 (Phillip Helbig, 1998, private communication)." However, the derivation in Peebles (1993, p. 336) seems to indicate that this result is independent of the sign of curvature.

There is no numeric problem because Ok0 * (dm/dh)**2 evaluates to -sin^2[sqrt(|Ok0| dc/dm], which has a lower bound of -1.
Using the expression above and comparing it to Dm(z1, z2)/(1+z2), where Dm(z1, z2) is computed by performing the integral over 1/E(z) from z1 to z2, gives numerically identical results.

@nhmc Am I missing something or is the claim in Hogg (1999) wrong and the code should allow negative Ok0?

Activity

  1. joergdietrich commented on Mar 2, 2016

    @joergdietrich
    ContributorAuthor

    PR in case I'm right.

  2. astrofrog commented on Mar 2, 2016

    @astrofrog
    Member

    @nhmc @aconley - can you comment on this?

  3. nhmc commented on Mar 2, 2016

    @nhmc
    Contributor

    I just copied the result from Hogg et al. So it sounds likely that the expression is indeed valid for Omega_k < 0, but it would be nice to have an independent check. I won't have time to check for the next couple of days though - maybe @aconley can take a look sooner?

  4. aconley commented on Mar 7, 2016

    @aconley
    Contributor

    Hmm.. looks like my email started swallowing everything from github about a month ago. So I just saw this (and a few other things!). I will try to investigate this tomorrow-ish.

  5. aconley commented on Mar 8, 2016

    @aconley
    Contributor

    As far as I can tell, @joergdietrich is right. Unfortunately, this seems to be very seldom discussed in cosmology books. But in the derivation I found in Weinberg 'Cosmology' (I threw my Peebles away a long time ago... ugh, what a horrible book), indeed, this doesn't depend on the sign of the curvature.

    I think the next thing to do is to ask Dave (Hogg) if he has any more details to back this up.

  6. joergdietrich commented on Mar 8, 2016

    @joergdietrich
    ContributorAuthor

    @davidwhogg We're discussing the origin and content of the priv. comm. (Helbig) in your 1999 arXiv manuscript on cosmic distances. Could you clarify the reasoning behind the claim that the expression from Peebles for the angular diameter distance does not hold for negative curvature? It seems to us that is valid in this case.

  7. davidwhogg commented on Mar 8, 2016

    @davidwhogg

    I don't think it is valid. But I have an old email chain on this, with a generalization that is valid. Keep pinging me.

  8. joergdietrich commented on Mar 9, 2016

    @joergdietrich
    ContributorAuthor

    I also contacted Phillip Helbig. He promised to look into this as well (by the weekend at the latest).

  9. joergdietrich commented on Mar 21, 2016

    @joergdietrich
    ContributorAuthor

    I exchanged quite a few lengthy emails with Phillip Helbig on this. The first point to note is that Omega_k < 0 as defined in astropy (and many textbooks) of course implies positive curvature, different from what the comment in the code and this issue title say. Nevertheless, the code does what it is supposed to do and reflects the discussion between Phillip and David in 1998.

    Phillip and I now come to the conclusion that Peebles equation, which is used in angular_diameter_distance_z1z2 is indeed correct for all values of the curvature parameter, just as Peebles (1993) says at the bottom of page 336. (@davidwhogg you may want to update your arXiv manuscript). However, while this expression is correct, it seems to have numerical issues for very closed universes:

    In [2]: cosmo = cosmology.LambdaCDM(H0=100, Om0=10, Ode0=0, Tcmb0=0)
    In [4]: cosmo.angular_diameter_distance_z1z2(2, 5) / cosmo.hubble_distance
    Out[4]: <Quantity 0.005551158453390785>
    

    while an analytic expression (Schneider, Ehlers, Falco, Gravitational Lenses, eq. 4.57) gives the value 1.9772409E-02. Numerical integration of a 2nd order differental equation as described in Kayser, Helbig, Schramm (https://ui.adsabs.harvard.edu/#abs/1997A%26A...318..680K/abstract) agrees with the SEF result to 9 decimals.

    Note that the expression currently used in astropy gives very good results for all cosmologies we would consider reasonable. However, there may be unforeseen use cases (educational maybe?) for extreme values of Omega_m. If this is a concern, #4462 should not be merged. Rather the problem can be mostly circumvented by computing the elliptical integral over E(z) to obtain the comoving los distance between two redshifts. Then use Hogg's eq. (16) to get the comoving transverse distance between z1 and z2 and divide by (1+z2) to arrive at the angular diameter distance between two redshifts. My old cosmology code does this and gets to within ~1% of the SEF result for this test case.

    I will prepare a PR for this over the next week or two unless the consensus if that #4662 is sufficient, in which case a warning should be added for very small values of Omega_k0.

  10. aconley commented on Mar 21, 2016

    @aconley
    Contributor

    My vote is for doing it 'right'.

    However -- are there significant performance implications? If so, it might
    be worth
    having a branch that uses the current method up to some |Omega_k|
    and then switches to the method you propose.

    Also: make sure whatever you are doing can handle arbitrary w(z). If not,
    more branching.

    On Mon, Mar 21, 2016 at 6:28 AM, Jörg Dietrich [email protected]
    wrote:

    I exchanged quite a few lengthy emails with Phillip Helbig on this. The
    first point to note is that Omega_k < 0 as defined in astropy (and many
    textbooks) of course implies positive curvature, different from what
    the comment in the code and this issue title say. Nevertheless, the code
    does what it is supposed to do and reflects the discussion between Phillip
    and David in 1998.

    Phillip and I now come to the conclusion that Peebles equation, which is
    used in angular_diameter_distance_z1z2 is indeed correct for all values
    of the curvature parameter, just as Peebles (1993) says at the bottom of
    page 336. (@davidwhogg https://github.com/davidwhogg you may want to
    update your arXiv manuscript). However, while this expression is correct,
    it seems to have numerical issues for very closed universes:

    In [2]: cosmo = cosmology.LambdaCDM(H0=100, Om0=10, Ode0=0, Tcmb0=0)
    In [4]: cosmo.angular_diameter_distance_z1z2(2, 5) / cosmo.hubble_distance
    Out[4]: <Quantity 0.005551158453390785>

    while an analytic expression (Schneider, Ehlers, Falco, Gravitational
    Lenses, eq. 4.57) gives the value 1.9772409E-02. Numerical integration of a
    2nd order differental equation as described in Kayser, Helbig, Schramm (
    https://ui.adsabs.harvard.edu/#abs/1997A%26A...318..680K/abstract) agrees
    with the SEF result to 9 decimals.

    Note that the expression currently used in astropy gives very good results
    for all cosmologies we would consider reasonable. However, there may be
    unforeseen use cases (educational maybe?) for extreme values of Omega_m. If
    this is a concern, #4462 #4462
    should not be merged. Rather the problem can be mostly circumvented by
    computing the elliptical integral over E(z) to obtain the comoving los
    distance between two redshifts. Then use Hogg's eq. (16) to get the
    comoving transverse distance between z1 and z2 and divide by (1+z2) to
    arrive at the angular diameter distance between two redshifts. My old
    cosmology code does this and gets to within ~1% of the SEF result for this
    test case.

    I will prepare a PR for this over the next week or two unless the
    consensus if that #4662 #4662 is
    sufficient, in which case a warning should be added for very small values
    of Omega_k0.

    —
    You are receiving this because you were mentioned.
    Reply to this email directly or view it on GitHub
    #4661 (comment)

  11. joergdietrich commented on Mar 23, 2016

    @joergdietrich
    ContributorAuthor

    The new PR #4723 should take care of all this. Arbitrary w(z) are naturally included because the new functions only integrate over 1/E(z) as did the old ones.

  12. aconley commented on Oct 29, 2016

    @aconley
    Contributor

    @joergdietritch -- I guess you must not have put closes #4661 in#4723.
    So would you mind closing this issue?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions