Observational constraints on varying fundamental constants in a minimal CPC model

The flat $\Lambda$CDM model of cosmology, known as cosmic concordance model, correctly predicts the overall evolution of the Universe, including the current phase of acceleration. It is also consistent with many cosmological and astrophysical observations, which makes general relativity (GR) the protagonist as the fundamental theory of gravity. However, even though it is successful at different scales and explains very diverse data sets, there exist some specific observational discrepancies and even theoretical shortcomings that open up the possibility for extensions of GR. As few examples, we mention the recent observation that the Universe is expanding faster than expected [116], the cosmic curvature problem [49, 72], small-scale problems [111] and the statistically significant tension between the predictions for the Hubble constant value by early-time probes assuming the $\Lambda$CDM model and the corresponding determinations of $H_{0}$ by a number of late-time model-independent fits from local measurements; this incongruity is the so called Hubble tension [102, 50]. For an interesting recent review on the $\Lambda$CDM model problems, we refer the reader to Ref. [27].

The above-mentioned problems of the concordance cosmological model expose some of GR’s shortcomings in the infrared regime; on top of that, the ultraviolet regime of gravity also hints limitation of its description by GR. Within the cosmological context, the early universe puzzles demanded the introduction of the inflation paradigm [25, 107, 18]. Moreover, very-high energy/curvature regimes would demand a quantum version of GR that is not currently available in a final consistent form—in spite of the loop quantum gravity [11] and string theory [83] proposals. The early works by Stelle and others were an attempt to alleviate the quantization issues of GR via the addition of $R^{2}$ [126, 127] and other higher-order curvature terms to the Einstein-Hilbert action [92]. Starobinsky model also makes use of the $R^{2}$ term to realise an inflationary scenario derived from an action principle beyond that of GR [71, 124, 125]. From a more theoretical standpoint, one could also argue that the scope of GR is somewhat limited since it is based on a Riemannian manifold that is metric compatible and torsion free; as such it does not include teleparallel gravity [6, 12], Riemann-Cartan gravity [48, 104] or metric-affine gravity [78]. Moreover, GR is not a typical gauge theory since it exhibits universality and soldering properties [5]. These limitations of GR and its applications to cosmology may serve as motivation for pursuing modified gravity theories.

Some examples of modified gravity theories are Horndeski [75, 32, 85], DHOST [89, 19, 88], and $f(R)$ theories of gravity [123, 47, 108]. There is also the class of $f\left(R,\nabla R,\dots,\nabla^{n}R\right)$ theories [40] which can be viewed as effective field theories stemming from a more complete and fundamental theory in certain energy/curvature regimes; the $f\left(R,\dots,\nabla^{n}R\right)$ proposal has applications in gauge theory of gravity [38], black-hole physics [41], gravitational waves [132], inflation [42, 43], and dark-energy [39]. Both $f\left(R\right)$ and $f\left(R,\dots,\nabla^{n}R\right)$ theories were shown to be equivalent to Brans-Dicke (BD) theory [26] for specific values of its characteristic dimensionless parameter $\omega$ [123, 40]. In fact, Brans-Dicke theory is the paradigm for a large class of alternative proposals to GR [57]. BD follows from other frameworks too, such as a result of Kaluza-Klein after compactification [80, 84, 13, 110]. BD can be derived as a limit of string theory: the lower-energy limit of bosonic string leads to a BD model with $\omega=-1$ [28, 59]. According to [57], BD theory is accommodated within Lyra’s geometry [97, 122]—see also [44].

