Alexandria  2.14.1
Please provide a description of the project.
CosmologicalDistances.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2012-2020 Euclid Science Ground Segment
3  *
4  * This library is free software; you can redistribute it and/or modify it under
5  * the terms of the GNU Lesser General Public License as published by the Free
6  * Software Foundation; either version 3.0 of the License, or (at your option)
7  * any later version.
8  *
9  * This library is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
12  * details.
13  *
14  * You should have received a copy of the GNU Lesser General Public License
15  * along with this library; if not, write to the Free Software Foundation, Inc.,
16  * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
17  */
18 
24 #include <cmath>
25 
27 #include "ElementsKernel/Real.h"
28 
33 
35 
36 namespace Euclid {
37 namespace PhysicsUtils {
38 
40  return Elements::Units::c_light * 1.0E3 / parameters.getHubbleConstant(); // in [pc]
41 }
42 
44  const CosmologicalParameters& parameters) const {
45  double square = (1. + z) * (1. + z);
46  double cube = square * (1. + z);
47  return std::sqrt(
48  parameters.getOmegaM() * cube + parameters.getOmegaK() * square
49  + parameters.getOmegaLambda());
50 }
51 
53  const CosmologicalParameters& parameters, double relative_precision) const {
54  if (Elements::isEqual(0., z)) {
55  return 0.;
56  }
57 
60  relative_precision, MathUtils::SimpsonsRule::minimal_order) };
61 
63  [&parameters,this](double x) {return 1./hubbleParameter(x,parameters);});
64 
65  return hubbleDistance(parameters)
66  * integrate(adpater, 0., z, std::move(integrationScheme));
67 
68 }
69 
71  const CosmologicalParameters& parameters) const {
72  double comoving = comovingDistance(z, parameters);
73  if (Elements::isEqual(0., parameters.getOmegaK())) {
74  return comoving;
75  }
76 
77  double dHOverSqrtOmegaK = hubbleDistance(parameters)
78  / std::sqrt(std::abs(parameters.getOmegaK()));
79 
80  if (parameters.getOmegaK() > 0) {
81  return dHOverSqrtOmegaK * std::sinh(comoving / dHOverSqrtOmegaK);
82  } else {
83  return dHOverSqrtOmegaK * std::sin(comoving / dHOverSqrtOmegaK);
84  }
85 }
86 
88  const CosmologicalParameters& parameters) const {
89  if (Elements::isEqual(0., z)) {
90  return 10.;
91  } else {
92  return (1. + z) * transverseComovingDistance(z, parameters);
93  }
94 }
95 
97  const CosmologicalParameters& parameters) const {
98  return 5. * std::log10(luminousDistance(z, parameters) / 10.);
99 }
100 
101 
103  const CosmologicalParameters& parameters) const{
104  double D_H = hubbleDistance(parameters);
105  double E = hubbleParameter(z,parameters);
106  double D_M = transverseComovingDistance(z,parameters);
107  return D_M*D_M/(E*D_H*D_H);
108 }
109 
110 }
111 }
double getOmegaLambda() const
Get Omega Lambda for the cosmology.
double distanceModulus(double z, const CosmologicalParameters &parameters) const
return the correction for the Magnitude due to the distance: DM =5*log_10(DL/10pc)
T sinh(T... args)
double transverseComovingDistance(double z, const CosmologicalParameters &parameters) const
return the transverse comoving distance in [pc]
Class implementing the NumericalIntegrationScheme interface.
T log10(T... args)
Adapt a std::function<double(double)> to the Function Interface.
T sin(T... args)
double getOmegaK() const
Get the Omega curvature (computed as 1 - Omega_m - Omega_L) for the cosmology.
double getOmegaM() const
Get Omega matter for the cosmology.
double getHubbleConstant() const
Get the Hubble constant H_0 in (km/s)/Mpc.
double hubbleDistance(const CosmologicalParameters &parameters) const
Get the computed Hubble distance for the cosmology.
T move(T... args)
ELEMENTS_API double integrate(const Function &function, const double min, const double max, std::unique_ptr< NumericalIntegrationScheme > numericalIntegrationScheme=nullptr)
double luminousDistance(double z, const CosmologicalParameters &parameters) const
return the luminous distance in [pc]. For z=0 the returned value is 10pc.
double comovingDistance(double z, const CosmologicalParameters &parameters, double relative_precision=0.0000001) const
return the comoving distance in [pc]. This value is obtained through a numerical integration....
STL class.
bool isEqual(const RawType &left, const RawType &right)
T sqrt(T... args)
static const int minimal_order
Definition: SimpsonsRule.h:51
double dimensionlessComovingVolumeElement(double z, const CosmologicalParameters &parameters) const
return the dimensionless comoving volume element.
double hubbleParameter(double z, const CosmologicalParameters &parameters) const
Returns the (unit-less) Hubble parameter E(z)
Model the cosmological parameters. Omega_m, Omega_lambda, Omega_k and hubble_constant....
constexpr double c_light