Direct tests of General Relativity under screening effect with galaxy-scale strong lensing systems

As a successful theory of a dynamical space time where general coordinate invariance acts an essential role, Einstein’s Theory of General Relativity (GR) has passed all observational test so far (Ashby, 2002; Bertotti et al., 2003), from the millimetre scale in the laboratory to the Solar System tests and consistency with gravitational wave emission by binary pulsars (Dyson et al., 1920; Pound & Rebka, 1960; Shapiro, 1964; Taylor et al., 1979; Shapiro, 2004; Williams et al., 2004). A recent review about the status of experimental tests of GR and of theoretical frameworks for analyzing them can be found in Will (2014). Extrapolating to the cosmological scales, GR seems to precisely characterize the history and evolution of the Universe, as well as the large scale structure of space-time and matter. On the other hand, a mysterious component with negative pressure, dubbed dark energy (DE) (Copeland et al., 2006) responsible for the accelerating expansion of the Universe should be invoked in the framework of GR, and the simplest candidate of DE is the cosmological constant (introduced for different reasons at the early years of the GR). However, the inconsistency between the observed value of this constant and the theoretical value of zero-point energy predicted by quantum field theory is considerable (Weinberg, 1989), which is a fairly challenging problem in theoretical physics and opens up a discussion of whether the GR could fail at larger, cosmological scales (Koyama, 2016). Therefore, besides presenting possible solutions to the cosmological constant problem (Sahni & Starobinsky, 2000; Padmanabhan, 2003; Peebles & Ratra, 2003; Mukohyama & Randall, 2004; Nojiri & Odintsov, 2006; Padilla, 2015), it would be equally important to formulate and quantitatively interpret the tests of GR on the extra-galactic scale with high precision, and the parameterized post-Newtonian (PPN) framework (Thorne & Will, 1971) offers an interesting approach to test the departure from GR. The post-Newtonian parameter denoted by $\gamma_{PN}$, which measures the amount of spatial curvature per unit mass, can be used to probe the deviations from GR (i.e. from the value $\gamma_{PN}=1$) on significantly larger length scales due to the scale independence.

Recently, applying time-delay measurements of strongly lensed quasars, Jyoti et al. (2019) proposed a new phenomenological model of gravitational screening as a step discontinuity in the measurement of $\gamma_{PN}$ at a cutoff scale $\Lambda$ to test the nonlinear departure from GR. In this model, Jyoti et al. (2019) made use of two key characteristics of most modified gravity (MG) theories: gravitational slip, meaning the difference between the gravitational potential created by temporal and spatial metric components (Daniel et al., 2008), and screening, which can restore GR on small scales but may lead to distinct signatures in the large scale structure of the universe (Vainshtein, 1972; Jain & Khoury, 2010; Joyce et al., 2015; Ferreira, 2019). It is worthy to note that the constraint given by Jyoti et al. (2019) is $|\gamma_{\rm PN}-1|\leq 0.2\times(\Lambda/100\,\rm kpc)$, with screening length $\Lambda=10-200$ kpc, which is limited by $10^{-1}\sim 10^{-2}$ precision of the strong lensing time delay measurements using quasars. More recently, (Abadi & Kovetz, 2021) proposed to use strongly lensed fast radio burst time-delay measurements to put constraints on $\gamma_{PN}$, and yield constraints as tight as $\left|\gamma_{PN}-1\right|\leqslant 0.04\times(\Lambda/100\;kpc)\times[N/10]^{-1/2}$ ($N$ denoting the sample size), which indicates that ten events alone could place constraints at a level of 10% in the range of $\Lambda=10-300\;kpc$. Therefore, extra cosmological probes are also required to make complementary and more precise studies of the screened MG model featuring a gravitational slip.