BD theory, its scalar-tensor generalizations—see e.g. [14, 24, 15]—, and other extensions of the current standard model naturally lead to the concept that fundamental couplings are indeed expected to vary. Some examples of varying physical couplings proposals are studies of variation of the fine-structure constant $\alpha$ [82, 61, 93, 129, 86, 35, 63, 96, 22], of the Newton’s gravitational constant $G$ [51, 79, 130, 90, 62, 109, 137, 131], and of the speed of light $c$ [1, 37, 119, 115, 135, 136, 29, 30, 118, 94, 95, 2, 103, 138]. The common feature of these works is to take only one of the couplings as a function of the cosmological time while disregarding possible variations of the others. Variation of more than one of the aforementioned couplings is a possibility that has been studied recently by [60, 36, 64, 65, 66, 67, 68, 69, 70, 21, 91, 45], which include the cosmological constant $\Lambda$—among others—to the set. More specifically, Refs. [60, 36, 64] consider the simultaneous variations of $\left\{c,G,\Lambda\right\}$ within the so called Co-varying Physical Couplings (CPC) framework. For example, Ref. [36] introduces a covariant $c$-flation mechanism to resolve the early-universe puzzles usually addressed by inflation. Ref. [64] applies the CPC scheme to the late-time universe for a particular ansatz for the time-varying trio $\left\{c,G,\Lambda\right\}$. The entangled variation of the couplings is a crucial feature of the CPC framework, one that distinguishes it from earlier approaches that account for varying physical constants, and a characteristic whose importance was already pointed out in Refs. [112, 16].

The Varying Speed of Light (VSL) proposals are examples of scenarios accounting for the variation of $c$ alone—see e.g. [105, 4], and [128] for a review. The use of Baryon Acoustic Oscillations (BAO) to constrain the value of $c$ was worked out recently by [119], based on a relation between the Hubble parameter function and the angular diameter distance at the same redshift. In [29, 118], the measurement of $c$ was done by utilising radio quasars calibrated as standard rulers. Ref. [115] used type Ia Supernovae (SNe Ia), BAO, $H(z)$, and CMB data to constrain the variation of $c$ as $\sim 10^{-2}$, a result confirmed by [95] through a combination of SNe Ia and strong gravitational lensing (SGL). In [103] a new method to test the eventual dependence of the speed of light on the redshift was implemented by combining the measurements of galaxy cluster gas mass fraction (GMF), $H(z)$ and SNe Ia; the analyses indicate negligible variation of the speed of light.

In this paper, we will advance in the proposal of studying the eventual variation of more than one fundamental constant simultaneously, namely, we shall investigate possible variations of $\{c,G\}$. Accordingly, we will perform a consistency check for the time-independence of the speed of light $c$ and of the gravitational coupling $G$ in the context of the Co-varying Physical Couplings framework; this will be done by constraining the minimal CPC model with the most recent GMF data obtained from galaxy clusters [101] and complementary data from SNe Ia [121]. The minimal CPC model is a particular realization of the CPC scheme that neglects possible variations of $\Lambda$ and entangles the eventual variations of $c(z)$ and $G(z)$ in a consistent way, viz. $G\left(z\right)=G_{0}\left(c/c_{0}\right)^{4}$. Three different parametrizations for $c(z)$ are analysed.

The paper is organized as follows: Section 2 contains the theoretical background focused on Brans-Dicke theory of gravity as an inspiration for the Co-varying Physical Couplings framework. The CPC scenario is built in Section 2.1 where the general constraint relating the simultaneous variation of $\{G,c,\Lambda\}$ is derived. In Section 2.2 we choose the minimal realization of the CPC scenario described in the end of the previous paragraph. Section 3 and the following deal with the fitting of the minimal CPC model to the galaxy cluster gas mass fraction data. Section 4 explains the data and how the related equations are adapted in the face of varying couplings; and, in Section 5, the plots and numerical values of the free parameters of the model are presented and discussed. Our conclusions and final comments are stated in Section 6.

