© The Authors 2026

1. Introduction

Magnetic fields are known to be ubiquitous in the Universe. They span a wide range of scales and strengths: from stars and magnetars (up to ∼10¹⁵ G; Shu et al. 2025) all the way to ∼10⁻⁹ G in cosmic web filaments (Carretti et al. 2025). In the overdensities of the cosmic web (i.e., superclusters and galaxy clusters) large-scale magnetic fields have also been found. The turbulent intracluster medium (ICM) is permeated by ∼μG magnetic fields at ∼1 − 100 kpc scales (see Govoni & Feretti 2004; Donnert et al. 2018 for reviews of cluster magnetic fields). These magnetic fields play important roles in the dynamics and understanding of the baryonic content of the structures they permeate.

The baryonic content of the local Universe is well constrained, both by the predictions of Big Bang nucleosynthesis through the analysis of the abundance of primordial light elements and observations of the Lyman-α forest at high redshift (Nicastro et al. 2008) and by cosmic microwave background (Planck Collaboration VI 2020). However, only ∼70% of these baryons have been found in observations of the local Universe (Shull et al. 2012). Nonetheless, the remaining ∼30%, the so-called missing baryons, are now claimed to be accounted for, although the uncertainties are yet to be lowered for a more definitive conclusion (Driver 2021). They are expected to lie in the warm hot intergalactic medium (WHIM; Nicastro et al. 2018; Macquart et al. 2020). The WHIM is made up of diffuse ionized gas that fills the filaments of the cosmic web. These missing baryons therefore lie in low-density regions, which makes them hard to detect. Cosmic web filaments and intercluster bridges are the type of low-density structures in the Universe where these baryons are expected be. Recently, evidence has been found for magnetic fields in these filaments (Carretti et al. 2025) and it is therefore natural to think that intercluster bridges, where diffuse radio emission and thermal Sunyaev-Zeldovich (tSZ) effect signals have been detected (Pignataro et al. 2024a,b; Tanimura et al. 2019; de Graaff et al. 2019; Radiconi et al. 2022), are magnetized as well. It is then interesting to determine their level of magnetization and their magnetic field strengths and properties.

The Shapley Supercluster core (SSC; z ≃ 0.048) contains one such intercluster bridge between the massive galaxy clusters Abell 3558 and Abell 3562 (hereafter A3558 and A3562), detected by Planck through the tSZ effect (Planck Collaboration XXII 2016). The tSZ emission of two massive groups of galaxies was also detected in this bridge (SC 1327 and SC 1329). The masses, M₅₀₀, and radii, r₅₀₀¹, of these four objects can be found in Table 1. The projected distance between the centers of the Abell clusters is ≃4.2 Mpc. Such a system provides a unique laboratory for the study of diffuse magnetized gas in intercluster regions to derive the properties of their magnetic field through Faraday rotation.

Linearly polarized light that travels through an ionized and magnetized medium experiences Faraday rotation, i.e., a rotation of the linear polarization vector χ, which is proportional to the square of its wavelength: Δχ = RM ⋅ λ², where the constant of proportionality is defined as the rotation measure (RM) (e.g., Ferrière et al. 2021). Since the RM is an integral along the line of sight, the quantity we measure is the sum of all the contributions from every Faraday screen, 𝒮_i

$$\begin{array}{c}\hfill {\mathrm{RM}}_{\mathrm{obs}}={\displaystyle \sum _i}{\mathrm{RM}}_i=\frac{e^3}{2\pi m_e^2{c}^4}{\displaystyle \sum _i}{\displaystyle {\int }_{{\mathcal{S}}_i}}{n}_e(l)B_{\Vert }(l)dl.\end{array}$$(1)

The integral goes from the source to the observing point. The constants are e, the absolute value of the charge of the electron, m_e, its mass, and c, the speed of light. The RM depends on the product of the number density of thermal free electrons (n_e) and the magnetic field component along the line of sight (B_∥)². Therefore, to properly quantify how much of the observed rotation is actually due to the Shapley Supercluster, we need to properly account for all other possible contributions, typically local contributions to the source itself, intervening extragalactic structures along the way, and most importantly the Galactic contribution, denoted by GRM. When analysing RM data we consider differences in the RM dispersion of different subsamples. Thus, any contribution to the dispersion that can be set as common to both subsamples, such as the intrinsic contribution from the sources, will cancel out. Given that the GRM contribution can be different at different scales throughout the Shapley region we aim to study, it is necessary to model it in detail.

Radio polarization data has been widely used to study magnetic fields in clusters through Faraday rotation (Murgia et al. 2004; Govoni et al. 2006; Guidetti et al. 2008; Bonafede et al. 2010; Govoni et al. 2010; Vacca et al. 2012; Govoni et al. 2017; Stuardi et al. 2021). However, one of the main sources of uncertainty in these RM grid studies is the limited number of polarized sources whose lines of sight go through the given object of study (Johnson et al. 2020). This limited sampling of the density and magnetic field information is one of the main drivers for the Australian Square Kilometer Array Pathfinder-Polarisation Sky Survey of the Universe’s Magnetism (ASKAP-POSSUM) collaboration to build the densest RM grid of the southern sky to date, with about 1 million expected extragalactic polarized sources, aiming for ∼30−50 RM/deg² (Gaensler et al. 2025). The previous densest wide-area catalog has ∼1 RM/deg² (Taylor et al. 2009). Recent studies with Square Kilometer Array (SKA) pathfinder telescopes data have shown its capabilities to map and detect ionized and magnetized gas in the outskirts of galaxy clusters (Anderson et al. 2021; Loi et al. 2025) as well as the RM signature of galaxy groups (Anderson et al. 2024). Other recent studies of magnetic fields in galaxy clusters have shown the potential of statistical stacking as an alternative to wide area dense RM grids (Osinga et al. 2022, 2025).

In this work, we study the magnetic field properties of the SSC with ASKAP-POSSUM³ data from the POSSUM Pilot II survey. We focus on the potential of these dense RM grids to estimate the magnetic properties of the gas in the intercluster bridge region between clusters A3558 and A3562. The underlying assumed cosmology throughout this paper is a ΛCDM with Ω_m = 0.3, Ω_Λ = 0.7, and H₀ = 70 km/s/Mpc. The linear scale at the redshift of Shapley is thus 0.941 kpc/″. Section 2 is dedicated to the data used in this work. Sections 3 and 4 show the results derived from the data and the different modeling approaches we followed, respectively. Section 5 contains the discussion of our results. In Sect. 6 we conclude our findings. Appendices A and B provide ancillary results and outline the methods used to extract information from the observations and the modeling, respectively. Throughout this paper, the symbol σ will be used as a measure of the scatter in the RM data, and it always refers to the median absolute deviation (MAD) standard deviation defined as σ = 1.4826 MAD, which is less sensitive to outliers than the traditional standard deviation.

2. Data

2.1. Radio polarization and Faraday rotation data from ASKAP-POSSUM

In this work we use observations from the ASKAP-POSSUM pilot survey by combining data from Band 1 (800 − 1088 MHz) and Band 2 (1296 − 1440 MHz) with a spectral resolution of 1 MHz (Gaensler et al. 2025). The target of these observations was the Shapley Supercluster in two different fields: “core” and “south” (as shown in Fig. 1). The tile center coordinates for “core” and “south” in Band 1 and Band 2 are listed in Table 2. ASKAP’s field of view is 30 deg², while for the typical integration time and rms noise consult Table 2.

Fig. 1. POSSUM RMs used in this work and the definition of the different RM subsamples. The background shows the Planck tSZ effect y map. The dash-dotted squares represent the (6 deg)2 Band 1 ASKAP fields, namely, the “core” and “south” fields, with their centers represented by red crosses. The black contour represents the ybdry = 4.24 × 10−6 value used to define the on-target (orange triangles) and off-target (purple dots) samples. These two are inside the 1.71 × 2.42 deg2 (5.8 × 8.2 Mpc2) y-map cutout used for the analysis (see Sect. 2.2). The zoomed-in region at the top left corner of the plot shows the bridge (cyan triangles) and cluster (blue triangles) subsamples, all of which belong in the on-target region. The bridge box used to define the bridge sources is also represented to ease the interpretation of the plot.

2.1.1. Deriving the RM catalog

For each field, the Stokes I, Q, and U image cubes for the fields were obtained from the CASDA repository⁴, which were produced by the ASKAP Observatory pipeline using the ASKAPsoft software package (Guzman et al. 2019). The on-axis polarization calibration was done for each of the 36 beams (closepack36 footprint with beam pitch of 0.9 deg) using the unpolarized bandpass-corrected primary calibrator (PKS B1934−638). The off-axis leakage in the Stokes Q and U cubes was corrected using beam models derived from holography observations (which was estimated to result in instrumental polarization levels of order 1% or less of Stokes I for these particular observations by the POSSUM data validation team Vanderwoude et al. 2024; Gaensler et al. 2025). A Stokes I catalog for each field is also provided by the ASKAP Observatory team using the Selavy software package (Whiting & Humphreys 2012; Whiting et al. 2017). To quantify the linear polarization and Faraday rotation properties of all Stokes I sources, a development version of the POSSUM Polarimetry Pipeline was used. Initially it was ensured that the image cubes were convolved both spatially and spectrally to a common angular resolution of 20″ (using the beamcon_3D routine from racs_tools⁵). The Q and U image cubes were then corrected for the average ionospheric Faraday rotation during the observation using the package FRion⁶. The 1D IQU spectra were then extracted at the position of each Stokes I catalog entry. To account for any potential contamination of the Stokes Q and U spectra by diffuse polarized emission from the Milky Way, the median Stokes Q and U emission in a 109″ annulus surrounding each source was subtracted (see e.g. Vanderwoude et al. 2024 for details). The optimal annulus size was determined through simulated source tests in an internal POSSUM memo (Oberhelman et al. 2024). The 109″ outer radius is large enough such that the error on the subtracted median value does not affect the noise in the subtracted spectra, while also being small enough such that the procedure does not subtract diffuse emission that is unrelated to that contaminating the target spectra. RM synthesis was then applied to the individual spectra using the RM-Tools software package rmsynth_1d⁷, which output an RM catalog in the standardised form (Van Eck et al. 2023).

2.1.2. Band 1 and Band 2

While both Band 1 and Band 2 data were processed with the development version of the POSSUM pipeline (Purcell et al. 2017), the Band 1 data were processed as single fields (individually for core and south fields), while the Band 2 data were first combined (interleaved for core and in the overlap region with the south) and then split into smaller tiles (∼3.5 deg²). The IQU spectra were extracted in the same manner from both the single fields and the smaller tiles, as described above.

To create the most reliable RM catalog from the Band 1 and Band 2 data, we initially applied a cut in signal to noise of S/N > 6 (in polarization) to Band 1 core data, and S/N > 8 to Band 2 data⁸. In both Band 1 and Band 2 the angular resolution is 20″, so for multiple sources whose separation is < 20″, we only kept the one with a higher S/N, thus avoiding statistical bias from duplicated RM components. Band 1 was split into core and south data, so we merged both catalogs dealing with the overlap region between them by removing the south counterparts and keeping the core ones, due to a more reliable instrumental polarization correction in the latter. In Band 1 and Band 2 we only kept sources with Ifitstat < 5 thus filtering out those that are too faint in total intensity or those with problems raised during the model fitting of Stokes I⁹. There were sources that are lying in the overlap region between Band 1 south and Band 2 observations. To deal with these, for those in Band 1 we kept sources with a distance to the nearest beam < 0.3 deg, and distance to the tile < 1 deg. This allowed us to remove sources potentially dominated by leakage and closer to the edges of the beams. For both Band 1 south and Band 2 sources we applied a cut in the degree of polarization > 1.5% (Band 1 core sources were applied a > 1% cut, due to a better widefield leakage correction in that field). Then, we only kept unique Band 2 sources, i.e., not detected in Band 1, since Band 2 has a smaller bandwidth in λ² thus higher uncertainty in the RMs.

2.1.3. Combined final catalog

In the end, 369 sources from Band 2 and 2106 out of Band 1 data were kept¹⁰. Nonetheless, a further selection process was made to end up with a total 149 rotation measures. These where chosen to be those in the relevant region of the SSC that allows us to study its clusters and groups, comprising the combination of the on and off-target regions (see Fig. 1). Out of the final total 149 sources, only 12 of them are from Band 2 and no “south” sources made it to this final version. The catalog has metrics to quantify the Faraday complexity of the sources, and we find that only 4% of the on-target sources have hints of Faraday complexity. Furthermore, the data used in this work is unsuitable for detailed complexity analyses (e.g. QU-fitting) due to the suboptimal instrumental polarization correction. Figure A.1 shows the Stokes I image of ASKAP’s fields with the positions of the final RMs overlaid. The median error of the sources in the final RM grid is δRM = 1.98 rad/m² and the density of sources is 36 RM/deg². The most relevant quantities for the polarized sources in the final catalog can be found summarized in Table A.3 for an example set of ten of them.

2.1.4. Embedded sources identification

We searched for sources embedded in Shapley, i.e., within ≲16 Mpc of A3558 (at z = 0.048). 4 sources were found out of the 46 sources in the detected subsample (see Sect. 3.4 for the definition of it) to be embedded. Their component names are: J132928−313924, J133048−314323 (in the bridge subsample), J133503−313918, and J132803−314527 (in the clusters subsample).

2.2. Thermal Sunyaev-Zeldovich map from Planck

When CMB photons propagate through a gas of hot thermal electrons (like a galaxy cluster) they are likely to scatter off them, and typically they gain energy from these collisions (inverse Compton scattering). The way the strength of this effect is parametrized is via the y-Compton parameter (Mroczkowski et al. 2019) defined as