In the last decade, strong gravitational lensing (SGL) of galaxies have become powerful and promising probes to study astrophysical issues (Treu, 2010), such as the structure and evolution of galaxies, the parameters that characterize the geometry, content, and expansion of the Universe. In recent years, great efforts have been made in estimating cosmological parameters through SGL (Zhu, 2000; Chae, 2003; Chae et al., 2004; Mitchell et al., 2005; Grillo et al., 2008; Oguri et al., 2008; Zhu & Sereno, 2008a; Zhu et al., 2008b; Cao & Zhu, 2012; Cao, Covone & Zhu, 2012; Cao et al., 2012; Cao, Zhu & Zhao, 2012; Biesiada, 2006; Biesiada, Piórkowska, & Malec, 2010; Collett & Auger, 2014; Cardone et al., 2016; Amante et al., 2020), in measuring the Hubble constant $H_{0}$ (Fassnacht et al., 2002; Bonvin et al., 2016; Suyu et al., 2017; Birrer et al., 2019; Wong et al., 2020), the cosmic curvature (Räsänen et al., 2015; Xia et al., 2017; Qi et al., 2019; Liu et al., 2020; Zhou & Li, 2020), and the distribution of matter in massive galaxies acting as lenses (Zhu & Wu, 1997; Mao & Schneider, 1998; Jin et al., 2000; Keeton, 2001; Kochanek & White, 2001; Ofek et al., 2003; Treu et al., 2006; Chen et al., 2019; Geng et al., 2021). It should be stressed that all of the studies mentioned above have been carried out under the assumption of GR. On account of some controversies in the framework of GR, for instance the aforementioned cosmological constant problem, and the Hubble tension between different cosmological probes, we are prompted to test GR in this paper. Moreover, SGL by galaxies provide a unique opportunity to test GR in kpc scales (Bolton et al., 2006, 2008; Smith, 2009; Schwab et al., 2010; Cao et al., 2017; Collett et al., 2018; Yang et al., 2020; Liu et al., 2022), with reasonable prior assumptions and independent measurements of background cosmology and appropriate descriptions of the internal structure of lensing galaxies. Inspired by Schwab et al. (2010) and Cao et al. (2017), who used galaxy-scale SGL systems with measured stellar velocity dispersions to test GR, we utilize the recently combined sample of 158 SGL systems, which is summarized in Cao et al. (2015); Shu et al. (2017). Based on this sample we investigate gravitational slip with screening effect in a phenomenological model. Considering the limitation of the sample size of available SGL systems, we also take advantage of the simulated future LSST SGL data to assess the expected improvement in precision and study the degeneracy between $\gamma_{PN}$ and $\Lambda$.

In this paper, we focus on demonstrating the possibilities of testing gravitational slip at super-galactic screening scales $\Lambda$ by the observed velocity dispersion of SGL systems. As already mentioned, our approach to test gravitational slip at super-galactic screening scales $\Lambda$, invokes a combination of lensing and stellar kinematics. Two different lens models, with a power law density profiles for the total mass density and luminous density (i.e. stellar mass density) assumed, are used to illustrate the influence of the lens mass distribution on the PPN parameter $\gamma_{PN}$. It would be interesting to figure out if the constraints from the current and future SGL systems reveal any specific cutoff scales, beyond which MG is relevant, or GR is still valid with high precision up to scales of $\Lambda\sim 300$ kpc. Following the simplifying assumptions presented by (Jyoti et al., 2019), we also assume that the screening radius $\Lambda$ is bigger than the Einstein radius of the lens galaxy, and thus the stellar dynamics within the lens are not changed, but the departure from GR at larger radius would impact the photon path. This assumption will allow us to investigate the regime of super-galactic screeing, complementary to the discussions in (Bolton et al., 2006; Smith, 2009; Collett et al., 2018), where the screening effects occur within the galaxy and the post-Newtonian parameter $\gamma_{PN}$ is constrained through two different mass measurements: the dynamical mass measured from the stellar kinematics of the deflector galaxy, and the lensing mass inferred from the lensing image.

This paper is organized as follows. In Sec. 2, we give a brief introduction of the model we used to evaluate the velocity dispersion of lensing galaxies. In Sec. 3, all the observational and simulated data sets are introduced. We perform a Markov chain Monte Carlo (MCMC) analysis using different data sets, and discuss our results in Sec.4. Finally, conclusions are summarized in Sec. 5.

