Dark Energy Survey Year 6 Results: Intra-Cluster Light from Redshift 0.2 to 0.5

Using the full six years of imaging data from the Dark Energy Survey, we study the surface brightness profiles of galaxy cluster central galaxies and intra-cluster light. We apply a ``stacking'' method to over four thousand galaxy clusters identified by the redMaPPer cluster finding algorithm in the redshift range of 0.2 to 0.5. This yields high signal-to-noise radial profile measurements of the central galaxy and intra-cluster light out to 1 Mpc from the cluster center. Using redMaPPer richness as a cluster mass indicator, we find that the intra-cluster light brightness has a strong mass dependence throughout the 0.2 to 0.5 redshift range, and the dependence grows stronger at a larger radius. In terms of redshift evolution, we find some evidence that the central galaxy, as well as the diffuse light within the transition region between the cluster central galaxy and intra-cluster light within 80 kpc from the center, may be growing over time. At larger radii, more than 80 kpc away from the cluster center, we do not find evidence of additional redshift evolution beyond the cluster mass dependence, which is consistent with the findings from the IllustrisTNG hydrodynamic simulation. We speculate that the major driver of intra-cluster light growth, especially at large radii, is associated with cluster mass growth. Finally, we find that the color of the cluster central galaxy and intra-cluster light displays a radial gradient that becomes bluer at a larger radius, which is consistent with a stellar stripping and disruption origin of intra-cluster light as suggested by simulation studies.


