Accurate computation of the screening of scalar fifth forces in galaxies

Screening mechanisms allow light scalar fields to dynamically avoid the constraints that come from our lack of observation of a long-range fifth force. Galactic scale tests are of particular interest when the light scalar is introduced to explain the dark matter or dark energy that dominates our cosmology. To date, much of the literature that has studied screening in galaxies has described screening using simplifying approximations. In this work, we calculate numerical solutions for scalar fields with screening mechanisms in galactic contexts, and use these to derive new, precise conditions governing where fifth forces are screened. We show that the commonly used binary screened/unscreened threshold can predict a fifth force signal in situations where a fuller treatment does not, leading us to conclude that existing constraints might be overestimated. We show that various other approximations of the screening radius provide a more accurate proxy to screening, although they fail to exactly reproduce the true screening surface in certain regions of parameter space. As a demonstration of our scheme, we apply it to an idealised Milky Way and thus identify the region of parameter space in which the solar system is screened.


Introduction
Attempts to resolve long-standing complications within the dark sector often introduce a new scalar field which naturally couples to the gravitational sector, causing modifications to standard gravity, General Relativity (GR), in the form of scalar fifth forces.However, such deviations from GR are tightly constrained by lab-based experiments testing the inversesquare law and tests of Solar System dynamics [1].To avoid such constraints, scalar-tensor theories with screening propose that the scalar field is environment-dependent, allowing for modifications to GR on certain scales whilst suppressing fifth forces on other scales [2][3][4][5][6][7][8].
However, investigating these scalar fifth forces on galactic scales is a challenging problem: it is difficult to isolate any potential fifth force signal from the numerous astrophysical effects that can obscure it, and forward-modelling the signal involves the solution of a non-linear equation of motion (EoM).As a consequence, much of the literature has relied on approximations to solve the equations of motion, assuming a simplistic setup and dropping terms to simplify the non-linear dynamics.The most common approximation is to simply assume that a galaxy, depending on its mass and environment, is either fully screened or fully unscreened.This disregards the possibility of partial screening, under which the central region of the galaxy is screened while its outskirts remain unscreened.An example of a work using this approximation is that of ref. [27], which uses it to derive the current strongest constraints on chameleon f (R) gravity, altogether ruling out the astrophysically-relevant part of parameter space 1 .The purpose of the present work is to test these approximations, and the possible impact that relaxing them might have on such constraints.
In this work, we investigate the non-linear nature of screening theories in the galactic regime by utilising a numerical solver, based on refs.[31,43,44], to determine the true scalar field profile inside our static model galaxies.Guided by these field solutions, we obtain empirical rules describing the locations of screening surfaces and compare these results with the results of the common approximations.We discuss the implications of the proper nonlinear treatment of these theories and how our findings can impact observables and constraints established using simplified methods.These results will form the foundation for a reevaluation of galaxy-based constraints in the near future, using methods that we have explored and tested in this paper.This reanalysis will enable us to provide a quantitative estimate of the impact of using various standard approximations to screening, compared to our precise numerical methods that consider partial screening.
The outline for this paper is as follows.In section 2 we present an overview of the general scalar-tensor theory, as well as our specific chameleon and symmetron models.In section 3 we describe our methods to solve the theories numerically, including the galactic model we use, the mechanics behind the numerical solver, and the criteria used to define the screened and unscreened regions.Section 4 outlines the results of our investigation.Our conclusions are laid out in section 5.

General Scalar-Tensor Model
In the present study, our investigation into modified theories of gravity is exclusively focused on scalar-tensor theories.These theories are described by the Einstein frame action, given by where g is the determinant of the Einstein frame metric, g µν , R is the Ricci scalar, G is the gravitational constant, ϕ is a canonically normalised (dimensionless) scalar field with a potential V (ϕ), and S m represents the action for the standard model fields, ψ SM i .The theories also include a conformal coupling to matter, meaning that matter fields move on geodesics of the Jordan frame metric gµν = A 2 (ϕ)g µν . (2.2) Deviations from GR arise from the non-minimal coupling, A(ϕ).Although both frames have identical classical observables, for ease of calculation, this analysis focuses solely on the Einstein frame, where there is a direct coupling between matter and the scalar field, and the pure gravitational action does not depend on the scalar field.Scalar-tensor theories of gravity modify gravitational forces by introducing fifth forces, which can be obtained by taking the Newtonian limit of the geodesic equations (see for example ref. [5]) to arrive at the following expression for the acceleration where β(ϕ) ≡ d ln A/dϕ.Meanwhile, the scalar field's EoM is derived by extremising the general scalar-tensor action, eq.(2.1), with an energy-momentum tensor for the matter fields that assumes matter is a non-relativistic, perfect fluid with mass density ρ where the effective potential, V eff (ϕ), is given by For certain choices of V (ϕ) and A(ϕ), the effective potential has a stable minimum and this density-dependent EoM admits screening solutions, i.e. solutions in which the fifth force, eq.(2.3), is rendered negligible in certain environments but unleashed elsewhere.Screening can arise via several mechanisms: a small matter coupling, β(ϕ), a large scalar mass resulting in short-ranged forces, or not all the mass sourcing (or coupling to) the scalar field [5].