The post-Newtonian variables are applied to quantify the behavior of gravity and deviations from GR, and we adopt the notation and conventions of Ma & Bertschinger (1995) to express the perturbed Robertson-Walker metric

to characterize the cosmological space-time, where $a(\eta)$ is the cosmological scale factor ($\eta$ being the conformal time), $\Phi$ and $\Psi$ are Newtonian and longitudinal gravitational potentials. In the weak-field limit GR predicts $\Phi=\Psi$ and a gravitational slip manifested as $\Phi\neq\Psi$ occurs in MG theories, for instance scalar-tensor theories (Schimd et al., 2005), $f(R)$ theories (Chiba et al., 2003; Sotiriou, 2010), Dvali-Gabadadze-Porrati (DGP) model (Dvali et al., 2000; Sollerman et al., 2009), and massive gravity (Dubovsky, 2004). We remark here that for a modified gravity theory where the Newtonian and longitudinal gravitational potentials are different, i.e. the gravitational slip is not zero, the Poisson equation $\nabla^{2}\Phi=4\pi Ga^{2}\rho$ usually is modified into a form $\nabla^{2}\Phi=4\pi G\mu a^{2}\rho$ which is different from that of GR, where the modified gravity parameters are included in the $\mu$ term (Koyama, 2016). Rigorous approach would thus require to focus on a specific theory of modified gravity. Similar to Jyoti et al. (2019) we are taking a heuristic approach assuming that the screening mechanism recovers GR exactly in the central parts of lensing galaxies, so that standard spherical Jeans equation is sufficient to recover the dynamical mass inside the Einstein radius.

The deviation from GR is quantified by the ratio $\gamma_{PN}=\Psi/\Phi$, and $\gamma_{PN}=1$ represents the GR. Gravitational screening mechanisms can restore GR in regions of high density, potential or curvature, such as galaxies, but allow modification of gravity at cosmological scales, which corresponds to suppressing the additional gravitational degrees of freedom within a certain region. Recent reviews about screening mechanisms and the experimental tests of such theories are available (Joyce et al., 2015; Ferreira, 2019). In order to model the effect of screening, we follow the notations in Jyoti et al. (2019) and suppose the gravitational slip to be a step-wise function, discontinuous at a screening radius $\Lambda$, which is assumed to be larger than the (physical) Einstein radius, $\Lambda>R_{E}=D_{L}\theta_{E}$.

Photon geodesics feel the sum of Newtonian and longitudinal gravitational potentials, which is defined as $\Sigma\equiv\Phi+\Psi$. Then, for a spherically-systematic mass distribution, the departure from GR can be expressed as

where r and $\Lambda$ are physical distances, and $\Theta$ is the Heaviside step function. For $r\leq\Lambda$, Eq. (2) reduces to $\Sigma=2\Phi$ as in GR. According to the definition in Narayan & Bartelmann (1996), the lensing potential is

where $D_{L}$, $D_{S}$, $D_{LS}$ are the angular diameter distances from observer to the lens, the source, and between the lens and source, and $\mathcal{Z}$ represents the distance along the line of sight. Then, the lensing potential can be decomposed as

In this model, $\psi_{\rm GR}$ and $\Delta\psi$, denote the lensing potential in GR and the correction due to screening, which can be expressed as

respectively. Assuming the power-law density profiles for the total mass density for elliptical galaxies,

where the constants: $\rho_{0}$ and $r_{0}$, set the total mass of the lens and solving the Poisson equation $\nabla^{2}\Phi=4\pi Ga^{2}\rho$, one can get the Newtonian potential:

This will lead to the lensing potential in GR (Suyu, 2012; Jyoti et al., 2019; Abadi & Kovetz, 2021)

where the $r_{0}$ and $\rho_{0}$ parameters have been subsumed into $\theta_{E,\rm GR}$, and $\theta_{E,\rm GR}$ in the Einstein radius of the lens inferred within GR. The deflection angle becomes

One can obtain the PPN correction to the lensing potential $\Psi$ by integrating the potential in Eq. (6) (Jyoti et al., 2019):