BD equations of motion are given by [57]:

	$\displaystyle R_{\mu\nu}-\frac{1}{2}g_{\mu\nu}R+\frac{V}{2\phi}g_{\mu\nu}$	$\displaystyle=\frac{8\pi}{\phi}T_{\mu\nu}+\frac{\omega}{\phi^{2}}\left(\nabla_{\mu}\phi\nabla_{\nu}\phi-\frac{1}{2}g_{\mu\nu}\nabla^{\rho}\phi\nabla_{\rho}\phi\right)$
		$\displaystyle\hphantom{=\frac{8\pi}{\phi}T_{\mu\nu}\;}+\frac{1}{\phi}\left(\nabla_{\mu}\nabla_{\nu}\phi-g_{\mu\nu}\square\phi\right),$		(1)

where $g_{\mu\nu}$ is the space-time metric, $R_{\mu\nu}$ is the Ricci tensor, $R$ is the curvature scalar, $T_{\mu\nu}$ is the matter energy momentum tensor, $\phi$ is BD scalar field and $V$ is the associated potential. Differently from GR, BD equations of motion fully incorporate Mach’s principle by realising the idea that the gravitational coupling $G$ is an effect due to the whole matter in the universe. In fact, comparison of (1) with Einstein’s equations of GR suggests that BD scalar field $\phi$ plays the hole of the inverse gravitational coupling:

$$G=\frac{1}{\phi}.$$

(2)

BD equation (1) is supplemented by the equation of motion for $\phi$, namely:

$$\frac{2\omega}{\phi}\square\phi+R-\frac{\omega}{\phi^{2}}\nabla^{\rho}\phi\nabla_{\rho}\phi-\frac{dV}{d\phi}=0.$$

(3)

Therefore, $\phi$ is a dynamical field. Due to (2), a dynamical scalar field $\phi=\phi\left(x^{\mu}\right)$ implies a $G$ that is a function of the space-time coordinates. In this way, BD theory admits a varying physical coupling, viz. Newton’s gravitational constant.

Notice that BD theory reduces to standard GR in the limit of $\omega\rightarrow\infty$ [133]. Solar system tests require large values of $\omega$, typically $\omega\gtrsim 40,000$ [134]. This means that pure BD should be very close to GR in local scales. This is a shortcoming that could be potentially eliminated by a screening mechanism [81, 58]. Otherwise, one could argue that an extension of BD is of interest in cosmological scales. Our proposal aligns with the latter possibility.

A quantity of interest for imposing limits on possible $G$ and $c$ variations is the so-called gas mass fraction of galaxy clusters, defined by $f_{gas}=M_{gas}/M_{tot}$, where $M_{gas}$ is the mass of the intracluster gas and $M_{tot}$ is the total mass, which includes baryonic gas mass and dark matter mass. As it is largely known, the physical process dominating in the galaxy cluster baryonic content is the intracluster gas X-ray emission, detected predominantly via thermal bremmsstrahlung [120].

Within a radius $R$, the gas mass $M_{gas}(<R)$ can be obtained by X-ray observations, being written as [120]:

	$\displaystyle M_{gas}(<R)$	$\displaystyle=$	$\displaystyle\left(\frac{3\pi\hbar m_{e}c^{2}}{2(1+X)e^{6}}\right)^{1/2}\left(\frac{3m_{e}c^{2}}{2\pi k_{B}T_{e}}\right)^{1/4}m_{H}$		(16)
			$\displaystyle\times\frac{1}{[\overline{g_{B}}(T_{e})]^{1/2}}{r_{c}}^{3/2}\left[\frac{I_{M}(R/r_{c},\beta)}{I_{L}^{1/2}(R/r_{c},\beta)}\right][L_{X}(<R)]^{1/2}\;,$

where $T_{e}$ is the gas temperature, $m_{e}$ and $m_{H}$ are the electron and hydrogen masses, respectively, $X$ is the hydrogen mass fraction, $\overline{g_{B}}(T_{e})$ is the Gaunt factor, $L_{X}(<R)$ is the bolometric luminosity, $r_{c}$ stands for the core radius and

$$I_{M}(y,\beta)\equiv\int_{0}^{y}(1+x^{2})^{-3\beta/2}x^{2}dx\;,$$