$$\begin{array}{c}\hfill y=\frac{{\sigma }_T}{m_e{c}^2}{\displaystyle \int P_e}(l)dl,\end{array}$$(2)

where again, the integral is along the line of sight, P_e = n_ek_BT_e is the pressure of the gas of electrons, and σ_T is the Thomson cross section.

We have used the MILCA (Modified Internal Linear Combination Algorithm) full mission ComptonSZ map from the 2015 Planck results (Planck Collaboration XXII 2016). MILCA is a component separation algorithm that was developed to reliably map spatially localized emission as opposed to the diffuse emission associated with the CMB (Hurier et al. 2013). We downloaded a cutout of the Shapley Supercluster from Planck’s Legacy Archive web page¹¹ with pixel size of 1.71′. The relevant SSC region is 1.71 × 2.42 deg² (see Fig. 1).

There are other y maps available in the literature that use the same initial Planck data and implement different algorithms for tSZ signal reconstruction, along with deprojection of the cosmic infrared background (CIB) dust contribution (McCarthy & Hill 2024). The root-mean-squared noise of this map in the same SSC region is 4% lower than for the MILCA Planck map. The median ratio of the map values from McCarthy & Hill (2024) and those from the MILCA Planck map is 1.1, thus indicating minimal impact on results derived from the y map.

3. Observational results

In this section we use the POSSUM RMs to study the statistical properties of the magnetized gas in the SSC, and how they vary as a function of distance from the cluster and group centers.

3.1. Galactic foreground subtraction and residual RMs

To study the RM properties that are only related to the SSC, we need to remove correlated RM structures on scales larger than the SSC. The dominant contribution on such scales is expected to come from the Milky Way insterstellar medium (ISM), which we refer more generally to as the Galactic RM (GRM), including the ISM, circum-galactic medium (CGM) and any magneto-ionized media associated with the Milky Way. We define the residual RM as

$$\begin{array}{c}\hfill \mathrm{RRM}\equiv {\mathrm{RM}}_{\mathrm{obs}}-\mathrm{GRM}.\end{array}$$(3)

To estimate the GRM we use the total 2475 POSSUM RMs derived in Sect. 2.1.3. This grid has an areal density of ρ_RM ≃ 41.5 deg⁻². As a consequence, we did not employ commonly used methods such as the all-sky map of the Milky Way contribution (Hutschenreuter et al. 2022), which effectively loses information on scales ≲1 deg². Instead, we have used the annulus method (Anderson et al. 2024), in which we can adapt the inner and outer annulus radii as relevant to the SSC system.

3.2. GRM from the annulus method

In this method, for each RM value, we estimate the GRM at the source’s position by defining an annulus and averaging¹² over the RM values that lie inside it. The annulus is defined by two parameters: the inner exclusion radius r_inner and the number of sources we want to average over. First we compute the distance from the source to every other source in the catalog and we keep the 50 nearest ones whose distances are > r_inner. Thus, the furthest source defines the outer radius of the annulus.

The parameters r_inner and the number of sources to average over, were chosen to avoid removing the actual SSC RM structure, while accounting for the foreground Galactic contribution. The chosen r_inner = 1.7 Mpc (0.5 deg) has a size comparable to the average r₂₀₀ of the four objects in the Shapley core (${\overline{r}}_{200}=1.58\mathrm{Mpc}\lesssim r_{\mathrm{inner}}$), thus we are not removing information inside the r₂₀₀ of the clusters and groups, while any contributions on larger scales are expected to be sub-dominant, as evidenced by cluster-derived magnetic field power spectra (Domínguez-Fernández et al. 2019). By choosing 50 sources to average inside the annulus, the average outer radius is ${\overline{r}}_{\mathrm{outer}}\approx 2.8\mathrm{Mpc}\gtrsim 1.5{\overline{r}}_{200}$. The RRM errors are computed as the quadratic sum of the measurement errors (δRM) and the standard error of the mean of the 50 sources averaged over in the annulus for that particular source i

$$\begin{array}{c}\hfill \delta {\mathrm{RRM}}_i=\sqrt{{(\delta {\mathrm{RM}}_i)}^2+{\left({\sigma }_{\mathrm{RM}}^{\mathrm{ann},i}\right)}^2/N}.\end{array}$$(4)

In Fig. 2 the RRM and GRM maps are shown. A smooth gradient over the SSC is seen in the GRM, while the RRM shows a more irregular pattern. We compared our GRM results with the all-sky reconstruction of the Milky Way contribution in Hutschenreuter et al. (2022). In contrast with the GRM gradient in Fig. 2, we found an almost constant GRM throughout the SSC compatible with the average value from the annulus method. This is consistent with the fact that this map loses information on scales lower than ∼1 degree (due to the sparsity of the underlying input RM data). The median δGRM in the Hutchenreuter et al. map is also approximately two times larger than for the annulus method. We also checked for a potential correlation between the GRM and the dust polarization map (Erceg et al. 2024) from Planck observations (Planck Collaboration XI 2020) finding no significant relation. This is supported by the fact that the thermal electrons, which are responsible for the GRM, do not necessarily reside in the dust-dominated region of the Galaxy.

Fig. 2. RRM and GRM maps made by interpolation to the nearest pixel of the estimated value at the particular position of a source from our catalog with the annulus method. The black contour corresponds to ybdry = 4.24 × 10−6 (same as in Fig. 1; see Sect. 3.3). The centers of the clusters and groups are represented by black crosses, while the dash-dotted circles represent their r500. Left: Residual rotation measure map. The random nature of the patches in the map with similar sizes between them indicates that we have removed the larger coherent RM structure of the foreground Galactic contribution, while retaining the information about the SSC. Right: GRM map. Contrary to the behavior of the patches in the RRM map, this map exhibits a large continuous gradient in the RMs, which is expected for the large-scale Galactic contribution.

A comparison of the statistical properties of the RM_obs and the RRM is shown in Table 3. The RRM distribution exhibits a mean and median closer to zero and a lower scatter than the observed RM distribution. Figure 3 shows the histograms of the RM_obs and RRM distributions. In Appendix A.2 we demonstrate the robustness of this method.

Fig. 3. Histograms of observed rotation measures (RMobs) and residual rotation measures (RRM) of the 149 sources used for our analysis.

3.3. Definition of subsamples

The 149 sources in the catalog we produced (see details in Sect. 2.1) sample the SSC, which corresponds to ${\rho }_{\mathrm{RM}}^{\mathrm{SSC}}=36$ RMs/deg² (3 RMs/Mpc²). We define the extent of the SSC region for study based on the Planck y-map, with the boundary between on-target and off-target regions defined as twice the y-map’s rms noise level y_bdry = 2y_rms = 4.24 × 10⁻⁶ (see Fig. 1 for context), where

• Off-target ≡ y < y_bdry,

• On-target ≡ y ≥ y_bdry.

The on-target region contains 46 sources. The off-target subsample contains 103 sources, which is a sufficiently large control sample that also avoids regions too far from the center. This avoids encountering contributions from other clusters and groups in the wider Shapley Supercluster area. The on-target region contains lines of sight through the Abell clusters and the SC groups. To study these contributions separately the on-target sample is further split into the “clusters” and “bridge” subsamples.

The bridge box (see Fig. 1) was defined to be parallel to the axis joining the centers of the Abell clusters. Despite the overlapping of the groups and clusters, this definition aims to isolate those lines of sight that mainly go through the SC groups (the bridge) inside their r₅₀₀, with 12 sources in this subsample. The remaining 34 sources that belong to the on-target region, are defined as the clusters subsample. These mainly represent lines of sight that cross the Abell clusters inside their r₅₀₀.

X-ray data from eROSITA could be used to trace thermal electrons and define the SSC boundaries. The main advantage of eROSITA compared to the XMM-Newton pointed observations is the larger and contiguous area coverage (Bulbul et al. 2024). While eROSITA covers the whole of Shapley, given that the deeper XMM-Newton observations cover the entirety of the SSC, which is the region of interest for this work, no significant improvement is expected from eROSITA data in this context. Furthermore, the linear dependence of the tSZ effect with the electron density makes it a more suitable and sensitive tracer to define the SSC boundary, compared to the squared dependence of the X-ray emission.

3.4. RM signature in the Shapley Supercluster core

Hereafter, we use “RM” for the results derived from the RRMs except where needed in case of a potential confusion. The higher the scatter in a distribution of RMs, the more likely it is that there is an excess contribution from magnetized gas. We calculate the excess RM scatter (MAD standard deviation, Sect. 1) of the SSC as

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}=\sqrt{{\sigma }_{\mathrm{RM},\mathrm{on}-\mathrm{target}}^2-{\sigma }_{\mathrm{RM},\mathrm{off}-\mathrm{target}}^2}=30.52\pm 4.55\mathrm{rad}/{\mathrm{m}}^2,\end{array}$$(5)

with a statistical significance of 6.7σ. We found no significant variations in the exact value of the scatter after exploring minor variations in the value of y_bdry. The scatter in the on-target region is σ_{RM, on − target} = 32.91 ± 4.18 rad/m², while the off-target sample presents a scatter of σ_{RM, off − target} = 12.32 ± 1.54 rad/m². This value is higher than the expected ∼6 − 8 rad/m² for a random extragalactic background (Schnitzeler 2010). It implies that there is a contribution from magnetized gas outside the SSC, or that there remains uncorrected GRM contributions in the background RRMs. However, ${\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}$ is not affected by this, since both subsamples are equally affected by these contributions. Removing the 19 sources with 6 < S/N < 8 leaves the statistics unchanged within the errors.

3.5. RM signature within the clusters and bridge

We investigate the RM substructure of the SSC by separately analysing the Abell clusters and the bridge region (where the SC two groups are). Figure 4 shows the y-map of the SSC, along with the |RRM| of the sources and the bridge box. It also highlights the r₅₀₀ of the clusters and the groups. The excess RM scatter detected from these regions follows Eq. (5), with the same off-target sample and replacing σ_{RM, on − target} with the RM scatter of the bridge and clusters subsamples (Table 4).

Fig. 4. tSZ Planck map of the A3558-A3562 clusters system as well as the two massive groups of galaxies SC 1327 and SC1329. The triangles represent the locations of the background ASKAP radiogalaxies and they are colored according to their |RRM| values. We also show the r500 of the four objects as reference for the size of the system at the plane of the sky, along with their centers. The black contour represents the threshold we used to define the boundary between the off-target and on-target regions: ybdry = 2yrms = 4.24 × 10−6. The rectangular box defines the bridge. The sources inside it sample the region between the Abell clusters and outside their r500, despite some overlapping effects. The counterpart on-target sources that lie outside the box mainly sample the inside the of r500 of the clusters, thus its name: the clusters subsample.

The excess RM scatter within the bridge implies that ${\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}$ is not completely dominated by the RM of the clusters. Furthermore, the fact that ${\sigma }_{\mathrm{RM}}^{\mathrm{Bridge}}$ is compatible within the uncertainties with ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$, implies that we are probing a very similar range of ${\sigma }_{\mathrm{RM}}\propto n_e{B}_{\Vert }\sqrt{\mathrm{\Lambda }L}$ (Eq. (8)) in the clusters and bridge regions. This is supported by the electron density range computed in Appendix A.5. We have also computed the excess RM scatter from each of the clusters, and found ${\sigma }_{\mathrm{RM}}^{\mathrm{A}3558}=30.16\pm 6.51$ and ${\sigma }_{\mathrm{RM}}^{\mathrm{A}3562}=27.69\pm 8.73$ at 4.6σ and 3.2σ respectively. These, further suggest that the range of $n_e{B}_{\Vert }\sqrt{\mathrm{\Lambda }L}$ values is similar in the entire SSC region and its substructures, i.e., clusters combined, bridge, and each cluster separately.

3.6. Observed RM scatter profile with distance to nearest cluster and/or group

To further study how the scatter in the RM varies with distance from the clusters and groups centers, we follow Anderson et al. (2024) and Osinga et al. (2025), and define the observable for the RM scatter profile with distance between independent sources as

$$\begin{array}{c}\hfill {\mathfrak{S}}_{\mathrm{RM}}(d_{\mathrm{nrst}})\equiv \sqrt{{\sigma }_{\mathrm{RM}}^2(d_{\mathrm{nrst}})-\frac{{\displaystyle \sum _{i=1}^N}\delta {\mathrm{RRM}}_i^2}{N-1}-{\sigma }_{\mathrm{RM},\mathrm{off}-\mathrm{target}}^2},\end{array}$$(6)

where we correct for the off-target contribution and for the RRM measurement errors (see Eq. (4)). The distance d_nrst ≡ r_nrst/r₅₀₀ is the (projected) distance from a given source to the nearest cluster/group normalized by its r₅₀₀. This definition deals with the fact that in the SSC, the complicated structure and geometry does not have spherical symmetry, as in single cluster studies. For this profile, we use the on-target region subsample (see Sect. 3.3, Figs. 1 and 4). The range of distances of our clusters and groups we are able to probe with these data is 0.3 ≲ d_nrst ≲ 1.8.

To bin our data, we used a sliding window approach, with fixed size of the window of N = 16 data points, which translates to a median window size of ∼0.3 r₅₀₀. We obtained the 68% confidence level (C.L.) error region through bootstraping with 10⁶ resamples. The resulting behavior (Fig. 5) is statistically consistent with a moderate anti-correlation between the scatter in the RM and d_nrst, as a Spearman rank test yields a correlation coefficient of −0.54 with a p-value of 2 × 10⁻³.