where ${}_{2}F_{1}$ is the hypergeometric function,

and $\Gamma$ is the Euler’s Gamma function. Then, the correction to the deflection angle from Eq. (11) can be expressed as $\Delta\alpha=\partial_{\theta}\Delta\psi$. It should be pointed out that the gravitational slip correction to the lensing potential $\Delta\psi$ influences the lens parameters, such as time delay, image positions, observational Einstein radius $\Theta_{E,\rm obs}$, and so on, which are deduced from the lensing observables (Jyoti et al., 2019). For each observed lensing event, one should carry out the Markov chain Monte Carlo (MCMC) method to analyze the entire set of lens parameters under the effect of $\gamma_{PN}$ at a screening radius $\Lambda$, which is a costly process for a sample of SGL systems. Therefore, we follow the same procedure as adopted in (Jyoti et al., 2019), relating the unobservable GR Einstein angle $\theta_{E,\rm GR}$ and the observed $\theta_{E,\rm obs}$, which would took into account the gravitational slip, through the lens equation (Narayan & Bartelmann, 1996)

where $\beta$, $\theta$, as well as $\alpha$ are the source position, the angular position of lens images, and the deflection angle respectively. Setting $\beta=0$ in the equation above, we obtain

According to $\Delta\alpha=\partial_{\theta}\Delta\psi$, the exact solution for $\theta_{E,\rm GR}$ is found to be

where $c^{\prime}$ is given by Eq. (12), and $\gamma_{\rm PN}=1$ will lead to $\theta_{E,\rm GR}=\theta_{E,\rm obs}$. After the slip-term is introduced into the GR Einstein angle $\theta_{E,\rm GR}$, we want to express the averaged observed velocity dispersion under the gravitational slip. According to the theory of gravitational lensing (Schneider et al., 1992; Schwab et al., 2010), the angular size of the Einstein radius corresponding to a point mass $M_{E}$ (or the spherically symmetric mass distribution within the Einstein radius) is expressed as

where $M(\theta_{E,\rm GR})$ is the mass enclosed within a radius of $\theta_{E,\rm GR}$ and should not be affected by the modification of GR. We can acquire a useful expression if we rearrange terms in Eq. (16) with $R_{E}=D_{L}\theta_{E,\rm GR}$ (Schneider et al., 1992; Abadi & Kovetz, 2021):

One of the simplest models to describe the mass density of an elliptical galaxy is a scale-free model based on power-law density profiles for the total mass density, $\rho$, and luminosity density, $\nu$, (Koopmans, 2006a):

Here $r$ is the spherical radial coordinate from the lens center: $r^{2}=R^{2}+\mathcal{Z}^{2}$. One can use the anisotropy parameter $\beta$ to characterize anisotropic distribution of three-dimensional velocity dispersion pattern, which can be written as:

where $\sigma^{2}_{t}$ and $\sigma^{2}_{r}$ represent the tangential and radial components of the velocity dispersion. In our analysis, $\beta$ is assumed to be independent of $r$, and two cases will be considered: isotropic dispersion $\beta=0$ and anisotropic dispersion $\beta=const.\neq 0$. Applying the spherical Jeans equation (Binney, 1980), the radial dispersion of luminous matter $\sigma_{r}^{2}(r)$ of the early-type galaxies can be written as

where $\beta$ is the anisotropy parameter. Applying the mass density profile in Eq. (18), one can obtain the mass contained within a spherical radius $r$ and the total mass $M_{E}$:

where $\lambda(x)=\Gamma\left(\frac{x-1}{2}\right)/\Gamma\left(\frac{x}{2}\right)$ denotes the ratio of Euler’s gamma functions. Following the same formulae with the notation used in (Koopmans, 2006a) : $\xi=\delta+\gamma-2$, we acquire the radial velocity dispersion by scaling the dynamical mass to the Einstein radius:

Additionally, (Schwab et al., 2010; Cao et al., 2016, 2017) have pointed out that the actually observed velocity dispersion is measured over the effective spectroscopic aperture $\theta_{ap}$ and effectively averaged over the line-of-sight effects. Given the entanglement of the aperture with atmospheric blurring and luminosity-weighted averaging, the averaged observed velocity dispersion takes the form:

