Observational Cosmology | Chapter 4

Weak Gravitational Lensing

How the subtle distortion of distant galaxy shapes by intervening mass reveals the distribution of dark matter, constrains cosmological parameters, and exposes tensions in our standard model

1. Introduction: Mass Bends Light

General relativity predicts that mass curves spacetime, and light follows geodesics through that curved geometry. When a massive object — a galaxy, a cluster, or even the large-scale cosmic web — lies between us and a distant source, the source's light is deflected. This phenomenon is gravitational lensing, and it provides one of the most powerful and direct probes of the total mass distribution in the universe, including the dominant dark matter component that emits no light.

In the strong lensing regime, massive clusters produce spectacular arcs, multiple images, and even complete Einstein rings. But these dramatic cases require nearly perfect alignment and probe only the densest regions. In the far more common weak lensing regime, the distortion of background galaxy shapes is tiny — typically a few percent — and can only be detected statistically by averaging over many galaxies. Despite its subtlety, weak lensing has become one of the most important tools in observational cosmology.

Why Weak Lensing Matters

Unlike galaxy surveys that trace the luminous matter distribution (and require assumptions about how galaxies trace mass), weak lensing responds directly to the total gravitational potential. It does not care whether the mass is baryonic or dark — both deflect light identically. This makes weak lensing the gold standard for:

Mapping dark matter in galaxy clusters and the cosmic web
Calibrating cluster masses for cosmology
Measuring the growth of cosmic structure over time
Constraining $\Omega_m$ and $\sigma_8$ independently of the CMB
Testing general relativity on cosmological scales

The first detection of cosmic shear — the weak lensing signal from large-scale structure — came in 2000, from four independent groups publishing nearly simultaneously. Since then, surveys of increasing scale and precision (CFHTLenS, DES, KiDS, HSC) have elevated weak lensing into a primary cosmological probe, and next-generation missions (Euclid, Vera C. Rubin Observatory) will survey billions of galaxies to push constraints to percent-level precision.

2. Derivation: Gravitational Lensing Basics

2.1 The Deflection Angle

Consider a photon passing a point mass $M$ with impact parameter $b$. In Newtonian gravity, treating the photon as a particle moving at speed $c$, one obtains a deflection angle $\hat{\alpha}_{\text{Newton}} = 2GM/(c^2 b)$. However, general relativity predicts exactly twice this value because of the contribution from spatial curvature (the $g_{rr}$ component of the Schwarzschild metric). The full GR calculation proceeds by integrating the geodesic equation for a null ray in the weak-field Schwarzschild metric:

The metric in isotropic form:

\[ds^2 = -\left(1 - \frac{2\Phi}{c^2}\right)c^2 dt^2 + \left(1 + \frac{2\Phi}{c^2}\right)(dx^2 + dy^2 + dz^2)\]

where $\Phi = -GM/r$ is the Newtonian potential. The effective index of refraction is $n = 1 - 2\Phi/c^2$, and by Fermat's principle the deflection is:

\[\hat{\alpha} = -\frac{2}{c^2}\int \nabla_\perp \Phi \, dl\]

For a point mass, evaluating the integral over the unperturbed straight-line path gives the celebrated Einstein deflection angle:

\[\boxed{\hat{\alpha} = \frac{4GM}{c^2 b}}\]

This is the result confirmed by Eddington's 1919 solar eclipse expedition, providing the first experimental test of general relativity. Note the factor of 2 enhancement over the Newtonian prediction.

2.2 The Lens Equation

For a thin gravitational lens, we define angular positions on the sky: the source has true position $\boldsymbol{\beta}$ and observed (image) position $\boldsymbol{\theta}$. The lens, observer, and source are separated by angular diameter distances $D_L$,$D_S$, and $D_{LS}$. The physical deflection angle $\hat{\alpha}$maps to a reduced deflection angle on the sky:

\[\boldsymbol{\alpha}(\boldsymbol{\theta}) = \frac{D_{LS}}{D_S}\hat{\boldsymbol{\alpha}}(D_L \boldsymbol{\theta})\]

The thin lens equation then relates source and image positions:

\[\boxed{\boldsymbol{\beta} = \boldsymbol{\theta} - \boldsymbol{\alpha}(\boldsymbol{\theta})}\]

This deceptively simple equation encodes all the geometry of gravitational lensing. For strong lenses, it can have multiple solutions (multiple images). For weak lensing, there is a single slightly displaced image.