For d_nrst > 0.7 r₅₀₀ we find sources that are inside the r₅₀₀ of two different structures, i.e., a cluster and a group, thus we define this as the boundary at which overlapping structures begins. We quantified the trend with d_nrst of sources unaffected by the overlapping of structures (d_nrst ≤ 0.7 r₅₀₀) with a Spearman rank test with correlation coefficient of −0.65 and p-value of 0.03. Beyond > 0.7 r₅₀₀ the profile appears to flatten (Fig. 5, right of the vertical pink line). A Spearman rank test yields a correlation coefficient of 0.02, supporting this argument, but not with enough statistical significance (p-value of 0.9). To assess the impact of the bridge and galaxy groups on the profile in Fig. 5, we examined a subsample of sources, focusing on 19 sources whose lines of sight pass through the outer side of the cluster halos (counter to the bridge direction). After dividing the sample into two bins based on d_nrst, we find that the scatter still closely follows the original profile for the entire sample. Therefore, the inner structures of the bridge and groups do not significantly affect the profile.

Fig. 5. Observed 𝔖RM(dnrst) profile (black line). The gray shaded region is the 68% C.L. error region. The blue crosses represent the amplitude of the RRMs of the on-target region, while the horizontal line is σoff − target. The arrow-shaped data point corresponds to an outlier value of ∼230 rad/m2. The vertical line is set at 0.7 r500, which is the scale at which overlapping effects start to take place. The horizontal segment at the lower left corner represents the median window size of ∼0.3 r500. The vertical purple error bar (lower right corner) indicates the median on-target δRRM = 4.6 rad/m2.

3.7. Impact of small-scale magnetic fields on RM statistics

Magnetic fields coherent on scales smaller than our beam resolution, i.e., ∼20 kpc, may play an important role contributing, for instance, to the total ${\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}$. They may also cause significant depolarization, such that within the SSC region we could be missing a population of polarized sources due to our selection criteria (Sect. 2.1.3). However, there is no clear evidence of a missing population given that the RM number densities are similar in the on-target and off-target regions (${\rho }_{\mathrm{RM}}^{\mathrm{off}-\mathrm{target}}=35$ RMs/deg² and ${\rho }_{\mathrm{RM}}^{\mathrm{on}-\mathrm{target}}=37$ RMs/deg²).

To further assess the RM scatter due to small-scale fields, which we denote by Σ_RM following the convention in Osinga et al. (2025), we can use Burn’s law (Burn 1966; Sokoloff et al. 1998)

$$\begin{array}{c}\hfill p=p_0{e}^{-2{\mathrm{\Sigma }}_{\mathrm{RM}}^2{\lambda }_0^4},\end{array}$$(7)

where λ₀ ≃ 32 cm corresponds to the weighted mean of the observed wavelength-squared channels of Band 1 observations (ν ≃ 932 MHz). We set p_{off − target} ≡ p₀ = 5.3%, which is the median degree of polarization of the RM sources in the off-target sample. The median degree of polarization of the on-target sample is lower at p_{on − target} ≡ p = 4.0%, and therefore on average these sources suffer additional depolarization due to a Σ_RM, on/off = 3.6 rad/m².

This small-scale RM contribution could contribute to a decrease in 𝔖_RM at low values of d_nrst. The median degree of polarization of on-target sources inside and outside 0.7 r₅₀₀ is p_on, ≤0.7 = 3.7% and p_{on, > 0.7} = 4.4%, respectively. Taking the ratio of these two, we find a depolarization of ∼0.84 between sources closer to the clusters’ centers with respect to those outside. Using Eq. (7), we find Σ_{RM, on − in/out} = 2.9 rad/m² (small in comparison to the total signal from the SSC of ≃30 rad/m²). These results indicate that beam depolarization appears to be a minor effect, likely playing a subdominant role in the overall RM statistics of our sample. Deeper observations and/or at higher frequencies, where depolarization effects are less pronounced, would be required to advance this analysis.

4. Modeling results

4.1. Magnetic field strength estimates with single-scale models

To estimate the required magnetic field strength to explain the RM observations in the SSC region, we can use a single-scale model such as that in Gaensler et al. (2001). This model assumes cells of fixed size Λ and constant electron density ${\overline{n}}_e$, where magnetic field strength is constant with random orientation with line of sight. The observed σ_RM and magnetic field strength are then related by

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{RM}}=\frac{e^3}{2\pi m_e^2{c}^4}\frac{{\overline{n}}_{e}B}{2\sqrt{3}}\sqrt{L\mathrm{\Lambda }},\end{array}$$(8)

where L = 1.8 Mpc is the path length (for details on their estimation we refer the reader to Appendix A.3). The coherence scale Λ was fixed to 50 kpc (Guidetti et al. 2008; Bonafede et al. 2010; Vacca et al. 2012; Govoni et al. 2017). The electron density estimates range from ∼10⁻⁴ to 10⁻³ cm⁻³ (Appendix A.5, Table A.2). From Eq. (8), we determined the magnetic field strength

$$\begin{array}{c}\hfill B\simeq 1\mathrm{\mu }\mathrm{G}\left(\frac{{\sigma }_{\mathrm{RM}}[\mathrm{rad}/{\mathrm{m}}^2]}{30}\right){\left(\frac{{\overline{n}}_e[{\mathrm{cm}}^{-3}]}{4\times {10}^{-4}}\right)}^{-1}{\left(\frac{L[\mathrm{Mpc}]}{1.8}\frac{\mathrm{\Lambda }[\mathrm{kpc}]}{50}\right)}^{-1/2}\end{array}$$(9)

in the bridge, clusters, and SSC regions, ranging from 0.2 μG to 1.7 μG in the bridge and cluster regions for a wide range of electron density assumptions (see Table A.2).

4.2. Gaussian random fields with MiRò

Here we investigate more realistic ICM models in order to obtain more robust magnetic field estimates, determine the field correlation with electron density (i.e., |B| ∝ n_e^η), and assess whether an additional gas component besides the clusters and groups in the SSC is required. Therefore, we use semi-analytic models for the ICM that go beyond constant density assumptions and single-scale magnetic field fluctuations by incorporating electron density profiles for clusters and power spectra motivated by magnetohydrodynamic (MHD) simulations onto a Gaussian random magnetic field model. In particular, we modeled the electron density for the clusters and groups with universal density profiles (Pratt et al. 2022; see our Appendix B.1.1). To implement multi-scale fluctuations in the magnetic field, we used the power spectrum from Domínguez-Fernández et al. (2019) that peaks at 230 kpc. The MiRò code (Bonafede et al. 2013; Stuardi et al. 2021) is capable of incorporating the aforementioned elements, thus we use MiRò to model the SSC. Additional details on how this is done, can be found in Appendix B.1.2.

4.2.1. RM signatures predicted by MiRò

We compute the excess RM scatter from the bridge and clusters regions for the MiRò models in the following parameter space: B_rms = 1 μG−3.5 μG, η = 0, 0.5, 1, where η parametrizes the correlation between magnetic field strength and electron density |B| ∝ n_e^η. Following a Monte-Carlo and bootstrap method (see Appendix B.1.3) we obtained ${\sigma }_{\mathrm{RM}}^{\mathrm{Bridge}}$, ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ shown in Fig. 6 in comparison to the observed values. The bottom line shows the observed excess RM scatter in the bridge and clusters with its error bars, and vertically we show the estimates from the MiRò models.

Fig. 6. Visual comparison of the excess RM scatter in the bridge and clusters regions between the observations (bottom line) and the Gaussian random field models from MiRò (computed following Appendix B.1.3). Horizontally, we represent the 68% C.L. error region for each model. Vertically we outline the statistical error estimated for the observations for visual comparison reasons. The values of Brms are in μG.

With this statistic only, despite the fact that there is a significant amount of overlap between models due to the large associated uncertainties, we can already rule out some models. The criterion to decide was whether the predicted ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ is compatible, within the uncertainties (at the 1σ level), with the observed one. Given that the power spectrum is derived from cluster simulations and there is no evidence that indicates it should also hold in the bridge region, MiRò is expected to perform best at the clusters. Table B.1 shows the predicted RM scatter only of these models.

Regarding the bridge region, note that only the inclusion of the two massive groups in the modeling is enough to explain the observations within the uncertainties for most models. Therefore, the contributions to the observed excess RM scatter from the clusters and groups dominates over other magnetized gas in the SSC. Hence, at the level of our data, there is no motivation to add an additional RM component (e.g. from filamentary gas) in the intercluster region.

4.2.2. Comparison between modeled and observed 𝔖_RM(d_nrst)

To better discriminate between models, we have computed scatter profiles of the RM with distance to the nearest cluster/group to compare the predictions from MiRò with the results in Sect. 3.6. To calculate the profiles from the mock RM maps, we binned the region with a sliding window of fixed step size Δd_nrst = 10⁻², and the 68% C.L. error region was computed through bootstrap with 10⁴ resamples. For more details see Appendix B.1.4. In Fig. 7 we plot the resulting profiles as well as the observed profile from Sect. 3.6.

Fig. 7. Rotation measure scatter profiles (𝔖RM). The color-shaded regions represent the 68% C.L. error regions for the MiRò profiles (made with a sliding window of fixed size 0.01 r500). We show the five MiRò models with highest Bayesian evidence. The black line represents the observed scatter profile and the gray shaded region represents its 68% C.L. error region. The vertical line is at 0.7 r500, where overlapping structures start to dominate.

To determine which model best reproduces the observed profile we defined a Gaussian likelihood function ℒ_k for each MiRò model k as follows:

$$\begin{array}{c}\hfill {\mathcal{L}}_k=\frac{1}{\sqrt{{(2\pi )}^N{\mathrm{\Pi }}_{i=1}^N{\sigma }_{\mathrm{tot},i,k}^2}}exp(-\frac{1}{2}{\displaystyle \sum _{i=1}^N}\frac{{({\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{MiR}\stackrel{`}{\mathrm{o}},\mathrm{k}}-{\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{Obs}})}^2}{{\sigma }_{\mathrm{tot},i,k}^2}),\end{array}$$(10)

where ${\sigma }_{\mathrm{tot},i,k}^2={(\mathrm{\Delta }{\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{MiR}\stackrel{`}{\mathrm{o}},k})}^2+{(\mathrm{\Delta }{\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{Obs}})}^2$ is the total uncertainty in 𝔖_RM from both the model and the data, as computed through bootstrapping. We computed the Bayesian evidence for each model which we can simply identify with the likelihood¹³. Among the models tested, the one that provides the best description of the observational data is the one with the largest value of lnℒ_k. To compare the performance between different models we compute Bayes factors ℬ_jk ≡ ℒ_j/ℒ_k and use the Kass & Raftery (1995) scale (K&R). For more details on the Bayesian model selection analysis see Appendix B.1.4.

To do this comparison, we used the value of the RM maps with a ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ compatible with the observed excess RM scatter in the clusters (see Fig. 6) only at the location of the on-target RMs and applied Eq. (6) with the same sliding window parameters applied to the observational profile. Even though MiRò is modeling overlapping structures, it is not accounting for the likely complex interactions between them. Therefore, we only compared the models to the data in the bins satisfying d_nrst ≤ 0.7, where the overlapping structures are not dominating the behavior of the profiles. For values > 0.7 r₅₀₀ all profiles exhibit a similar flattening and most of them lie within the uncertainties of the observed profile. Thus, in order not to bias our conclusions from the Bayesian model selection with the performance of the profiles in that range of distances, we do not include them in this analysis.

In Table 5 we show the results of the Bayesian model selection analysis. The model with the largest Bayesian evidence corresponds to the model with the best predictive performance: B_rms = 2.5 μG, η = 0. We show the Bayes factors of all other models with respect to it sorted by their performance up to the first one with strong evidence against it. We find barely worth mentioning evidence against the B_rms = 3 μG, η = 0 model, and positive evidence against the best performing η = 0.5 model.

Another approach is to use the scatter of the RM as a function of electron density, 𝔖_RM(n_e). This quantity could, in principle, provide a more direct mapping of the profile of the RM scatter from the center of the cluster toward its outskirts. However, the overlapping of structures in the SSC geometry, gives 𝔖_RM(d_nrst) the advantage of clearly showing at what scale these become important. Moreover, to estimate the electron density from the observations some assumptions need to be made, such as constant density and temperature over the path length (see Appendix A.5). Nonetheless, we show the profiles with respect to electron density and the comparison between models and data with the same Bayesian model selection analysis in Appendix B.4. From these profiles the data favors the model with B_rms = 2.5 μG, η = 0, which is in agreement with the above results from the comparison with the profiles with d_nrst. The five models with highest Bayesian evidence are the same as those from the profiles with d_nrst. However in this case, the second best performing model is B_rms = 1 μG, η = 0.5, with barely worth mentioning evidence against it.

4.2.3. Estimated magnetic field strength from MiRò

The average magnetic field strength can be estimated in the bridge and clusters regions from the magnetic field 3D cubes generated by MiRò, but retains significant systematic uncertainty due to the model selection ambiguity. Applying Eq. (B.6) to the model with the best performance according to the results from Sect. 4.2.2, i.e., B_rms = 2.5 μG, η = 0, we obtained ⟨|B|⟩≃1.7 μG both in the bridge and clusters regions with a statistical error of ∼0.1 μG.

4.3. Comparison to cosmological MHD simulations: SLOW

MiRò provides a significant improvement in the ICM and magnetic field modeling in clusters with respect to single-scale models. Nevertheless, it lacks a realistic description of the complex interactions between the clusters and groups in the SSC (it only assumes noninteracting spheres with a given density profile), or variable magnetic field coherence lengths (the same power spectrum was assumed to hold at all scales and distances in the entire SSC). However, cosmological MHD simluations provide a solution to these problems.

Therefore, we extend our analysis to compare the RM data from Sect. 3 to modeling obtained from the SLOW (Simulating the Local Web)¹⁴ simulation set (Dolag et al. 2023; Böss et al. 2024). SLOW is a cosmological MHD simulation constrained to match the local Universe large-scale structure. We refer the reader to Appendix B.5 for the simulation setup details. From the simulation, we identified analogs of the SSC groups and clusters with the following properties (see discussion in Seidel et al. 2025):

• SLOW-A3562: M₅₀₀ = 5.3 × 10¹⁴ M_⊙, r₅₀₀ = 1.2 Mpc

• SLOW-A3558: M₅₀₀ = 1.46 × 10¹⁵ M_⊙, r₅₀₀ = 1.7 Mpc

• SLOW-SC1329: M₅₀₀ = 2.25 × 10¹⁴ M_⊙, r₅₀₀ = 0.91 Mpc

• SLOW-SC1327: M₅₀₀ = 2.04 × 10¹⁴ M_⊙, r₅₀₀ = 0.88 Mpc.

The radii were obtained from the M₅₀₀ following the usual self-similarity relation M₅₀₀ ∝ r₅₀₀³. The masses from the analog clusters and groups differ from the observed ones by a factor ranging from 0.02 (SC 1327) to 3.5 (SC 1329) larger. The r₅₀₀ is underpredicted by a factor 0.02 for SC 1327 and overpredicted by a factor 0.1 − 0.5 for A3558 and SC 1329, respectively, while SLOW-A3562 has a radius equal to the observed one.

Within the simulation adiabatic compression and dynamo processes amplify a seed magnetic field. This simulated magnetic field (B_sim) produces realistic values in galaxy clusters but the radial decline is too steep due to resolution effects (see Steinwandel et al. 2022), hence filament magnetic fields are under-estimated. To study the RM in filamentary regions between the clusters, several magnetic field models are applied to rescale the simulated magnetic field. Namely: a constant plasma-β ratio between magnetic and thermal pressures (B_β); turbulent and magnetic pressure equilibrium (B_ℱ); a magnetic flux conservation during amplification by 3D compression mechanism (B_ff); and two turbulent dynamo scenarios, where one of them operates for a range of densities n_e ≳ 10⁻⁴ cm⁻³ (B_dyn, ↓) with an attached decline at n_e < 10⁻⁴ cm⁻³ to match observations for filament magnetic fields (Carretti et al. 2023), and the other is a pure saturated dynamo (B_dyn, ↑); see Böss et al. (2024) for details.

To understand how these models predict different correlations between magnetic field strength and electron density (in analogy to the η parameter in the modeling with MiRò), we have fit these two quantities in the electron density range 4 × 10⁻⁵ ≤ n_e [cm⁻³] ≤ 8 × 10⁻⁴, where the comparison between the SLOW models and the data is performed (see Appendix B.5.2). We found that all amplification scenarios are well fitted by a power law with: η_sim = 0.44, η_β = 0.53, η_ℱ = 0.50, η_ff = 0.67, η_dyn, ↓ = 0.51 and η_dyn, ↑ = 0.50.

Using SLOW it is possible to obtain maps of different physical effects and observables, e.g., X-ray surface brightness, tSZ effect, temperature, electron density. As an example, Fig. 8 shows the electron number density and y-Compton maps, along with the RM map for one of the dynamo mechanisms of the Shapley analog found in the SLOW simulation. While generally well reproducing the clusters and groups’ masses and sizes, it is worth mentioning that there are two main differences with respect to the observations: the projected distance between the Abell clusters is a factor 0.8 lower, and the positions of the groups do not exactly match the observations.

Fig. 8. Zoom-in simulation of the Shapley analog found in the SLOW simulation. From left to right: electron number density, y-Compton parameter, RM map for the dynamo ↓ amplification mechanism. We show the r500 of the two clusters and the two groups, namely: 1.7, 1.2, 0.91, and 0.88 Mpc in decreasing order of mass.

Using mass weighted magnetic field maps for the different amplification mechanisms, we have computed the average magnetic field in the bridge and clusters regions following Eq. (B.6)¹⁵. All models provide magnetic field strengths of the order of ∼1 to 3 μG, both in the bridge and clusters regions, except the B_sim model that yields ∼0.1 μG. The statistical errors are ≲10⁻³ μG.

4.3.1. RM signatures predicted by SLOW

We determine ${\sigma }_{\mathrm{RM},\mathrm{SLOW}}^{\mathrm{Bridge}}$ and ${\sigma }_{\mathrm{RM},\mathrm{SLOW}}^{\mathrm{Clusters}}$ from the SLOW SSC analog in a similar manner to the observations (see Appendix B.5.2). We have estimated the excess RM scatter in these regions using the RM maps derived from all SLOW magnetic field amplification mechanisms.

Based on ${\sigma }_{\mathrm{RM},\mathrm{SLOW}}^{\mathrm{Clusters}}$ we find that the B_sim model is inconsistent with observations at 4.6σ, and therefore can be ruled out (see Fig. 9). This implies that we can rule out the possibility that in the gas between the Abell clusters no further amplification of the fields is taking place.

Fig. 9. Visual comparison of the excess RM scatter in the bridge and clusters regions between the observations (bottom line) and all magnetic field amplification mechanisms implemented in SLOW. Horizontally, we represent the 68% C.L. error region for each model. Vertically we outline the statistical error estimated for the observations for visual comparison reasons.

All other mechanisms are able to reproduce the observed scatter in the clusters. The B_β model, while showing a larger RM scatter in the bridge than in the clusters region, it also predicts a ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ that is consistent with observations at the 1σ level. Moreover, similarly to what we found from the Gaussian random field modeling with MiRò, with this statistic we find overlap between these models within the uncertainties, thus being unable to discriminate between them. Table B.4 lists the values of the excess RM scatter shown in Fig. 9.

4.3.2. Comparison between simulated and observed 𝔖_RM(d_nrst)

To better discriminate between models, using the RM maps derived from the simulation, we explore how the RM scatter behaves with distance to the center of the clusters for all amplification mechanisms considered. Details of how these profiles were constructed are provided in Appendix B.5.3. We cannot directly compare with the observational results in Sect. 3.6, because of the differences in projected distance between the clusters and the location of the groups. Therefore, we defined d_nrst this time as the distance to the nearest cluster only, without including the groups. We computed a profile of the RM scatter from the observations with this new definition of distance, using Eq. (6).

Figure 10 shows the new observational 𝔖_RM(d_nrst) profile and the equivalent SLOW profiles of all amplification mechanisms. We have performed a Spearman rank test to the observed profile to asses its flatness, obtaining a correlation coefficient of −0.36 with a p-value of 0.05. This indicates that no statistically significant anti-correlation is found, in contrast with the result from Sect. 3.6 when we included the groups into the analysis.

Fig. 10. Rotation measure scatter profiles (𝔖RM) with distance to the nearest cluster only. The black line represents the observed scatter profile, which was computed using a fixed sliding window of N = 16 points. The gray shaded region represents its 68% confidence intervals obtained with bootstrapping resampling. The observational dnrst values have been rescaled by a factor of ×0.8 to account for the difference in projected distance between the SLOW analog clusters and the real A3558-A3562. The colored regions are the SLOW 68% C.L. error regions (made with a sliding window of fixed size 0.01 r500) for the RM maps that correspond to all amplification mechanisms of the magnetic field in the intercluster region.

To determine the performance of these models we implemented a Bayesian analysis similar to that of Sect. 4.2.2. The main difference from the previous analysis is that this time we are not able to directly map the positions of the observed RMs in the SLOW SSC analog maps. Therefore, we used the value of the RM maps inside the entire analog on-target region¹⁶ to produce the profiles. However, we binned the SLOW profiles in the same way as the observational profile to make the comparison with the data by evaluating ${\mathfrak{S}}_{\mathrm{RM}}^{\mathrm{SLOW}}$ at the same d_nrst, i bins (see Appendix B.5.3).

In Table 6 we show the results from the Bayesian analysis. The turbulent scenario for amplification of the magnetic field in the intercluster region, i.e., B_ℱ is the one that provides the best description of the data. We have very strong evidence against all other mechanisms according to the K&R scale.

5. Discussion

The main objective of this work is to characterize the properties of the detected magnetized gas in the SSC (Sects. 3.4 and 3.5). To do this, we first modeled the SSC with Gaussian random magnetic fields with the MiRò code (Sect. 4.2). These models provide insights into the correlation of magnetic field strength and electron density in the ICM through the η parameter defined as |B| ∝ n_e^η. However, after estimating σ_RM in the bridge and clusters regions predicted by MiRò (Sect. 4.2.1), we found a model degeneracy, i.e., with this information alone, we are not able to distinguish between models with different values of η. Moreover, we have also estimated how the uncertainties in σ_RM behave when increasing the RM grid density, and found that even for SKA-expected grid densities ∼100 RMs/deg², this statistic suffers from large uncertainties (see Appendix B.3). These statistical difficulties are not model-dependent, i.e., the local-constrained cosmological MHD simulation SLOW (Sect. 4.3), despite providing a much more detailed prescription of ICM physics, suffers from them as well (Fig. 9). Therefore, this motivates the need to study the trend in the RM scatter as a function of distance from the cluster and group centers (𝔖_RM(d_nrst) profiles) to gain further insight into the physical properties of the magnetized gas in the SSC.

5.1. Flatter than expected 𝔖_RM profiles: The η < 0.5 puzzle

As mentioned in the previous subsection, while Gaussian random field models (MiRò) with η = 0.5, corresponding to a constant thermal-to-magnetic energy scaling, are able to match the observations in terms of ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ (Fig. 6), there are degeneracies with other models with η < 0.5. The analysis of the RM scatter as a function of distance to the nearest cluster/group, 𝔖_RM(d_nrst) (Sect. 4.2.2), revealed that the model with highest Bayesian evidence, hence, the one with best predictive performance, is B_rms = 2.5 μG, η = 0, which has an unphysical scaling of η = 0 (Table 5). Furthermore, we found positive evidence against the best-performing η ≥ 0.5 model, i.e., B_rms = 1 μG, η = 0.5. In an attempt to overcome the difficulties in the analysis introduced by the complicated geometry and overlapping structures, we analyzed profiles with respect to electron density 𝔖_RM(n_e) and followed the same Bayesian analysis for all MiRò models Appendix B.4. While these further indicate that B_rms = 2.5 μG, η = 0 is still the best-performing model, which has the unphysical η = 0 scaling, there is only barely worth mentioning evidence against the B_rms = 1 μG, η = 0.5 model (Table B.3).

The observed profile with respect to distance to the nearest cluster and/or group (Fig. 5), while showing an anti-correlation as suggested by a Spearman rank test, is still considerably flatter than the profiles derived from MiRò, which causes these η < 0.5 models to be consistent with the data. However, the observed profile is consistent with a flat dependence when considering the distance to the nearest cluster only (Sect. 4.3.2), as well as with electron density (Fig. B.1).

Considering that MiRò was designed mainly for single-cluster studies, it does not account for interactions between the clusters and groups of the SSC or a spatially variable magnetic field coherence length. However, the SLOW simulations overcome these limitations, providing detailed ICM cluster physics with a full interaction and merger history¹⁷, amplification of a seed magnetic field by adiabatic compression and dynamo processes throughout the dynamical evolution of the system, and prescriptions for further amplification mechanisms outside the cluster cores. As a consequence, the analog SSC provided by SLOW (Fig. 8) is the best possible description of the system to date.

The 𝔖_RM(d_nrst) profile comparison between all SLOW amplification mechanisms and the observed profile (Fig. 10) indicates, through the Bayesian model selection analysis, that the turbulent velocity scenario (B_ℱ) has the best predictive performance (Table 6). There is very strong evidence against all other amplification mechanisms. Nonetheless, there is a large discrepancy at scales d_nrst ≳ 0.8. It is possible that this could be solved by more directly comparable positions of the analog SC groups with respect to the observed ones (Fig. 8 vs. Fig. 4).

In the B_ℱ scenario, the magnetic field strength is rescaled as a fraction ℱ of the turbulent pressure, which was assumed to be in equilibrium with magnetic pressure (ℱ = 1). This model defines turbulent velocity as the deviation from the gas bulk velocity (similar to e.g., Simonte et al. 2022) and therefore preferentially amplifies magnetic fields in regions of shocks. This opens an avenue for future investigations with more elaborate sub-grid models for turbulence (e.g., Iapichino et al. 2017) of magnetic field amplification at shocks via the Bell instability. Recent results by Zhou et al. (2024) show magnetic field generation and amplification in low-density accretion shock regimes which could provide stronger magnetic fields than can be captured in the SLOW simulation. It will also be worth investigating the impact of astrophysical seeds for magnetic fields such as supernovae and active galactic nuclei (AGN), and how feedback effects contribute to distributing, mixing and amplifying magnetic fields (see e.g., Donnert et al. 2018, for a review). This could have a further impact on the η < 0.5 puzzle and will need to be investigated with the next generation of SLOW simulations.

Other recent works like Osinga et al. (2025), Khadir et al. (2025), have found similar results, where the data suggest models with η < 0.5. Understanding what is driving the flattening of the observed RM profile is key to solving this problem. As pointed out by Osinga et al. (2025), there seems to be an inherent difference between studies that use resolved radio galaxies, such as Murgia et al. (2004), Bonafede et al. (2010), Vacca et al. (2012), Govoni et al. (2017), and Stuardi et al. (2021), which find η ≥ 0.5, and those using unresolved radio sources. The latter have found this flatter-than-expected 𝔖_RM profile as a function of distance from the cluster center (Böhringer et al. 2016; Stasyszyn & de los Rios 2019; Osinga et al. 2022, 2025; Khadir et al. 2025), which favor the η < 0.5 models. The size of the samples is considerably larger in studies where the η < 0.5 models are favored, e.g., Osinga et al. (2025) stacks a sample of order ∼100 clusters with ∼600 RMs, and Khadir et al. (2025) uses 111 RM. While studies with resolved radio galaxies derive their conclusions from a limited number of sources, typically ≲10, their focus on individual clusters allows a more precise control of some of the possible confounding factors in stacking studies and interacting systems (e.g., Bonafede et al. 2010; Pagliotta et al. 2025).

Note that the RM scatter profile derived from observations in this work (Fig. 5), matches the behavior of the profile derived by Osinga et al. (2025) for sources behind the clusters (see their Fig. 6), i.e., they have the same amplitude of the profile at similar scales and a similar flattening at ≳0.7 r₅₀₀. This implies that the behavior found for the SSC profile is not a special case, but rather likely representative of the general population.

5.2. The possibility of a missing RM population due to depolarization

One major difference between our observations and the SLOW simulations is that the observations suffer from depolarization effects. One expects that sources that probe lower projected distances from the clusters’ centers have a large scatter in the RM, and it is at these scales precisely where we could be missing sources due to depolarization. Given our beam’s linear size of ∼20 kpc, we cannot resolve spatial fluctuations below this scale. This means that there could be sources with an intrinsic scatter in their RMs that we are not accounting for. We studied this in Sect. 3.7, to find missing additional scatter of only Σ_RM ≲ 4 rad/m², on average, due to this effect. This is clearly not enough to boost the scatter in the first bins of the profiles in Figs. 7 and 10.

In fact, we can estimate what is the Σ_RM scatter we would need to reconcile the observed 𝔖_RM(d_nrst) value in the first bin (much flatter than expected; see Fig. 10) with the value predicted by the B_ℱ magnetic field amplification scenario (the best performing model according to the Bayesian analysis). We would need 16 sources (the size of the sliding window used to compute the profile) with Σ_RM ≃ 50 rad m⁻² to yield a ${\mathfrak{S}}_{\mathrm{RM}}^{\mathrm{Obs}}(d_1)\simeq {\mathfrak{S}}_{\mathrm{RM}}^{B_{\mathcal{F}}}(d_1)$, where d₁ corresponds to the first bin in d_nrst. This implies we would be losing a scatter as high as ${\mathrm{\Sigma }}_{\mathrm{RM}}\simeq 1.65{\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}$. Moreover, we do not find any significant difference in the density of the RM grid in regions closer to the clusters centers than at larger distances. Note though that in this work, we set a threshold of p ≥ 1% to produce the catalog, while future ASKAP observations of the Shapley Supercluster may allow us to use a lower threshold of 0.1% (Gaensler et al. 2025). This could reveal more sources potentially enabling an analysis close to the cluster and group centers.

6. Conclusions

In this work we study the properties of the magnetized gas in the Shapley Supercluster core (SSC), with a special emphasis on the ≃4.2 Mpc (projected) intercluster region (or bridge) between the massive Abell clusters (A3558-A3562). We use radio polarization data from the ASKAP-POSSUM Pilot II survey in combination with thermal Sunyaev-Zeldovich effect data (from Planck).

Upon the detection of this magnetized gas, we followed different modeling approaches to characterize it, namely: semi-analytic galaxy cluster modeling with universal density profiles and Gaussian random fields following a power spectrum (MiRò), and constrained local-Universe cosmological MHD simulations (SLOW). Both of these models generate mock RM maps that we used to compare with the RMs we derived from ASKAP-POSSUM observations. The main conclusions derived from this work are the following:

1. We detected a highly significant (6.7σ) excess RM scatter from the magnetized gas present in the SSC, ${\sigma }_{\mathrm{RM}}^{\mathrm{SSC}}=30.52\pm 4.55$ rad/m². Furthermore, we determined that this excess is not dominated by the Abell clusters’ contribution alone, given that the bridge region where the massive groups SC 1327 and SC 1329 are has a comparable excess RM scatter to the clusters ${\sigma }_{\mathrm{RM}}^{\mathrm{Bridge}}=25.32\pm 8.48$ rad/m².

2. We also studied the trend of the RM scatter as a function of distance (normalized by r₅₀₀) to the nearest cluster and/or group 𝔖_RM(d_nrst), in the range sampled by our data, 0.3 ≲ d_nrst ≲ 1.8. A Spearman rank test revealed a moderate anti-correlation (correlation coefficient −0.54) with a p-value of 2 × 10⁻³. Overlapping structures in the SSC appear in our data for d_nrst > 0.7 r₅₀₀.

3. We estimated the average magnetic field strength in the bridge and clusters regions (Fig. 1). Single-scale model estimates, for a variety of assumptions, yield 0.2−1.7 μG (Table A.2); MiRò estimates are ≃1.7 μG for the best performing model (Sect. 4.2.3), while SLOW predicts 1−3 μG for all amplification mechanisms, except for B_sim, which yields ∼0.1 μG.

4. When comparing MiRò with observations, models with η < 0.5 are preferred by the data, as revealed by a comparison of RM scatter profiles both with distance (𝔖_RM(d_nrst), Sect. 4.2.2) and electron density (𝔖_RM(n_e), Appendix B.4) through the Bayesian model selection analysis.

5. The SLOW simulations provide the most detailed and complete description of the ICM physics and, despite minor differences in projected distance between the clusters and the position of the groups, an analog SSC system was identified (Fig. 8). Out of the six different amplification mechanisms for the magnetic field, the RM scatter profile analysis (Sect. 4.3.2) indicates that the best performing model is the one with amplification by turbulent velocity B_ℱ ∝ n_e^1/2v_turb in the intercluster region, with equilibrium between turbulent and magnetic energies.

6. The observational RM profiles are flatter than expected compared to all models considered in this work. We refer to this as the η < 0.5 puzzle, which cannot be easily solved by accounting for depolarization effects (considering the instrumental polarization of the ASKAP pilot survey observations of ≲1%). We find no evidence of a large, missing population of sources close to the cluster centers that would contribute to the RM scatter sufficiently (see Sect. 5.2).

Given that the RM profile of the SSC is not a special case and is likely representative of the RM behavior of the general population (Sect. 5), it provides the perfect laboratory for enhancing our understanding of cluster physics and magnetic field amplification scenarios (e.g., does magnetic reconnection play a significant role?). In the near future, the full POSSUM survey data, with an order of magnitude better threshold in the degree of polarization (Gaensler et al. 2025), will be able to detect more RMs toward the core of the clusters and groups. Deeper MeerKAT and eventually SKA-Mid observations can further improve RM grid densities, while their higher angular resolution (Braun et al. 2019) should allow for more detailed analysis on scales closer to the cluster centers. Higher resolution tSZ maps from the Simons Observatory (Ade et al. 2019) will also significantly improve the electron density mapping of complex systems like the SSC.

Acknowledgments

D. Alonso-López would like to thank the GUAIX group at the UCM for useful help and comments throughout the development of this work. D. Alonso-López also acknowledges support from the Universidad Complutense de Madrid and Banco Santander through the predoctoral grant CT25/24 and PID2022-138621NB-I00, funded by MCIN/AEI/10.13039/501100011033/FEDER, EU. D. Alonso-López and SPO acknowledge support from the Comunidad de Madrid Atracción de Talento program via grant 2022-T1/TIC-23797, and grant PID2023-146372OB-I00 funded by MICIU/AEI/10.13039/501100011033 and by ERDF, EU. This scientific work uses data obtained from Inyarrimanha Ilgari Bundara/the Murchison Radio-astronomy Observatory. We acknowledge the Wajarri Yamaji People as the Traditional Owners and native title holders of the Observatory site. The Australian SKA Pathfinder is part of the Australia Telescope National Facility (https://ror.org/05qajvd42) which is managed by CSIRO. Operation of ASKAP is funded by the Australian Government with support from the National Collaborative Research Infrastructure Strategy. ASKAP uses the resources of the Pawsey Supercomputing center. Establishment of ASKAP, the Murchison Radio-astronomy Observatory and the Pawsey Supercomputing Centre are initiatives of the Australian Government, with support from the Government of Western Australia and the Science and Industry Endowment Fund. The POSSUM project (https://possum-survey.org) has been made possible through funding from the Australian Research Council, the Natural Sciences and Engineering Research Council of Canada, the Canada Research Chairs Program, and the Canada Foundation for Innovation. AB acknowledges support from the ERC-CoG project BELOVED n. 101169773. KD acknowledges support by the COMPLEX project from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program grant agreement ERC-2019-AdG 882679. BAS acknowledges support by the grant agreements ANR-21-CE31-0019/490702358 from the French Agence Nationale de la Recherche/Deutsche Forschungsgemeinschaft, DFG through the LOCALIZATION project. BAS acknowledges support by the COMPLEX project from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program grant agreement ERC-2019-AdG 882679. LMB would like to thank Frederik Groth for help with the vortex-p code. LMB is supported by NASA through grant 80NSSC24K0173. C.S.A. acknowledges funding from the Australian Research Council in the form of FT240100498. POSSUM is partially funded by the Australian Government through an Australian Research Council Australian Laureate Fellowship (project number FL210100039 awarded to NM-G). The Dunlap Institute is funded through an endowment established by the David Dunlap family and the University of Toronto. Basic research in radio astronomy at the Naval Research Laboratory is supported by 6.1 Base funding. SM acknowledges support from grant PID2023-146372OB-I00, funded by MICIU/AEI/10.13039/501100011033 and by ERDF, EU.

References

Appendix A: Observational methods and ancillary results

A.1. Stokes I emission of the ASKAP-POSSUM radio galaxies

Figure A.1 shows the Band 1 (800-1088 MHz) Stokes I image from the ASKAP-POSSUM observations used to derive the RM grid of this work.

Fig. A.1. Stokes I image of the ASKAP fields with the location of the 149 final RMs in the catalog and the centers of the clusters and groups of the Shapley Supercluster core relevant for this work. The circles represent their r500 = 1.5, 1.2, 0.9, 0.6 Mpc. The colorbar is in units of Jy/beam.

A.2. Robustness of GRM estimation and structure function analysis

We have performed a series of tests on mock RM data for the full 2475 sources used during the GRM estimation to determine the robustness of the annulus method. The mock RMs were divided into the on-target region, and the rest of the core and south POSSUM fields (which we denote by not-on). The mock RMs of this not-on distribution, were drawn from a Gaussian distribution with constant standard deviation ${\sigma }_{\mathrm{RM}}^{\mathrm{not}-\mathrm{on}}=12$ rad/m² over the whole field The mock on-target distribution, was further split into two bins based on their d_nrst values. Those sources with d_nrst ≤ than the median distance, were drawn from a Gaussian distribution with ${\sigma }_{\mathrm{RM}}^{\mathrm{on},\mathrm{below}}=40$ rad/m² while those at d_nrst > than the median from another Gaussian with ${\sigma }_{\mathrm{RM}}^{\mathrm{on},\mathrm{above}}=20$ rad/m². The mean is common for all three input distributions ⟨RM⟩ = 20 rad/m².

To test the performance of the annulus method, we analyzed to what level the standard deviations of the input mock RMs are recovered by the final RRMs. We parametrize this with a deviation parameter that reflects the change in the standard deviation before and after applying the annulus method:

$$\begin{array}{c}\hfill ϵ\equiv 1-\frac{{\sigma }_{\mathrm{RM}}}{{\sigma }_{\mathrm{RRM}}}·\end{array}$$(A.1)

Thus, |ϵ| ≃ 0 would indicate a good performance of the annulus method. However, we need to take into account the statistical errors computed as the standard error of the standard deviation (δσ_RRM). The relative error Δσ ≡ δσ_RRM/σ_RRM is then completely determined by the number N of mock RMs of each subsample Δσ = (2(N − 1))^−1/2. These errors go from 7% to 15% for N = 103 (off-target) and N = 23 (on-target split in two halves), respectively. Therefore, to test, within the statistical errors the performance of the annulus method on recovering the input statistics we need to evaluate the quotient |ϵ|/Δσ.

We ran the simulation n = 1000 times, and obtained a median deviation parameter of |ϵ| = 0.05, thus indicating a great performance of the annulus method. The median quotient |ϵ|/Δσ goes from 1.4-3.2 for the different subsamples, which implies that the relative errors are typically a factor ≲3 lower than the deviation parameter. Therefore, the results of this robustness test are not dominated by the statistical uncertainties.

Then, we tested by what percentage the annulus GRM method can flatten the 𝔖_RM(d_nrst) profile by checking the relative slope between the sources in the two d_nrst bins before (input mock RMs), and after applying the annulus (RRMs). After running the simulation n = 1000 times, we found the annulus GRM method to flatten the RM scatter profile by ∼0.3%. This result illustrates the robustness of our GRM modeling as well as proving that it does not play any significant role in flattening the observed RM profiles of Sect. 3.6.

We also carried out a structure function analysis¹⁸ (Thomson 2024) on the residual rotation measures. It revealed a flat behavior with respect to the angular size on scales between 0.2-2 degrees. This means that for pairs of sources whose projected separation is different, no significant variation is found in the difference of their RMs. This is interpreted as a successful foreground removal. A Galactic-type contribution in the data, is associated with a non-flat structure function, since large gradients of coherent patches expected from it would make the structure function peak at angular scales across them (Haverkorn et al. 2004).

A.3. Path length estimation

One of the key parameters in the single-scale model approach is the path length L. This represents the extent of the Shapley Supercluster that Faraday rotates the angle of polarization of the light coming from our background sources. Thus, to properly determine its value, one needs to assume certain aspects of Shapley’s geometry. We have followed a simple approach, where we have assumed that both A3558 and A3562 as well as the two galaxy groups between them are at the same redshift, i.e., they lie on the same plane. We have further assumed that the bridge is axisymmetric, and the axis of symmetry is parallel to the plane of the sky. Then, we can assume that the depth is equal to the width of the projection of the bridge in the plane of the sky.

Under these assumptions, we have estimated L by studying transverse profiles through the bridge in the y-map. In particular, the axis of symmetry is parallel to the line joining the centers of A3558 and A3562, so we have examined profiles in lines perpendicular to that. Then, the value for y_bdry defines the width of the bridge (Fig. A.2). The median of all the estimated path lengths yields the value L = (1.8 ± 0.1) Mpc used throughout this paper.

This approach gives a result for the path length of lines of sight that go through the bridge. However, we have also used the value derived in this way to estimate the magnetic field strength using the single-scale model (Sect. 4.1) for the clusters subsample. Note that, using Eq. (A.3) we can see that the ratio between electron densities estimated with different path lengths, say, L and $\stackrel{\sim }{L}$, is $n_e(L)/n_e(\stackrel{\sim }{L})=\stackrel{\sim }{L}/L$. Then, from Eq. (9) it follows that the magnetic field strength ratio is

$$\begin{array}{c}\hfill \frac{B(L)}{B(\stackrel{\sim }{L})}=\sqrt{\frac{L}{\stackrel{\sim }{L}}}\equiv {\mathcal{R}}_B(L,\stackrel{\sim }{L}).\end{array}$$(A.2)

Fig. A.2. Transverse y-map profiles for the bridge. 1 pix = 1.17’. The horizontal line represents ybdry. For each profile, the difference between where the two values that cut this line represent the estimation for the path length. The resulting value of L ≃ 1.8 Mpc is the median of all transverse profiles.

Then, once we estimate what is the path length through the clusters $\stackrel{\sim }{L}$, we can determine the robustness of our assumption by evaluating Eq. (A.2). Following the same method as the one outlined before for the bridge, we have found that the (median) path lengths for the clusters are: ${\stackrel{\sim }{L}}_{\mathrm{A}3558}=2.6$ Mpc and ${\stackrel{\sim }{L}}_{\mathrm{A}3562}=2.3$ Mpc, which gives ratios ℛ_B(1.8, 2.6) = 0.86 and ℛ_B(1.8, 2.3) = 0.88, respectively. It should be clear then, that for a given polarized source in the clusters subsample, we expect a factor ≈0.9 ∼ 1 difference in the rms magnetic field strength estimate using a single-scale model, which is not a significant effect and will not affect any conclusions derived from our results.

We have explored yet another way of path length estimation via the universal density profile cube we created (see B.1.3). In this case, we have made use of the convenience of having a 3D cube and we have directly extracted the density profile along the lines of sight of our sources: n_e(x_s, y_s, z), where x_s, y_s are the plane-of-the-sky coordinates of a given source in the detected sample.

Fig. A.3. Line of sight universal density profiles through our sources in the on-target sample. The horizontal line represents ne, thresh = 1.33 × 10−4 cm−3. The mean path length for the clusters is Lclust ≃ 2.6 Mpc, while for the bridge Lbridge ≃ 2.1 Mpc, just like with the y-map transverse profiles.

To estimate the path length from these, we proceeded in the same way as before, by setting a threshold¹⁹, and estimating the depth of each line of sight as the width of the profile above the threshold. We found that the median bridge path length is L_bridge = (2.1 ± 0.3) Mpc, while for the clusters L_clust = (2.6 ± 0.2) Mpc. These values are similar to what we found with the y-map method. These results indicate that our path length estimations are consistent with each other.

A.4. X-ray counts and temperature maps from XMM-Newton

We have used a temperature map and an X-ray counts map, both derived from XMM-Newton observations (proposal IDs 060198 & 065159, PI Herve Bourdin). The counts map has a pixel size of 0.11′. The temperature map is the same as in Venturi et al. (2022), made with a multi-scale spectral mapping algorithm by means of a B2-spline wavelet analysis (Bourdin & Mazzotta 2008), but smoothed using the tophat function with radius of 20 pixels (132”) in DS9. The pixel size of this map is then 0.23′. The smoothing was done to fill the holes in the map that come from removing emission from AGNs that do not belong to Shapley. The particular scale of 132″ corresponds to ∼135 kpc, that ensures we minimize small-scale fluctuations in the map. This temperature map was used in Appendix A.5 to estimate the electron density of the SSC.

A.5. Breaking the n_eB_∥ degeneracy using tSZ and X-ray

From the RM information alone, we can only probe the product of density of free thermal electrons with the magnetic field strength parallel to the line of sight (Eq. (1)); thus, further information is needed to estimate the density in Shapley before handling the magnetic field strength. From Eq. (2), we see that we can determine n_e after assuming a dependence with line of sight distance of density and temperature. The most simple approach is to assume that density and temperature are constant throughout a given path length L across Shapley. We estimated the path length using transverse profiles of the y-map through the bridge and obtained L = 1.8 Mpc (see Appendix A.3). Then, for a given polarized source we can estimate its electron density simply as

$$\begin{array}{c}\hfill n_e=\frac{m_e{c}^2}{{\sigma }_T{k}_B}\frac{y}{T_{e}L}·\end{array}$$(A.3)

In particular, using Eq. (A.3) and the path length L, we have employed four different approaches to estimate the density, three of them based on the XMM temperature map (Sect. A.4), and one of them using the universal density cube we created (see B.1.1):

1. T = 10⁷ K: we assume a constant temperature throughout the path length L of 10⁷ K. It corresponds to the average of the temperature map. This is a lower limit of the average temperature in the ICM, thus giving an upper limit for the electron densities.

2. T = 5 × 10⁷ K: same method as above but with a temperature five times higher, which is the average of the temperature map at the location of the on-target sources.

3. XMM T-map: we use the exact value of the temperature map at the location of each of the on-target sources.

4. Universal n_e profiles: we evaluate the projected density map generated from the universal density profiles for each object (see Appendix B.1.1) at the position of the on-target sources.

All of these methods were used to estimate the magnetic field strength from a single-scale model (see Sect. 4.1, Appendix A.6 and Table A.2). To define the on-target region in MiRò we related y_bdry with an electron density value to define the 𝒩_e regions using method (c) (see Appendix B.1.3). The resulting n_e range that we can probe in all three observed subsamples as obtained with method (c) can be found in Table A.1.

A.6. Single-scale magnetic field model details

One of the most used models to estimate the magnetic field strength in clusters assumes an idealized scenario where the ICM is split into cells of equal size, with a constant electron density ${\overline{n}}_e$ throughout the path length L. In each of these cells, of size Λ, the magnetic field and the line of sight vector form a random angle (cos θ) drawn from a uniform distribution between −1 and 1 (see e.g., Gaensler et al. (2001) and references therein). For a given line of sight, provided that Λ ≪ L, i.e., for sufficiently many number of cells N, summing over the random angles yields a expected value for the rotation measure of ⟨RM⟩ = 0. Therefore, the RM grid will be Gaussian distributed as a consequence of the Central Limit Theorem. The standard deviation of this distribution is found to be Eq. (8). The factor ${(2\sqrt{3})}^{-1}$ has two origins: $1/\sqrt{3}$ comes from the Gaussianity of the magnetic field distribution, so that, on average, $B_i=|\mathit{B}|/\sqrt{3}$, while 1/2 is introduced in Gaensler et al. (2001) as a result of considering the possibility of lines of sight with randomly aligned cells.

We have fixed the coherence scale to Λ = 50 kpc (Guidetti et al. 2008; Bonafede et al. 2010; Vacca et al. 2012; Govoni et al. 2017), that satisfies the Λ ≪ L assumption. Finally, we used the electron density estimates in Sect. A.5 to use Eq. (9) to compute the magnetic field strength. In particular, using all four methods in that section, we computed the thermal electron density for each observed RM, and used the median value of each subsample as our estimate of ${\overline{n}}_e$.

In Table A.2 we present the magnetic field estimate results from all of the mentioned methods and assumptions. The lower and upper bounds are: 0.34 μG ≲ B_SSC ≲ 1.71 μG, 0.23 μG ≲ B_Bridge ≲ 1.60 μG, 0.33 μG ≲ B_Clusters ≲ 1.67 μG. Methods (a) and (b) for the n_e estimation provide the lower and upper bound for the SSC and clusters subsamples, respectively. While in the bridge, method (a) provides the lower bound and method (d) provides the upper bound.

The assumed path length L = 1.8 Mpc was estimated by studying transverse profiles through the bridge in the y-map (Appendix A.3). It is worth exploring how variations in the path length can affect our estimates of the magnetic field. For instance, we can think of a source in the clusters subsample and the possibility that the path length through the clusters may be larger than the one through the bridge. Then, the magnetic field estimate would be affected by a factor $B(L)\approx 0.9B(\stackrel{\sim }{L})$ (Appendix A.3), which is order 1, where L is the bridge path length and $\stackrel{\sim }{L}\simeq 2.45$ Mpc is the average path length of the clusters. Thus, we have used the estimated path length for the bridge in all cases.

A.7. Polarized source catalog

In Table A.3 we present a sample of three polarized sources of our final RM catalog with some of the most relevant columns.

Appendix B: Modeling methods and ancillary results

B.1. Modeling Faraday Rotation from clusters with MiRò

In galaxy clusters, turbulence is mainly driven by gas flows caused during mergers. Given that the ICM is ionized, the magnetic field gets twisted by the turbulent gas motions (small-scale dynamo), hence, its orientation has an inherent random nature (Donnert et al. 2018). It is natural then to model the magnetic field as a random field. A single line of sight will encounter then many reversals of the field’s orientation. For a sufficiently large path length L through the cluster, with enough reversals on scales much smaller ≪L, it is equally likely to encounter magnetic fields pointing toward us and away from us, averaging out to zero. Therefore, because of the Central Limit Theorem, the distribution of RMs tends to a Gaussian with zero mean, thus Gaussian random fields are the most common way to model the magnetic field distribution in galaxy clusters e.g., (Murgia et al. 2004; Bonafede et al. 2010, 2013; Stuardi et al. 2021; Osinga et al. 2022, 2025; De Rubeis et al. 2024). A widely used way to generate mock RM maps that can be compared with observations is with the MiRò code (Bonafede et al. 2013). This code implements these Gaussian magnetic fields on top of a three dimensional density cube.

B.1.1. Universal density profiles

As previously mentioned, it is important to map the electron density distribution for a robust estimation of the magnetic field. Therefore, we created a (9 Mpc)³ density cube, the size of which was chosen to be comparable to that covered by the on and off-target regions colectively, which we denote by 𝒜 (see Fig. 1). The pixel size is 10 kpc. We did a 1:1 match of the relative geometry (at the plane of the sky²⁰) of A3558, A3562 and the two massive groups SC 1329, SC 1327, and assigned each of these objects a universal density profile following Pratt et al. (2022). The main advantage of these profiles is their universality, which accounts for departures from self-similarity by introducing a mass and redshift dependence. These profiles are scaled with the r₅₀₀ of each cluster/group, and since we work at a fixed redshift of z = 0.048, the only physical parameter related to the clusters that remains free is their mass (see Sect. 1 to find their values). Defining x ≡ r/r₅₀₀, then the universal density profile follows the equation (Pratt et al. 2022)

$$\begin{array}{c}\hfill n_e(x,z,M_{500})=\frac{{\rho }_{\mathrm{gas}}(x,z,M_{500})}{{\mu }_e{m}_p}=\frac{{\rho }_{500}(z)}{{\mu }_e{m}_p}A(z,M_{500})f(x),\end{array}$$(B.1)

where μ_e = 1.148 is the mean molecular weight per free electron, m_p is the proton’s mass and ρ₅₀₀(z) = 500ρ_c(z)∝E(z)² = Ω_m(1 + z)³ + Ω_Λ. The departure from self-similarity function A(z, M₅₀₀) encodes how the self-similarity of these profiles breaks down at larger radii Pratt et al. (2022). The function f(x) follows a generalized Navarro-Frenk-White (GNFW) shape

$$\begin{array}{c}\hfill f(x)=\frac{f_0}{{(x/x_s)}^{\alpha }{(1+{(x/x_s)}^{\gamma })}^{(3\beta -\alpha )/\gamma }},\end{array}$$(B.2)

and the scaling parameters where fixed to the results obtained in (Pratt et al. 2022). Note that the cube we generated implicitly assumes noninteracting clusters and groups, i.e., $n_{e,\mathrm{tot}}={\displaystyle \sum _{i=1}^4}{n}_{e,i}$.

B.1.2. RM maps with MiRò: magnetic fields

To model the magnetic field distribution we used Gaussian random fields as implemented by the MiRò code. A brief summary of how the code works is:

• To model a Gaussian random magnetic field we work in Fourier space, where the relation between Fourier modes (k, in pixels) and spatial scales (Λ) is given by k = ℓ_cube/Λ, where ℓ_cube is the full 1D size of the cube computed as the product of the number of pixels times the resolution.

• To make sure that the final magnetic field is divergenceless, we start from the vector potential A(k), whose modulus is drawn from a Rayleigh distribution, while the phase is drawn from a uniform distribution between 0 and 2π. Then, the magnetic field is B(k) = ik × A(k). The code makes sure that the output magnetic field follows a given (input) power-spectrum.

• Afterwards, everything is transformed back to real space.

The total number of pixels of our simulation is 900 and the pixel size is 10 kpc, chosen such that the Nyquist theorem ensures that the minimum fluctuating scale will be 20 kpc, to match the resolution of our ASKAP-POSSUM observations.

The power spectrum we used was from Domínguez-Fernández et al. (2019), which is motivated by MHD simulations of clusters. The analytical expression for the magnetic energy spectrum is²¹

$$\begin{array}{c}\hfill E(k)=A{k}^{3/2}[1-\mathrm{erf}(Bln\left(\frac{k}{C}\right))],\end{array}$$(B.3)

where A is a normalization parameter, B sets the width of the spectrum and C is a characteristic wavenumber that corresponds to the inverse outer scale. We have fixed these parameters to B = 1.054, C = 4.354 Mpc⁻¹, which correspond to the case of a merging cluster at z = 0, whose power spectrum peaks at ∼230 kpc²² (see Domínguez-Fernández et al. 2019). There are two free parameters in our modeling:

• B_rms: this parameter sets the overall normalization of the magnetic field cube (fixing the value of A in Eq. (B.3)). Since the magnetic field is Gaussian, we can completely characterize it by its mean, and its standard deviation. For each spatial component, the mean is set to zero ⟨B_i⟩ = 0, that follows from the process outlined before by which the field was produced, and the standard deviation is related to the root-mean-square, which is fixed as an input in the code satisfying: ${\sigma }_{\mathit{B}}^2={\displaystyle \sum _{i=1}^3}⟨B_i^2⟩=3{\sigma }_{B_i}^2=B_{\mathrm{rms}}^2$.

• η: it is widely assumed in RM studies of clusters that the magnetic field scales somehow with density, i.e., |B| ∝ n_e^η. Some typical values for this parameter are η = 0.5, corresponding to a thermal-to-magnetic constant ratio scenario (predicted by dynamo models); the conservation of magnetic flux (flux-freezing, see e.g., Klein & Fletcher 2015) and adiabatic compression lead to η = 2/3 in the case of 3D compression, or η = 1 in the case of 1D compression (e.g., merger shocks), which leads to higher magnetic field strengths in the central regions of the cluster with a steeper radial decline. See Table 6 in Stuardi et al. (2021) for a summary of the η values found in the literature in RM studies of clusters.

With these inputs and these two free parameters, we have explored the following parameter space: B_rms = 1 μG − 3.5 μG (in steps of 0.5 μG), η = 0, 0.5, 1. This particular choice was made to test the aforementioned scenarios of magnetic field correlation with electron density and to find models that reproduce the observed scatter in the clusters region (Sect. 3.5).

The code outputs three magnetic field cubes (B_x, B_y, B_z) for each realization and a mock RM map computed by integrating the density and magnetic field cubes along the z axis.

B.1.3. Measuring the excess RM scatter with MiRò

From the observations we can only sample the magnetic field distribution from an RM grid, i.e., a statistical sample of background radio sources with individual RMs that probe a limited number of lines of sight through the extended Faraday screen of the Shapley supercluster.

We want to compare the observed ${\sigma }_{\mathrm{RM}}^{\mathrm{\Omega }}$ found in Sect. 3.5 with the predictions from MiRò. Here, Ω ⊂ 𝒜 refers to a 2D region of the maps, either the bridge or clusters region. First, we need to define these regions in the MiRò maps. This is carried out by means of two different criteria: matching the n_e range constrained by observations (Sect. A.5, Tab. A.1), and defining the same bridge box as in Fig. 4.

• To map the observed electron density range into the MiRò analysis, we first projected the density cube by integrating along the z-axis over the path length we determined observationally: L = 1.8 Mpc (Sect. A.3):

  $$\begin{array}{c}\hfill n_{e,\mathrm{proj}}(x,y)=\frac{1}{L}{\displaystyle \sum _{z=z_0-L/2}^{z=z_0+L/2}}{n}_e(x,y,z),\end{array}$$(B.4)

  where z₀ is the z-axis pixel where the clusters and groups centers are set to in the cube. Then, the projected density boundaries are defined in the following way:

  • 𝒩_e, Bridge ≡ 1.33 × 10⁻⁴ ≤ n_e, proj [cm⁻³] ≤ 5.24 × 10⁻⁴

  • 𝒩_{e, Clusters} ≡ 1.33 × 10⁻⁴ ≤ n_e, proj [cm⁻³] ≤ 9.24 × 10⁻⁴.

• Next, we define the bridge box to distinguish between pixels that sample the bridge and the clusters regions, just like in Fig. 4. Given that we matched the exact same relative geometry between the clusters and groups of the system in the density cube, we can easily define where the vertices of the bridge box are in terms of pixels. We denote the region delimited by the bridge box by ℬ.

• Combining these two we have:

  • Bridge: pixels that satisfy n_e, proj ∈ 𝒩_e, Bridge and inside the bridge box (Ω_Bridge ≡ 𝒩_e, Bridge ∩ ℬ)

  • Clusters: pixels that satisfy n_e, proj ∈ 𝒩_{e, Clusters} and outside the box (Ω_Clusters ≡ 𝒩_{e, Clusters} \ ℬ).

Then, to compute the excess RM scatter in these regions we perform a Monte-Carlo and bootstrapping method:

• First, we draw randomly generated pixels from the mock RM map. The number of these pixels is set to match the observed areal density of the observed RMs in the on and off-target regions (𝒜), i.e., ∼36 RM deg⁻².

• These pixel locations represent lines of sight through the full system. Thus, to extract ${\sigma }_{\mathrm{RM}}^{\mathrm{\Omega }}$ in the desired Ω region we select only the random pixels that satisfy the corresponding definitions: (x, y)∈Ω, and evaluate the rotation measure map at them RM(x, y) and compute its MAD standard deviation (σ_RM).

• This process is done 10⁶ times, which yields a distribution for ${\sigma }_{\mathrm{RM}}^{\mathrm{\Omega }}$. We quote the median of this distribution as the estimate of the scatter and the 68% confidence level interval sets the lower/upper error bars from the 16%/84% quantiles, respectively.

This bootstrapping was made to fully sample the SSC region we have modeled, restricted to the same areal density of RMs from our observations. The ASKAP-POSSUM RMs used in this work only represent a sample through this Shapley regions, while with our mock maps we can learn more about them through its statistical properties. The results are in Sect. 4.2.1.

B.1.4. 𝔖_RM(d_nrst) profiles and Bayesian model selection

To constrain magnetic field models obtained from MiRò we computed 𝔖_RM(d_nrst) profiles to directly compare with the observed profile from Sect. 3.6. In contrast with Eq. (6), no measurement error or off-target contribution corrections are needed in this case. Thus, computing these profiles reduces to computing the MAD standard deviation of the mock RM maps as a function of d_nrst. We remind the reader that we created a density cube with the exact same relative distances between the clusters and groups at the plane of the sky. As a consequence, we can directly compare the d_nrst values from the mock maps with those computed from the observations.

For the observational profile we used the data of the on-target subsample. This motivates the decision to use the pixels inside the 𝒩_{e, Clusters} region of the maps, which was defined to explicitly reproduce the region inside y_bdry, to compute 𝔖_RM(d_nrst). We used a sliding window with a fixed step size of Δd_nrst = 10⁻² that corresponds to a physical size of ∼10 kpc for an average r₅₀₀ of ≃1 Mpc. The 68% C.L. error region was obtained from a bootstrap resampling with 10⁴ resamples. The resulting profiles can be found in Sect. 4.2.2.

To statistically quantify which model best matches the data from the 𝔖_RM profiles we have defined a Gaussian likelihood function (ℒ_k) for each model k (see Eq. (10)). N = 11 is the number of bins used for the model selection based on these profiles. For this comparison, we evaluated the mock RM maps at the exact locations of the polarized sources, and derived 𝔖_RM(d_nrst) profiles with the same sliding window used for the observed profile (Sect. 3.6). Then, we only used bins in distance that satisfy d_nrst ≤ 0.7 to compare only in the region where the profile is not dominated by the overlapping structures. Given the fact that the error bars are asymmetric, we took into account the upper (lower) error bars whenever ${\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{MiR}\stackrel{`}{\mathrm{o}}}>{\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{Obs}}$, (${\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{MiR}\stackrel{`}{\mathrm{o}}}<{\mathfrak{S}}_{\mathrm{RM},i}^{\mathrm{Obs}}$).

To assess the performance of the models in this Bayesian framework, we computed the Bayesian evidence for each model k, given the data d which we can determine as

$$\begin{array}{c}\hfill p(d|k)={\displaystyle \int {\mathcal{L}}_k}(\mathit{\theta })p(\mathit{\theta },k)d\mathit{\theta }={\displaystyle \int {\mathcal{L}}_k}(\mathit{\theta })\delta (\mathit{\theta }-{\mathit{\theta }}_k)d\mathit{\theta }={\mathcal{L}}_k({\mathit{\theta }}_k),\end{array}$$(B.5)

where we have assumed a fixed parameter space θ_k (as it is the case for us). Therefore, the Bayesian evidence is just the likelihood for a given model (defined by Eq. (10)). Then, the model with the largest lnℒ_k is the model with the highest probability of being true given the data. To quantify the performance of the model with highest evidence (labeled by j) over the rest of the models (k), we compute the Bayes factor ℬ_jk = ℒ_j/ℒ_k with respect to it. Following the Kass & Raftery (1995) scale we say that if lnℬ_jk ≤ 1 there is “barely worth mentioning” evidence in favor of model j with respect to model k; if lnℬ_jk ≤ 3 there is “positive” evidence; if lnℬ_jk ≤ 5 there is “strong” evidence and if lnℬ_jk > 5 there is “very strong” evidence.

B.1.5. Average magnetic field strength

We estimated the average magnetic field strength over 2D regions of interest Ω from the 3D cubes generated by MiRò by

$$\begin{array}{c}\hfill {⟨|\mathit{B}|⟩}_{\mathrm{\Omega }}=\frac{\sqrt{3}}{N_z{N}_{\mathrm{\Omega }}}{\displaystyle \sum _{k=1}^{N_z}}{\displaystyle \sum _{j=1}^{N_{\mathrm{\Omega }}}}|B_i(x_j,y_j,z_k)|,\end{array}$$(B.6)

where we have made use of the Gaussianity of the distribution: B_x² + B_y² + B_z² = 3B_i², where this relation holds for any i = x, y, z. Here, N_z = 900 is the number of pixels of the cube along the line of sight, and N_Ω is the number of pixels in the Ω region.

B.2. MiRò σ_RM estimates in the clusters and bridge regions

Using the Monte-Carlo and bootstrap method in Appendix B.1.3, we estimated the RM signature in the bridge and clusters regions as predicted by MiRò. Table B.1 shows these values for the models whose ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ is compatible, within the uncertainties at the 1σ level, with the observed value (see Fig. 6).

B.3. RM grid density effect on uncertainties. Prospects with the SKA

The error bars derived with the Monte-Carlo and bootstrap sampling method (Appendix B.1.3) to compute σ_RM in the bridge and clusters regions are large. As an example of this situation, Fig. 6 illustrates how there is a degeneracy in the predicted RM scatter in the clusters region between models with very different physical meaning, e.g., B_rms = 1.5 μG, η = 0.5 and B_rms = 3.5 μG, η = 0.

Therefore, we have studied how the relative error on ${\sigma }_{\mathrm{RM}}^{\mathrm{Bridge}}$ and ${\sigma }_{\mathrm{RM}}^{\mathrm{Clusters}}$ would be affected by increasing the density of RMs on the sky (Table B.2), given that this is one of the main sources of uncertainty in magnetic field estimates we can derive from Faraday Rotation studies (Johnson et al. 2020).

We found that with the current density derived from the ASKAP-POSSUM data in this study (ρ_RM ≃ 36 RMs/deg²), the error bars are of the order of ∼25-37%, for the clusters and bridge regions, respectively. Even for a density of 100 RM/deg², about ∼3 times higher, which is close to the expected density for the SKA (Heald et al. 2020), we would still obtain a relative error of ∼15-21%. To get errors lower than 10%, we would have to increase the density of sources to an unrealistic 500 RM/deg². Note that, given the large angular size and the redshift of the SSC, these estimates serve as a lower limit in the uncertainties for the expectations of similar studies of objects of smaller angular size and at higher redshifts. This result indicates that the uncertainties we can expect for single object studies with RM grids, will remain of the order of 10% or greater in the next decades as far their density is concerned.

These limitations serve as motivation to do further analysis with the RM scatter profiles with distance and using other statistics to distinguish between models, such as the Bayesian model selection analysis performed in Sects. 4.2.2 and 4.3.2.

B.4. 𝔖_RM(n_e) profiles from MiRò

We computed RM scatter profiles as a function of electron density. For the observational profile, the electron density was estimated using Eq. (A.3), with L = 1.8 Mpc (Appendix A.3), temperature T_e = 5 × 10⁷ K (Appendix A.5), and the values of the y-map at the location of the on-target sources. The sliding window and the bootstrapping used was the same as in Sect. 3.6. For the profiles derived from the mock RM MiRò maps, we used their values evaluated at the pixels inside the 𝒩_{e, Clusters} region. The electron density values were obtained from this same region of the projected density map n_e, proj (Eq. (B.4)). We used a sliding window of fixed size Δn_e = 10⁻⁵ cm⁻³ and bootstrap with 10⁴ resamples for the errors.

We also carried out a Bayesian model selection analysis to determine the predictive performance of the models, and compare between them following the same procedure as in B.1.4. Analogously, we used Eq. (10) but the profiles are now binned in electron density values instead. We compared the data with the models only using the values of the mock RM and n_e, proj maps at the position of the on-target sources and the same sliding window and bootstrap parameters used for the observational profile. We restricted the electron density bins to those that satisfy n_e ≥ 2.5 × 10⁻⁴ cm⁻³ to avoid overlapping structures dominating the comparison.

Table B.3 shows the Bayesian model selection analysis results. The model with highest Bayesian evidence which is an indicator for best agreement between the model and the data has parameters B_rms = 2.5 μG, η = 0, which is the same one found from the analysis carried out in Sect. 4.2.2. The model with B_rms = 1 μG, η = 0.5 is the one with the second best performance, and there is only barely worth mentioning evidence against it.

Fig. B.1. 𝔖RM(ne) profiles. The color-shaded regions represent the 68% C.L. error regions for the five MiRò profiles with the best performance according to the Bayesian analysis. The black line represents the observed scatter profile and the gray shaded region represents its 68% C.L. error region. The vertical dash-dotted line marks ne ≥ 2.5 × 10−4 cm−3, which is the region where overlapping structures do not dominate.

B.5. SLOW simulations

SLOW (Dolag et al. 2023) uses initial conditions constructed from the CosmicFlows-2 dataset (Tully et al. 2013), to closely resemble the observed local Universe at z = 0, using the realization number 8 of CLONES (see Sorce 2018, and references therein for details on the construction process). The SLOW simulation produces a close match to observed local galaxy clusters (Hernández-Martínez et al. 2024) as well as local super clusters, like the Shapley Supercluster (see discussion in Seidel et al. 2025). For this work we performed a high-resolution zoom-in re-simulation of the Shapley Supercluster region based on initial conditions constructed from the collapse volumes of these structures (Seidel et al., in prep). We used the same physical setup as the COMA zoom simulation presented in Böss et al. (2024), performed with the cosmological TREE-SPH code OPENGADGET3 code (Groth et al. 2023). This simulation is a non-radiative MHD simulation, using the MHD solver described in Dolag & Stasyszyn (2009), the nonideal MHD extension by Bonafede et al. (2011), constrained hyperbolic divergence cleaning (Tricco et al. 2016, Steinwandel & Price; in prep) as well as higher order smoothing kernels and improvements to the hydro solver summarized in Beck et al. (2016).

Fig. B.2. Rotation Measure maps of all SLOW magnetic field amplification mechanisms used in this work.

B.5.1. Magnetic field amplification mechanisms

In this simulation, a seed magnetic field of strength B = 10⁻¹⁴ G at redshift z = 120 was set up. This field evolves and is amplified by adiabatic compression and dynamo processes. On the one hand, these mechanisms are enough to provide magnetic field strengths of ∼μG in the cluster centers, on the other hand, the magnetic field strengths in filaments are negligible, at the current resolution. To further amplify the magnetic field in the filaments, several mechanisms are implemented: two of them scaling with pressure, and three with density:

• Plasma-β: magnetic field pressure and thermal pressure are related by a constant factor $B_{\beta }=\sqrt{8\pi P_{\mathit{th}}/\beta }$, with β = 20.

• Turbulent velocity: the magnetic field strength is a fraction ℱ of the turbulent pressure $B_{\mathcal{F}}=\mathcal{F}\sqrt{4\pi \rho v_{\mathrm{turb}}^2}$²³. The fraction was set to ℱ = 1 assuming pressure equilibrium between turbulent and magnetic pressures²⁴.

• Flux-freezing: from the conservation of the magnetic flux, amplification by 3D compression implies B_ff ∝ n_e^2/3.

• Turbulent dynamo ↓: this model assumes that turbulent dynamo amplifies the field down to densities of n_e ∼ 10⁻⁴ cm⁻³. For lower densities than that, a fitting was made to match Carretti et al. (2023) observations of ∼30 nG fields in filaments down to densities of n_e ∼ 10⁻⁵ cm⁻³. The scaling behavior is B_dyn, ↓ ∝ n_e^1/2.

• Turbulent dynamo ↑: pure saturated dynamo with the same B_dyn, ↑ ∝ n_e^1/2 scaling.

For more details see Böss et al. (2024). Please note that the models B_ff, B_dyn, ↓ and B_dyn, ↓ are up-scaled by a factor of two in respect to the previous work for a better match with observations. To model the RM signal we rescale the parallel component of the simulated magnetic field according to the different magnetic field models.

To compare the results from the SLOW analog with the observations, we applied very similar techniques to those used to compare with the outputs from MiRò (Appendix B.1). The RM maps for the six different magnetic field amplification mechanisms implemented by SLOW can be see in Fig. B.2.

B.5.2. Excess RM scatter in the Shapley analog

To estimate the excess RM scatter in 2D regions (${\sigma }_{\mathrm{RM}}^{\mathrm{\Omega }}$), we based the analysis on B.1.3. This is, we defined the region Ω in terms of an appropriate electron density range and in terms of a bridge box. Regarding density, the boundaries are

• 𝒩_e, Bridge ≡ 4 × 10⁻⁵ ≤ n_e [cm⁻³] ≤ 4.5 × 10⁻⁴

• 𝒩_{e, Clusters} ≡ 4 × 10⁻⁵ ≤ n_e [cm⁻³] ≤ 8 × 10⁻⁴.

This choices were made trying to probe a similar range of densities as in the observations, avoiding lines of sight that go directly through the clusters and groups centers. The region Ω is also defined in terms of a bridge box, whose orientation is parallel to the line that joins the centers of the massive clusters. However, the geometry of the system is not directly comparable to that in the real Shapley, since the positions of the analog SC groups are not exactly the same. In fact, the projected distance between the clusters is a factor 0.8 lower. Despite this fact, the definition of the bridge box ensures that the groups lie inside it, so that the scatter of the RM in the bridge region accounts for their contribution.

In summary, the bridge and clusters regions defined in the SLOW maps are: Bridge: pixels that satisfy n_e ∈ 𝒩_e, Bridge and inside the bridge box; Clusters: pixels that satisfy n_e ∈ 𝒩_{e, Clusters} and outside the box. Once these regions are defined, to obtain ${\sigma }_{\mathrm{RM}}^{\mathrm{\Omega }}$ the same Monte-Carlo and bootstrap methods in Appendix B.1.3 were used.

B.5.3. 𝔖_RM(d_nrst) profiles with SLOW

We also computed profiles of the RM scatter with respect to distance. Given that the SLOW SSC analog does not have the groups in the same locations as in the real Shapley, the distance d_nrst was defined as the distance to the nearest cluster only, without considering the distances to the groups. Therefore, we also made a profile from the observations considering the distance to the clusters only. To compare the SLOW profiles with the observational one, we have re-scaled the latter by a factor ×0.8 that corresponds to the ratio between the projected distances between the A3558 and A3562 clusters of the SLOW analog of Shapley, and the real one (∼3.5 Mpc/4.2 Mpc)

To produce these profiles, we have proceeded in a similar way to Appendix B.1.4. Therefore, we have employed a sliding window approach with a fixed step size of Δd_nrst = 0.01 and bootstrap resampling with 10⁴ resamples. We derived profiles for the RMs derived from all amplification mechanisms in the region 𝒩_{e, Clusters}. Nevertheless, to define the likelihood function for the Bayesian model selection analysis, we have binned these profiles according to the binning of the observational profile. This was done by finding the values of ${\mathfrak{S}}_{\mathrm{RM}}^{\mathrm{SLOW}}$ at the binned values of d_nrst, i from ${\mathfrak{S}}_{\mathrm{RM}}^{\mathrm{Obs}}$. Then, the likelihood function and the model selection analysis is completely analogous to Appendix B.1.4. The results can be found in Sect. 4.3.2.

B.6. SLOW σ_RM estimates in the clusters and bridge regions

Table B.4 shows the excess RM scatters in the bridge and clusters regions for all SLOW amplification mechanisms for the magnetic field computed using the method in Appendix B.5.2.

All Tables

All Figures

Fig. 1. POSSUM RMs used in this work and the definition of the different RM subsamples. The background shows the Planck tSZ effect y map. The dash-dotted squares represent the (6 deg)2 Band 1 ASKAP fields, namely, the “core” and “south” fields, with their centers represented by red crosses. The black contour represents the ybdry = 4.24 × 10−6 value used to define the on-target (orange triangles) and off-target (purple dots) samples. These two are inside the 1.71 × 2.42 deg2 (5.8 × 8.2 Mpc2) y-map cutout used for the analysis (see Sect. 2.2). The zoomed-in region at the top left corner of the plot shows the bridge (cyan triangles) and cluster (blue triangles) subsamples, all of which belong in the on-target region. The bridge box used to define the bridge sources is also represented to ease the interpretation of the plot.

In the text

Fig. 2. RRM and GRM maps made by interpolation to the nearest pixel of the estimated value at the particular position of a source from our catalog with the annulus method. The black contour corresponds to ybdry = 4.24 × 10−6 (same as in Fig. 1; see Sect. 3.3). The centers of the clusters and groups are represented by black crosses, while the dash-dotted circles represent their r500. Left: Residual rotation measure map. The random nature of the patches in the map with similar sizes between them indicates that we have removed the larger coherent RM structure of the foreground Galactic contribution, while retaining the information about the SSC. Right: GRM map. Contrary to the behavior of the patches in the RRM map, this map exhibits a large continuous gradient in the RMs, which is expected for the large-scale Galactic contribution.

In the text

	Fig. 3. Histograms of observed rotation measures (RMobs) and residual rotation measures (RRM) of the 149 sources used for our analysis.
In the text

Fig. 4. tSZ Planck map of the A3558-A3562 clusters system as well as the two massive groups of galaxies SC 1327 and SC1329. The triangles represent the locations of the background ASKAP radiogalaxies and they are colored according to their |RRM| values. We also show the r500 of the four objects as reference for the size of the system at the plane of the sky, along with their centers. The black contour represents the threshold we used to define the boundary between the off-target and on-target regions: ybdry = 2yrms = 4.24 × 10−6. The rectangular box defines the bridge. The sources inside it sample the region between the Abell clusters and outside their r500, despite some overlapping effects. The counterpart on-target sources that lie outside the box mainly sample the inside the of r500 of the clusters, thus its name: the clusters subsample.

In the text

Fig. 5. Observed 𝔖RM(dnrst) profile (black line). The gray shaded region is the 68% C.L. error region. The blue crosses represent the amplitude of the RRMs of the on-target region, while the horizontal line is σoff − target. The arrow-shaped data point corresponds to an outlier value of ∼230 rad/m2. The vertical line is set at 0.7 r500, which is the scale at which overlapping effects start to take place. The horizontal segment at the lower left corner represents the median window size of ∼0.3 r500. The vertical purple error bar (lower right corner) indicates the median on-target δRRM = 4.6 rad/m2.

In the text

Fig. 6. Visual comparison of the excess RM scatter in the bridge and clusters regions between the observations (bottom line) and the Gaussian random field models from MiRò (computed following Appendix B.1.3). Horizontally, we represent the 68% C.L. error region for each model. Vertically we outline the statistical error estimated for the observations for visual comparison reasons. The values of Brms are in μG.

In the text

Fig. 7. Rotation measure scatter profiles (𝔖RM). The color-shaded regions represent the 68% C.L. error regions for the MiRò profiles (made with a sliding window of fixed size 0.01 r500). We show the five MiRò models with highest Bayesian evidence. The black line represents the observed scatter profile and the gray shaded region represents its 68% C.L. error region. The vertical line is at 0.7 r500, where overlapping structures start to dominate.

In the text

	Fig. 8. Zoom-in simulation of the Shapley analog found in the SLOW simulation. From left to right: electron number density, y-Compton parameter, RM map for the dynamo ↓ amplification mechanism. We show the r500 of the two clusters and the two groups, namely: 1.7, 1.2, 0.91, and 0.88 Mpc in decreasing order of mass.
In the text

	Fig. 9. Visual comparison of the excess RM scatter in the bridge and clusters regions between the observations (bottom line) and all magnetic field amplification mechanisms implemented in SLOW. Horizontally, we represent the 68% C.L. error region for each model. Vertically we outline the statistical error estimated for the observations for visual comparison reasons.
In the text

Fig. 10. Rotation measure scatter profiles (𝔖RM) with distance to the nearest cluster only. The black line represents the observed scatter profile, which was computed using a fixed sliding window of N = 16 points. The gray shaded region represents its 68% confidence intervals obtained with bootstrapping resampling. The observational dnrst values have been rescaled by a factor of ×0.8 to account for the difference in projected distance between the SLOW analog clusters and the real A3558-A3562. The colored regions are the SLOW 68% C.L. error regions (made with a sliding window of fixed size 0.01 r500) for the RM maps that correspond to all amplification mechanisms of the magnetic field in the intercluster region.

In the text

	Fig. A.1. Stokes I image of the ASKAP fields with the location of the 149 final RMs in the catalog and the centers of the clusters and groups of the Shapley Supercluster core relevant for this work. The circles represent their r500 = 1.5, 1.2, 0.9, 0.6 Mpc. The colorbar is in units of Jy/beam.
In the text

	Fig. A.2. Transverse y-map profiles for the bridge. 1 pix = 1.17’. The horizontal line represents ybdry. For each profile, the difference between where the two values that cut this line represent the estimation for the path length. The resulting value of L ≃ 1.8 Mpc is the median of all transverse profiles.
In the text

	Fig. A.3. Line of sight universal density profiles through our sources in the on-target sample. The horizontal line represents ne, thresh = 1.33 × 10−4 cm−3. The mean path length for the clusters is Lclust ≃ 2.6 Mpc, while for the bridge Lbridge ≃ 2.1 Mpc, just like with the y-map transverse profiles.
In the text

Fig. B.1. 𝔖RM(ne) profiles. The color-shaded regions represent the 68% C.L. error regions for the five MiRò profiles with the best performance according to the Bayesian analysis. The black line represents the observed scatter profile and the gray shaded region represents its 68% C.L. error region. The vertical dash-dotted line marks ne ≥ 2.5 × 10−4 cm−3, which is the region where overlapping structures do not dominate.

In the text

	Fig. B.2. Rotation Measure maps of all SLOW magnetic field amplification mechanisms used in this work.
In the text