where $\tilde{\sigma}_{\rm atm}\approx\sigma_{\rm atm}\sqrt{1+\chi^{2}/4+\chi^{4}/40}$ and $\chi=\theta_{\rm ap}/\sigma_{\rm atm}$ (Schwab et al., 2010), the detailed expression of $\theta_{E,\rm GR}$ is given by Eq. (15), $\sigma_{\rm atm}$ is the seeing recorded by the spectroscopic guide cameras during observing sessions, and specific values of seeing for different surveys have been summarized in (Cao et al., 2016). According to Eq. (24), one can probe the PPN parameter $\gamma_{PN}$ and the screening radius $\Lambda$ on a sample of lenses with available information about measured redshifts of the lens and the source, velocity dispersion, as well as the Einstein radius $\theta_{E,obs}$. In this work, a flat $\Lambda CDM$ model is assumed as a fiducial cosmology with $H_{0}=67.36\mathrm{~{}km~{}s^{-1}~{}Mpc^{-1}}$ and $\Omega_{m}=0.315$ (Planck Collaboration et al., 2018). In light of some evidence (Koopmans et al., 2006b; Ruff et al., 2011; Bolton et al., 2012; Sonnenfeld et al., 2013; Cao et al., 2016) supporting that the total density profile $\gamma$ for massive galaxies shows a trend of the cosmic evolution, we consider the lens models where the total mass density evolves as a function of redshift for two cases: $\gamma=\delta$ and $\gamma\neq\delta$. Meanwhile, we can examine the possible influence of different lens models on the constraints of $\gamma_{PN}$ and $\Lambda$. As for the anisotropy parameter, $\beta=0$ is assumed in the case $\gamma=\delta$, and $\beta$ is being marginalized using a Gaussian prior with $\beta=0.18\pm 0.13$ for the case $\gamma\neq\delta$ (Gerhard et al., 2001; Schwab et al., 2010; Cao et al., 2016).

In this section, we present the details of observational and simulated SGL systems which are used to put constraints on the post-Newtonian slip parameter $\gamma_{PN}$ and the screening radius $\Lambda$.

In this paper, we used a mass-selected sample including 99 SGL systems, 5000 simulated SGL systems expected from future LSST survey and four sub-samples extracted from the forecasted samples with different selection criteria to probe the gravitational slip with screening effects and the parameters $\gamma_{0},\gamma_{1},\delta$ characterizing the structure of elliptical galaxies. In order to determine these parameters, which are assumed to be the same for all lensing galaxies, we used MCMC method to sample their probability density distributions based on the likelihood ${\cal L}\sim\exp{(-\chi^{2}/2)}$, where

was derived from the measured values of velocity dispersion $\sigma_{ap}$, and theoretical prediction Eq. (24) using the observed Einstein radius $\theta_{E,obs}$, and aperture radius $\theta_{ap}$. For the term $\sigma_{atm}$, which requires the information of seeing recorded by the spectroscopic guide cameras during observing sessions, we have used the seeing summarized in Cao et al. (2016) for the current intermediate-mass sample of SGL systems. For the simulated lensing systems from LSST, the median seeing in i-band, which is 0.75 arcsec reported in Collett (2015) was used. This value was also adopted to calculate the likely yield of observable gravitationally lensed quasars and supernovae based on the properties of LSST Oguri & Marshall (2010). The uncertainties of $\sigma_{ap}$ and $\theta_{E,\rm obs}$ were propagated to the final uncertainty ${\Delta\bar{\sigma}_{*,i}}$. Following the SLACS team Bolton et al. (2008) and Cao et al. (2016), the fractional uncertainty of $\theta_{E,\rm obs}$ is taken as 5%, and the specific strategy for estimating the uncertainty of $\sigma_{ap}$ and $\theta_{E,\rm obs}$ for the simulated LSST lenses can be found in section 3.1. As for the sub-sample of forecasted LSST lenses using the logarithmic polynomial cosmographic reconstruction through SNe Ia, the uncertainties of $D_{L}(z)$ reconstructed from SNe Ia are also propagated to the final uncertainty ${\Delta\bar{\sigma}_{*,i}}$. The numerical results for $(\Lambda,\gamma_{PN})$ and the lens model parameters with 68.3% confidence level are summarized in Table 1, and the 1D probability distributions and 2D contours with $1\sigma$ and $2\sigma$ confidence levels are displayed in Figs. 5-9.

