Numerical Strategies of Computing the Luminosity Distance

We propose two efficient numerical methods of evaluating the luminosity distance in the spatially flat {\Lambda}CDM universe. The first method is based on the Carlson symmetric form of elliptic integrals, which is highly accurate and can replace numerical quadratures. The second method, using a modified version of Hermite interpolation, is less accurate but involves only basic numerical operations and can be easily implemented. We compare our methods with other numerical approximation schemes and explore their respective features and limitations. Possible extensions of these methods to other cosmological models are also discussed.