$$I_{L}(y,\beta)\equiv\int_{0}^{y}(1+x^{2})^{-3\beta}x^{2}dx\;.$$

Then, from Eq. (16),

$$M_{gas}(<R)\propto c^{3/2}.$$

(17)

On the other hand, the galaxy cluster total mass within a given radius $R$ is calculated by considering that the intracluster gas is in hydrostatic equilibrium [120, 7, 8], i.e.

$$M_{tot}(<R)=-\left.\frac{k_{B}T_{e}R}{G\mu m_{H}}\frac{d\ln n_{e}(r)}{d\ln r}\right|_{r=R},$$

(18)

where $n_{e}$ is the electron number density and $\mu$ is the mean molecular weight [120]. It can be seen from the above equation that

$$M_{tot}\propto G^{-1}.$$

(19)

From (17) and (19), one concludes that the gas mass fraction $f_{gas}=M_{gas}/M_{tot}$ scales as:

$$f_{gas}\propto c^{3/2}(z)\,G(z).$$

(20)

Finally, by using Eqs. (14) and (15), the expression of $f_{gas}$ is cast in an appropriate form for the minimal CPC model:

$$f_{gas}\equiv\phi_{c}^{11/2}(z)\frac{M_{gas}}{M_{tot}}.$$

(21)

It is important to emphasize here that, due to Eqs. (13) and (15), the function $\phi_{c}(z)$ is dimensionless. This is consistent with the main arguments found in Refs. [54, 53] conditioning the physical meaningfulness of theories with varying constants to the use of their corresponding dimensionless ratios.

We shall put limits on the parameter $c_{1}$ in Eq. (15) by using the most recent X-ray gas mass fraction measurements of hot ($kT\geq 5$ keV), massive and relaxed galaxy clusters from Chandra archive. The sample consists of 44 galaxy clusters in redshift range $0.018\leq z\leq 1.160$ [101] (see Fig. 1–left). The bias in the mass measurements from X-ray data arising by assuming hydrostatic equilibrium was calibrated using robust mass estimates for the target clusters from weak gravitational lensing [9]. Moreover, the authors of Ref. [101] measured the gas mass fraction in spherical shells at radii near $r_{2500}$, instead of using the cumulative fraction integrated over all radii ($<r_{2500}$) as in previous works. This procedure reduces significantly the systematic uncertainties, cf. detailed discussion by [100].

On the other hand, the expected constancy of the $f_{gas}$ within massive, hot and relaxed galaxy clusters with redshifts in the interval of interest has been used to constrain cosmological parameters via the following equation:

$$f_{gas}(z)=\gamma(z)K(z)A(z)\left[\frac{\Omega_{b}}{\Omega_{M}}\right]\left(\frac{D_{L}^{*}}{D_{L}}\right)^{3/2}\,.$$

(22)

(Refs. [7, 8, 56, 100, 74, 101] give the details on this technique.) In Eq. (22), $K(z)$ stands for the calibration constant which is equal to $0.96\pm 0.12$ [100] and $A(z)$ represents the angular correction factor which is close to unity. The asterisk denotes the corresponding quantities for the fiducial model used in the observations to obtain $f_{gas}$ (usually a flat $\Lambda$CDM model with Hubble constant $H_{0}=70$ km s^-1 Mpc^-1 and the present-day total matter density parameter $\Omega_{M}=0.3$), finally, $\gamma(z)$ is the depletion parameter, which indicates the amount of cosmic baryons that are thermalized within the cluster potential—see details in Refs.  [7, 8, 17, 113, 100]. It is very important to comment that the cosmological analyses with gas mass fraction measurements are model-independent due to the ratio in the parenthesis of Eq. (22). This factor takes into account the expected variation in the gas mass fraction measurement when the underlying cosmology is modified. We perform our statistical analysis using $\gamma(z)$ as a constant factor, whose value was obtained from hydrodynamical simulations:⁵⁵5Recent observational analyses have supported this value, cf. Refs.  [73, 23]. $\gamma_{0}=0.85\pm 0.08$ [17, 113].

