Strong gravitational lensing, which can make a background source galaxy appears multiple times due to its light rays being deflected by the mass of one or more foreground lens galaxies, provides astronomers with a powerful tool to study dark matter, cosmology and the most distant Universe. PyAutoLens is an open-source Python 3.6+ package for strong gravitational lensing, with core features including fully automated strong lens modeling of galaxies and galaxy clusters, support for direct imaging and interferometer datasets and comprehensive tools for simulating samples of strong lenses. The API allows users to perform ray-tracing by using analytic light and mass profiles to build strong lens systems. Accompanying PyAutoLens is the autolens workspace, which includes example scripts, lens datasets and the HowToLens lectures in Jupyter notebook format which introduce non-experts to strong lensing using PyAutoLens. Readers can try PyAutoLens right now by going to the introduction Jupyter notebook on Binder or checkout the readthedocs for a complete overview of PyAutoLens's features.

We investigate how strong gravitational lensing can test contemporary models of massive elliptical (ME) galaxy formation, by combining a traditional decomposition of their visible stellar distribution with a lensing analysis of their mass distribution. As a proof of concept, we study a sample of three ME lenses, observing that all are composed of two distinct baryonic structures, a 'red' central bulge surrounded by an extended envelope of stellar material. Whilst these two components look photometrically similar, their distinct lensing effects permit a clean decomposition of their mass structure. This allows us to infer two key pieces of information about each lens galaxy: (i) the stellar mass distribution (without invoking stellar populations models) and (ii) the inner dark matter halo mass. We argue that these two measurements are crucial to testing models of ME formation, as the stellar mass profile provides a diagnostic of baryonic accretion and feedback whilst the dark matter mass places each galaxy in the context of LCDM large scale structure formation. We also detect large rotational offsets between the two stellar components and a lopsidedness in their outer mass distributions, which hold further information on the evolution of each ME. Finally, we discuss how this approach can be extended to galaxies of all Hubble types and what implication our results have for studies of strong gravitational lensing.

A defining prediction of the cold dark matter (CDM) cosmological model is the existence of a very large population of low-mass haloes. This population is absent in models in which the dark matter particle is warm (WDM). These alternatives can, in principle, be distinguished observationally because halos along the line-of-sight can perturb galaxy-galaxy strong gravitational lenses. Furthermore, the WDM particle mass could be deduced because the cut-off in their halo mass function depends on the mass of the particle. We systematically explore the detectability of low-mass haloes in WDM models by simulating and fitting mock lensed images. Contrary to previous studies, we find that halos are harder to detect when they are either behind or in front of the lens. Furthermore, we find that the perturbing effect of haloes increases with their concentration: detectable haloes are systematically high-concentration haloes, and accounting for the scatter in the mass-concentration relation boosts the expected number of detections by as much as an order of magnitude. Haloes have lower concentration for lower particle masses and this further suppresses the number of detectable haloes beyond the reduction arising from the lower halo abundances alone. Taking these effects into account can make lensing constraints on the value of the mass function cut-off at least an order of magnitude more stringent than previously appreciated.

The distribution of dark and luminous matter can be mapped around galaxies that gravitationally lens background objects into arcs or Einstein rings. New surveys will soon observe hundreds of thousands of galaxy lenses, and current, labour-intensive analysis methods will not scale up to this challenge. We develop an automatic, Bayesian method which we use to fit a sample of 59 lenses imaged by the Hubble Space Telescope. We set out to leave no lens behind and focus on ways in which automated fits fail in a small handful of lenses, describing adjustments to the pipeline that ultimately allows us to infer accurate lens models for all 59 lenses. A high success rate is key to avoid catastrophic outliers that would bias large samples with small statistical errors. We establish the two most difficult steps to be subtracting foreground lens light and initialising a first, approximate lens model. After that, increasing model complexity is straightforward. We put forward a likelihood cap method to avoid the underestimation of errors due to pixel discretization noise inherent to pixel-based methods. With this new approach to error estimation, we find a mean $\sim 1{{\%}}$ fractional uncertainty on the Einstein radius measurement which does not degrade with redshift up to at least z = 0.7. This is in stark contrast to measurables from other techniques, like stellar dynamics, and demonstrates the power of lensing for studies of galaxy evolution. Our PyAutoLens software is open source, and is installed in the Science Data Centres of the ESA Euclid mission.

The elliptical power-law (EPL) mass model of the mass in a galaxy is widely used in strong gravitational lensing analyses. However, the distribution of mass in real galaxies is more complex. We quantify the biases due to this model mismatch by simulating and then analysing mock {\it Hubble Space Telescope} imaging of lenses with mass distributions inferred from SDSS-MaNGA stellar dynamics data. We find accurate recovery of source galaxy morphology, except for a slight tendency to infer sources to be more compact than their true size. The Einstein radius of the lens is also robustly recovered with 0.1\% accuracy, as is the global density slope, with 2.5\% relative systematic error, compared to the 3.4\% intrinsic dispersion. However, asymmetry in real lenses also leads to a spurious fitted `external shear' with typical strength, $\gamma_{\rm ext}=0.015$. Furthermore, time delays inferred from lens modelling without measurements of stellar dynamics are typically underestimated by $\sim$5\%. Using such measurements from a sub-sample of 37 lenses would bias measurements of the Hubble constant $H_0$ by $\sim$9\%. The next generation cosmography must use more complex lens mass models.