In this work, we have studied the gravitational slip under a phenomenological model of gravitational screening, where GR is maintained for small radii and the departures from GR take the form of a gravitational slip beyond the screening scale $\Lambda$. Based on mass-selected galaxy-scale strong gravitational lenses from the SLACS, BELLS, LSD, as well as SL2S surveys and simulated future measurements of 5000 galaxy-scale SGL from the forthcoming Large Synoptic Survey Telescope (LSST) survey, we were able to evaluate this screened MG model with screening length $\Lambda=10-300\;kpc$, which is broad enough to cover one typical massive galaxy. This is also the first attempt to use the stellar kinematics of SGL systems to assess constraints on the PPN parameter $\gamma_{PN}$, screening radius $\Lambda$ within this screening range under two different lens models simultaneously. Here we summarize our main conclusions in more details.

Considering two different lens models where the total mass density and luminosity density of lensing galaxies are modeled as power-law density profiles $\rho(r)=\rho_{0}(r/r_{0})^{-\gamma}$ and $\nu(r)=\nu_{0}(r/r_{0})^{-\delta}$ respectively and $\gamma$ is assumed to evolve with redshift, our results indicate that the current intermediate mass early-type galaxies are not able to provide tight constraints on the PPN parameter in this new theories of modified gravity, but GR ($\gamma_{PN}=1$) is still valid with screening length $\Lambda=10-300\;kpc$ in both lens models. On the other hand, our work studies the complementary range $\Lambda\geq R_{E}$ compared with previous researches (Bolton et al., 2006; Smith, 2009; Schwab et al., 2010; Collett et al., 2018; Liu et al., 2022) where the screening is assumed to take place within the galaxy. Then, one interesting thing is to figure out if there exists any specific cutoff scale for galaxy-scale SGL systems beyond which MG is relevant. The constraint achieved in the case of $\gamma=\delta$ shows a significant value of $\Lambda\sim 100$ kpc, beyond which the departure from $\gamma_{PN}$ is noticeable within $1\sigma$ confidence level. Nevertheless, in the case of $\gamma\neq\delta$, the current sample shows no specific cutoff value in the range $10\;kpc<\Lambda<300\;kpc$. It supports the claim that the intermediate mass lenses may shed new light on testing the validity of GR under this screened MG model. In addition, the 68 % confidence level constraints on $\Lambda$ and $\gamma_{PN}$ are quite different depending on the lens model. Setting the luminosity profile of elliptical galaxies as a free parameter, we obtained larger best fit values of $\Lambda=142.92_{-106.94}^{+97.51}\;kpc$ and $\gamma_{PN}=0.937_{-0.767}^{+1.384}$ compared with $\Lambda=99.56_{-57.55}^{+97.70}\;kpc$ and $\gamma_{PN}=0.378_{-0.269}^{+0.522}$ in the case of $\gamma=\delta$. Moreover, in this paper, we assessed the constraints for the total mass-profile and light-profile shapes of elliptical galaxies. Allowing for the cosmic evolution of the total mass density profile exponent in the form of $\gamma=\gamma_{0}+\gamma_{1}z$, there is no obvious evidence suggesting that the total density profile of intermediate mass early-type galaxies has become steeper over cosmic time (up to $z\sim 1$), and the singular isothermal sphere model is well favored by the current mass-selected sample in both lens models.