Astrophysical Screening Approximations
The crux of this paper is to determine the location of the screened region inside galaxies, defined by what we term the 'screening surface': the interface between the screened central region and the unscreened outskirts of a partially screened object.To determine this, we use numerical solutions of scalar EoM.These results will then be compared with the approximate methods that are commonly used in the literature.
One prevalent method to approximate screening surfaces considers a static, sphericallysymmetric source, δρ(r), in a constant density background, ρ, e.g. a star embedded in the galactic medium or a dark matter halo embedded in the cosmological background, as in refs.[45][46][47][48].For such an object, the screening surface will necessarily be a sphere with radius r s , the 'screening radius'.The general EoM, eq.(2.4), is then solved in a piecewise manner, described in full in appendix A, to obtain the expression where φ = ϕ(ρ) is the background field solution and Φ N (r) is the Newtonian potential.To obtain this χ parameter, hereby referred to as the critical potential, several assumptions are required, specifically that the field reaches its minimum value inside the object, that this minimum value can be neglected, that the theory can be linearised in the region containing the screening radius, that the scalar mass is negligible in this region, and that the coupling function β(ϕ) is constant throughout.If eq.(2.6), or equally eq.(2.7), has no solution for r s this corresponds to a fully unscreened object, i.e. none of the object is screened.Since Φ N < 0 and Φ ′ N > 0, an object must be fully unscreened when χ > |Φ N |/c 2 , where the Newtonian potential is measured at the surface of the object.In reverse, this indicates only objects with χ < |Φ N |/c 2 are partially screened (or fully screened, i.e. r s is larger than the extent of the object, in the case The critical potential is a starting point for two approaches commonly seen in the literature for approximating the screening of astrophysical objects.In the first approach, the implicit equation (2.6) is solved for the screening radius, which appears as the lower limit of the integral.As we shall see, this can readily be done numerically, but analytic solutions can also be found for certain choices of the density distribution (e.g., for idealised dark matter haloes; as in ref. [49]).
The second approach, which we refer to as the 'binary screening condition', disregards partial screening altogether and classifies objects as either fully screened or fully unscreened by comparing χ to an estimate of Φ N at the surface of the object.Neglecting the contribution of its environment, an object is classed as screened if |Φ N |/c 2 > χ, and unscreened otherwise.For example, at the level of current constraints in f (R) gravity, refs.[24,25,27,32], objects with deep potentials, e.g.main sequence stars with |Φ N |/c 2 ∼ 10 −6 are likely screened, while objects with shallower potentials e.g.dwarf galaxies with |Φ N |/c 2 ∼ 10 −8 can still be unscreened.In the context of galaxies, this approach was suggested by ref. [50], who proposed using Φ vir = GM vir /R vir , M vir and R vir will be defined in sec.3.1.The same approach was then subsequently used in a string of recent papers deriving modified gravity constraints on galactic scales (refs.[24,25,27,51]).However, it is worth noting that these works also included environmental contributions, classifying objects as screened when where |Φ env.| is the environmental contribution to the local gravitational potential.
One of the main aims of the present work is to compare the results obtained from these two approximate approaches with more accurate results obtained from exact numerical solutions of the scalar equations of motion.

Chameleon f (R)
Chameleon models are a class of scalar-tensor theories that screen primarily via a densitydependent scalar mass, leading to the presence of short-range fifth forces.Screened regions in chameleon models are characterised by flat scalar field profiles, which result in the suppression of the fifth force, eq.(2.3).We choose to work with f (R) chameleon screening because of its simplicity as a model and its ubiquity in the literature.
f (R) theories modify GR by replacing the Ricci scalar, R, in the Einstein-Hilbert action with the generalised R + f (R), allowing for higher-order curvature terms in the action (see ref. [52] for a general review), where R = R(g).Here, f R ≡ df (R) dR acts as the scalar field of the theory.This action can be recast to the general scalar-tensor action eq.(2.1) with the field redefinition [53] (2.9) For specific choices of the functional form of f (R), the theory can exhibit the chameleon screening mechanism.Several forms of f (R) can achieve this, refs.[54][55][56], we adopt the commonly used Hu-Sawicki model, [57].In this model where a and b are positive dimensionless parameters, and m is a dimensionful parameter.For simplicity, we choose b = 1, in line with the bulk of the current literature.The two remaining free parameters a and m can be related by requiring that gravity returns to GR + ΛCDM in the high curvature limit (i.e.f (R) ≈ −2Λ when R ≫ m 2 ), resulting in where Λ is the cosmological constant.This allows us to describe the theory in terms of one free model parameter, a or m, which we instead express as the present-day background field value f R0 : where Ω m and Ω Λ are the cosmological density parameters for matter and dark energy respectively.Throughout this article we assume values of Ω m = 0.3 and Ω Λ = 0.7.
Extremising the action eq.(2.8) with respect to the metric results in a set of modified Einstein field equations.Taking the trace and then the Newtonian limit, with the assumption |f R | ≪ 1 and the quasi-static approximation |∇f R | ≫ ∂ t f R , produces the f (R) EoM, Here, δρ is the density perturbation and δR is the curvature perturbation, given by under the Hu-Sawicki model, with zero subscripts denoting background values today.The fifth force acceleration, eq. ( 2.3), caused by our scalar field is given by where we have assumed f R ≪ 1, and the background Compton wavelength is given by as in ref. [50], based on ref. [57].Lastly, using our field redefinition eq.(2.9) we can define the critical potential eq.(2.6) in terms of our background scalar field value (2.17)

Symmetron
In contrast with chameleon models, where screening arises due to a variable scalar mass, symmetron models screen by suppressing the coupling to matter in dense regions.This is achieved via spontaneous symmetry breaking (SSB), given the correct form of potential and coupling in the scalar-tensor action, eq.(2.1).In our case, we choose to work with a quartic potential and quadratic coupling, where M 2 pl ≡ ℏc/8πG is the reduced Planck mass and L 2 pl ≡ 8πGℏ/c 3 is the reduced Planck length.In the coupling function above, and throughout the rest of this section, we have assumed that ϕ ≪ M sym /M pl to suppress higher order terms in A(ϕ).This choice of potential and coupling results in an effective potential eq.(2.5) with the necessary SSB properties.This model has three free parameters: the dimensionless λ governing the scalar self-interactions, the mass scale M sym controlling the strength of the coupling to matter, and µ, representing the bare tachyonic mass of the scalar field.
Extremising the action, eq.(2.1), with respect to the scalar field gives us the static symmetron EoM, for non-relativistic matter, where L c and ρ SSB are the background Compton wavelength and SSB density scale, expressed as It is possible to reduce this three-parameter EoM to a two-parameter one by rescaling the field by its vacuum expectation value.We find the vacuum field value by considering solutions to eq. (2.20) in an infinite box of density ρ.When ρ > ρ SSB the field is screened and ϕ = 0.If ρ < ρ SSB , then our scalar field falls into one of two solutions where ϕ vev is the vacuum expectation value given by Now, rescaling the scalar field, φ = ϕ/ϕ vev , in eq.(2.20) returns the simpler two-parameter EoM, The acceleration a test particle experiences due to the scalar fifth force is given by (2.26) Note that despite reducing the EoM to two parameters, the theory remains a three-parameter one as the scaling of the fifth force with ϕ2 vev reintroduces a dependence on λ.The critical potential eq.(2.7) for this symmetron model is (2.27) Note that χ depends only on M sym .However, as we shall see later, screening is a function of both M sym and L c .This discrepancy is a result of the approximations used to derive the general form of χ, eq.(2.6); specifically a combination of assuming the field can be linearised at the screening surface, the background scalar mass is light and thus negligible, and that β(ϕ) is constant throughout.
The symmetry-breaking nature of the symmetron provides a third possible approximation for the screening surface -in addition to the two described in section 2.2 for scalartensor theories in general.One can assume that the screening surface universally corresponds to the surface where ρ(r = r s ) = ρ SSB , neglecting other relevant factors such as the coupling strength and Compton wavelength of the theory.We will also examine the efficacy of this approximation in section 4.

Numerical Screening Solutions
This section describes the numerical methods we use to obtain isolated galaxy solutions for the scalar field equations of motion under the chameleon f (R), eq.(2.13), and symmetron theories, eq.(2.25).Section 3.1 describes our choice of galaxy model, section 3.2 then describes our numerical solver for the equations of motion, and section 3.3 outlines the conditions we use to identify the screening surface within our scalar field solutions.
Our methods are summarised in figure 1, the individual panels of which show a typical density profile under our galaxy model, the coordinate grid of our scalar field solver, the scalar field solutions for the same galaxy, and the screening surfaces identified from the scalar field solutions.These various panels of the figure will be examined in more detail in the corresponding parts of this section.

Galactic Model
Given the f (R) and symmetron EoMs, equations (2.13) and (2.25) respectively, we can solve the scalar field profile given a suitable galactic density model.For this, we use a twocomponent axisymmetric model comprising a dark matter halo and a galactic disc (note that we are limiting our study to late-type spiral disc galaxies, not elliptical or dwarf spheroidal galaxies).We model the dark matter halo with a Navarro-Frenk-White (NFW) profile [58]  where r is the spherical radial coordinate, and ρ NFW and r NFW are the characteristic density and length scale.The galactic disc is represented by a simple double exponential profile2 [59] where R represents the cylindrical radial coordinate and z represents the coordinate perpendicular to the galactic disc, and Σ disc , R disc and z disc are the characteristic density and length scales.In actuality, a more sophisticated galaxy model (e.g.[60]) could also include various other components, such as central bulges, stellar haloes, and multiple discs; however, we have chosen to exclude these for the sake of keeping our model as simple as possible while still capturing the key characteristics.Including these density components will introduce five extra input parameters used to scale the density profiles: ρ NFW , r NFW , Σ disc , R disc and z disc .Nonetheless, for the remainder of this article, we shall quote our galactic model in terms of one input parameter: the virial mass, M vir .This is defined as the mass enclosed within the virial radius, R vir , the radius within which the mean density inside is 200 times the cosmological critical density.We reduce the five-dimensional parameter space into a one-dimensional space by using empirical relations between the various parameters.In doing so, we typically adopt median values on the relations, ignoring any scatter.Thus, for a given M vir , the galaxy we model is in some sense a maximally typical one.We describe this process in further detail in appendix B. The top left panel of figure 1 shows an example density profile for a galaxy with M vir = 10 11.5 M ⊙ , corresponding to the following values for the parameters appearing in our density model: The adoption of the NFW profile for the dark matter halo presents two numerical challenges: the central density singularity and the logarithmically divergent mass.To resolve these issues, we set a constant density within the central 50 pc, and introduce an outer cutoff at 2.2R vir .A sharp cutoff at the halo's edge is well motivated astrophysically, e.g.ref. [61].We tested our analysis for various inner and outer cutoffs.In changing the inner cutoff neither the derived position of the screening surfaces, nor the conclusions we make are altered.We found that increasing the outer cutoff moved the screening surface farther outward, primarily due to the increase in the total mass of the system.However, we find this shift is mirrored in true and approximate solutions, and thus our conclusions remain robust.
In this work, we assume that the galaxy is embedded in a constant-density background equal to the cosmological critical density.Therefore, we do not consider environmental effects, such as inhomogeneities in the matter distribution throughout the Universe, which can manifest as regions of lower density (voids) and higher density (clusters, walls, and filaments), as well as neighbouring galaxies and satellites.This simplification enables us to better investigate other variables that impact screening in galaxies.Overall, higher-density environments tend to enhance screening effects, resulting in more fully screened galaxies and larger screened regions in partially screened galaxies, while the opposite is true for lower-density environments.We comment on where this will impact our results throughout section 4, and plan to account for these effects in follow-up work.

Numerical Solver
Using the axisymmetric galaxy density profiles from the previous section, we solve the f (R) and symmetron equations of motion, given by eqs.(2.13) and (2.25) respectively.We do this using a numerical solver based on that employed by ref. [31], updated to incorporate some optimisations suggested by [44] and to include the symmetron.The solver is ultimately a two-dimensional version of the f (R) solver used in the cosmological simulation code mggadget [43]: it constructs a coordinate grid, discretises the EoM on this grid, then solves this discretised EoM using an iterative Gauss-Seidel method.
Although a cylindrical polar coordinate system might seem the natural choice for an axisymmetric system, we use a two-dimensional spherical polar system (neglecting the azimuthal coordinate).This is because the radial range of our grid should stretch from the galactic centre to several Compton wavelengths (eqs.(2.16) and (2.21)) beyond the galaxy's edge; the dynamic range is such that a spherical coordinate system better enables us to use a variable grid resolution with smaller cells in the central region and larger cells further out.Specifically, we use a logarithmic radial discretisation, i.e. cells are equally spaced in the coordinate u ≡ ln(r/[chosen length scale]), producing higher resolution at smaller radii, as required.Similarly, in the polar direction, θ, because the galaxy has a disc, we need a higher resolution near the disc-plane θ ≈ π/2 than at the poles.We can achieve this by discretising in the coordinate v ≡ cos θ.In principle, one could assume mirror symmetry about the disc plane and restrict v to the range [0, 1], but we are often interested in the behaviour of the scalar field in the disc plane, and so wish to avoid pushing the disc plane to the boundary of the grid.We thus discretise v in the range [−1, 1].A visualisation of the coordinate grid, out to a maximum radius of 5R disc , is given in the top-right panel of figure 1.
In u−v coordinates, the Laplace operator acting on a function f (dropping the azimuthal term) is For later convenience, we have only expressed the derivatives in terms of u, while keeping other factors that contain r, as the conversion to u is straightforward.We construct a regular grid in u − v space.We use 512 radial (u) cells in the range [0.05, 10 4 ] kpc, and 201 angular (v) cells.The maximum radial extent of 10 4 kpc is several times larger than any galaxy scale considered in this work, so the field can reach the background value outside of the galaxy.By using an odd number of angular cells, we ensure that one row of the grid is centred on the disc plane.Each grid cell is labelled with indices (i, j), where i represents the u-coordinate and j represents the v-coordinate, i.e. u = u i and v = v j at cell (i, j).All other quantities (density, scalar field etc.) can be discretised on this grid.For example, we can write ρ ij to represent the density at the centre of cell (i, j).
The Laplace operator acting on a function f , eq. ( 3.3), can then be discretised as where h u and h v are the (constant) grid spacings in u and v respectively, and s ≡ sin θ, so that s 2 j = 1 − v 2 j .Note that half-integer subscripts (e.g.r i+1/2 ) indicate quantities evaluated at the cell edges rather than the cell centres.
We can now proceed to discretising the EoMs, eqs.(2.13) and (2.25).In the f (R) case, we follow reference [44] in redefining the scalar field as φ ≡ (f R /f R0 ) 1/2 .Consequently, the discretised forms of the symmetron and f (R) EoM both turn out to be cubic equations of the form where the coefficients p ij and q ij are given in the f (R) case by and in the symmetron case by and the objects A ij and B ij are given in turn by where ψ represents the argument of B i,j , i.e. φ 2 in the f (R) case and φ in the symmetron case.
Before proceeding to solve eq.(3.5) we must specify our boundary conditions and initial field profile guesses, which we will iteratively converge to the true solution.Our grid has four boundaries: v = −1, 1 (i.e.θ = π, 0) and u = ln r min , ln r max .At the first three, we impose ∂φ/∂v, ∂φ/∂u = 0.At the outer radial boundary, we instead impose φ = φ ∞ , i.e. the cosmic background value of the scalar field.In the f (R) case, this is simply f R0 ; in the symmetron case, this is the field value that minimises the effective potential eq.(2.23) evaluated at the cosmic matter density.For initial guesses, in the f (R) case we simply set φ = 1 everywhere, and allow the solver to converge on the correct profile.For the symmetron case, we speed up the convergence by setting the field to 1 everywhere, except in the region where ρ > ρ SSB which is assumed to be screened and thus set to 0.
We solve eq.(3.5) using a Gauss-Seidel iteration method.Here, our 2D grid is decomposed into a chequerboard pattern, i.e. we alternately assign cells with 'red' and 'black' labels.Each iteration n is then comprised of two halves: a red sweep and a black sweep.During the 'red sweep', the scalar field value is changed in all of the red cells using a Newton-Raphson update, where L ′ ≡ ∂L/∂φ = 3φ 2 + p.Note that L i,j here is evaluated using φ n−1 , i.e. the values of φ calculated at the previous Gauss-Seidel iteration (or the initial guess in the case n = 1).
Next, we perform a 'black sweep' in which the scalar field is updated in all of the black cells using the same expression.There is however one key difference in the black sweep: rather than using the leftover values from the previous Gauss-Seidel iteration (n − 1) as was done during the red sweep, one uses the neighbouring red cell values that have just been updated during the red sweep, i.e. in the first half of iteration n.
This iterative process repeats, updating the field values at each grid point until a converged solution is obtained.The iteration stops when the difference between the field in the previous step and the current one is less than a given threshold.In the f (R) case we use a threshold of 10 −7 , and in the symmetron case we use 10 −14 .We find these values give robust, converged results.
The left semicircles in the bottom row of figure 1 show outputs of the solver for the same example galaxy shown in the upper left panel of the figure, with a mass M vir = 10 11.5 M ⊙ , theory parameters of f R0 = −10 −6.4 in the f (R) case, M sym = 10 −4.5 M pl and L c = 10 −1 kpc in the symmetron case.In both cases, the field solutions are plotted to a maximum radius of 5R disc .

Screening Conditions
This section aims to establish a universal and reliable screening condition to locate the (generally non-spherical) screening surface, r s (θ), in a galaxy's scalar field profile.We will find that our conditions to calculate r s (θ) are given by c 2 δR(r, θ) 8πGδρ(r, θ) r=rs(θ) = 0.9 , (3.10a) for the f (R) and symmetron models respectively.In other words, for a given angle, the screening surface is the surface at which the f (R) curvature-density ratio or the symmetron field value matches a given threshold.The right sides of the bottom row panels in figure 1 display the f (R) curvature-density ratio and the symmetron field with a colour map truncated at the threshold.The sharp transition in both of these cases vindicates our chosen threshold.We arrived at the screening conditions presented in this section, after considering a variety of different properties of the scalar field solutions, for example the fifth force profiles.The screening conditions discussed here are the ones that could be most consistently applied across the range of model and galaxy parameters explored in this work.
In the f (R) case, the screening condition, eq.(3.10a), is a threshold on the curvaturedensity ratio, motivated by the two terms on the right-hand side of the f (R) EoM eq.(2.13).Specifically, inside the screened region the field is flat and the Laplacian is zero, thus the two terms on the right-hand side of the EoM are approximately equal in magnitude, and the curvature-density ratio is unity.Conversely, in the unscreened regime, the field adopts its background value, causing δR, eq.(2.14), and thus the curvature-density ratio to become zero.We find that the curvature-density ratio provides a clearer demarcation between the screened and unscreened region than the scalar field, as can be seen in the bottom left plot of figure 1.Nevertheless, there are still instances where the transition is more gradual, particularly for smaller values of |f R0 |.We find that a threshold on the curvature-density of 0.9 correctly separates the screened/unscreened regions and produces universal behaviour across the whole parameter space.
For the symmetron condition, we use a threshold value on the field itself.This is motivated by the coupling function, eq.(2.19), having a first derivative proportional to the scalar field; therefore the density-dependent term in the EoM is then also proportional to the scalar field.Setting this threshold as a ratio of the background scalar field allows for a screening condition that can be universally applied over the whole parameter space.We determine that a threshold value of 10% of the background scalar field value provides consistent results, with the scalar field profile being visually flat, and thus suppressing fifth forces, within a substantial region of r s .
Figure 2 shows the value of the f (R) curvature-density parameter and the symmetron field along the plane of the galactic disc, as well as the magnitude of the fifth forces of both theories.We rescale the radial coordinates by the derived r s along the disc plane to better visualise the behaviour over the whole parameter space.Taking figures 1 and 2 together, we are able to justify our choice of thresholds, namely 0.9 for the curvature-density parameter and 0.1 for the symmetron field, finding that this choice allows for consistent results, that suppress the fifth forces inside of the screened region, across our full range of parameters.
We classify a galaxy as fully screened if the screening surface determined according to eqs.(3.10) is greater than the extent of the galaxy, or if SSB has not yet occurred in the symmetron case, i.e. ρ > ρ SSB .Conversely, a galaxy is considered fully unscreened if there is no screened region within the galaxy, i.e., if the conditions given by eq.(3.10) are not satisfied anywhere.However, we also find a small subset of galaxies for which the screening conditions are satisfied somewhere but nonetheless have non-negligible fifth forces in their central regions.To catch these cases, we consider the innermost field value.In the f (R) case, solutions satisfying f R (r ≈ 0) > 10 −3 f R0 are fully unscreened.We find that an analogous threshold on the field central field value can not be used in the symmetron case, as some field solutions yield small central field values but quickly grow away from the centre, giving an unscreened central region.So, for the symmetron, we instead use the Laplacian of the field: solutions satisfying the condition ∇ 2 φ(r ≈ 0) > 10 −1 max ∇ 2 φ are fully unscreened.
The symmetron solutions that avoid the central field threshold condition (i.e.satisfying ϕ(r ≈ 0) < 10 −3 ϕ ∞ ) but are correctly classified as fully unscreened by the central Laplacian threshold condition are highlighted in figure 2.Here we see that these edge case field profiles grow to be large well within the screening radius, causing a non-zero fifth force in this region, as shown in the bottom panel.4 The Screening of Galaxies

Position and Shape of Screening Surfaces
Figure 3 shows an array of solutions for the screening surface (dashed lines), over a range of f (R) and symmetron parameters, for a galaxy with M vir = 10 11.5 M ⊙ .Generally, decreasing the adopted value of |f R0 | or M sym , or increasing L c , pushes the screening surface further out.We separate these solutions into those for which the screening surfaces lie within the discdominated region (left-hand semi-circles) and the halo-dominated region (right-hand semicircles).The figure shows that screening surfaces that occur in the disc-dominated region are elongated along the galactic disc, whereas larger screening surfaces are spherical and trace the dark matter component.
We also include approximate screening surfaces (dashed-dotted lines) for comparison.In the case of f (R), the approximate screening surface is computed by numerically solving eq.(2.6), where χ is given by the f (R) critical potential, eq.(2.17), and δρ(r) is the galactic density averaged over spherical shells.The approximate symmetron screening surfaces are isodensity contours with a value equal to the symmetry breaking scale, i.e. ρ(r = r s (θ)) = ρ SSB .
In the f (R) case (top left panel), the true screening surfaces (dashed lines) in the discdominated region follow the isodensity contours away from the disc plane, but deviate nearer the disc, where the isodensity contours are flattened to a larger extent.We find that the screening surfaces lie somewhere between the highly flattened isodensity contours of our galactic model and the isopotential contours (not shown), which are visually more spherical than the screening surfaces.
The non-spherical nature of the screening surface in the disc-dominated region of the galaxy (left-hand side of the f (R) plot) indicates that using approximations (dashed-dotted lines) that result in spherical screening radii, such as eq.(2.6), can fail to capture the true nature of the screened region.This is especially relevant when searching for observables along the plane of the disc, where r s will be underestimated as the spherical screening surface meets the disc plane at a smaller radius than the true, flattened screening surface.As a result, for a given value of the theory parameter, one may expect the presence of fifth force observables at a certain location in the disc plane because it lies outside the assumed spherical screening surface, whereas, in reality, the location lies within the true, flattened, screening surface, and so the fifth force is suppressed.This in turn could lead to exaggerated constraints on the theory parameters.
In the halo-dominated region of the galaxy (right-hand side of the f (R) plot) we see the approximation overestimates the screening radius for smaller values of |f R0 | -corresponding to larger values of r s , and smaller Compton wavelengths.This leads to an overestimation of the amount of dark matter that is screened by this approximation.This behaviour is likely due to a breakdown of the approximations utilised to derive the critical potential at large r s values, namely that the field can be linearised around the screening surface, and that the scalar mass is negligible in this region.Specifically, when the Compton wavelength is smaller than the extent of the galaxy, not all the mass within the galaxy impacts the field profile, whereas the approximation, which assumes negligible scalar mass and thus infinite Compton wavelength, responds to mass in the extremities of the galaxy.The effect that this overestimation of r s has on observables is discussed further in section 4.2.
By contrast with the f (R) screening surfaces, the screening surfaces under the symmetron theory (bottom plots) follow the isodensity contours very closely.However, this behaviour breaks down in the large L c , large M sym regime, as seen with the larger screening radii in the left-sided symmetron plots, where the screening surface fails to trace the finer details of the isodensity contour, especially near the disc plane.The reasoning for this breakdown is due to the physics associated with each parameter.The Compton wavelength, L c , is the characteristic scale on which the field can vary, i.e. the field will not respond to changes on scales much smaller than the Compton wavelength, but in the large L c regime the scale height of the disc profile is much smaller than L c .In other words, a larger Compton wavelength responds to mass averaged over a larger distance, causing the field profile to smooth out the finer details of the disc density profile.Meanwhile, for the coupling constant, M sym , the coupling function, and thus the density-dependent term in the EoM, eq.(2.25), have an inverse proportionality to M 2 sym , meaning a large change in density has less impact on the field profile.
The approximation of symmetron screening surfaces as isodensity contours matching ρ SSB is a natural consequence of the highly sensitive relationship between the screened region and the symmetry-breaking scale.These approximate solutions are shown as grey dasheddotted lines on the symmetron plots of figure 3.For the most part, this approximation does a good job of calculating the true value of r s , especially in the small M sym , small L c regime.Outside of this region of parameter space, the approximation always overestimates the value of r s when compared to the true value.This overestimation is especially prevalent in the large L c , large M sym regime.

Amounts of Galactic Mass Screened
One of the main results of our analysis is contained within figure 4, where we show the percentage of mass contained within the screening surface over a range of galactic and theory parameters.The left column represents the disc mass and the right column is the halo mass.In the symmetron plots, we vary one parameter with the other parameter fixed.This figure allows us to study how our system transitions from screened to unscreened across the parameter space, which we can compare with the transition from screened to unscreened given by the binary screening condition, eq.(2.17) in the f (R) case and eq.(2.27) in the symmetron case.As a proxy for the 'true' answer, we include a dashed contour line to show the parameter boundary where 50% of the mass is screened in each case 3 .This line can then be directly compared with the binary screening condition boundary, which we represent with a dotted line in the f (R) and symmetron M sym plots respectively (not included in the L c plot since the symmetron critical potential is a function of M sym only).
All the sub-figures contained in figure 4 exhibit some similar behaviours.The galactic disc component transitions sharply from being fully screened to fully unscreened, with the two regimes separated in the parameter space of galaxy and model parameters by only a narrow band in which the disc is partially screened.In contrast, the fraction of screened dark matter varies gradually, due to the significant dynamical range of the halo, which spans several orders of magnitude in density.Another notable feature is the discrepancy between the regions considered fully screened and unscreened by the binary condition and the region where half of the mass is screened.This is explored further throughout this section.In all cases, the r s approximations (dashed-dotted lines) are closer to the true 50% screened solution than the binary screening condition.
Focusing first on the f (R) case, we see that the f (R) binary screening condition, eq.(2.17), does not correctly characterise the screening of the galactic disc.Specifically, the region of parameter space considered 'fully unscreened' by the binary condition (above the dotted line in the top left panel) includes roughly an order of magnitude in which the galactic disc is actually entirely screened.This ultimately means that searching for signatures of screening theories in 'fully unscreened' galaxies may lead to false negative results since the visible stellar and gas components are fully screened.For example, when looking at stellar-gas offsets, ref. [24,26], in this region of parameter space the disc is fully screened and there is therefore no displacement between the stellar and gas components.However, using eq.(2.17) as a binary screening threshold would lead one to assume that the disc is unscreened and thus incorrectly predict the presence of an offset.In practice, this could entail relaxing the current tightest constraints.
By contrast with the disc, for the halo the binary screening condition kicks in at a point where the halo is only partially screened, so there will be haloes (below the dotted line in the top right panel) which will be classed as screened but actually carry a substantial unscreened component in their outskirts, as much as 60-70% of the total halo mass.Similarly, there will be haloes (above the dotted line) which will be classed as unscreened but actually have screened central regions, comprising as much as 30-40% of the total halo mass.This is the fundamental flaw in the binary screening condition: it assumes a step change in a quantity which actually changes gradually over the course of several orders of magnitude in parameter space.In summary, the f (R) binary screening condition underestimates unscreened mass when it classifies an object as screened, and it underestimates screened mass when it classifies an object as unscreened.
The various panels in figure 4 also include dashed-dotted lines to indicate when 50% of the mass is screened when r s is derived by approximate means.In the f (R) case (top row), this is when r s is derived by numerically integrating eq.(2.6) for a spherically averaged version of our galaxy density profile.This approximation solution yields results closer to the exact solution than the binary screening condition.For the galactic disc (top left panel), this approximation contour is almost identical to the equivalent numerical solution contour, and the difference in the proportion of mass screened between the two calculations is no more than 5% across the range of parameters we consider.However, as we saw in figure 3, the true and approximate solutions for r s can have vastly different shapes and the spherical nature of the approximate screening surface results in an underestimation of r s along the galactic disc and an overestimation along the perpendicular axis.
For the dark matter component (top right panel), the approximate value of r s , and therefore the percentage of mass screened, is very close to the actual value for low values of r s .However, as r s becomes larger, the approximate solution begins to overestimate r s compared to the actual solution, as previously seen in figure 3.This overestimation occurs in the limit of large M vir and small f R0 , and emerges because, in this regime, the Compton wavelength becomes smaller than the extent of the galaxy.Given that the approximation assumes a negligible scalar mass and thus an infinite Compton wavelength, it considers mass beyond a Compton wavelength (which does not significantly impact the field).The biggest discrepancy between the approximate and true screening surfaces, located in the lower right of the f (R) parameter space of figure 4, results in the estimated proportion of screened dark matter being as much as 20% larger than the actual solution.In other words, if the 50% threshold were increased to, e.g., 80%, then the two lines would be further separated.Ultimately this will result in the predicted fifth forces sourced by the galaxy being larger than in the true case, but not to the same extent as the binary screening condition would predict.
For the symmetron model (middle and bottom rows), the binary screening condition is much less accurate.As discussed in section 2.4, the condition depends only on M sym , but even in the middle row, with varying M sym and fixed L c , this estimate of screening disagrees with our solution by several orders of magnitude.Looking at the row with varying L c and fixed M sym , it is clear that the location of the screening surface is a strong function of L c .Therefore, this simplified version of the symmetron binary screening condition should be avoided.
Instead, a more useful approximation for screening is obtained by assuming the screening surface is coincident with the surface where the total density equals ρ SSB .The parameter space that results in 50% of each mass component being enclosed by this approximate screening surface is shown in the symmetron panels of figure 4 as a dashed-dotted contour.As discussed in section 4.1, this approximation will never underestimate the value of r s , and will overestimate r s for large values of M sym or L c .Looking first at the disc mass, in both cases (varying M sym and varying L c ) the approximation contour follows a similar behaviour to the true contour, but somewhat overestimates the amount of disc mass screened.Notably, the distance in parameter space between the true and approximate 50% screened contours is reduced in the case where L c is smaller (towards the right of the bottom-left panel with varying L c and fixed M sym ).For the dark matter component, we see that the approximate and true 50% screening lines are indistinguishable in the middle right panel with varying M sym and fixed L c .This is due to the fact that L c is fixed at a low value of 10 −1 kpc, and so to get a large screening radius that screens 50% of the halo we need a small value of M sym , thus falling into the regime where this approximation works best (small M sym , small L c ).In the bottom right panel with varying L c and fixed M sym , the 50% line is a good match for larger galaxies, but for smaller galaxies the approximation overestimates the amount of mass screened.As length scales get smaller they become comparable with L c and the field can not respond to the density changes quickly enough to screen exactly along the isodensity contour.In summary, the approximation does not account for the full effects of varying M sym and L c , so it falters in the regimes in which these effects are strongest.
One noteworthy caveat to this discussion is that we do not consider any environmental effects, e.g., the presence of over-/underdensities within a few Compton wavelengths of our target galaxies.In general, these would serve to vary the location of the screening surface of a given galaxy from that calculated in isolation.We have done some preliminary tests in this regard, and find the conclusions to this section are robust.We experimented with mimicking overdense environments by embedding our galaxies in large overdense spheres of constant density.As expected, the additional overdensity leads to the galaxies being screened more easily.However, we find that either the true and approximate solutions for the screening surface shift equally or the gulf widens.However, a full reanalysis with proper treatment of environmental effects is required to make a clear, quantitative conclusion.

Milky Way Observables
Figure 5 shows the derived screening radius along the plane of the galactic disc, over the entire model parameter spaces, for a Milky Way-like galaxy with M vir = 1.5 × 10 12 M ⊙ [60].The approximate distance from the Milky Way centre to the Solar System (8 kpc) is marked by the dashed lines.The dashed-dotted lines show the screening radius derived via approximate methods, namely, numerically integrating eq.(2.6) given the f (R) critical potential, eq.(2.17), or defining the screening surface as the isodensity contour equal to ρ SSB in the symmetron case.The dotted lines represent the binary screening boundary, separating the parameters that cause solutions to be considered 'fully screened' and 'fully unscreened'.Specifically, all |f R0 | values less than 10 −6.7 and all M sym values less than 10 −3.1 M pl are considered 'fully screened' by the binary screening conditions.In the f (R) plot, we again see the approximate r s begins to overestimate the true solutions when the Compton wavelength becomes comparable to the extent of the galaxy, at f R0 ≲ 10 −7 .
Stringent constraints based on tests of celestial gravitational dynamics require that the Milky Way screening radius extends to at least 8 kpc, encompassing the entire Solar System [62]. Figure 5 allows us to determine the parameter space permitted by such constraints, namely |f R0 | ≲ 10 −6.0 and the portion of the symmetron parameter space satisfying M sym ≲ 10 −3.6 M pl and to the right of the dashed line on the lower panel.There are a few caveats to these constraints.Firstly, they only apply to this particular galaxy model and choice of mass.We do not expect the Milky Way to be a maximally-typical galaxy but for consistency with our previous analysis, we have modelled it as such.In addition, there remains significant empirical uncertainty regarding the total mass of the Milky Way halo [63], which would have a non-negligible effect on the bounds we obtain.We have also neglected any environmental effects, such as the influence of the Local Group and the galaxy's wider cosmographic context.
Given the strict constraint on the minimum value of the Milky Way's screening radius, the question arises as to whether there can still be observable effects in the periphery of the galactic disc.For this galaxy, the disc component becomes subdominant along the disc plane after roughly 14 kpc.Looking at either panel of figure 5, between the r s = 8 kpc contour and the parameter space where r s = 14 kpc, there is only a small slice of parameter space in which the Solar System is screened but the bulk of the galactic disc is not.This means it would require a great deal of fortune to be in a universe adopting these parameter values and allowing us to explain observable effects in the disc via fifth forces.
Nevertheless, there remains several orders of magnitude of unconstrained parameter space in which the screening radius lies between the edge of the galactic disc and the virial radius of ∼ 300 kpc.As such, observable consequences of scalar fields with screening mechanisms located in the Milky Way's dark matter halo, are entirely plausible.For example, the kinematics of halo tracers or the detection of equivalence-principle-violating asymmetries in the leading and trailing tails of stellar streams, from tidally disrupted dwarf galaxies in the dark matter halo, as in refs.[28,[64][65][66].
An additional feature in figure 5 worth mentioning is the turn-off in the top right of the symmetron plot.Away from this region, lines of equal screening radius are diagonals across parameter space.These diagonals correspond with lines of equal ρ SSB (ρ SSB ∝ M sym /L 2 c ; cf.eq.(2.22)), and as we saw in figure .3, screening surfaces closely match the ρ = ρ SSB surfaces in the low M sym , low L c regime.The turn-off is another symptom of the breakdown of this behaviour in the large M sym , large L c limit.

Conclusions
In this work, we have investigated the effects of scalar modifications of gravity on galactic scales, for both chameleon f (R) and symmetron screening theories, with a particular focus on how screening differs between full numerical solutions and common approximations used to derive constraints.To determine the screened region, we solved the quasi-static scalar field EoM numerically, given an input density described by a maximally-typical galaxy consisting of a double-exponential disc profile and a dark matter NFW halo.We developed new and robust screening conditions to compute the position of the screening surface based on the scalar field values, eqs.(3.10).
We examined the shape of our calculated screening surfaces.We find that screening surfaces in the halo-dominated region follow the spherical isodensity contours.In the discdominated region, f (R) screening surfaces almost follow the isodensity contours, but are overall more spherical than the disc profile.Symmetron screening surfaces closely follow isodensity contours in the disc-dominated region, except in the large M sym , large L c limit where the individual microphysics of each parameter becomes important.
The simplest approximation we explore is the binary screening condition, ref. [50], where if a galaxy has a virial potential, Φ vir = GM vir /R vir , greater/smaller than the critical potential then it is classified as fully screened/unscreened, ignoring the effects of partial screening.We note that the symmetron binary screening condition, eq.(2.27), is not suitable as a proxy for screening.This is because the approximations required to derive the critical potential, eq.(2.6), are incompatible with the symmetron theory, namely the assumptions that the field can be linearised, the scalar mass is negligible and β(ϕ) is constant.
We explored the effectiveness of the binary screening condition as a proxy to f (R) screening, for observables in the visible disc and dark matter halo.We found that the binary screening condition drastically underestimates the amount of parameter space in which the galactic disc is fully screened.In this fully screened region of parameter space, using the binary screening condition would erroneously define the galaxy as 'fully unscreened', i.e. it would lead one to predict an observable signature that would not occur.It follows that constraints set using the binary screening condition in this contested region of parameter space would be tighter than if one numerically computed the true screening surface.
We also compare our true solutions against approximate screening surfaces.In the f (R) case, r s can be approximated by numerically solving the critical potential, eq.(2.6), for a spherically averaged galactic density.In this case, the resulting approximate screening surfaces are spherical.For the symmetron model, it is possible to approximate r s as the location where the density is equal to the symmetron breaking scale, i.e. ρ(r = r s ) = ρ SSB .
For the halo, we found the binary screening condition sits in the region of parameter space where the haloes are partially screened.This means objects classified as "screened" actually have 60-70% unscreened mass, while objects classified as "unscreened" actually have 30-40% screened mass.This is the fundamental flaw of the binary screening condition, the fact that the binary screening condition assumes a discrete step change between fully screened and unscreened, but actually the degree of screening varies gradually over many orders of magnitude in parameter space.
Direct approximation of r s allows for partial screening solutions and is overall a closer approximation to the true screening solutions than the binary screening condition.In regions where the disc component is not negligible, the spherical approximate screening surface will underestimate r s along the disc and overestimate r s perpendicular to the disc.As a result of this, we might expect to see an observable feature at a location along the plane of the disc, but in fact, the location is within the true screening surface and therefore there is no observable feature.In addition, this approximation overestimates r s for small values of |f R0 |, corresponding to smaller Compton wavelengths.This is a result of the approximation assuming a negligible scalar mass, and thus being influenced by additional mass that the field does not respond as strongly to.The effect is amplified in the large M vir , small f R0 limit, where the extent of the galaxy is largest and the Compton wavelength is smallest.
We have shown that defining the symmetron r s as the radius at which the total density equals ρ SSB is a useful approximation for a large portion of the parameter space.However, in the large M sym , large L c limit, varying the two parameters separately has consequences for the scalar field profile that are not captured purely through a dependence on ρ SSB and the approximation breaks down.In general, the approximate screening surface will never underestimate r s and the amount of overestimation depends on how large M sym and L c are.In smaller galaxies, the approximation may break down if the galactic scale length becomes comparable to L c .
Finally, we explored the possibility of finding observables in a Milky Way-like galaxy.Assuming that the Solar System is screened by the Milky Way, we found that (for our isolated test Milky-Way galaxy) |f R0 | must be less than 10 −6.0 and M sym must be less than 10 −3.6 M pl with L c being to the right of the dashed line in the bottom panel of figure 5. We saw that there is only a small slice of model parameter space in which there could be observable effects in the galactic disc.This suggests we would have to be rather fortunate to exist in a universe where we could detect the modifications to GR considered in this work within our own galactic disc -however, searches for observables in the halo remain promising, thanks in part to the large dynamical range of the halo density.
As cosmological surveys continue to improve the precision of their measurements it is important that theoretical predictions continue to improve in accuracy alongside them.It is in this spirit that we have presented, in this paper, new definitions of screening for both the f (R) chameleon model and the symmetron model of screening.We have shown that it is important to differentiate between situations where the galactic disc is screened, but much of the dark matter halo is unscreened, and situations where the majority of the total mass of the galaxy is screened.We have shown that some approximations currently used in the literature could lead to order-of-magnitude errors in the constraints placed on model parameters from galactic observations.We intend to present a reanalysis of current constraints in the near future.

A χ Parameter Derivation
This appendix provides a complete derivation of the critical potential, eq.(2.6), introduced in section 2.2.This prevalent approximation to screening comes from solving the general EoM eq.(2.4) for a static spherically-symmetric source, δρ(r), of radius R, in a constant density background, ρ.This method can be applied to scenarios such as a star embedded in the galactic medium or a dark matter halo embedded in the cosmological background.
Far from the object, the field settles into the minimum of the effective potential, i.e. ϕ → φ ≡ ϕ min (ρ).Near the object, if the field is able to reach the minimum of the potential inside the object then ∇ 2 ϕ = 0, meaning ϕ = ϕ min (ρ), and the fifth force is screened due to the vanishing field gradients.Otherwise, the EoM needs to be solved, which can be done if we assume the object is static and spherically symmetric, and therefore □ → ∇ 2 (r).If the field remains close to the background value, φ, we linearise as ϕ = φ + δϕ, and assume constant β(ϕ) throughout, resulting in the EoM

B.2 Galactic Disc Parameters
To derive the three parameters (Σ disc , R disc and z disc ) needed to construct the galactic disc profile, eq.(3.2), we must pass from the virial mass to the stellar mass (contributing to the bulk of the galactic disc).To go from halo to disc mass, we employ the empirical stellar-halo mass relation (SHMR), namely eq. 2 in ref. [68] where the four fitting parameters and related errors are in table 3 of ref. [68].Now that we have the disc mass, we pass via the half-light radius to the disc scale length, R disc .To achieve this, we use the observational relation defined by eq.18 in ref. [69].Specifically, where R hl denotes the median half-light radius and fitting parameters (α, β, γ, M 0 , σ 1 and σ 2 ) take the values from the bottom row of table 1 in ref. [69].Then, by integrating the density, eq.(3.2), to calculate the mass, and assuming a constant mass-to-light ratio so that we can utilise M disc (R = R hl ) = 1 2 M disc (R → ∞), we arrive at the relation where W −1 (x) is the n = −1 solution to the Lambert W n function.
The disc mass can be determined by integrating the density profile, eq.(3.2), at which point it is easily reversed and combined with our value of R disc to calculate Σ disc .
Finally, we employ one last observational relation to determine the disc scale height, z disc , in terms of the disc scale radius R disc .For this we utilise eq. 1 of ref. [70], R disc z disc = 0.367 log 10 R disc kpc + 0.708 ± 0.095 .(B.9) With this, we have defined all five density and length scales (ρ NFW , r NFW , Σ disc , R disc and z disc ) in terms of a single input, M vir .

Figure 1 .
Figure 1.Various panels illustrating the different stages of our numerical procedure, as described throughout section 3.All plots are zoomed in to 5R disc , showing only the disc-dominated central region.Top left: input density profile, for an axisymmetric galaxy with M vir = 10 11.5 M ⊙ , rescaled by the cosmic mean density.Top right: coordinate grid, showing every fifth radial and azimuthal angle gridline.Bottom left, left half : f R field profile rescaled by the background field value, |f R0 | = 10 −6.4 .Bottom left, right half : the curvature-density ratio, eq.(3.10a).Bottom right, left half : symmetron field profile rescaled by its background field value, for the input parameters M sym = 10 −4.5 M pl and L c = 10 −1 kpc.Bottom right, right half : truncated scalar field profile to illustrate the symmetron screening condition, eq.(3.10b).

Figure 2 .
Figure 2. Values of the f (R) curvature-density ratio (top left), the f (R) fifth force (top right, the symmetron field (bottom left) and the symmetron fifth force (bottom right), in the plane of the galactic disc, plotted against radial coordinates rescaled by r s (the radius at which the screening surface meets the disc plane), across a range of theory parameters and galaxy masses.In the f (R) case, log 10 (|f R0 |) ∈ [−8, −5] and log 10 (M vir /M ⊙ ) ∈ [10, 13.5], both in increments of 0.2.In the symmetron case, log 10 (M sym /M pl ) ∈ [−6.5, −3], log 10 (L c /kpc) ∈ [−3, 3] and log 10 (M vir /M ⊙ ) ∈ [10, 13.5], all in increments of 0.5.Profiles are cut off at the extent of the galaxy, 2.2R vir .Dashed vertical lines represent the location of the screening radius.Dashed horizontal lines correspond to the screening condition threshold.Purple lines on the symmetron plot indicate solutions that evade the fully unscreened central scalar field definition but are properly considered fully unscreened by the central Laplacian definition.

Figure 3 .
Figure 3. Screening surfaces across a range of theory parameters, in a galaxy with mass M vir = 10 11.5 M ⊙ , for f (R) (top) and symmetron (bottom) models.In the symmetron case, we vary M sym holding L c fixed (bottom left) and vary L c holding M sym fixed (bottom right).Sky blue dashed lines are the calculated screening surfaces with theory parameters labelled.Grey dashed-dotted lines are approximate solutions: numerically integrating the critical potential (2.6) in the f (R) case; isodensity contours equal to ρ SSB in the symmetron case.In each case, the right-hand semi-circle plots the radial range of the galaxy out to 10r NFW (approximately 150 kpc), while the left-hand semi-circle zooms into the central region, out to 5R disc (approximately 8 kpc).The coloured background regions represent the adopted matter density profile, plotted as a fraction of the cosmic mean density.

Figure 4 .
Figure 4. Fraction of disc mass (left) and dark matter (right) screened, for a range of f (R) (top) and symmetron parameters, with fixed L c = 10 −1 kpc (middle) and with fixed M sym = 10 −4.5 M pl (bottom).Sky blue dashed contours show where 50% of mass is screened.Grey dashed-dotted contours indicate where 50% of mass is screened when r s is approximated: by numerically integrating eq.(2.6) for the f (R) critical potential, or by finding isodensity contours equal to ρ SSB for the symmetron model.Grey dotted lines correspond to the binary screening condition; solutions below this line would be considered 'fully screened'.Solid white lines separate the fully unscreened and partially screened solutions.

Figure 5 .
Figure 5.Positions where non-spherical screening surface intersects the disc plane of a Milky Way-like galaxy over f (R) (top) and symmetron (bottom) parameter space.In the f (R) plot, the dashed line shows the true r s value, whereas the dashed-dotted line shows the approximate r s values, derived by numerically integrating eq.(2.6) for the f (R) critical potential.In the symmetron plots, the colourmap shows the true r s value, the dashed line indicates where r s is equal to the Milky Way centre-Solar System separation of 8 kpc, and the dashed-dotted line shows where the approximate r s , derived by finding isodensity contours equal to ρ SSB , is equal to the Milky Way centre-Solar System separation.Dotted lines indicate the binary screening condition.In the symmetron plot, a solid white line separates the fully unscreened and partially screened solutions, and the white region shows where SSB does not occur, i.e. the cosmic background density, ρ, is greater than ρ SSB .