Therefore, for a redshift-dependent speed of light given by $c(z)=c_{0}\phi_{c}(z)$ and $z$-dependent gravitational coupling $G(z)\propto c^{4}(z)$, Eq. (22) must be modified according to Eq. (21), i.e.:

$$f_{gas}(z)=\phi_{c}^{11/2}(z)\gamma(z)K(z)A(z)\left[\frac{\Omega_{b}}{\Omega_{M}}\right]\left(\frac{D_{L}^{*}}{D_{L}}\right)^{3/2}\,.$$

(23)

It is possible to constrain $\phi_{c}(z)$ if one knows the luminosity distance to the galaxy cluster. Here, the luminosity distance for each galaxy cluster of the sample is obtained by using SNe Ia data in [121]. Our sample of $D_{L}$ is constructed from the apparent magnitude ($m_{b}$) of the Pantheon catalog (Fig. 1–right) by considering $M_{b}=-19.23\pm 0.04$ (the absolute magnitude) through the relation

$$D_{L}=10^{(m_{b}-M_{b}-25)/5}\mathrm{Mpc}.$$

(24)

However, in order to perform the appropriate tests, we must use SNe Ia at (approximately) the same redshift of the each galaxy cluster. Thus, we make a selection of SNe Ia according to the criterion: $|z_{GC}-z_{SNe}|\leq 0.005$. Then, we perform the weighted average for each measurement utilising:

$$\bar{D}_{L}=\frac{\sum_{i}{D_{L}}_{i}/\sigma_{{D_{L}}_{i}}^{2}}{\sum_{i}1/\sigma_{{D_{L}}_{i}}^{2}},$$

(25)

$$\sigma_{\bar{D}_{L}}^{2}=\frac{1}{\sum_{i}1/\sigma_{{D_{L}}_{i}}^{2}}.$$

(26)

It is worth mentioning that Pantheon [121] is the most recent wide refined sample of SNe Ia observations found in the literature; it consists of 1049 spectroscopically confirmed SNe Ia and covers a redshift range of $0.01\leq z\leq 2.3$.

At this point, it should be emphasized that Eq. (25) does not use a luminosity distance where $c=c(z)$. Therein, $\bar{D}_{L}$ is the the average value weighted by the uncertainties associated with the measurements of the SNe Ia distance modulus $\mu\equiv\left(m_{b}-M_{b}\right)$, cf. Eq. (24). Our working hypotheses are that the SNe Ia apparent magnitude data directly determine their luminosity distance and that any eventual correction to $D_{L}$ coming from $c(z)$ is subdominant with respect to the pre-factor $\phi_{c}^{11/2}(z)$ already present in Eq. (23). ⁶⁶6In actuality, an eventual dependence of $D_{L}$ on $c(z)$—such as that considered by [46]—would contribute to reduce the power of $\phi_{c}(z)$ in Eq. (23), thus decreasing the effect produced by a hypothesized varying $c$ on $f_{gas}$. This would make even more difficult for us to assess if $c$ is in fact $z$-dependent via $f_{gas}$ data. Our current analysis (disregarding the contribution of $c(z)$ in $D_{L}$) therefore provides a prediction of the maximum effect of the varying $c$ according to the $f_{gas}$ observational window.

Eq. (23) requires the value of the ratio $\Omega_{b}/\Omega_{M}$. In order to make our analyses independent of the Planck satellite results, we consider the value of this ratio as obtained by [101]: $0.150\pm 0.021$. This bound is obtained from six gas mass fraction measurements of galaxy clusters with $z<0.16$ and by using the SH0ES results [117] for the Hubble constant. The advantage in doing so is that this bound is insensitive to the particular model of dark energy used, since that its equation of state does not evolve strongly over this redshift range. Then, for further statistical analyses, we exclude these six gas mass fraction data with $z<0.16$, and thereby use only the remaining 38 $f_{gas}$ data (with $z>0.16$ ).