2.3 The Einstein Radius

For a point mass with perfect source-lens-observer alignment ($\beta = 0$), the lens equation becomes $\theta = \alpha(\theta)$. Using the point-mass deflection angle $\hat{\alpha} = 4GM/(c^2 D_L \theta)$ and the reduced angle relation:

\[\theta = \frac{D_{LS}}{D_S} \cdot \frac{4GM}{c^2 D_L \theta}\]

Solving for $\theta$:

\[\boxed{\theta_E = \sqrt{\frac{4GM}{c^2}\frac{D_{LS}}{D_L D_S}}}\]

This is the Einstein radius, the characteristic angular scale of gravitational lensing. For a galaxy cluster of mass $M \sim 10^{15} M_\odot$ at$z \sim 0.3$ lensing sources at $z \sim 1$:

\[\theta_E \sim 30'' \quad \text{(galaxy cluster)}\]

For a single galaxy of mass $M \sim 10^{12} M_\odot$:

\[\theta_E \sim 1'' \quad \text{(galaxy-scale lens)}\]

2.4 Strong vs. Weak Lensing Regimes

The nature of lensing depends on whether the source lies inside or outside the Einstein radius:

Strong lensing ($\beta \lesssim \theta_E$): Multiple images, giant arcs, or complete Einstein rings form. The magnification can be enormous ($\mu \gg 1$). This regime probes the inner mass profile of clusters.
Weak lensing ($\beta \gg \theta_E$): A single slightly distorted image. The distortion (shear) is typically $|\gamma| \sim 0.01\text{--}0.05$, far smaller than the intrinsic ellipticity of galaxies ($\sigma_\epsilon \sim 0.3$). The signal must be extracted statistically from many galaxies.
Microlensing ($\theta_E$ unresolved): The Einstein radius is too small to resolve, but the magnification is detectable. Used for detecting exoplanets and compact dark matter (MACHOs).

3. Derivation: Shear and Convergence

3.1 The Lensing Potential

For a general mass distribution, we define the lensing potential as the scaled, projected Newtonian potential:

\[\psi(\boldsymbol{\theta}) = \frac{2}{c^2}\frac{D_{LS}}{D_L D_S}\int \Phi(D_L\boldsymbol{\theta}, z)\, dz\]

The reduced deflection angle is the gradient of this potential:

\[\boldsymbol{\alpha}(\boldsymbol{\theta}) = \nabla_\theta \psi(\boldsymbol{\theta})\]

3.2 Convergence and Critical Surface Density

The convergence $\kappa$ is the dimensionless surface mass density, defined as:

\[\boxed{\kappa(\boldsymbol{\theta}) = \frac{\Sigma(\boldsymbol{\theta})}{\Sigma_{\text{cr}}}}\]

where $\Sigma(\boldsymbol{\theta})$ is the projected (surface) mass density of the lens, and $\Sigma_{\text{cr}}$ is the critical surface density. To derive $\Sigma_{\text{cr}}$, note that for a uniform sheet of surface density$\Sigma$, the deflection angle at angular position $\theta$ is:

\[\hat{\alpha} = \frac{4G}{c^2}\pi \Sigma D_L \theta\]

The reduced deflection becomes $\alpha = (D_{LS}/D_S)\hat{\alpha}/(D_L) \cdot D_L = (4\pi G \Sigma D_L D_{LS})/(c^2 D_S) \cdot \theta$. Setting $\kappa = \alpha/\theta = 1$ (the condition for a critical lens) gives:

\[\boxed{\Sigma_{\text{cr}} = \frac{c^2}{4\pi G}\frac{D_S}{D_L D_{LS}}}\]

For typical cosmological distances, $\Sigma_{\text{cr}} \sim 0.3 \text{--} 0.5 \; \text{g/cm}^2$. When $\kappa > 1$, we are in the strong lensing regime. Weak lensing corresponds to $\kappa \ll 1$.

Note that $\kappa$ is also related to the lensing potential via the 2D Poisson equation:

\[\kappa = \frac{1}{2}\nabla^2_\theta \psi = \frac{1}{2}\left(\frac{\partial^2 \psi}{\partial \theta_1^2} + \frac{\partial^2 \psi}{\partial \theta_2^2}\right)\]

3.3 The Distortion Matrix and Shear

Weak lensing distorts the shapes of background galaxies. The mapping from source to image plane is described by the distortion matrix (or amplification matrix), obtained from the Jacobian of the lens mapping:

\[\mathcal{A}_{ij} = \frac{\partial \beta_i}{\partial \theta_j} = \delta_{ij} - \frac{\partial^2 \psi}{\partial \theta_i \partial \theta_j}\]

Writing out the second derivatives matrix explicitly:

\[\mathcal{A} = \begin{pmatrix} 1 - \kappa - \gamma_1 & -\gamma_2 \\ -\gamma_2 & 1 - \kappa + \gamma_1 \end{pmatrix}\]

where we have decomposed the second derivatives into the convergence and the two components of the complex shear:

\[\boxed{\gamma = \gamma_1 + i\gamma_2}\]

with the shear components defined as:

\[\gamma_1 = \frac{1}{2}\left(\frac{\partial^2 \psi}{\partial \theta_1^2} - \frac{\partial^2 \psi}{\partial \theta_2^2}\right), \quad \gamma_2 = \frac{\partial^2 \psi}{\partial \theta_1 \partial \theta_2}\]

The physical interpretation is clear: the convergence $\kappa$ produces an isotropic magnification (all dimensions scaled equally), while the shear $\gamma$ produces an anisotropic distortion that turns circles into ellipses. The magnification is:

\[\mu = \frac{1}{\det \mathcal{A}} = \frac{1}{(1-\kappa)^2 - |\gamma|^2}\]

In the weak lensing regime ($\kappa, |\gamma| \ll 1$), the magnification simplifies to $\mu \approx 1 + 2\kappa$ and the observed ellipticity of a galaxy is shifted by $\delta \epsilon \approx \gamma$.

3.4 Measuring Shear from Galaxy Shapes

In practice, each galaxy has an intrinsic ellipticity $\epsilon^s$ drawn from a distribution with zero mean and dispersion $\sigma_\epsilon \approx 0.3$. The observed ellipticity is:

\[\epsilon^{\text{obs}} \approx \epsilon^s + \gamma\]

Since $\langle \epsilon^s \rangle = 0$, averaging over $N$ galaxies gives:

\[\langle \epsilon^{\text{obs}} \rangle = \gamma \pm \frac{\sigma_\epsilon}{\sqrt{N}}\]

For a typical shear of $|\gamma| \sim 0.03$, one needs $N \sim 100$ galaxies per resolution element to detect the signal at $1\sigma$. This is why weak lensing surveys require enormous sky areas with high galaxy number density.

4. Derivation: Cosmic Shear Power Spectrum

4.1 From 3D to 2D: The Limber Approximation

Cosmic shear measures the cumulative lensing effect of all the large-scale structure between us and distant galaxies. To connect the observed shear field to the underlying matter power spectrum $P(k)$, we use the Limber approximation, which is valid for small-angle (high-$\ell$) modes.

The convergence field $\kappa$ at angular position $\boldsymbol{\theta}$ is a weighted projection of the 3D matter overdensity $\delta$ along the line of sight:

\[\kappa(\boldsymbol{\theta}) = \int_0^{\chi_H} W(\chi)\, \delta\!\left(\frac{\ell}{\chi}, \chi\right) d\chi\]

where $\chi$ is the comoving distance and $W(\chi)$ is the lensing efficiency kernel:

\[W(\chi) = \frac{3H_0^2 \Omega_m}{2c^2}\frac{\chi}{a(\chi)}\int_\chi^{\chi_H} n(z)\frac{dz}{d\chi'}\frac{\chi' - \chi}{\chi'}\, d\chi'\]

Here $n(z)$ is the normalized redshift distribution of source galaxies, and$a(\chi)$ is the scale factor. The kernel peaks roughly midway between observer and sources.

4.2 The Shear Power Spectrum

Since the shear and convergence are related by second derivatives of the lensing potential, their power spectra are equal: $C_\ell^{\gamma\gamma} = C_\ell^{\kappa\kappa}$. Applying the Limber approximation to the convergence projection gives:

\[\boxed{C_\ell^{\gamma\gamma} = \int_0^{\chi_H} \frac{W^2(\chi)}{\chi^2}\, P\!\left(\frac{\ell}{\chi},\, \chi\right) d\chi}\]

This is the central equation of cosmic shear cosmology. It connects the observable angular power spectrum of galaxy shape correlations to the 3D matter power spectrum $P(k, \chi)$, weighted by the lensing kernel $W(\chi)$.

The Limber approximation replaces $k = \ell/\chi$ (the flat-sky approximation), which is accurate to $\sim 1\%$ for $\ell \gtrsim 10$. For next-generation surveys, corrections to the Limber approximation may be needed at low $\ell$.

4.3 Sensitivity to Cosmological Parameters

The cosmic shear power spectrum depends on cosmology through two channels:

The matter power spectrum $P(k, z)$: The amplitude is set by $\sigma_8$ (the rms fluctuation in spheres of$8\, h^{-1}\text{Mpc}$), and the shape depends on $\Omega_m$,$\Omega_b$, $n_s$, and the transfer function.
The lensing kernel $W(\chi)$: This depends on the geometry ($\Omega_m$, $\Omega_\Lambda$,$w$) through the distance-redshift relation, and on the source redshift distribution $n(z)$.

The dominant degeneracy is between $\Omega_m$ and $\sigma_8$: increasing$\Omega_m$ while decreasing $\sigma_8$ can leave the amplitude of $C_\ell^{\gamma\gamma}$ approximately unchanged. This leads to the famous “banana-shaped” contours in the $\Omega_m$–$\sigma_8$ plane.

5. Derivation: The S₈ Parameter and Tensions

5.1 Defining S₈

To capture the well-constrained combination of $\sigma_8$ and $\Omega_m$that weak lensing measures, we define:

\[\boxed{S_8 \equiv \sigma_8 \sqrt{\frac{\Omega_m}{0.3}}}\]

The exponent arises because the cosmic shear amplitude scales approximately as $C_\ell^{\gamma\gamma} \propto \sigma_8^2 \Omega_m^\alpha$ with$\alpha \approx 0.5$ for typical survey configurations. The parameter $S_8$is then the direction of minimal uncertainty in the $\sigma_8$–$\Omega_m$ plane.

5.2 Weak Lensing Constraints on S₈

Modern weak lensing surveys have converged on a consistent picture:

Current S₈ Measurements

Survey	S₈ Value	Year	Area (deg²)
DES Y3	$0.776 \pm 0.017$	2022	4143
KiDS-1000	$0.759^{+0.024}_{-0.021}$	2021	1006
HSC Y3	$0.776^{+0.032}_{-0.033}$	2023	416
Planck CMB	$0.834 \pm 0.016$	2020	Full sky

5.3 The S₈ Tension

The most striking feature is the systematic offset between weak lensing and CMB measurements. Weak lensing surveys consistently find $S_8 \approx 0.76 \pm 0.02$, while the Planck CMB predicts $S_8 \approx 0.83 \pm 0.02$ assuming $\Lambda$CDM. The discrepancy is at the $\sim 2\text{--}3\sigma$ level — not decisive on its own, but persistent across independent surveys with different telescopes, pipelines, and systematics.

The tension can be quantified as:

\[\frac{S_8^{\text{Planck}} - S_8^{\text{WL}}}{\sqrt{\sigma_{\text{Planck}}^2 + \sigma_{\text{WL}}^2}} \approx \frac{0.834 - 0.770}{\sqrt{0.016^2 + 0.020^2}} \approx 2.5\sigma\]

While not yet at the $5\sigma$ discovery threshold, the consistency of the offset across multiple independent surveys makes it difficult to dismiss as a statistical fluke.

Possible explanations for the tension include:

Unaccounted systematics: Shear calibration errors, photo-z biases, intrinsic alignments, or baryonic feedback suppressing the matter power spectrum at small scales.
New physics: Decaying dark matter, dark energy with $w \neq -1$, massive neutrinos ($\sum m_\nu > 0.06$ eV), or modified gravity could reduce structure growth relative to the CMB prediction.
CMB systematics: Possible unaccounted-for effects in the Planck analysis, though this is considered less likely given extensive cross-checks.

6. Applications and Surveys

6.1 Current Surveys

Dark Energy Survey (DES): A 5000 deg² survey using the 570-megapixel DECam on the Blanco 4m telescope in Chile. DES Year 3 results used ~100 million galaxies to produce the most precise weak lensing constraints to date. The “3x2pt” analysis combines cosmic shear, galaxy-galaxy lensing, and galaxy clustering.

Kilo-Degree Survey (KiDS): A 1350 deg² survey using the OmegaCAM instrument on the ESO VLT Survey Telescope. Known for excellent image quality and rigorous shear calibration. KiDS-1000 results are among the most precise weak lensing measurements.

Hyper Suprime-Cam (HSC): A deep survey using the 870-megapixel HSC on the Subaru 8.2m telescope. Though smaller in area (~1100 deg²), its depth probes higher redshifts and complements DES and KiDS.

6.2 Cluster Mass Calibration

Beyond cosmic shear, weak lensing plays a critical role in cluster mass calibration. Galaxy clusters are the most massive gravitationally bound objects in the universe, and their abundance as a function of mass and redshift is highly sensitive to $\Omega_m$ and $\sigma_8$. However, cluster masses cannot be observed directly — they must be inferred from observables like X-ray luminosity, the Sunyaev-Zel'dovich (SZ) effect, or richness (galaxy count). Weak lensing provides the most unbiased mass calibration because it directly probes the gravitational potential.

6.3 CMB Lensing

The CMB photons themselves are weakly lensed by the intervening large-scale structure. CMB lensing probes the matter distribution at higher redshifts ($z \sim 0.5\text{--}5$) than galaxy lensing ($z \sim 0.1\text{--}1.5$), providing a complementary view. Cross-correlating CMB lensing maps (from Planck, ACT, SPT) with galaxy surveys enables powerful joint analyses that break degeneracies in cosmological parameters.

6.4 Next-Generation Surveys

Euclid (ESA, launched 2023): A space-based survey of 15,000 deg² with exquisite image quality from above the atmosphere. Euclid will measure the shapes of ~1.5 billion galaxies, achieving percent-level precision on $S_8$ and sub-percent constraints on dark energy parameters.

Vera C. Rubin Observatory LSST: A ground-based survey of 18,000 deg² using a 3.2-gigapixel camera, reaching unprecedented depth with 10 years of repeated imaging. The co-added data will measure ~4 billion galaxy shapes, making it the definitive weak lensing survey of the 2030s.

Roman Space Telescope (NASA): A 2.4m space telescope surveying ~2000 deg² with very high resolution. Though smaller in area, its depth and resolution complement Euclid and Rubin for high-redshift lensing and intrinsic alignment studies.

7. Historical Context

Timeline of Gravitational Lensing

1801 — Johann Georg von Soldner calculates the Newtonian deflection of light by the Sun, obtaining half the correct (GR) value.
1915 — Einstein completes general relativity and calculates the full deflection: $\hat{\alpha} = 4GM/(c^2 b) = 1.75''$ for the Sun.
1919 — Eddington's solar eclipse expedition confirms Einstein's prediction, making GR front-page news worldwide.
1936 — Einstein publishes the lens equation and the concept of an “Einstein ring,” but considers it practically unobservable.
1937 — Fritz Zwicky proposes that galaxy clusters could act as gravitational lenses, correctly predicting strong lensing decades before its discovery.
1979 — First observation of a multiply-imaged quasar (Q0957+561), confirming strong gravitational lensing.
1990s — Theoretical framework for weak lensing developed (Kaiser, Squires, Broadhurst 1995; Bartelmann & Schneider).
2000 — First detection of cosmic shear by four independent groups: Bacon, Refregier & Ellis; Kaiser, Wilson & Luppino; Van Waerbeke et al.; Wittman et al. A landmark year for observational cosmology.
2006 — The Bullet Cluster (1E 0657-558): weak lensing mass maps show dark matter separated from baryonic gas after a cluster merger — considered the most direct evidence for particle dark matter.
2013–present — DES, KiDS, and HSC produce increasingly precise cosmic shear measurements, revealing the S₈ tension with Planck.
2023 — Euclid launches, beginning the next era of space-based weak lensing cosmology.

8. Python Simulation: Weak Lensing Shear and S₈ Constraints

The following simulation computes the tangential shear profile around a galaxy cluster using an NFW (Navarro-Frenk-White) density profile, and also visualizes the S₈ constraints from different surveys. All computations use numpy only — no scipy.

Python

script.py247 lines

#!/usr/bin/env python3
"""
weak_lensing_simulation.py
Weak gravitational lensing: NFW tangential shear profile
and S8 constraints from major surveys.
"""
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# ========================================
# Constants
# ========================================
G = 6.674e-11       # m^3 kg^-1 s^-2
c = 2.998e8          # m/s
Msun = 1.989e30      # kg
Mpc = 3.086e22       # m
pc = 3.086e16        # m
arcsec = np.pi / (180.0 * 3600.0)  # radians per arcsec

# ========================================
# Cosmological parameters (flat LCDM)
# ========================================
H0 = 70.0           # km/s/Mpc
h = H0 / 100.0
Omega_m = 0.3
Omega_L = 0.7
rho_crit = 3.0 * (H0 * 1e3 / Mpc)**2 / (8.0 * np.pi * G)  # kg/m^3

# ========================================
# Angular diameter distance (flat LCDM)
# ========================================
def comoving_distance(z, nsteps=1000):
    """Comoving distance via trapezoidal integration."""
    zz = np.linspace(0, z, nsteps)
    dz = zz[1] - zz[0]
    Ez = np.sqrt(Omega_m * (1 + zz)**3 + Omega_L)
    integrand = 1.0 / Ez
    chi = (c / (H0 * 1e3 / Mpc)) * np.trapezoid(integrand, dx=dz)
    return chi  # meters

def angular_diameter_distance(z):
    return comoving_distance(z) / (1 + z)

def angular_diameter_distance_12(z1, z2):
    chi1 = comoving_distance(z1)
    chi2 = comoving_distance(z2)
    return (chi2 - chi1) / (1 + z2)

# ========================================
# NFW profile
# ========================================
z_lens = 0.3
z_source = 1.0
M200 = 1.0e15 * Msun   # kg
c200 = 5.0              # concentration

D_L = angular_diameter_distance(z_lens)
D_S = angular_diameter_distance(z_source)
D_LS = angular_diameter_distance_12(z_lens, z_source)

# Critical surface density
Sigma_cr = (c**2 / (4.0 * np.pi * G)) * (D_S / (D_L * D_LS))  # kg/m^2

# NFW scale radius
rho_c_z = rho_crit * (Omega_m * (1 + z_lens)**3 + Omega_L)
R200 = (3.0 * M200 / (4.0 * np.pi * 200.0 * rho_c_z))**(1.0/3.0)
r_s = R200 / c200

# NFW characteristic density
delta_c = (200.0 / 3.0) * c200**3 / (np.log(1 + c200) - c200 / (1 + c200))
rho_s = delta_c * rho_c_z

# ========================================
# NFW projected surface density Sigma(R)
# ========================================
def nfw_sigma(R):
    """NFW projected surface density (Bartelmann 1996)."""
    x = R / r_s
    result = np.zeros_like(x)

# x < 1
    mask1 = x < 1.0
    x1 = x[mask1]
    arg1 = np.sqrt((1 - x1) / (1 + x1))
    f1 = 1.0 / (x1**2 - 1) * (1 - np.log((1 + arg1) / (1 - arg1)) / np.sqrt(1 - x1**2))
    result[mask1] = 2 * rho_s * r_s * f1

# x = 1
    mask2 = np.abs(x - 1.0) < 1e-6
    result[mask2] = 2 * rho_s * r_s / 3.0

# x > 1
    mask3 = x > 1.0
    x3 = x[mask3]
    f3 = 1.0 / (x3**2 - 1) * (1 - 2.0 * np.arctan(np.sqrt((x3 - 1) / (x3 + 1))) / np.sqrt(x3**2 - 1))
    result[mask3] = 2 * rho_s * r_s * f3

return result

# ========================================
# NFW mean surface density within R
# ========================================
def nfw_sigma_mean(R):
    """Mean surface density within projected radius R."""
    x = R / r_s
    result = np.zeros_like(x)

# x < 1
    mask1 = x < 1.0
    x1 = x[mask1]
    arg1 = np.sqrt((1 - x1) / (1 + x1))
    g1 = np.log((1 + arg1) / (1 - arg1)) / np.sqrt(1 - x1**2) + np.log(x1 / 2.0)
    result[mask1] = 4 * rho_s * r_s * g1 / x1**2

# x = 1
    mask2 = np.abs(x - 1.0) < 1e-6
    result[mask2] = 4 * rho_s * r_s * (1 + np.log(0.5))

# x > 1
    mask3 = x > 1.0
    x3 = x[mask3]
    g3 = np.arctan(np.sqrt((x3 - 1) / (x3 + 1))) / np.sqrt(x3**2 - 1) + np.log(x3 / 2.0)
    result[mask3] = 4 * rho_s * r_s * g3 / x3**2

return result

# ========================================
# Tangential shear and convergence
# ========================================
R_min = 0.1 * Mpc   # 100 kpc
R_max = 5.0 * Mpc   # 5 Mpc
R = np.logspace(np.log10(R_min), np.log10(R_max), 200)

Sigma_R = nfw_sigma(R)
Sigma_mean_R = nfw_sigma_mean(R)

kappa = Sigma_R / Sigma_cr
gamma_t = (Sigma_mean_R - Sigma_R) / Sigma_cr  # tangential shear
delta_sigma = Sigma_mean_R - Sigma_R  # excess surface density

# Convert R to arcminutes
theta = (R / D_L) * (180.0 * 60.0 / np.pi)  # arcmin

# ========================================
# Plot 1: Tangential shear profile
# ========================================
fig, axes = plt.subplots(1, 2, figsize=(14, 6))

ax1 = axes[0]
ax1.loglog(theta, gamma_t, 'purple', linewidth=2.5, label=r'$\gamma_t(\theta)$')
ax1.loglog(theta, kappa, 'magenta', linewidth=2, linestyle='--', label=r'$\kappa(\theta)$')
ax1.set_xlabel('Angular radius (arcmin)', fontsize=13)
ax1.set_ylabel('Amplitude', fontsize=13)
ax1.set_title('NFW Weak Lensing Profile\n' +
              r'$M_{200} = 10^{15} M_\odot$, $c = 5$, $z_L = 0.3$, $z_S = 1.0$',
              fontsize=12)
ax1.legend(fontsize=12)
ax1.set_xlim(0.3, 80)
ax1.set_ylim(1e-3, 1)
ax1.grid(True, alpha=0.3)
ax1.axhline(y=1, color='gray', linestyle=':', alpha=0.5, label=r'$\kappa = 1$ (strong lensing)')

# Add R scale on top axis
ax1_top = ax1.twiny()
R_ticks_mpc = np.array([0.2, 0.5, 1.0, 2.0, 5.0])
theta_ticks = (R_ticks_mpc * Mpc / D_L) * (180.0 * 60.0 / np.pi)
ax1_top.set_xscale('log')
ax1_top.set_xlim(ax1.get_xlim())
ax1_top.set_xticks(theta_ticks)
ax1_top.set_xticklabels([f'{r:.1f}' for r in R_ticks_mpc])
ax1_top.set_xlabel('Physical radius (Mpc)', fontsize=11)

# ========================================
# Plot 2: S8 constraints from surveys
# ========================================
ax2 = axes[1]

# Survey data: (name, S8, sigma_low, sigma_high, color)
surveys = [
    ('DES Y3',       0.776, 0.017, 0.017, '#9333ea'),
    ('KiDS-1000',    0.759, 0.021, 0.024, '#a855f7'),
    ('HSC Y3',       0.776, 0.033, 0.032, '#c084fc'),
    ('DES+KiDS+HSC', 0.770, 0.017, 0.017, '#7c3aed'),
    ('ACT CMB lens', 0.840, 0.028, 0.028, '#f472b6'),
    ('Planck CMB',   0.834, 0.016, 0.016, '#ec4899'),
]

y_positions = np.arange(len(surveys))
names = [s[0] for s in surveys]
s8_vals = np.array([s[1] for s in surveys])
err_low = np.array([s[2] for s in surveys])
err_high = np.array([s[3] for s in surveys])
colors = [s[4] for s in surveys]

for i, (name, s8, el, eh, col) in enumerate(surveys):
    ax2.errorbar(s8, i, xerr=[[el], [eh]], fmt='o', color=col,
                 markersize=10, capsize=5, capthick=2, linewidth=2)

# Weak lensing band
ax2.axvspan(0.770 - 0.020, 0.770 + 0.020, alpha=0.15, color='purple',
            label='WL combined band')
# Planck band
ax2.axvspan(0.834 - 0.016, 0.834 + 0.016, alpha=0.15, color='pink',
            label='Planck band')
# Vertical reference lines
ax2.axvline(x=0.770, color='purple', linestyle='--', alpha=0.4)
ax2.axvline(x=0.834, color='deeppink', linestyle='--', alpha=0.4)

ax2.set_yticks(y_positions)
ax2.set_yticklabels(names, fontsize=11)
ax2.set_xlabel(r'$S_8 = \sigma_8 \sqrt{\Omega_m / 0.3}$', fontsize=13)
ax2.set_title(r'$S_8$ Constraints: The Lensing Tension', fontsize=13)
ax2.set_xlim(0.68, 0.92)
ax2.legend(fontsize=10, loc='lower right')
ax2.grid(True, axis='x', alpha=0.3)

# Annotate tension
ax2.annotate('', xy=(0.770, -0.7), xytext=(0.834, -0.7),
             arrowprops=dict(arrowstyle='<->', color='white', lw=1.5))
ax2.text(0.802, -0.9, r'$\sim 2.5\sigma$ tension', ha='center',
         fontsize=10, color='white', style='italic')

plt.tight_layout()
plt.savefig('output.png', dpi=150, bbox_inches='tight',
            facecolor='#0a0a0a', edgecolor='none')
plt.close()
print("Figure saved: NFW tangential shear profile and S8 constraints")
print(f"\nCluster parameters:")
print(f"  M200 = {M200/Msun:.1e} Msun")
print(f"  c200 = {c200}")
print(f"  R200 = {R200/Mpc:.2f} Mpc")
print(f"  r_s  = {r_s/Mpc:.3f} Mpc")
print(f"\nDistances (z_L={z_lens}, z_S={z_source}):")
print(f"  D_L  = {D_L/Mpc:.1f} Mpc")
print(f"  D_S  = {D_S/Mpc:.1f} Mpc")
print(f"  D_LS = {D_LS/Mpc:.1f} Mpc")
print(f"\nCritical surface density:")
print(f"  Sigma_cr = {Sigma_cr:.3e} kg/m^2")
print(f"  Sigma_cr = {Sigma_cr * 1e-4:.2f} g/cm^2")
print(f"\nEinstein radius:")
theta_E = np.sqrt(4*G*M200/(c**2) * D_LS/(D_L*D_S))
print(f"  theta_E = {theta_E/arcsec:.1f} arcsec = {theta_E/arcsec/60:.2f} arcmin")
print(f"\nPeak tangential shear: gamma_t_max = {np.max(gamma_t):.3f}")
print(f"Peak convergence: kappa_max = {np.max(kappa):.3f}")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

References

Bartelmann, M. & Schneider, P. (2001). “Weak gravitational lensing,” Physics Reports 340, 291–472. arXiv:astro-ph/9912508. — The definitive review of weak lensing theory.
Kilbinger, M. (2015). “Cosmology with cosmic shear observations: a review,” Reports on Progress in Physics 78, 086901. arXiv:1411.0115. — Modern review of cosmic shear methodology and results.
DES Collaboration (2022). “Dark Energy Survey Year 3 results: Cosmological constraints from galaxy clustering and weak lensing,” Physical Review D 105, 023520. arXiv:2105.13549.
Asgari, M. et al. (KiDS Collaboration) (2021). “KiDS-1000 Cosmology: Cosmic shear constraints on the amplitude of matter fluctuations,” Astronomy & Astrophysics 645, A104. arXiv:2007.15633.
Hikage, C. et al. (HSC Collaboration) (2019). “Cosmology from cosmic shear power spectra with Subaru Hyper Suprime-Cam first-year data,” Publications of the Astronomical Society of Japan 71, 43. arXiv:1809.09148.
Planck Collaboration (2020). “Planck 2018 results. VI. Cosmological parameters,” Astronomy & Astrophysics 641, A6. arXiv:1807.06209.
Kaiser, N., Squires, G. & Broadhurst, T. (1995). “A Method for Weak Lensing Observations,” The Astrophysical Journal 449, 460. — The foundational KSB method for shear measurement.
Clowe, D. et al. (2006). “A Direct Empirical Proof of the Existence of Dark Matter,” The Astrophysical Journal Letters 648, L109. arXiv:astro-ph/0608407. — The Bullet Cluster weak lensing analysis.
Navarro, J. F., Frenk, C. S. & White, S. D. M. (1997). “A Universal Density Profile from Hierarchical Clustering,” The Astrophysical Journal 490, 493. arXiv:astro-ph/9611107. — The NFW profile used in cluster lensing.
Euclid Collaboration (2024). “Euclid. I. Overview of the Euclid mission,” Astronomy & Astrophysics. arXiv:2405.13491. — The next generation of space-based weak lensing.

← Baryon Acoustic Oscillations Course Overview →

Survey	S₈ Value	Year	Area (deg²)
DES Y3	\(0.776 \pm 0.017\)	2022	4143
KiDS-1000	\(0.759^{+0.024}_{-0.021}\)	2021	1006
HSC Y3	\(0.776^{+0.032}_{-0.033}\)	2023	416
Planck CMB	\(0.834 \pm 0.016\)	2020	Full sky