INTRODUCTION
Galaxy clusters contain a diffuse stellar component of intra-cluster light (ICL).First discovered more than half a century ago (Zwicky 1951(Zwicky , 1952)), ICL is abundant around the cluster central galaxies (CGs) or the brightest cluster galaxies (BCGs), and contains stars dispersed into the intra-cluster space.It has been studied using optical or infrared imaging, and spectroscopic observations, which have been reviewed in Abraham et al. (2017); DeMaio (2017); Contini (2021); Montes (2022); Arnaboldi & Gerhard (2022).Because of the ICL's faint brightness and the overall difficulties of studying low surface ★ E-mail: yuanyuan.zhang@noirlab.edubrightness features (Abraham et al. 2017;Mihos 2019), ICL has remained a poorly understood subject until recently, when the number of studies jumped with refreshed interests due to new data (Montes & Trujillo 2022), simulations (Shin et al. 2022), and techniques (Marini et al. 2022).
Simulation and semi-analytical studies have investigated ICL formation in many different channels (e.g., Sommer-Larsen 2006; Rudick et al. 2006;Barai et al. 2009;Rudick et al. 2009;Puchwein et al. 2010;Contini et al. 2014;Martel et al. 2012) including galaxy disruption, stellar stripping, merging (Murante et al. 2004(Murante et al. , 2007;;Contini et al. 2018) and preprocessing (Chun et al. 2022).The ICL's formation is often studied together with the evolution of cluster CGs or even the cluster's overall galaxy distribution due to difficulties in separating the two (e.g., Conroy et al. 2007; Monaco et al. 2006;Cooper et al. 2015;Pillepich et al. 2018;Cañas et al. 2020).Different channels of ICL formation carry implications for the ICL's observational properties and their redshift evolution in terms of age, color, and metallicity (Harris et al. 2017;Contini et al. 2019), fraction of ICL in the cluster stellar light (Purcell et al. 2007;Murante et al. 2007;Cui et al. 2014;Tang et al. 2018), morphology (Rudick et al. 2006), or scaling relation to cluster mass or mass distribution (Alonso Asensio et al. 2020;Contini & Gu 2021).For example, Contini et al. (2019) analyzed ICL color and metallicity using semianalytical models which contain ICL formed through tidal stripping of cluster satellite galaxies as well as through merging relaxation; they found a negative radial color and metallicity gradient.From hydrodynamic simulations, Pillepich et al. (2018) found that ICL stellar mass strongly correlates with the host halo mass, but this correlation appears to evolve little from redshift 1 to 0.
In observational studies, the formation and evolution of ICL has been studied using its color (Mackie 1992;Krick et al. 2006;Krick & Bernstein 2007;Montes & Trujillo 2014;DeMaio et al. 2015;Iodice et al. 2017;Morishita et al. 2017;DeMaio et al. 2018;Montes & Trujillo 2018;Ko & Jee 2018;Montes et al. 2021;Ragusa et al. 2021;Yoo et al. 2021;Golden-Marx et al. 2023a;Martínez-Lombilla et al. 2023), stellar mass (Krick et al. 2011;Burke et al. 2012;DeMaio et al. 2020;Spavone et al. 2020;Barfety et al. 2022), stellar population spectroscopic studies, (e.g., Coccato et al. 2010bCoccato et al. ,a, 2011;;Ventimiglia et al. 2011;Arnaboldi et al. 2012;Longobardi et al. 2015a;Edwards et al. 2016;Barbosa et al. 2016) and is often investigated together with BCG evolution (e.g., Gonzalez et al. 2005;Zhang et al. 2016;Golden-Marx et al. 2022).For example, DeMaio et al. (2020) studied the stellar mass of BCG and ICL between redshift 0.05 to 1.75 and found its growth rate to be greater than that of the cluster by a factor of two.They also found that the core of the BCG formed early while the BCG outskirt and ICL were built at later times.On the other hand, detailed analysis of local BCG and ICL stellar populations by Edwards et al. (2020) indicate that while the stellar population in the ICL is old, it is still younger (≈ 9 Gyr) than the BCG (≈ 13 Gyr), pointing towards a late and continuous formation of ICL through minor merging.
In this work, we continue the observational study of ICL evolution by examining its properties in the redshift range of 0.2 to 0.5.Our work is based on thousands of galaxy clusters and the full six years of observations from the Dark Energy Survey (DES) (Abbott et al. 2021), a wide-field imaging survey (DES Collaboration 2005) designed to probe cosmic structures in the late Universe (e.g., Abbott et al. 2020Abbott et al. , 2022b,a;,a;DES Collaboration et al. 2022).We use a "stacking" method (Zibetti et al. 2004(Zibetti et al. , 2005;;Tal & van Dokkum 2011;Zhang et al. 2019c;Sampaio-Santos et al. 2021;Chen et al. 2022;Ahad et al. 2023) with the DES galaxy cluster sample to reduce measurement noise.Our goal is to acquire high signal-to-noise measurements of the ICL surface brightness profile, color, and luminosity and quantify their evolution between redshift 0.2 and 0.5.This paper presents one of the largest ICL redshift evolution studies, based on a cluster sample a few times larger than that in Golden-Marx et al. (2023a) which used a Cosmic Microwave Background (CMB) selected cluster sample from the Atacama Cosmology Telescope.
One challenge to ICL studies is the difficulty in disentangling ICL from cluster CGs (see discussions in Dolag et al. 2010;Contini et al. 2022).Although stars in the intra-cluster space may have different stellar composition or dispersion dynamics (e.g., Longobardi et al. 2015bLongobardi et al. , 2018a,b;,b;Hilker et al. 2018;Gu et al. 2020;Pérez-Hernández et al. 2022) than cluster CGs or BCGs -from imaging data alone, it is often difficult to separate the ICL from the low-surface brightness outskirts of those galaxies.Different separation methods have been suggested (Rudick et al. 2011), including analytical decomposition, using a machine learning algorithm (Marini et al. 2022), surface brightness limits (Presotto et al. 2014), or using physical distance apertures to separate those components.In this paper, we follow the practice of Pillepich et al. (2018) who analyzed CG and ICL as the "diffuse light" of galaxy clusters.We use the phrase diffuse light interchangeably with CG+ICL in this paper.When needed, we use a physical aperture to separate CG and ICL, with ICL defined as the diffuse light beyond 30 kpc from the CG center (an outer radius limit is defined according to the context), while CG is defined as the diffuse light component within 30 kpc.
The remainder of this paper is organized as the following.In section 2 we review our data sets and the methods.Section 3 presents our measurements of the diffuse light surface brightness, while Section 4 quantifies the cluster mass and redshift dependence of the diffuse light luminosities.Section 5 discusses observational effects that may impact the interpretation of our results, and Section 6 discusses our results in the context of simulations and other observational studies.Section 7 summarizes our findings.Throughout this paper, we assume a flat ΛCDM cosmology with Ω  = 0.3, and ℎ = 0.70.

The redMaPPer Galaxy Cluster Catalog
The red sequence Matched-filter Probabilistic Percolation cluster finder algorithm (redMaPPer ) (Rykoff et al. 2014) has been used by the DES Collaboration to derive galaxy cluster catalogs from the Science Verification data, (Rykoff et al. 2016), the Year 1 observations (McClintock et al. 2019), and the Year 1 to Year 3 observations (O'Donnell et al. 2021).redMaPPer is a red-sequence based algorithm that provides excellent cluster richness () and photometric redshift estimates (Rykoff et al. 2014).It also provides a random point catalog that tracks the sky footprint and depth covered by the cluster-finding algorithm.
This paper is based on the redMaPPer cluster catalog, version 6.4.22+2, derived from the DES Year 3 Gold data sets (Sevilla-Noarbe et al. 2021).A relevant difference between this catalog and the DES Year 1 version (McClintock et al. 2019) is the much larger sky coverage.As a result, this catalog contains more than 21,000 galaxy clusters with richness above 20, which approximately corresponds to a halo mass threshold of 10 14.1 M ⊙ (McClintock et al. 2019;Farahi et al. 2019).These galaxy clusters are detected from the DES single-object fit (SOF) catalog (Drlica-Wagner et al. 2018;Sevilla-Noarbe et al. 2021) which contains objects detected and deblended by Source Extractor (Bertin & Arnouts 1996), while the photometry was derived from single-object fitting using the ngmix algorithm on multi-epoch image stamps of each object, and the deblended nearby objects are masked on each single-epoch image.For the DES Year 3 data processing campaign, the SOF photometry measurements are preferred in many applications because of the tighter photometry constraints compared to the Source Extractor measurements derived using the coadded images.
Of particular importance to this analysis, redMaPPer provides CG (Central Galaxy) candidates for each cluster.Unlike algorithms that aim to select the BCG (Brightest Cluster Galaxies), redMaP-Per aims to select a relatively luminous cluster galaxy that is nearest to the cluster's gravitational center, and the goal of this selection is to find the central galaxy of the massive dark matter halo as defined in simulation modeling studies (e.g.De Lucia & Blaizot 2007;Yang et al. 2008).redMaPPer provides five CG candidates for each cluster.We use the most likely CG candidate; multi-wavelength studies have shown that this candidate is the correct one with an ∼ 80% frequency (Saro et al. 2015;Zhang et al. 2019a;Bleem et al. 2020).
For this diffuse light analysis, we apply a few additional selection criteria to the redMaPPer clusters as well as to the redMaPPer random points.(1) Around the cluster center (or a random point), in a circular region with a radius of 0.15 deg, we require at least one DES exposure image in each of the , , , and  filters.(2) Around the cluster center (or a random point), in a circular region with a radius of 0.15 deg, we require the 10  depth magnitude of the DES Auto measurements to be deeper than a redshift dependent "masking" magnitude (see next section for details).This selection criterion ensures that our diffuse light measurements are comparable between different redshift slices.The cut has minimal effect on the clusters/randoms below redshift 0.4 but excludes a small fraction of the clusters between redshift 0.4 and 0.5 and most of the clusters above redshift 0.5.(3) Around the cluster center (or the random point), in a circular region with a radius of 0.2 deg, we exclude areas containing famous or bright stars (the Yale bright stars or 2MASS stars of  < 8) 1 , nearby galaxies including the Large Magellanic Cloud, and globular clusters to reduce scattered light in the images.
After applying these selection criteria, we are left with a sample of over 4000 clusters in the redshift range of 0.2 to 0.5.The number of clusters in each richness/redshift bin is listed in Table 1.We also include clusters in the redshift range of 0.5 to 0.6 in some of the analyses in this paper.However, because of the DES depth limit, we are concerned that our galaxy masking procedure (see Section 2.3) may be incomplete in this redshift range 2 and above.Therefore, the measurements of the 0.5 to 0.6 clusters are presented only for illustrative purposes and are not included in our quantification of diffuse light evolution.less than 2. 2 The redshift 0.5 to 0.6 clusters will be masked to 22.32 mag in z-band in mag_auto.Only ∼ 1% of the redMaPPer clusters reach 22.3 mag in the DES 10 z-band depth map continuously in the whole 0.15 deg 2 region around them and a depth-based cut would eliminate most of the clusters.For comparison, The redshift 0.4 to 0.5 clusters will be masked to 21.9 mag in z-band in mag_auto, while 67% of the redMaPPer clusters reach 21.9 mag in the DES 10 z-band depth map continuously in the whole 0.15 deg 2 region around them.Note that the DES coadd catalogs are generally over 95% complete above 23.7 mag in z-band, and the decision of a depth cut is made out of an abundance of caution.
We note that a redMaPPer galaxy cluster catalog based on the full six years of DES observations is also internally available to the DES collaboration.However, we opt to use the Year 3 version described here because its richness definition has better consistency with the Year 1 version in McClintock et al. (2019), which provides the richness-mass relation used in our estimations.The redMaPPer catalog based on Year 1 to Year 6 observations goes to higher redshift ( ∼ 0.9) than the Year 3 version in this paper (which is based on Year 1 to Year 3 observations), but both versions have excellent redshift coverage in the 0.2 to 0.5 redshift range studied in this paper.

DES Object Catalogs and Images
In this paper, we use the DES images and catalogs produced by the Dark Energy Survey Data Management (DESDM) project (Sevilla et al. 2011;Morganson et al. 2018).A detailed description of the DESDM pipeline can be found in Abbott et al. (2018).To summarize, the DESDM pipeline takes raw images from the Dark Energy Camera (DECam) (Flaugher et al. 2015), performs instrumental signature removals and corrections (Plazas et al. 2014;Gruen et al. 2015), flatfield corrections (Bernstein et al. 2018), full-focal-plane background subtractions (Bernstein et al. 2017b), as well as photometric (Burke et al. 2018) and astrometric (Bernstein et al. 2017a) calibrations to produce science-ready single exposure images.Those images are coadded (Bertin et al. 2002) into multi-epoch coadd images, which are used to produce object catalogs and photometry measurements by the Source Extractor software (Bertin & Arnouts 1996).The science-ready single exposure images are also used as input for photometry measurements such as the ngmix photometry measurements mentioned in the previous section.
For this work, we benefit from the full 6 years of DES operations (Diehl et al. 2018), and the DES Data Release 2 (DR2) processing campaign (Abbott et al. 2021) which includes not only more data, but also improved processing since the previous data release.Changes and improvements relevant to our analysis include: coadded images based on single-epoch full-focal-plane background subtraction, which do not include local background subtraction as applied to previous versions of DES coadd images (we use the "_nobkg" version of the coadd images, which do not have the local short-scale sky background subtracted as mentioned in Section 2.3); combining DES ,  and  band images into detection images as an average to create more robust faint objects detection; finally, changing the source detection threshold from 10 to 5 to produce more com-plete object catalogs.The DR2 coadd images have a combined sky coverage of 4913 deg 2 in DES , , ,  and  bands, and the 95% completeness of the coadd catalogs reaches 23.7 magnitude in the DES -band, with a 10  magnitude limit of 23.1.
We use both DES images and catalogs in this paper.Other than the redMaPPer cluster catalog, the coadd object catalogs used in this paper are constructed from DES coadd images using the Source Extractor software to detect and deblend objects.Moreover, we make use of the object's Auto measurements from Source Extractor to determine each object's masking area, which is based on Kron apertures and magnitudes (Kron 1980).The images used in this paper include both the science-ready single exposure images and the multi-epoch coadd images.In the next section, we describe how we use these images and catalogs in our workflow.

The Averaging/Stacking method
In this paper, we again use a "stacking" method described in Zhang et al. (2019c), which has also been adopted in Leung et al. (2020);Sampaio-Santos et al. (2021);Golden-Marx et al. (2023a).We present ICL and CG properties averaged over a large cluster sample.The "stacking" method proceeds as the following: (i) Coadd images of each galaxy cluster, are downloaded from the DESDM database3 .For each image, we mask all objects above a -band magnitude limit determined by the cluster's redshift.We exclude the redMaPPer CG from masking to preserve the CG light.The masking magnitude is chosen to be 0.2 * with  * being the characteristic luminosity of a cluster red galaxy luminosity function measurement (Zhang et al. 2019b).Assuming a faint end slope of -1, this masking limit would remove 82% of the light from cluster galaxies.
(ii) After masking, the radial surface brightness (SB) profile is derived from the masked images, as the mean pixel value in radial annuli on the images.
(iii) Similarly, radial SB profiles are derived for a sample of random points that cover the same sky area as the cluster catalog.In later steps, those random profiles are subtracted from the cluster profiles to eliminate residual backgrounds in the cluster profiles.
(iv) Using a Jackknife resampling method (see examples of applications in Norberg et al. 2009;Melchior et al. 2017), we divide the cluster samples, and the random points into 40 subsamples according to their sky coordinates 4 .For each sky coordinate subsample, we derive the CG+ICL SB profile by averaging the profiles of clusters and randoms, and then subtracting the random profiles from that of the clusters.We then apply the Jackknife resampling method (Efron 1982) to those sky coordinate subsamples to derive the means and uncertainties of the final CG+ICL measurements.
(v) Additional quantities, such as the CG+ICL color and luminosities are further computed from the CG+ICL SB profiles.In this paper, we also analyze the radial SB profiles of the cluster total light including CG, ICL and cluster satellite galaxies.Those measurements are derived using the same procedures listed here, but without the objects masking described in step (i).
We highlight one difference between this paper and Zhang et al. (2019c); Leung et al. (2020); Sampaio-Santos et al. (2021); Golden-Marx et al. (2023a).In Step (i), for each cluster, the previous analyses processed and coadded single exposure images from DESDM.For this work, we make use of the already coadded images from the DESDM database.Those DESDM images are based on the single exposure images, but coadded without applying local background subtraction steps in the swarp and Source Extractor software.In section 5.3, we compare the profiles derived from those coadded images and those based on the single exposure images Zhang et al. (2019c); Leung et al. (2020);Sampaio-Santos et al. (2021); Golden-Marx et al. (2023a).The results are highly consistent.
The computational resources needed for this method are not trivial.The masking of a 0.15 × 0.15 deg 2 region centered on one galaxy cluster, depending on the masking magnitude, can take a few minutes to hours with a single CPU processor.For this work, the masking and profile measurements of each cluster/random point is performed on the Open Science Grid5 , a High Throughput Computing Consortium.The processing of tens of thousands of clusters or random locations is distributed to thousands of parallel processes on the Open Science Grid in an "opportunistic" mode.Given the need to test and validate the analyses with different set-ups, we estimate that up to hundreds of thousands of CPU hours have been used in this work.

Richness and Redshift Dependence
Our first analysis is the surface brightness (SB) radial profiles of the diffuse light (CG+ICL), as well as the SB of the total cluster stellar content including the rest of the cluster galaxies.The goal of this analysis is to visually examine the shapes of those profiles and their general redshift/richness trends.The galaxy clusters are split into redshift and richness sub-samples, and those surface brightness profiles are presented in Figure 1.
We use three redshift bins from 0.2 to 0.5 to analyze the clusters.In each redshift bin, the clusters are further divided into four richness bins, with the richness binning defined in previous DES clusterlensing studies (McClintock et al. 2019) and listed in Table 1.As mentioned in Section 2, we apply the "stacking" procedure to each redshift/richness binned cluster subsample.The residual background for each subsample is derived from random points, but the masking magnitude limit for the random points is adjusted according to the cluster subsample's redshift range.For each redshift bin, we apply distance corrections so that the measurements are shifted to be in the observer frame of redshift 0.25.Those measurements are presented in Figure 1.In each of those richness and redshift bins, we measure the diffuse light profiles up to 1 Mpc from the center.Our measurements agree with previous studies which show the radial extension of ICL up to several hundreds of kpc, or even one Mpc from the cluster center (Krick & Bernstein 2007;Zibetti et al. 2005;Kluge et al. 2021;Li et al. 2022;Chen et al. 2022).
The upper panels of Figure 1 show the SB profiles of the galaxy clusters first split by redshift and then by richness.In each subpanel, the redshift range of the clusters is fixed to be the same and each line represents a different richness range.Those radial SB profiles show a clear richness dependence: richer galaxy clusters generally are brighter in SB, while less rich clusters are fainter.The trends are observed for both the diffuse light, and for the total light including cluster satellite galaxies.Moreover, the distinctions between the different richness subsamples are present throughout the three redshift bins, indicating robust richness dependence across the 0.2 We show both the surface brightness of the diffuse light (CG+ICL, red hues) and also the surface brightness of total cluster light (blue hues).Both profiles display strong richness dependence across the four redshift panels.Middle Row: Clusters in the same richness ranges in each panel, with different lines representing different redshift subsamples.Again, we show both the surface brightness of the diffuse light (CG+ICL, red hues) and also the surface brightness of total cluster light (blue hues).We do not observe consistent redshift trends across the richness panels, indicating weak, if any, signs for redshift evolution.We quantify the significance of the redshift/richness trends in the next Section.Lower Row: Surface brightness profiles after the cluster's radius has been scaled by a percolation radius (corresponding to the cluster subsample's average  200 ).The radial profiles of the diffuse light as well as the clusters' total stellar contents are "self-similar" after radial scaling.to 0.5 redshift range.Our result supports previous findings that detect strong ICL correlations with cluster mass (e.g., Gonzalez et al. 2005;Montes & Trujillo 2019;Huang et al. 2020;Sampaio-Santos et al. 2021;Kluge et al. 2021;Huang et al. 2022;Golden-Marx et al. 2023a;Ragusa et al. 2022;Chen et al. 2022).Further, the SB richness dependence in Figure 1 grows more prominent with enlarging radius.
In the middle panel of Figure 1, we present the SB profiles, first split by richness and then by redshift.In each subpanel, the richness range of the clusters is fixed to be the same and each line represents a different redshift range.Interestingly, those redshift-divided profiles appear similar within their uncertainty ranges; while fixing the cluster's richness range, we do not observe a consistent trend of the SB profiles being either brighter or fainter at a lower redshift.The lack of a consistent trend does not necessarily indicate that there is no redshift evolution for a richness-fixed sample, but potentially an evolution that is too small to be noticeable in those SB figures.In Section 4, we further quantify the redshift-related trends using their luminosity measurements.

ICL "Self-Similarity"
Previously we have noted the remarkable similarity of the ICL SB profiles after scaling by the cluster's radius in Zhang et al. (2019c); Sampaio-Santos et al. (2021).Similarly, we investigate this effect with a much larger sample in this analysis.
Given the relatively small richness range of each bin, we use a similar procedure as described in Sampaio-Santos et al. ( 2021) to scale the radial profiles.For each richness bin, we rescale the clusters' SB profiles by one radius determined by the average richness of each richness bin.Because there are no weak-lensing mass measurements for the galaxy cluster samples studied in this paper, we scale their radial SB profiles by the redMapper percolation radius which is a function of richness,  200 = 1.95 × (/100) 0.45 Mpc/ℎ.This radius relation is based on the  200m to richness relation derived using the DES Year 1 richness-mass relation (McClintock et al. 2019), and meant to be an approximation of  200m derived from cluster richness.We note that the richnesses are not necessarily consistent between different versions of the redMaPPer catalogs based on different DES data releases (a small difference has been found in preliminary comparisons).The percolation radii used here is a close but not necessarily accurate estimation of the clusters average  200 .
The last row of Figure 1 shows the SB radial profiles after scaling by the percolation radius.This row is meant to be compared to the top row of the same figure (without radial rescaling).These scaled profiles, both of the diffuse light and the cluster light, indeed appear to be much more similar across the richness bins, especially outside a transition radial range around 0.04  200 .In the central regions, the profiles do not appear to be similar after rescaling, suggesting that the CG SB profiles can not be well-described by scaled cluster radii.This phenomenon can be explained by an inside-out growth scenario (Oser et al. 2010;van Dokkum et al. 2010), such that CG stellar cores formed early at  > 2, while the accretion of CG outskirts and the ICL profiles happen later and are more influenced by the galaxy cluster's mass accretion process.

Volume-Limited Cluster Sample
In the previous subsection, we have shown that when fixing the cluster's richness, clusters in different redshift ranges have similar SB profiles.However, this does not answer the question of how ICL and CGs evolve with time in a single galaxy cluster whose richness will also evolve with time -more likely, their richnesses/masses will increase over time because of ongoing merging events.
In this sub-section, we account for cluster richness evolution by constructing a volume-limited cluster sample in different redshift bins.Specifically, in the highest redshift bin of 0.4 to 0.5, we compute the cosmic volume contained in this redshift bin, and select the clusters above richness of 50.For the lower redshift bins, 0.3 to 0.4 and 0.2 to 0.3, we again compute their respective cosmic volumes, and then adjust the richness thresholds for the cluster selections, so that each redshift sub-samples have the same cluster density given their different cosmic volumes.For the redshift slice of 0.3 to 0.4, the richness threshold becomes 55 and for the redshift slice of 0.2 to 0.3, the richness threshold becomes 61.Only clusters above those richness thresholds are selected for comparison in this sub-section.These The radial profiles of these redshift subsamples appear to be consistent within the measurement uncertainties, indicating a redshift evolution that is below a detection level.Again, we show both the surface brightness of the diffuse light (CG+ICL, red hues) and the total cluster light (blue hues).We will further quantify the redshift evolution of the cluster luminosities in the next section.
selections ensure that the cluster samples have the same spatial densities in each redshift bin.A similar volume-limited selection method is also use in Golden-Marx et al. (2023a), with the distinction that Golden-Marx et al. (2023a) select clusters based on SZ-computed masses, while this analysis is using optical richness.
The SB-profiles of those volume-limited samples are presented in Figure 2. We again do not observe significant differences between the redshift sub-samples, as the previous subsection has already noted no visible differences between redshift subsamples within fixed richness ranges.In this exercise, we limit the analysis to a high richness threshold which tends to have lower richness-to-mass scatter (Farahi et al. 2019;Anbajagane et al. 2020) and is less subjective to potential systematic effects such as line-of-sight projections (Costanzi et al. 2019;Abbott et al. 2020;Wetzell et al. 2022;Wu et al. 2022).For the same reason, we also do not sub-divide the clusters according to richness.The consequently smaller sample size lowers the significance of a possible redshift trend.We further quantify the redshift-related trend in Section 4.

Color Radial Profiles
We derive the  −  color of the ICL radial profiles using measurements of the DES  and  band SB profiles.These color measurements are shown in Fig. 3. Given that such measurements require highly-significant ICL SB profiles in two bands, we only show color measurements out to a radius slightly beyond 400 kpc.
The color profiles are presented in redshift subpanels with clusters further divided in richness bins.As previously noted in the literature, the color of the diffuse light displays a radial gradient, becoming bluer at a larger radius.Interestingly, we also notice a consistent, although not significant, richness trend in those colors, which appear to be redder in richer clusters.In addition to the diffuse light color profiles, we have acquired the color measurements of the cluster's total stellar content, but the measurements are much more uncertain because of Poissonian noise.However, we note that the average color of the cluster's total stellar content is generally consistent with the color of the diffuse light.
To further quantify the colors, we fit those measurements as a The diffuse light color is consistent with the average color of the cluster's total stellar content in the center.In addition, a radial gradient can be seen in all of the redshift and richness ranges, such that the diffuse light becomes bluer at larger radii.This is consistent with previous studies that ICL consists of more metal-poor and younger stars (Edwards et al. 2020) and suggests an ICL origin from galaxy disruption and tidal stripping (Contini et al. 2019).
function of radius: Here,  is the radial slope of the colors, and  is the intercept of the profile at  = 1.0 kpc.The fitted parameters are shown in Figure 4.
In each redshift/richness bin,  appears to be negative, indicating a robust detection of a radial gradient.This is consistent with previous measurements of ICL radial color gradient (e,g., Zibetti et al. 2005;Chen et al. 2022;DeMaio et al. 2018;Yoo et al. 2021) and that the ICL consists of more metal-poor and younger stars (e.g., Edwards et al. 2020) compared to the CG.As mentioned and analyzed in Contini et al. (2019), the color radial gradient suggests that the ICL's origin is from galaxy disruption and stripping: if clusters acquire ICL mainly through merging, ICL would have a relatively uniform color because of stellar population mixing.On the other hand, the disruption and stripping of cluster satellite galaxies would produce a radial gradient because of the radial dependence of those processes.
In addition to the radial gradient, we also detect a possible richness dependence in color, as the intercept () of the fitting results appears to be redder (higher positive value) in richer clusters.This is possibly related to richer and thus more massive clusters containing a higher fraction of red sequence galaxies than the less massive clusters (Hansen et

Richness and Redshift Dependence
To further quantify the ICL's richness and redshift dependences, we examine the luminosities of the diffuse light and the cluster's total stellar content.Those luminosities are derived by integrating the SB profile in radial ranges as the following: Here, (,  = 0.25) is the SB measurements presented in the previous section, which have been distance-corrected as if it was observed at redshift 0.25.Thus,  ( < ,  = 0.25) is the luminosity measurements enclosed within radius , but shown as the apparent magnitude in the observer frame of redshift 0.25.We choose this pivot redshift because it is close to the median redshift value in the lowest redshift subsample.Figure 5 shows the radial profiles of the integrated luminosity as a function of the outbounding .
There are two significant trends in those luminosity profiles.First, richer and thus more massive clusters contain more diffuse and total light.The richness dependence is observed across the four redshift ranges.Second, at small radii, within ∼ 50 kpc, diffuse light contributes to the bulk of the cluster total stellar light.Outside ∼ 50 kpc, the cluster's total stellar light increases significantly because of the contribution from cluster satellite galaxies.As a result, diffuse light appears to grow less significantly than the total stellar content with radius.We also show the luminosity radial profiles as a function of radius scaled by  200 .The richness/mass dependence of the luminosity becomes even more pronounced in those scaled radius plots.Luminosities (distance-corrected to be the apparent magnitude in the observer frame of  = 0.25) enclosed within 4 radial bin (0 to 30 kpc, 30 to 80 kpc, 80 to 300 kpc and 300 to 600 kpc).We examine how the luminosities change with redshift and richness.The luminosities of the total cluster stellar content (blue lines) and the cluster diffuse light (CG+ICL, red lines), increase significantly as the cluster's richness increases.However, those luminosities do not appear to change significantly with redshift, except in the lowest richness range and within 30 kpc.See Section 4 for quantitative analyses.We further investigate how these luminosities change with redshift and radial apertures.For clusters in a fixed richness range, we derive their luminosities enclosed within 30 kpc and in 30 to 80, 80 to 300 and 300 to 600 kpc annuli.The innermost 30 kpc radial bin is chosen to match the CG size, while the second radial bin, out to 80 kpc (we have experimented with 50 kpc, 75 kpc and 100 kpc, and found 80 kpc to be most representative in terms of the redshift trends), is chosen to probe the CG to ICL transition range.Finally, the 80 to 300 and 300 to 600 kpc annuli are chosen to probe the extended components of ICL.Those luminosity measurements in apertures/annuli are shown in Figure 6.
Interestingly, with the lowest richness sample -the clusters in the richness range of 20 to 30 -we notice that both their diffuse light and total light appear to get brighter towards lower redshift, indicating redshift evolution.In some of the higher richness bins, the luminosities of both the cluster's total and stellar light appear to be getting brighter or unchanged towards lower redshift within 30 kpc.However, outside of 30 kpc, there is no sign of consistent redshift evolution.
To further quantify the richness dependence and redshift evolution in those different apertures, we fit the measurements to the following relation: In this relation, the total amount of light ( ( 0 <  <  1 ), in the unit of magnitudes) contained within an aperture or annulus, is fitted with a linear relation to the logarithmic values of the cluster subsample's average richness and redshift.The richness and redshift dependences are described by parameters  and  respectively, and a non-zero value would indicate detection of dependence.In the relation, the pivot richness is chosen as 33, which is the mean richness value of the sample, while the pivot redshift is chosen as 0.25, which is close to the median redshift value in the lowest redshift subsample.Thus, the intercept of the relation, , can be interpreted as the apparent magnitude of a richness 33 and redshift 0.25 cluster.The fitting of those parameters , , and , is performed with a Markov Chain Monte Carlo (MCMC) method, and the likelihood is constructed from the  2 values between the measurements and the relations (using the uncertainties of the measurements as the weighting).Table 2 shows the derived posterior values of the , , and  parameters.
The fitted results confirm our fore-mentioned observations.First, the values of  deviate from 0 at very significant levels, between ∼ 7 to 16  levels in all of the analyzed apertures for both the diffuse and the total stellar content.This result confirms the significant richness, and thus cluster mass dependence for the diffuse and total stellar content.Further, the value of  is increasingly negative at large radii, suggesting stronger richness, and thus mass, dependence for the diffuse light as well as the cluster total stellar content.
Also consistent with our fore-mentioned observations, for the diffuse light, the value of parameter , which quantifies redshift evolution, is consistent with 0 outside 30kpc.This indicates that when fixing cluster richness, the amount of diffuse light is not obviously increasing (or decreasing) towards lower redshift.Cluster mass may be the main driver for diffuse light evolution in a large radial bin.However, within 30 kpc, the amount of diffuse light appears to be increasing towards lower redshift, suggesting that the amount of stellar light associated with the CG is building up over time.According to the fitting results, the luminosity of the CG becomes brighter by 0.113 magnitude (calculated from  × log 10 (1.45/1.25))from redshift 0.45 to 0.25 within 30 kpc.This brightening corresponds to a flux increase of 11%.The  value of the diffuse light between 30 to 80 kpc is consistent with 0, but we have also adjusted the 80 kpc outbound range between 50 and 100 kpc, and found that the  value is larger in apertures closer to 30 kpc, indicating possible additional redshift evolution closer to the CG.
For the cluster's total light,  is positive within 30 kpc, also indicating some growth with time (brightening towards lower redshift).Because CG is the dominating component of cluster total light within 30 kpc, the growth in this radial range mostly reflects CG growth.In annuli outside 30 kpc,  is generally consistent with 0, indicating no evidence of significant redshift evolution.

Volume-Limited Cluster Sample
We further investigate the luminosity evolution of the volume-limited cluster sample discussed in Section 3.3, which helps answer the question of ICL growth when tracking the same cluster's evolution over time.Similarly, we calculate the luminosities enclosed within radial bins for both the cluster diffuse light and total light and show how they change with redshift in Figure 7.A sign of redshift evolution can be seen within 80 kpc.We do not find evidence of redshift evolution in the rest of the radial bins.
To quantify the redshift evolution of this cluster sample, we also fit their measurements to the following relation: In this relation,  quantifies the redshift evolution of this cluster sample.A positive value would indicate brightening luminosity over time (towards lower redshift), while a negative value would indicate the opposite.The fitting procedure is performed with the curve_fit function of Scipy, and the derived values and uncertainties of  are noted in Figure 7.
With both the diffuse light and the cluster total light, we detect positive values of  above the significance level of 1 within 30 kpc in the CG as well as in the CG to ICL transition region, between 30 to 80 kpc range.In both of the radial ranges, the redshift evolution, indicated by the positive value of , appears to be larger in the diffuse light than in the cluster total light.This result, together with the results from the richness/redshift subsamples, indicates that the increase in stellar content with time associated with diffuse light is partially driven through the deposition of new material onto the CG and at the CG's outskirt.For the larger radial bins, we again find  to be consistent with 0, and thus, do not find evidence of redshift evolution.It is possible though, that the redshift evolution in those radial bins falls below our measurement limit as the uncertainties of  are large.We do not observe consistent redshift or richness-dependent trends (except in the 0.021 to 0.056  200 bin) in the measurements.However, the diffuse fractions appear to be dropping at larger radii.See Section 4.3 for more detailed discussions.

Diffuse fraction
The fraction of CG and ICL in the cluster's total stellar content is an important quantity (e.g., Gonzalez et al. 2013;Behroozi et al. 2013;Ragusa et al. 2022;Joo & Jee 2023).The build-up of the ICL, CG, and cluster total light is not necessarily aligned over time.For example, the ICL and cluster total light may have gone through more or less significant growth in the recent era compared to the CG.Thus, one may observe a change over time in the CG/ICL to cluster total light ratio.
Based on the luminosity measurements, we quantify the fractions of CG and ICL in the total cluster stellar content, with the following equation, In this equation, Lum total ( 0 <  <  1 ) enclosed within radius  0 and  1 is the luminosity of the cluster's total stellar content derived from its brightness, and Lum diffuse ( 0 <  <  1 ) is the luminosity of the diffuse light derived from its surface brightness.We refer to the ratio of these measurements as the diffuse fraction in this analysis.
Figure 8 shows those fractions derived for clusters in different redshift/richness ranges.Other than the 0.021 to 0.056  200 radial range, the diffuse fraction appears to stay unchanged with richness and redshift, indicating similar richness and redshift dependence between diffuse light and cluster total stellar content.However, the diffuse fraction does appear to decrease at a large radius, which is close to 40% in the 30 to 80 kpc range, but decreases to ∼ 20% in the 300 to 600 kpc range.Between the 0.021 to 0.056 200 radial range, the diffuse fraction does appear to decrease with richness, but this is likely because of the scaling of  200 with richness; a richer cluster would have a higher  200 value, which excludes more of the BCG outskirt with a 0.021 200 cut.
In addition to those trends, Figure 8 highlights the importance of selecting a radial and cluster mass/richness range when studying diffuse fractions.The fractions appear to drop with an increasing radius, but also is dependent on whether or not the measurements are made in physical radii or radial units scaled by cluster radius.Given the discrepancies in literature reports on diffuse fractions, fair comparisons will need to be made between cluster samples of comparable masses in similar radial scales.

PSF effect
Because the point spread function is known to have extended wings (Moffat 1969;King 1971;Racine 1996;Bernstein 2007), it would contribute to the extended low surface brightness features of galaxies or galaxy clusters.The radial scales we probe in this paper are significantly larger than the PSF FWHM of the DES images, therefore we expect minimal PSF contributions to the ICL detection (see discussion in Zhang et al. (2019c)).On the other hand, those contributions may change with redshift given the change of angular distance scale with redshift.Thus, we perform image simulations to probe the possible effect of PSF on the results presented in this paper.
To do so, we convolve a PSF model with an analytical diffuse light profile model and examine the differences before and after PSF convolution.We generate mock 2D images of diffuse light using an analytical model, setting the angular scale of the analytical diffuse light profile models at four redshifts, 0.25, 0.35, 0.45 and z = 0.55 (but without adjusting their surface brightness level as we are only looking at before-and-after PSF convolution differences).These 2D images are then convolved with a 2D PSF image model.Both the analytical models and the PSF models are based on the DES-Year1 measurements in Zhang et al. (2019c) in r-band, as the PSF models have similar large radial behaviors outside 2 arcseconds.We then derive the SB measurements in radial bins before and after PSF convolution.The results are shown in Figure 9.
The top panel shows the flux changes of the profiles before and after the PSF convolution for the three profiles at different redshifts.PSF convolution flattens the central regions of those profiles limited by the pixel scale of the images (0.263 arcsecond pixel scale).The middle panel of Figure 9 shows the relative changes in those profiles before and after convolution.Outside of 10 kpc, PSF convolution has a minor effect on SB measurements which change by less than 10%, but the change depends on redshift.Outside of 100 kpc, PSF effects appear to be negligible for all of the four redshifts, which is less than 1% at 100 kpc for  = 0.55.
As with the integrated (within radius) brightness measurements, similarly, the PSF effect appears to be negligible if integrating to 20 kpc, affecting less than 5% of the flux measurement, or around 2% if integrating to 30 kpc.Within 10 to 20 kpc, the PSF may affect the CG flux measurements by up to 12%, depending on the redshift.Within 10 kpc, the integrated luminosity needs to be carefully interpreted due to the PSF effect.
We conclude that PSF effect alone can not account for the redshift evolution in the diffuse light luminosity measurement within 30 kpc, which shows a change of ∼ 0.2 mag, or ∼ 20% in flux from redshift 0.45 to 0.25 (Section 4).With a carefully designed CG aperture (30 kpc in this analysis), our luminosity redshift evolution results should be minimally affected by the PSF effect.

Masking Magnitude Limit
The masking magnitude limit we use for this work varies with redshift.This may affect the results of this paper when cluster galaxies below the masking limit contribute a noticeable amount of light to the diffuse light measurements.We acknowledge this issue as a limitation in our analysis, as we do not explicitly account for the contributions from the fainter cluster member galaxies below the masking limit.
We test how much our results may have been affected by these magnitude limits.In this test, we redo the measurements of the diffuse light using a masking magnitude that is fainter by 0.7526 mag (or Figure 10.Diffuse light profiles derived when using a deeper masking magnitude limit (upper panels) and the relative differences to the fiducial measurements presented in previous sections (lower panels).The surface brightness measurements outside 100 kpc can be lowered by 10% when using a deeper masking magnitude, which removes more contamination from the faint cluster satellite galaxies.masking to 0.1 * of the cluster luminosity function), and compare the results to those from the fiducial analyses presented earlier.We have not applied this deeper magnitude limit in our fiducial analysis because of the increasingly incomplete galaxy detection associated with this magnitude limit, which would render the results in the redshift 0.4 to 0.5 bin less reliable.Nevertheless, we show the SB measurements with this deeper magnitude and the comparison to the fiducial analysis in Section 4.
Indeed, using a deeper masking limit notably reduces the surface brightness measurements of the diffuse light throughout the redshift and richness bins.Outside of 100 kpc, the reduction in flux consistently reaches a ∼ 10% level, although there are significant fluctuations as indicated by the uncertainties.Given that satellite galaxies 2.5 to 5 times brighter than the ICL in this radial range (Section 4) are excluded, a reduction of 10% in flux means that the deeper magnitude is only further removing 2% to 4% of the faint cluster satellite galaxy contribution.A deeper masking magnitude is unlikely to significantly further reduce ICL brightness unless there is a noticeable upturn in the cluster galaxy luminosity function at the faint end (e.g., Lan et al. 2016).
Other than the masking limit as well as the PSF effect, there are other additional effects that influence our results.Another issue related to masking is that the masking aperture does not enclose all of the light from cluster satellite galaxies.A galaxy's light can reach tens or even hundreds of kpcs.In Zhang et al. (2019c), we found that the aperture of masking only affects diffuse light measurements at a percentage level.In addition, in this analysis, we have adjusted the masking radius to be 3.5 Kron radii rather than 2.5 Kron radii which will further reduce the effect.Moreover, the cluster galaxy luminosity function may evolve with redshift.However, recent literature studies find that the redshift evolution of the cluster galaxy luminosity func-tion is very mild at most (e.g., Hansen et al. 2009;Sarron et al. 2018;Zhang et al. 2019b;Puddu et al. 2021).

Sky Background
Accurate diffuse light measurements require accurate evaluation and removal of the sky background in optical images.Similar to Zhang et al. (2019c), in this paper, the images we use have removed sky background that is estimated over the whole field-of-view (FOV) of DECam, approximately 3 deg 2 , using a PCA method (Bernstein et al. 2017b).Given that one galaxy cluster, even at redshift 0.2, only makes up a very small area in the DECam FOV, the sky background estimation is not sensitive to the presence of galaxy clusters, thus avoiding a background over-estimation issue that often plagues ICL measurements.Zhang et al. (2019c) tested the DECam FOV PCA background evaluations for ICL measurements, and it was shown that the PCA sky estimations at the cluster centers and at a large cluster radius (1.36 arcmin from the cluster center) are highly consistent.
After removing the full FOV sky background level, the images still possess a residual background.Since we average the measurements for several hundreds and sometimes several thousands of clusters, we estimate a residual background for those averaged measurements, through processing "sky randoms" that track the area coverage of the cluster sample.A surface brightness profile of the sky randoms is acquired using the same procedure of the cluster measurements.Those "random" profiles are subtracted from the "raw" cluster measurements to acquire the final cluster-related measurements.The top panel of Figure 11 illustrates the procedure.
In Figure 11, because of the residual background, the "raw" cluster measurements still have a SB level of ∼ 2 /kpc 2 in flux measured at large radii (∼ 2 Mpc), but this residual is also present in the Panel: Differences in these profiles.These two approaches yield consistent surface brightness measurements to an accuracy level of over 30 mag/kpc 2 in terms of raw diffuse light and random profile measurements.After random profile subtractions, the differences vanish at a surface brightness level of 40.5 mag/kpc 2 .
Figure 12.Surface brightness of the randoms when using different distance cuts to bright objects in the sky.Our analysis requires the bright objects identified in a DES masking file to be 0.2 deg away from the cluster center.Using further cuts would lower the surface brightness measurements of randoms because of less contamination from the bright objects.
"random" measurements.After subtracting the randoms, the final cluster measurements fluctuate around 0 at very large radii (∼ 2 Mpc).Note that in Figure 11, we are showing the measurements in DES "flux" (10 −12 of a maggy), where the "flux" used here is a linear measure of an object's brightness, as opposed to the logarithmic "magnitude" unit of brightness with the following relation mag = 30 − 2.5 × log 10 (flux).
We note the importance of using random catalogs that faithfully trace the sky coverage of the redMaPPer cluster catalog.The raw profile measurements of randoms in Figure 12 are sensitive to the selection of the random catalogs (and thus the redMaPPer cluster catalogs).These two catalogs are selected to avoid sky regions that contain bright foreground galaxies and stars -at least 0.2 deg away from their centers.If we adjust the distance cuts to 0.3 deg or 0.4 deg, the random's profile value would become lower, indicating different "residual" background levels.
Finally, a crucial difference between this paper and Zhang et al. ( 2019c) is that we use the coadded images from the Dark Energy Survey directly, which is based on coadding single epoch images after the PCA sky background subtraction.The DES coadd images (the "no-bkg" coadd images in the DES data release, which did not subtract local background) are based on the procedure in Zhang et al. (2019c) to better preserve low-surface brightness features.
To test that the DES coadds are indeed suitable for detecting intracluster light, we separately process the redshift 0.2 to 0.35 clusters by coadding single epoch images using the same procedures in Zhang et al. (2019c) and compare the measurements to the DES special coadd-based measurements.Their differences are shown in Figure 11.
The raw SB measurements of those clusters and the randoms from both sets of images are offset at a flux level of 0.2, corresponding to a surface brightness level of 31.7 mag/kpc 2 .Those raw measurement differences between the two coadding procedures are likely caused by pixel weighting differences.After the random subtraction, the measurements are similar at a surface flux level of 0.015, which means that the two methods are similar to a surface brightness level of 40.5 mag/kpc 2 , and thus highly consistent.

Comparison to Simulation
To gain theoretical insights into the evolution of ICL, we turn to the IllustrisTNG simulation suite (Nelson et al. 2019;Pillepich et al. 2018;Springel et al. 2018;Nelson et al. 2018;Naiman et al. 2018;Marinacci et al. 2018) to examine how the diffuse stellar components of galaxy clusters change with cluster mass and redshift.This has already been a subject of investigation in Pillepich et al. (2018).
Our analysis here is based on the TNG300-1 simulation, which has the largest volume (300 Mpc 3 ) in the IllustrisTNG simulation suite, and also the highest simulating resolution among the 300 Mpc 3 volume series.The TNG300-1 simulation contains 263 dark matter halos above the mass threshold of 6 × 10 13 M ⊙ /ℎ at redshift 0.27.It has an advantage over the smaller-volume simulations (for example the TNG 100 Mpc 3 and 50 Mpc 3 series) which contain much smaller samples of cluster-sized dark matter halos despite their higher simulation resolutions.
We select the redshift snapshots at 0.27, 0.42 and 0.58 for this analysis, to represent the redshift range studied in this paper.In each redshift snapshot, we select dark matter halos with  200m above 6 × 10 13 M ⊙ /ℎ as "galaxy clusters".After cutting dark matter halos that are too close to the simulation box boundaries (within 20 cMpc/h6 ), we are left with 205, 155 and 115 dark matter halos respectively in the three redshift snapshots.Those dark matter halos will be referred to as galaxy clusters in the rest of this analysis.For each simulated cluster, centered on its weighted mass center, we select the diffuse stellar particles contained in 3D distance apertures and compute their total stellar masses.Those stellar masses are shown together with the host halo mass in each redshift snapshot in Figure 13.In addition to the diffuse stellar component, we also included the dark matter halo's total stellar content (subhalo+diffuse) within those radial apertures for comparison.
In this simulation, both the diffuse and total stellar components of galaxy clusters steadily increase as the galaxy cluster mass increases.This mass dependence grows steeper in larger radial apertures.On the other hand, examining the mean of those stellar masses ( * ) as a function of halo mass ( 200m ), there does not appear to be tangible differences in different redshift snapshots, indicating no redshift evolution.
To further quantify the mass dependence and redshift evolution in the simulation, we fit the halos's stellar masses  * , halo masses  200m and redshifts  0 to the following stellar-mass and halo-mass relation.
This relation is similar to the one adopted in Section 4, substituting richness-dependence for mass-dependence.
The fitted constraints on the relation are listed in Table 3.In the relation, parameter  quantifies the mass-dependence of the stellar masses.Its value is positive in all of the radial bins, and becomes even more positive at larger radii, in agreement with the steeper mass dependence we have seen in Figure 13.On the other hand, parameter  quantifies the redshift evolution.Its value is consistent with 0 in all ) that did not find redshift evolution in the stellar-mass to halo-mass relation of galaxy clusters, in either the diffuse component or the subhalo component.However, both stellar mass components scale strongly with halo mass.We note that in small radial ranges, the properties of halo central galaxies or massive galaxies in the simulation do not always match observations (e.g., Pillepich et al. 2018;DeMaio et al. 2020;Li et al. 2019;Cannarozzo et al. 2023).Nevertheless, those simulation results qualitatively agree with our measurements in the large radial ranges outside of 80 kpc.

Comparison to Literature
Perhaps the most surprising result from this paper is the relative lack of ICL evolution at a radius larger than 80 kpc.Many analyses characterizing CG and ICL growth, including the work of the coauthors of this paper, have predicted significant growth of the ICL (e.g., Behroozi et al. 2013;Zhang et al. 2016;Contini et al. 2018;Golden-Marx et al. 2022) as a mechanism to explain the relatively slow CG growth observed below redshift 1.5 (e.g., Stott et al. 2010;Lidman et al. 2012;Lin et al. 2013).However, we do find signs of diffuse light redshift evolution in the CG as well as in the CG to ICL transition within 80 kpc.
Prior to this analysis, there have been few works that analyze large samples of ICL profiles over a broad range of redshift to directly quantify their redshift evolution.One of the most comparable literature studies to our work is presented in DeMaio et al. (2020), which analyzed 42 clusters in the redshift range of 0.05 to 1.75. DeMaio et al. (2020) measured the BCG and ICL growth out to about 100 kpc from the cluster center, and found that the stellar masses of BCG and ICL increase more rapidly than the cluster's total mass from redshift 1.5 to the present.They conclude that BCG+ICL growth is not solely driven by cluster mass growth.In this analysis, we indeed observe that the CG and ICL luminosity increases mildly within 30 kpc.
There are some differences between our work and that of DeMaio et al. (2020).DeMaio et al. (2020) find that the ICL grows by a factor of 1.08 ± 0.21 from redshift 1.55 to 0.4 when examining the 10 kpc to 100 kpc range.While in this paper, we find evidence of ICL growth by 11% within the 30 kpcfrom redshift 0.45 to 0.25.In our analysis, the results are derived for clusters in a time span of roughly 1.7 Gyr (redshift 0.45 to 0.25).The ICL growth observed in DeMaio et al. ( 2020) occurs over an extended period of 4.97 Gyr from redshift 1.55 to 0.4.Interpolating from their measurements, the ICL measured between 10kpc to 100 kpc can grow by 37% in 1.7 Gyr, significantly higher than the 11% measured in our work in the 30 kpc range where we see the most growth.On the other hand, DeMaio et al. ( 2020) have noted a slow-down in the BCG and ICL growth after redshift 0.4, that there's no change in the diffuse light stellar mass (between 10 and 100 kpc) to halo mass relation from redshift 0.4 to 0.1.
Furthermore, in our work, we do not find signs of ICL growth outside 80 kpc.On the other hand, Golden-Marx et al. (2023a) studied ICL growth from redshift 0.8 to 0.2, but also do not find much evidence for ICL growth.Golden-Marx et al. (2023a) defines ICL with a large radial aperture of between 50 to 300 kpc and those results are based on the same imaging data set and processing method as in this paper.In both Golden-Marx et al. (2023a) and this work, we are limited by the PSF resolution (as discussed in Section 5) to probe a smaller radial range such as 10 to 30 kpc.Combining the findings from Golden-Marx et al. (2023a);DeMaio et al. (2020) and this work, we speculate that the CG, as well as the region close to the CG within 100 kpc, rather than the ICL at a very large cluster radius, holds the key for explaining CG and ICL growth.However, the growth may not be very noticeable below redshift 0.45.
Another comparable analysis is from Furnell et al. (2021), which studied ICL growth over the redshift range of 0.1 to 0.5, using 18 X-Ray selected clusters with Hyper Suprime Cam Subaru Strategic Program observations.Using a radial aperture of  500 and a surface brightness limit of 25mag/arcsec 2 , Furnell et al. (2021) find that the ICL fraction increases by a factor of 2 − 4 over the 0.1 to 0.5 redshift range with no obvious mass dependence.However given that the ICL definition in Furnell et al. (2021) is based on a surface brightness limit, a radius aperture of  500 that scales with cluster mass, as well as a "divot" correction due to background subtraction in the images, it is possible that the ICL definitions in their and our analyses are not directly comparable.

SUMMARY AND PROSPECTS
In this paper, we present measurements of the CG and ICL radial profiles using the full 6 years of DES data.The major findings from those measurements can be summarized as the following: (1) the diffuse light (CG+ICL) extends to 1 Mpc in the redshift range of 0.2 to 0.5 investigated in this analysis.Prior to this analysis, multiple studies have already detected ICL in the several hundreds of kpc to Mpc radial range, which includes both "stacking" based analysis like this paper (e.g., Zibetti et al. 2005;Chen et al. 2022), and deep imaging studies of individual galaxy clusters (e.g., Krick & Bernstein 2007;Kluge et al. 2021;Golden-Marx et al. 2022).Our finding again showcases the wide radial reach of ICL.There may be much to study with the radial properties of ICL.
(2) We find that the diffuse light surface brightness and luminosity strongly depend on richness -a galaxy cluster mass proxy.This dependence is stronger at large radii outside of 50 kpc from the cluster center.The richness and thus cluster mass dependence appears to be the major factor behind the differences between diffuse light observations in different subsamples, as their radial profiles scale well with the cluster's radius ( 200 ) and their fractions in the cluster's total stellar luminosity appears to be richness-independent.The results agree with previous studies that find a strong mass correlation between ICL luminosity or stellar mass, or a possible correlation between the cluster mass distribution and ICL surface brightness (e.g., Montes & Trujillo 2019;Huang et al. 2020;Sampaio-Santos et al. 2021;Kluge et al. 2021;Ragusa et al. 2022).
Perhaps most interesting of all to cluster cosmology studies, this finding, again, suggests the potential of ICL as a cluster mass proxy (Golden-Marx et al. 2023a), or to help improving cluster finding algorithms (Huang et al. 2022).Cosmology studies based on galaxy cluster abundance measurements have long emphasized the importance of developing accurate and precise cluster mass proxies (i.e., galaxy cluster observables that scale well with masses), because a mass proxy with low scatter to the cluster's true mass can significantly reduce the requirement for follow-up observations, and thus reduce the derived uncertainties on cosmological parameters such as Ω  and  8 (Rozo et al. 2010).Further, the precision of those cosmology studies also depends on having an accurate mass proxy that is not affected by the cluster's large-scale structure environment (Wu et al. 2022).It will be particularly interesting to incorporate diffuse light quantities in developing cluster mass proxies or cluster finding algorithms (Huang et al. 2022).
(3) Perhaps with a bit of a surprise, we find that the diffuse light at large cluster radii (outside of 80 kpc from the cluster center) is not evolving much with redshift in the 0.2 to 0.5 range.Closer to the cluster center, within 80 kpc, we have found some evidence that the diffuse light's luminosity increases with time (towards lower redshift).We speculate that ICL build-up may be more pronounced closer to the CG, while at large radii, on the scale of hundreds of kpc, ICL build-up is more in tune with the cluster mass build-up, which also explains the stronger cluster mass dependence at large radii.
In the context of CG and ICL co-evolution studies, many (including the authors of this paper) have speculated a more rapid growth of ICL than the BCGs below redshift 1.Given that ICL and CG is often vaguely defined in those studies, our findings suggest that ICL growth happens at a much smaller radius (i.e., in the CG or at the CG outskirt) than we previously expected.On the other hand, our finding of little redshift evolution at large cluster radius is in excellent agreement with the hydrodynamic simulation study of IllustrisTNG (Pillepich et al. 2018), which finds little redshift evolution in diffuse light stellar mass once the cluster's halo mass is fixed.
(4) We have measured additional properties of ICL: the color profile of diffuse light has a radial gradient, which becomes bluer at a larger radius, and also bluer in less rich/massive clusters.In addition, the diffuse light SB profiles appear to be "self-similar" after scaling by the cluster radius, and that ICL fraction in the total cluster stellar light appears to be dropping at a larger radius.
Moving forward, there are multiple follow-up opportunities from our measurements.For example, in this paper, we have only studied the average properties of galaxy clusters in richness-redshift subsamples using a "stacking" method.As demonstrated in Golden-Marx et al. (2023a), it is possible to acquire diffuse light measurements of individual galaxy clusters, especially within the 300 kpcs radial range.This would allow us to study how diffuse light properties may change with cluster ellipticity, dynamical state, or with CG properties.It may also be interesting to compare the diffuse light to other galaxy cluster measurements, such as their weak lensing signals as done in Sampaio-Santos et al. (2021).
That said, there are also limitations in this study, especially related to the masking depth as discussed in Section 5.The redshift evolution results reported here are limited by the masking depth of cluster galaxies detected by DES.Faint or undetected cluster galaxies below the masking magnitude limit would have blended into our diffuse light measurements.In this analysis, we use the luminosity function and a test with a deeper magnitude limit to argue that the contribution from those faint galaxies does not affect our redshift evolution conclusions.However, this masking issue can be largely avoided by using a much deeper photometric catalog to more thoroughly mask the contribution of cluster galaxies.Future cosmic surveys like the Legacy Survey of Space and Time (LSST) from the Vera C. Rubin Observatory will be able to provide such a photometric catalog.
On a different note, those future surveys will also provide many more photons, and a much larger cluster sample for this "stacking" (averaging) method, significantly improving the accuracy of diffuse light measurements.Meanwhile, space-based cosmic survey programs like Euclid and the Nancy Grace Roman Telescope can provide imaging data that are less affected by sky background.We look forward to using those data in the coming years.

DATA AVAILABILITY AND ACKNOWLEDGEMENTS
This paper is dedicated to the memory of the pioneering Egyptian/American astronomer, Sahar Allam, a woman whose wisdom, bravery, care, sensitivity, and sense of humor have guided and supported us in the past decade and during the most difficult times.We miss you dearly.
Our analyses are performed with a few software packages, including Astropy (Astropy Collaboration et al. 2013, 2018, 2022), Numpy (Harris et al. 2020), Scipy (Virtanen et al. 2020), andEmcee (Foreman-Mackey et al. 2013).The data underlying this article were accessed from the DES and Illustris-TNG database.The derived data generated in this research will be shared on reasonable request to the corresponding author.
The IllustrisTNG simulations were undertaken with compute time awarded by the Gauss Centre for Supercomputing (GCS) under GCS Large-Scale Projects GCS-ILLU and GCS-DWAR on the GCS share of the supercomputer Hazel Hen at the High Performance Computing Center Stuttgart (HLRS), as well as on the machines of the Max Planck Computing and Data Facility (MPCDF) in Garching, Germany.
Funding for the DES Projects has been provided by the U.S. Department of Energy, the U.S. National Science Foundation, the Ministry of Science and Education of Spain, the Science and Technology Facilities Council of the United Kingdom, the Higher Education Funding Council for England, the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign,

Figure 1 .
Figure 1.Surface brightness of the clusters in richness/redshift ranges.Upper Row: Clusters in the same redshift ranges in each panel, with different lines representing different richness subsamples.We show both the surface brightness of the diffuse light (CG+ICL, red hues) and also the surface brightness of total cluster light (blue hues).Both profiles display strong richness dependence across the four redshift panels.Middle Row: Clusters in the same richness ranges in each panel, with different lines representing different redshift subsamples.Again, we show both the surface brightness of the diffuse light (CG+ICL, red hues) and also the surface brightness of total cluster light (blue hues).We do not observe consistent redshift trends across the richness panels, indicating weak, if any, signs for redshift evolution.We quantify the significance of the redshift/richness trends in the next Section.Lower Row: Surface brightness profiles after the cluster's radius has been scaled by a percolation radius (corresponding to the cluster subsample's average  200 ).The radial profiles of the diffuse light as well as the clusters' total stellar contents are "self-similar" after radial scaling.

Figure 2 .
Figure2.Surface brightness of a volume-limited cluster sample.The radial profiles of these redshift subsamples appear to be consistent within the measurement uncertainties, indicating a redshift evolution that is below a detection level.Again, we show both the surface brightness of the diffuse light (CG+ICL, red hues) and the total cluster light (blue hues).We will further quantify the redshift evolution of the cluster luminosities in the next section.

Figure 3 .
Figure3.The DES  −  color profile of the diffuse light in cluster subsamples of different redshift and richness ranges.The diffuse light color is consistent with the average color of the cluster's total stellar content in the center.In addition, a radial gradient can be seen in all of the redshift and richness ranges, such that the diffuse light becomes bluer at larger radii.This is consistent with previous studies that ICL consists of more metal-poor and younger stars(Edwards et al. 2020) and suggests an ICL origin from galaxy disruption and tidal stripping(Contini et al. 2019).

Figure 4 .
Figure 4. We fit the color profiles in Figure 3 to a linear model with radius Color() = × () + to further quantify those profiles.The radial slope parameter  (upper panel) is negative for all of the richness and redshift bins, which is a robust detection of the radial gradient of the diffuse light's color profile.The  parameter (lower panel), which is the diffuse light's color at  = 1 kpc, appears to be slightly larger for richer clusters, possibly reflecting a more massive and redder satellite population in richer clusters.

Figure 5 .
Figure5.Integrated luminosity (K-corrected to be the apparent magnitude in the observer frame of  = 0.25) as a function of radius (upper panels), or radius scaled by  200 (lower panels).This figure illustrates the luminosity measurements used in our analyses and shows that diffuse luminosity (red lines) and total cluster stellar luminosity (blue lines) increase as the radial range increases.See Section 4 for discussions on the trends manifested in this figure.Each panel shows a different redshift range, while the last redshift panel (redshift 0.5 to 0.6) is coded in a gray color, indicating that this redshift slice may be less reliable because of potential incomplete masking.We note that the last redshift range is not used in quantitative analyses.

Figure 6 .
Figure6.Luminosities (distance-corrected to be the apparent magnitude in the observer frame of  = 0.25) enclosed within 4 radial bin (0 to 30 kpc, 30 to 80 kpc, 80 to 300 kpc and 300 to 600 kpc).We examine how the luminosities change with redshift and richness.The luminosities of the total cluster stellar content (blue lines) and the cluster diffuse light (CG+ICL, red lines), increase significantly as the cluster's richness increases.However, those luminosities do not appear to change significantly with redshift, except in the lowest richness range and within 30 kpc.See Section 4 for quantitative analyses.

Figure 7 .
Figure 7. Luminosity (distance-corrected to be the apparent magnitude in the observer frame of  = 0.25) as a function of redshift in a volume-limited cluster sample.Again we analyze the brightness enclosed within 4 radial bins (0 to 30 kpc, 30 to 80 kpc, 80 to 300 kpc, and 300 to 600 kpc, top to bottom panels) and examine how they change with redshift.The luminosities of the total cluster stellar content (blue lines) and the cluster diffuse light (CG+ICL, red lines) both show some signs of becoming brighter over time within 30 kpc, and between 30 to 80 kpc.Those trends indicate growth in the CG and in the CG to ICL transition region.In section 4, we include quantitative analyses of those trends.

Figure 8 .
Figure 8. Diffuse fractions in the total cluster stellar content, calculated within physical radii (right column) and scaled radii by  200 (left column).We do not observe consistent redshift or richness-dependent trends (except in the 0.021 to 0.056  200 bin) in the measurements.However, the diffuse fractions appear to be dropping at larger radii.See Section 4.3 for more detailed discussions.

Figure 9 .
Figure 9. Testing the effect of PSF on SB and luminosity measurements.Upper panel: diffuse light SB models are convolved with a PSF model at different redshifts.The PSF flattens out the SB distribution in the center.Middle Panel: Relative change in SB after the diffuse light profile models are convolved with a PSF model.The changes are most significant within 10 kpc.Lower panel: Relative changes in luminosity (derived by integrating the SB profiles radially as described in Section 4) after a PSF model convolution.The integrated luminosity is most affected within 20 kpc.

Figure 11 .
Figure 11.Upper Panel: Diffuse profiles derived from DES special coadded images (red lines, Y6ACoadd, fiducial results in this paper) vs. those derived from single epoch images as in Zhang et al. 2019 (Y6SE, blue lines) .Lower Panel:Differences in these profiles.These two approaches yield consistent surface brightness measurements to an accuracy level of over 30 mag/kpc 2 in terms of raw diffuse light and random profile measurements.After random profile subtractions, the differences vanish at a surface brightness level of 40.5 mag/kpc 2 .

Figure 13 .
Figure 13.The stellar mass of diffuse light and of the cluster's total stellar content in the IllustrisTNG 300-1 simulation, as a function of the host halo's mass.The different lines indicate the running means in different redshift snapshots.These stellar mass-halo mass relations do not seem to vary with redshift.However, those relations depend on the apertures used to calculate the stellar masses, and are steeper in large radial ranges.

the
Kavli Institute of Cosmological Physics at the University of Chicago, the Center for Cosmology and Astro-Particle Physics at the Ohio State University, the Mitchell Institute for Fundamental Physics and Astronomy at Texas A&M University, Financiadora de Estudos e Projetos, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvimento Científico e Tecnológico and the Ministério da Ciência, Tecnologia e Inovação, the Deutsche Forschungsgemeinschaft and the Collaborating Institutions in the Dark Energy Survey.The Collaborating Institutions are Argonne National Laboratory, the University of California at Santa Cruz, the University of Cambridge, Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas-Madrid, the University of Chicago, University College London, the DES-Brazil Consortium, the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH) Zürich, Fermi National Accelerator Laboratory, the University of Illinois at Urbana-Champaign, the Institut de Ciències de l'Espai (IEEC/CSIC), the Institut de Física d'Altes Energies, Lawrence Berkeley National Laboratory, the Ludwig-Maximilians Universität München and the associated Excellence Cluster Universe, the University of Michigan, NSF's NOIRLab, the University of Nottingham, The Ohio State University, the University of Pennsylvania, the University of Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the University of Sussex, Texas A&M University, and the OzDES Membership Consortium.Based in part on observations at Cerro Tololo Inter-American Observatory at NSF's NOIRLab (NOIRLab Prop.ID 2012B-0001; PI: J. Frieman), which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation.The DES data management system is supported by the National Science Foundation under Grant Numbers AST-1138766 and AST-1536171.The DES participants from Spanish institutions are partially supported by MICINN under grants ESP2017-89838, PGC2018-094773, PGC2018-102021, SEV-2016-0588, SEV-2016-0597, and MDM-2015-0509, some of which include ERDF funds from the European Union.IFAE is partially funded by the CERCA program of the Generalitat de Catalunya.Research leading to these results has received funding from the European Research Council under the European Union's Seventh Framework Program (FP7/2007-2013) including ERC grant agreements 240672, 291329, and 306478.We acknowledge support from the Brazilian Instituto Nacional de Ciência e Tecnologia (INCT) do e-Universo (CNPq grant 465376/2014-2).This manuscript has been authored by Fermi Research Alliance, LLC under Contract No. DE-AC02-07CH11359 with the U.S. Department of Energy, Office of Science, Office of High Energy Physics.

Table 1 .
The Galaxy Cluster Sample in this Analysis

Table 2 .
constraints on parameters in  ( 0 <  <  1 ) =  × log 10  , which is used to quantify the relation between stellar luminosity and cluster richness and redshift.