The constraints on the model parameter $c_{1}$ can be obtained from Eq. (23) by maximizing the likelihood function, ${\cal{L}}\propto e^{-\chi^{2}/2}$, where $\chi^{2}$ is given by:

$$\chi^{2}=\sum_{i=1}^{38}\frac{\left(\phi_{c}^{11/2}(z)-\left[\frac{f_{gas}}{K(z)A(z)\gamma_{0}}\right]\left[\frac{\Omega_{M}}{\Omega_{b}}\right]\left[\frac{D_{L}}{D^{*}_{L}}\right]^{3/2}\right)^{2}}{\sigma_{i}^{2}}.$$

(27)

Here, $\sigma_{i}$ denotes statistical errors associated with the gas mass fraction measurements and the Pantheon sample; they are obtained by using standard error propagation methods.

We consider the three cases for $\phi_{c}(z)$ in Eq. (15):

• (I)

  $\quad\phi_{c}(z)=1+c_{1}z$

• (II)

  $\quad\phi_{c}(z)=1+c_{1}z/(1+z)$

• (III)

  $\quad\phi_{c}(z)=1+c_{1}\log\left(1+z\right)$

We plot the results from our analyses in Fig. 2. The blue-dashed, red-solid, and black-solid lines correspond to results from the forms (I), (II), and (III) of the function $\phi_{c}(z)$. We obtain (at 1 and 2$\sigma$ c.l.): $c_{1}=-0.025\pm 0.012\pm 0.027$, $c_{1}=-0.037\pm 0.025\pm 0.038$, and $c_{1}=-0.003\pm 0.005$ for the parameterizations (I), (II), and (III), respectively. Hence, the standard result ($c_{1}=0$) is recovered within 2$\sigma$ c.l. for all the three cases. It is worth comparing our results with those in Ref. [103]. Therein, the authors considered that only $c$ was varying, with a $\phi_{c}(z)$ function such as that in our case (I); they obtained $c_{1}=-0.069\pm 0.056$ and $c_{1}=-0.198\pm 0.052$ from analyses involving gas mass fraction plus $H(z)$ data and from gas mass fraction plus SNe Ia observations, respectively. Clearly, the results are in total mutual agreement.

The negative results regarding a possible variation of $c$—and, consequently, of $G$—occur for the linear parameterizations of $c$ in terms of $z$ and $a$—P (I) and P (II). It is interesting to note that the non-linear parameterization in terms of $\log(z)$—P (III)—also favors a non-varying $c$, and it does so within an even more stringent uncertainly interval than that of P (I) and P (II). Of course, the constancy of the pair $\{c,G\}$ was checked for the particular parameterizations in Eq. (15) with the particular data set from galaxy cluster gas mass fraction. Other parameterizations for $c(z)$ and other observational windows could lead to different conclusions.

The $f_{gas}$ data is related to relatively low redshifts. For these redshift values one expects that slowly varying functions for $c(z)$ would fit the data well. In fact, our results support this conjecture. This is true even for parameterization (III) of $c(z)\propto\log(z)$. However, for higher redshifts, one might need more rapidly varying functions for $c$ in terms of the redshift. Ref. [67] suggests a parameterization of this kind by adopting $c(a)\propto\exp(a^{\alpha}-1)$ with $\alpha=1.8$. [67] is then able to to resolve the primordial lithium abundance problem. For doing so, the numerical code AlterBBN [10] utilized to predict the Big-Bang Nucleosynthesis (BBN) of light elements was carefully modified for including the effect of varying constants compliant with the exponential type of ansatz adopted for $c(a)$. It is unclear if our slowly changing parameterizations would be as suited for the enterprise of studying BBN as the form in Ref. [67]: recall that BBN occurred in the very early universe, in the first few minutes after the Big Bang [55]. This suggests another line of investigation: to tailor the BBN-analysis code for including our proposals P (I)–(III) and check for the predictions. We deem this task too extensive to be included in the present paper; we intend to pursue the BBN analysis of our model in a separate paper—just like [67] dedicated an entire article to this important problem. The BBN analysis was also the subject of relevant works by [76, 77, 87, 33, 34, 52, 20, 106]. Therein, the lithium abundance problem is addressed by allowing variations of some constants (such as the fine structure constant) and/or parameters (e.g. deuterium binding energy, quark masses, etc.). This shows the subject is highly non-trivial.