Furthermore, we elaborated what kind of results one could acquire making use of the future measurements of a well-selected sample containing 5000 LSST lenses. The final results imply that much more severe constraints on the $\gamma_{PN}$ will be achieved with $10^{-2}\sim 10^{-3}$ precision in the regime of screening radii $\Lambda=10-300\;kpc$. Benefiting from LSST’s wide field of view and sensitivity, the SGL systems detectable in the future can be very helpful for testing GR on kpc-Mpc scales in modified theories of gravity. Interestingly, the degeneracy between the screening scale $\Lambda$ and the PPN parameter $\gamma_{PN}$ derived from the full simulated SGL samples and four sub-samples is very similar to the results presented in Jyoti et al. (2019); Abadi & Kovetz (2021), where the authors do not carry out full MCMC analysis in parameter space. Meanwhile, we still do not find any particular cutoff scale for these simulated SGL samples, which is consistent with the assumption that GR is valid and no screening effects are involved in simulating the LSST SGL systems. With the increasing number of available galaxy-scale lenses, our results imply that it would be advantageous to use velocity dispersion measurements of the intermediate-mass elliptical galaxies to probe departures from GR under the screened MG model, where the gravitational slip is modeled as a step-wise discontinuous phenomenon with the screening radius $\Lambda$. Other cosmological probes, such as strongly lensed fast radio bursts (FRBs), lensed gravitational-wave signals and so on would be beneficial complementary probes in the case $\Lambda\geq R_{E}$. Besides, the SIS model ($\gamma_{0}=2,\gamma_{1}=0$) is included within $1\sigma$ for all simulated samples selected with different criteria, which is consistent with the prior $\gamma=2.01\pm 1.24$ used to model the total mass density of the forecasted LSST lenses.

In this work, we also considered four sub-samples derived from the well-defined sample of 5000 LSST lenses and applied the logarithmic polynomial cosmographic reconstruction of distances based on the SNe Ia pantheon sample. It should be noted that the slight deviation from $\gamma_{PN}=1$ ($\gamma_{PN}=0.981_{-0.020}^{+0.021}$ for $\gamma=\delta$, and $\gamma_{PN}=0.935_{-0.131}^{+0.051}$ for $\gamma\neq\delta$) is a bit more noticeable than in the case of the full simulated SGL sample ($\gamma_{PN}=0.998_{-0.007}^{+0.003}$ for $\gamma=\delta$, and $\gamma_{PN}=0.973_{-0.071}^{+0.027}$ for $\gamma\neq\delta$) in both lens models. This indicates the possibility that over the redshift range $0<z<2.5$, the SNe Ia pantheon sample serving as standard candles may provide a valuable supplement to the assumed fiducial cosmology in testing departure from GR. For the sub-samples defined by different redshift ranges, the constraints on $\gamma_{PN}$ became more diverse in the lens model where $\gamma\neq\delta$ is assumed. Significant departures from $\gamma_{PN}=1$ are present, which are $\gamma_{PN}=0.547_{-0.369}^{+0.346}$ and $\gamma_{PN}=0.546_{-0.375}^{+0.394}$ corresponding to the sub-samples $0.3<z<0.65$ and $z\geq 0.65$, respectively. However, in the case of $\gamma=\delta$, the same sub-samples gave $\gamma_{PN}=0.999_{-0.013}^{+0.007}$ and $\gamma_{PN}=0.999_{-0.097}^{+0.054}$ respectively. Therefore, the lens mass model seems to have a great influence on the limits on the PPN parameter with screening length $\Lambda=10-300\;kpc$, which also can be concluded from the current intermediate mass SGL systems. Additionally, we have re-simulated the LSST samples with modified gravity effect present in the fiducial model and our results demonstrate the effectiveness of our methodology. In this paper, we just adopted a power-law profile to characterize the distribution of the luminous component. On the other hand, there are other more complicated but also more realistic descriptions of the luminosity density profiles in the literature (Hernquist, 1990; Navarro et al., 1997). It should be emphasized that more appropriate and accurate modeling of the structure of lensing galaxies will contribute to the precision and accuracy of testing the validity of GR using the SGL systems, and future systematic surveys such as LSST (Abell et al., 2009), DES (Frieman et al., 2004), and Euclid survey (Pocino et al., 2021) will greatly conduce to such studies.