In this paper we have recalled the Co-varying Physical Coupling framework which admits the Einstein Field Equations for gravitational interaction while promoting the couplings $\left\{c,G,\Lambda\right\}$ to space-time varying functions. Requiring the validity of Bianchi identity and conservation of the stress-energy tensor leads to the General Constraint entangling the simultaneous variation of $\left\{c,G,\Lambda\right\}$. The simplest model within CPC framework is the one where $\Lambda=$ constant and $G\propto c^{4}$. This model was dubbed the minimal CPC model. It was constrained by using the most recent $f_{\text{gas}}$ data set: 42 X-ray gas mass fraction measurements of hot ($kT\geq 5$ keV), massive and relaxed galaxy clusters from Chandra archive. The data set consists of galaxy cluster in the redshift range $0.078\leq z\leq 1.063$.

The data analysis constrained the couplings $c=c(z)$ and $G=G(z)$ to be consistent with no variation of these constants, at least for the particular reasonable parameterizations we have used. This conforms with the assumptions in the currently accepted physical theories and it is an important consistency check to their underlying hypotheses. However, it must be emphasised that this conclusion proceeds from a single observational window and should be cross checked with other data sets, such as $H(z)$ cosmic chronometers and type Ia supernovae data. This analysis is relegated to a different paper [46] because it is model dependent and relies on the use of the cosmology field equations in the CPC framework. On the other hand, we did not require the CPC analogue of Friedmann equations in this paper. In fact, we did not need to device the luminosity distance expression in the context of CPC framework because any eventual modification to $D_{L}$ due to the evolution of the co-varying couplings cancels out in the expression $f_{gas}(D_{L})$—see Section 4.

Although our observational constraints favour the standard interpretation of $c$ and $G$ as constant couplings, it should be emphasised that more general models within the CPC framework might indicate otherwise. In this direction, we point the reader to the works by [64, 65, 66, 67, 68, 69, 70], which not only demand $\Lambda$ to vary and adopts more sophisticated ansatze for the time-dependent couplings, but also relax the constancy of other fundamental physical constants, such as $\{\hbar,k_{B}\}$. The resulting phenomenology is rich and should not be disregarded.

Future perspectives include the study of CPC models including possible variation of $\Lambda$. The solution of the General Constraint for $\Lambda=\Lambda(z)$ is more evolved then the one considered in this paper since it necessarily involves specifying the matter-energy content. This modelling is currently under development. A first step in this direction is Ref. [46].

From a observational point of view, it is important to mention that the method discussed here is ready to be used in the near future with the upcoming data set from the X-ray survey eROSITA [114]; this survey is expected to detect $\approx$ 100,000 galaxy clusters.

RRC is grateful to Prof. Rajendra P. Gupta (University of Ottawa) for the hospitality and to CNPq-Brazil (Grant 309984/2020-3) for partial financial support. SHP acknowledges financial support from Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) (No. 303583/2018-5 and 308469/2021-6.). The authors thank Prof. Rajendra P. Gupta for useful discussions and comments and an anonymous referee whose questions and suggestions lead us to produce a scientifically stronger version of our paper.

All data used in this paper are from the references cited.

Authors have no conflict of interest.