Simba - C : the evolution of the thermal and chemical properties in the intragroup medium

The newly updated GIZMO and Simba based simulation, Simba-C , with its new stellar feedback, chemical enrichment

From X-ray observations, the following three group properties have attracted the most attention: (i) Entropy of the hot diffuse gas (S = k B T e /n 2/3 e ) within R 500 -This quantity contains the time-integrated history of the heating and cooling that the gas has experienced, including non-gravitational heating induced by stellar-powered galactic outflows and/or AGNs (Lewis et al. 2000;Sijacki et al. 2008;Le Brun et al. 2014).Entropy offers a more precise time-integrated representation of the energy flow within these groups compared to temperature or density measurements alone (Balogh et al. 1999;Babul et al. 2002).
(ii) Hot gas fractions in central regions -The hot gas fraction in the central region of groups, i.e. within R 2500 , is considerably lower than that found in massive clusters (Balogh et al. 1999;Vikhlinin et al. 2006;Gastaldello et al. 2007;Sun et al. 2009).This discrepancy could arise from various processes, such as the nongravitational heating of the gas described in (i) driving its expulsion (e.g.Balogh et al. 1999;Davé et al. 2008;Liang et al. 2016), or the depletion of gas due to efficient cooling (e.g.Voit & Bryan 2001a).Solely the latter would lead to significantly higher stellar mass fractions than those observed (Davé et al. 2002), but when combined with feedback processes, it could become an important contributing factor.
(iii) Metal content of the hot diffuse gas -Heavier elements originate from stellar nucleosynthesis, primarily through core-collapse supernova explosions of massive stars (SNe II) and thermonuclear detonations of accreting white dwarf stars (SNe Ia), and are dispersed into the broader environment via galactic outflows (Davé et al. 2008) and/or ram pressure stripping (Domainko et al. 2006).It is believed that SNe II seeded the Universe with metals during the early stages of galaxy formation, whereas SNe Ia are associated with the late stages of stellar evolution of smaller stars, dominating metal production over longer periods (McNamara & Nulsen 2007).Consequently, X-ray observations of the relative abundances of the byproducts of SNe II vs. SNIa enable the constraint of galaxy groups' star formation histories (McNamara & Nulsen 2007).Furthermore, larger groups and clusters are observed to have an iron abundance of approximately [Fe/H] ∼ 0.3 (Edge & Stewart 1991;Peterson et al. 2003;De Grandi et al. 2004;de Plaa et al. 2007), indicating that a significant fraction of the metals produced escape the interstellar medium (ISM), potentially with a considerable fraction of the metals having been ejected prior to the formation of the groups themselves (Oppenheimer et al. 2012).Recent studies (e.g.Appleby et al. 2021;Saeedzadeh et al. 2023) consider strong outflows/winds as the most plausible candidate for enrichment, although ram pressure stripping also contributes (Saeedzadeh et al. 2023).Once the metals leave the ISM, they are further redistributed by outflow-driven turbulence (Rennehan et al. 2019;Bennett & Sijacki 2020;Lochhaas et al. 2020;Rennehan 2021;Li et al. 2023).
From these three properties, it becomes evident that feedbackdriven outflows play a pivotal role in the formation and evolution of galaxies within group environments.These winds are more ubiquitous at high redshifts but are now also locally observable (Martin 2005;Sturm et al. 2011;Veilleux et al. 2013;Turner et al. 2014).It is suggested that outflows arise from either stellar processes (e.g.SNe explosions) or AGNs and possess the capacity to alter the physical properties of their environments, since the feedback energy is comparable to the binding energies of groups (Lovisari et al. 2021).Specifically, AGN-powered winds originate as high-velocity outflows on parsec scales at the centre of galaxies and have a profound impact on the gas in the central (∼ 1 kpc) region of their host galaxies (Sturm et al. 2011;Veilleux et al. 2013;Villar Martín et al. 2014).Consequently, it is interesting to ask questions about their impact on the large-scale group environment (Yang et al. 2024).
Stellar-powered outflows, driven by SNe energy and momentum injection, along with photoheating and radiation pressure from massive stars (Murray et al. 2005(Murray et al. , 2010;;Krumholz & Thompson 2013), exhibit remarkable efficacy in the removal of metal-enriched ISM gas (Taylor & Kobayashi 2015).This effectiveness stems from the launch sites of the winds being embedded in star-forming regions throughout the ISM, with metal-enriched winds observed to reach velocities exceeding the escape velocity with mass outflow rates comparable to the galaxy's star formation rate (SFR).Simulations such as the Feedback In Realistic Environments (FIRE) simulations (Hopkins et al. 2014) and the Numerical Investigation of a Hundred Astrophysical Objects (NIHAO) simulations (Wang et al. 2015) have demonstrated that stellar feedback is capable of launching powerful galaxy-wide winds with substantial mass loading, transporting metals into the intragroup medium (IGrM) and beyond.The diverse underlying galaxy processes combined with the galactic outflows generate dynamic changes observable in group halo properties, particularly in the hot X-ray gas.Such observations provide an opportunity to constrain the as-yet poorly understood stellar and AGN feedback mechanisms.
In this paper, we study intragroup gas properties using our latest Simba-C galaxy evolution simulation (Hough et al. 2023), which is based on the Simba galaxy formation model (Davé et al. 2019).This updated version incorporates enhancements such as an improved chemical enrichment and stellar feedback model, a reintegrated dust model, and a modified AGN feedback model.Our investigation aims to understand the influence of metal-enriched outflows on various characteristics of hot diffuse gas X-ray properties, including entropy, temperature, luminosity, mass, and the hot gas mass fractions within central regions.Using the newly updated chemical enrichment model (Hough et al. 2023), we also examine the metal content present in this hot diffuse gas, contrasting it with the approximation of instantaneous recycling of the metals model (Davé et al. 2019).
Sec. 2, discusses the input physics of the Simba-C simulation, highlighting its significantly modified chemical enrichment and feedback modules.Additionally, this section also provides a brief description of how we compile our catalogue of simulated galaxy groups and outlines the distinctions between the various versions of the Simba/Simba-C simulations.In Sec. 3, our focus lies on discussing the global X-ray properties of our simulated galaxy groups at z = 0, emphasizing the three commonly discussed group X-ray scaling relations: (i) Luminosity-temperature, (ii) mass-temperature, and (iii) entropy-temperature.In Sec.3.4 and 3.5, we discuss the gas mass fractions and the metal content of the IGrM, respectively.In Sec. 4, we discuss the evolution (z = 2 to 0) of the 1 keV temperature groups for each of the various scaling relations and physical properties discussed in Sec. 3. Finally, we summarize our findings and conclusions in Sec. 5.

Simba-C
In this paper, we use Simba-C as our main galaxy evolution model.This simulation (Hough et al. 2023) is a forked version of the original Simba simulation (Davé et al. 2019), a comprehensive large-volume cosmological simulation utilizing the hydrodynam-ics+gravity solver GIZMO (Springel 2005;Hopkins 2015Hopkins , 2017)).This section provides an overview of the Simba-C simulation, with further detailed information available in Davé et al. (2019) and Hough et al. (2023) for interested readers.GIZMO evolves the hydrodynamic equations using the mesh-free finite-mass (MFM) method (Lanson & Vila 2008a,b;Gaburov & Nitadori 2011).GIZMO also handles shocks using a Riemann solver without artificial viscosity, while preserving the mass within each fluid element at simulation time, facilitating the simple tracking of gas flows (Hopkins 2015;Davé et al. 2019;Alonso Asensio et al. 2023).
Both the Simba and Simba-C simulations utilize the GRACKLE-3.1 library (Smith et al. 2017), to handle the radiative cooling and photoionization heating of the gas.This also includes metal cooling and the non-equilibrium evolution of the primordial elements.Within this framework, adiabatic and radiative terms evolve simultaneously during a cooling sub-time-step, ensuring a stable thermal evolution.The model also accounts for self-shielding self-consistently based on a local attenuation approximation (Smith et al. 2017;Davé et al. 2019).Star formation is modelled using an H 2 -based SFR, where the fraction of H 2 is based on the sub-grid model that considers the metallicity and local column density of H 2 , based on the approach outlined in Krumholz & Gnedin (2011).
Simba-C differentiates itself from Simba by adopting the Chem5 cosmic chemical enrichment model, developed and refined in various studies, including Kobayashi (2004), Kobayashi et al. (2007), Taylor & Kobayashi (2014), Kobayashi & Nakasato (2011), Kobayashi et al. (2020a), andKobayashi et al. (2020b).The Chem5 model is the 'version-5' of a self-consistent1 chemodynamical enrichment model, tracking all elements from Hydrogen (H) to Germanium (Ge).It accounts for a range of physical stellar feedback processes: Stellar winds, Asymptotic Giant Branch (AGB) stars, super AGB stars, Type Ia SNe, Type II SNe -including Hypernovae (HNe), and 'failed' SNe for the most massive stars.Importantly, Chem5 does not use the instantaneous metal recycling approximation or a simplified delayed feedback model for Type Ia SNe and AGB stars.Instead, it treats each star particle as an evolving stellar population that ejects energy, mass, and metals into its nearby environment.This fundamental difference between the original Simba simulation and the Simba-C simulation leads to a more time-resolved modelling of stellar enrichment, relying on updated yields and stellar evolution models.This approach offers a more detailed and comprehensive representation of stellar enrichment dynamics.We refer interested readers to Hough et al. (2023) for a more in-depth discussion of the implementation of the Chem5 model into the Simba-C simulation, as well as the comparative tests to the Simba simulation.
A key benefit of Simba-C over its predecessor is its improved ability to predict observed abundance ratios.
Simba-C uses Simba's model for star formation-driven galactic winds, which employs decoupled two-phase winds with a mass load factor.It also uses the on-the-fly approximated friend-of-friends (FOF) finder, as described in Davé et al. (2016), to compute various galaxy properties to which the mass loading and wind velocity are scaled.While the scaling of mass loading remains unchanged from Simba, there is a reduced normalization modification in the velocity scaling, with Simba-C utilizing 0.85 instead of the previous value of 1.6.For a detailed understanding of the rationale behind this recalibration, refer to Hough et al. (2023).
Simba-C mostly follows Simba's black hole physics with some updates.Black holes are seeded when the galaxy is initially resolved (M * ≳ 6 × 10 8 M ⊙ ).However, the black hole accretion is suppressed exponentially by a factor of 1 − exp(−M BH /10 6 M ⊙ ), aimed at mimicking the effect of star formation suppressing BH growth in dwarf galaxies as seen in simulations as described in Anglés-Alcázar et al. (2017a) and Hopkins et al. (2022).At simulation time, both the BHs' dynamical mass (inherited from the parent star particle) and their physical black hole mass (set to M BH,seed = 10 4 M ⊙ h −1 and allowed to grow via accretion) are tracked.
Simba-C continues to use Simba's two-mode black hole accretion model.For the cool gas (T < 10 5 K), we compute the accretion rate based on the torque-limited accretion model (limited to three times the Eddington accretion rate) presented in Hopkins & Quataert (2011) and Anglés-Alcázar et al. (2017b) ṀTorque (1) Hot gas, meanwhile, is accreted following the standard Bondi-Hoyle-Lyttleton accretion rate (limited to the Eddington rate) presented in Bondi (1952 (2) As in Simba, the total large-scale accretion rate onto each black hole is given as the sum of the two modes, taking into account the conversion of matter into radiation where η = 0.1 is the radiative efficiency (Yu & Tremaine 2002).
Simba-C uses the accretion energy to drive the black hole feedback that quenches galaxies.This is achieved through a kinetic subgrid model for black hole feedback combined with X-ray energy feedback, as described in Davé et al. (2019).Furthermore, Simba-C's AGN feedback model mimics the energy injection into the large-scale surrounding gas by using purely bipolar feedback in the angular momentum direction of the black hole accretion radius.Regarding the mass scale at which the jets are permitted, Simba-C sets a range between 7×10 7 M ⊙ −1×10 8 M ⊙ .Additionally, each black hole particle is effectively assigned its own jet onset mass which it retains throughout the simulation run.Another distinction from Simba in Simba-C is that the black hole maximum jet velocity boost is allowed to reach a value that scales with the halo escape velocity, rather than a constant value of 7000 km/s (see Hough et al. (2023) for more details).

Dust integration and results
The sole difference between the Simba-C simulation results presented in Hough et al. (2023) and the simulation utilized in this article is the inclusion of Simba's dust model.The model was intentionally omitted in Hough et al. (2023), due to its possible influence on the metal content in a simulation dedicated to testing the new metal enrichment model.
For this paper, we re-introduce Simba's dust model to strengthen the precision of our representations of a physical system.As described in Davé et al. (2019), the dust is passively advected following the gas particles.Dust is produced by the condensation of metals from the ejecta of SNe and AGB stars.Simba uses fixed dust condensation efficiencies of δ AGB i,dust = 0.2 and δ SNII i,dust = 0.15, based on theoretical models from Ferrarotti & Gail (2006) and Rémy-Ruyer et al. (2014) to match the low metallicity end of the observed dustto-gas mass ratios.It should be noted that Type Ia SNe condensation is omitted due to its low impact on dust production (Nozawa et al. 2011;Dwek 2016;Gioannini et al. 2017).Dust can grow, be destroyed, and undergo thermal sputtering (Li & Mattsson 2020).However, it is essential to note that the dust model is still applied only to the original 11 elements that were tracked in Simba.The additional elements introduced in Simba-C are assumed to not participate in the formation of dust, as they constitute a very small fraction of the total metal mass.
In Fig. 1, we show the comparison between the Simba-C and Simba simulations' dust mass functions at z = 0.For comparison with observational results, we show results from Dunne et al. (2011), who used data from Herschel Astrophysical Terahertz Large Area Survey (Herschel-ATLAS), and observational results from Clemens et al. (2013), who used data from the Wide field Infrared Survey Explorer (WISE), the Spitzer space telescope, Infrared Astronomical Satellite (IRAS), and Herschel.The general trend shows that the dust mass function in the Simba-C simulation obtained lower values compared to the Simba simulation.However, it remains in line with the trend observed in Dunne et al. (2011).From subsequent comparisons conducted, we can also confirm that there were no other discernible changes resulting from the inclusion of dust in the simulations.This is likely due to the dust sputtering back into metals in the hot diffused gas, therefore minimizing its effect in the group regime.However, this could be different for cold gas, which is not explored here.Due to this, the Simba-C NoDust simulation will not be shown in the subsequent comparison tests.

Runs and analysis
We use four different Simba/Simba-C simulations, each incorporating a different physics module, until we reach the complete Simba-C model as described in Sec.2.1, with dust included.The complete Simba-C and Simba simulations (iii and iv) have volumes of side length 100 Mpc h −1 , consisting of 10243 gas particles and 1024 3 dark-matter particles, while the other two simulations (i and ii) have volumes of side-length 50 Mpc h −1 with 512 3 gas particles and 512 3 dark-matter particles and are used only for reference.All simulations ran from an initial redshift of z = 249, down to z = 0 and follow a Planck  We analyze the simulation outputs using a friends-of-friends galaxy finder, assuming a spatial linking length of 0.0056× the mean interparticle spacing, as well as the AMIGA Halo Finder (AHF), a tool developed by Knebe et al. (2008) and Knollmann & Knebe (2009) specifically designed for halo identification.Postprocessing involves cross-matching galaxies and haloes using two distinct Python packages: Caesar and XIGrM, each serving specific purposes in computing group properties.
(i) Caesar:2 This yt-based Python package performs galaxy finding which is applied to all stars, black holes, and cool gas elements above a specified minimum SF threshold density of n H > 0.13 H atoms cm −3 .Black holes are associated with galaxies to which they are most gravitationally bound, and the most massive black hole within a galaxy is designated as the central black hole.(ii) XIGrM: 3 This Python package specializes in computing the X-ray properties of the IntraGroup Medium, to compute the group properties as detailed below.
Previous studies exploring the effects of AGN feedback has been done by comparing the Simba NoAGN simulation to the Simba simulation (see Robson & Davé (2020, 2021) and Kar Chowdhury et al. (2022), and references therein).Despite these previous investigations, we include the Simba NoAGN simulation for reference, since we recalibrated the strength of the AGN feedback in Hough et al. (2023).Discussions based on the impacted properties are given in Sec.3.3.

Finding the galaxy group haloes
To identify galaxies and galaxy groups within each simulation, analysis is performed for outputs corresponding to redshifts z = 2, 1, 0.5 and z = 0. We follow the method used in Liang et al. (2016), Jung et al. (2022), andSaeedzadeh et al. (2023).Similar studies using different cosmological simulations, such as the Feedback Acting on Baryons in Large-scale Environments (FABLE) simulation (Henden et al. 2018) have also been done.We utilized the AHF software to find hierarchy structures nested for haloes and sub-haloes, by locating peaks in the adaptively smoothed density field, through the identification of all particles (gas, stars, dark matter, and black holes) that are gravitationally bound to each other (Jung et al. 2022;Saeedzadeh et al. 2023).We then proceed up in the hierarchy to find the larger structures, while the centres of these haloes are located by applying the shrinking-sphere approach (Power et al. 2003).The determination of halo masses (M ∆ ) involves the construction of a sphere and the expansion of its radius until the total mass density (enclosed) equals the viral density for a specific cosmology at a certain redshift.
with E(z) ≡ H(z)/H 0 being the dimensionless Hubble parameter given by: and ρ crit being the critical cosmology density. 4We also use the frequently used ∆ values, that is, 200, 500, and 2500 in addition to the virial radius/mass5 (Babul et al. 2002;Lovisari et al. 2021;Saeedzadeh et al. 2023).The following conversion equations are also used: We utilize the Python package XIGrM to determine the various halo quantities, including the radius, mass, luminosity function, entropy of the system, and the temperature of the host halo at each ∆-value.

Group mass functions
In Fig. 2, we present the z = 0 halo mass function for the full Simba-C simulation (the halo mass functions of the other three simulations are very similar).The black curve represents all haloes, whereas the red, blue, and magenta curves depict haloes with at least three, two, and one 'luminous' galaxies, respectively.A 'luminous' galaxy is defined as having a stellar mass of at least M * ∼ 1.16 × 10 9 M ⊙ , equivalent to at least 64 star particles, with each particle having a mass resolution of approximately ∼ 1.8 × 10 7 M ⊙ .Per definition, groups and clusters are identified as haloes with three or more luminous galaxies (Liang et al. 2016).Consistent with Liang et al. (2016), we observed that on mass scales of ≥ 10 13 M ⊙ , nearly all haloes have at least three luminous galaxies.
In Fig. 3, the average galaxy stellar mass function (GSMF) for all luminous galaxies in their host groups are shown for the complete Simba-C (solid lines) and the original Simba simulations (dotted lines).The groups are categorized into mass bins corresponding to 12.5 < log M vir ≤ 13.0 M ⊙ (magenta), 13.0 < log M vir ≤ 13.2 M ⊙ (blue), and 13.4 < log M vir ≤ 14.0 M ⊙ (red).These mass bins align with those used in Davé et al. (2008) and Liang et al. (2016).In addition, a comparison is made with the observed GSMF for low-mass X-ray detected groups in the Cosmic Evolution Survey (COSMOS) of Giodini et al. (2012).It is important to note that, for comparison purposes, the galaxy's stellar mass function in all three mass ranges has been artificially reduced by a factor of three for both simulations.This adjustment is made because of the unknown volume of the galaxy groups within the catalogues relative to the observations, which is treated as a free parameter.This artificially reduced value has no physical meaning.Hence, the only predictive power for comparison with the data is in the shape of the group stellar mass function.
The shape of the binned galaxy stellar mass functions for both simulations is shown in Fig. 3, matches the overall observational trend fairly well, particularly for the highest mass bin that is closest to the observed sample.Notably, unlike in Liang et al. (2016), there is no excess of galaxies with very large stellar masses in Simba-C.Simba does have an excess, although much smaller than in Liang et al. (2016), which did not include AGN feedback.Liang et al. (2016) suggested that the introduction of AGN feedback would quench the 'supersized' central galaxies.This assertion is confirmed when using the Simba NoAGN simulation, which did indeed show an excess of 'supersized' central galaxies.Therefore, similarly to the overall galaxy mass function, the group GSMF also requires AGN feedback to be accurately reproduced in models.

Global X-ray galaxy group properties
Once we identify the haloes, the next step is to compute the various observed X-ray properties.To simplify our analysis, we focus on the hot diffused IGrM gas particles.These are defined as gas particles with temperatures above the threshold of T > 5 × 10 5 K with a hydrogen number density (n h < 0.13 H atoms cm −3 ), which is below the density threshold to allow star formation.
The first X-ray property computed is the rest-frame 0.5 − 2.0 keV luminosity within R 500 , denoted as L X,0.5−2.0 .To calculate this property, the luminosity of each individual IGrM gas particle within a distance r ≤ R 500 from the halo centre is summed.The emission characteristics of the gas particles are computed using the Python package PyAtomDB,6 which in turn uses the Atomic Data for Astrophysicist database (AtomDB), which itself consists of two components: (i) Astrophysical Plasma Emission Database (APED) and (ii) Astrophysical Plasma Emission Code (APEC; Smith et al. 2001).The gas is assumed to be optically thin and to be in collisional ionisation equilibrium.APEC uses the particle's mass, SPH-weighted Figure 3.The galaxy stellar mass function (GSMF) for all luminous galaxies in the simulated groups, sorted into three mass bins: 12.5 < log M vir ≤ 13.0 M ⊙ (magenta), 13.0 < log M vir ≤ 13.2 M ⊙ (blue), 13.4 < log M vir ≤ 14.0 M ⊙ (red).Simulated galaxies from Simba-C is shown with solid lines, while Simba is shown with dotted lines.We compare this to observations of the GSMF for low-mass X-ray detected groups that span the same mass range as our simulated groups (Giodini et al. 2012).We artificially reduced the galaxy's stellar mass function in the simulation by a factor of three, as described in the text.The vertically dashed black line shows our luminous galaxy's stellar mass resolution of M * ∼ 1.16 × 10 9 M ⊙ .density, temperature, and metallicity as input while returning an Xray spectrum.By summing the photon energy intensities in the specified energy range, the luminosity is obtained.Contributions to line and continuum emissions associated with each tracked element in Simba-C are computed separately, except for Gallium and Germanium due to the 30-element limit of AtomDB.
We use mass-weighted (mw) and emission-weighted (ew) abundances to study the metal content of the IGrM.They are defined as where q is the metal species, Z q,i is the SPH kernel-weighted abundance of the i th particle, m i and L i its mass and X-ray luminosity, respectively.The sum runs over all IGrM particles within the halo's volume.Furthermore, all metal abundance estimates are in terms of solar 'photosphere abundance' values from Anders & Grevesse (1989). 7n the calculation of X-ray properties, the spectroscopic temperature T spec is used rather than the frequently used T X temperature, where T X represents an average of the temperatures of the individual components weighted by their contributions to radiative emission.The choice of T spec over T X is motivated by the tendency of T X to be biased high by approximately 25% for clusters (Mazzotta et al. 2004).In contrast, T spec is designed as a weighting scheme to produce a temperature that is comparable to the temperatures of a hot gas determined by observations of groups and clusters (Vikhlinin et al. 2006).Specifically, T spec is determined by identifying a singletemperature thermal model whose spectrum best matches the observed spectrum.
In the X-ray analysis, the decoupled hot-wind particles, representing a fraction of the stellar feedback-driven wind particles that are heated, are excluded.The fraction of these particles in Simba-C follows the trends of the FIRE simulations (Pandya et al. 2023), while in Simba, it was a fixed fraction of 30%.Although these particles constitute a very small portion of the total mass (Appleby et al. 2021), due to their high density as they emerge from the ISM, they have extreme X-ray luminosities.Two different X-ray calculations are performed: (i) All the hot diffuse gas particles within R 500 are used, or (ii) only the hot diffuse gas particles within the radial range 0.15R 500 ≤ r ≤ R 500 are used.The latter is referred to as 'coreexcised' (e.g., Liang et al. 2016), which is often more robust, as the cores of clusters and groups can have widely varying temperatures.Hence, observations often present 'core-excised' temperatures to which we will compare where appropriate.
We first study direct comparisons to observations at present-day (z ≈ 0) for each of our four simulations, to obtain an understanding of which physics models/properties play the largest role in improving the simulation for each scaling relation in Sec. 3. The expected largest impact physics model is the addition of AGN as shown in Robson & Davé (2020); however, other effects may also have a significant contribution, depending on the tested property; therefore, all combinations must be tested.Then, using this knowledge, we look at how each simulation evolves with redshift in Sec. 4. This will give us the necessary insight to understand how the introduction of each physics module played its role in improving the evolution of the galaxy groups within simulations.Lastly, given the proven successes of the Simba simulation with certain X-ray properties, e.g. the mass fractions in Robson & Davé (2020), we also show the trends of the Simba simulation to confirm that the improvements made in Simba-C did not adversely affect these predictions and to identify any improvements that Simba-C has compared to its predecessor.
We follow the convention in the literature and plot the quantities motivated by the self-similar model for group and cluster haloes (Kaiser 1986).In this model, the scaling relations are preserved when using the quantities: (i) L X (z)E(z) −1 , (ii) M ∆ (z)E(z), and (iii) S ∆ (z)E(z) 4/3 , with E(z) ≡ H(z)/H 0 being the dimensionless Hubble parameter.

GLOBAL SCALING RELATIONS AND HALO STRUCTURE COMPARISONS TO PRESENT-DAY OBSERVATIONS
We begin by presenting comparisons to observations of the X-ray scaling relations and baryon mass fractions at z = 0 for each of our four simulations.This comparison will help assess the capability of the Simba-C simulation to capture the underlying physics that drives the X-ray properties of galaxy groups.In this section, the spectroscopic weighted temperature is used, as discussed in Section 2.4.3.

Luminosity-Temperature scaling relation
We first consider the scaling relation between the X-ray luminosity L X and the X-ray spectroscopic temperature T spec .Observations generally indicate a steeper scaling relation (L X ∝ T 3−5 ) for the lowest mass groups, whereas massive clusters generally align with the predicted slope of the self-similar model (Balogh et al. 1999;Rob-son & Davé 2021;Lovisari et al. 2021), that is, L X ∝ T 1−2 ,8 depending on the X-ray passband under consideration.Two primary physical effects impact this scaling relation.First, there is the radiation mechanism.In clusters (T X > ∼ 1 keV), the self-similar model predicts the bolometric luminosity scaling as L ∼ T 2 due to the bremsstrahlung being the dominant mechanism.In groups, line radiation begins to dominate emission.However, this predicts a flatter, not a steeper L − T relationship (Balogh et al. 1999).Second, feedback affects the IGrM, influencing both the amount of gas within the IGrM and its radial profile.The overall impact is more pronounced in the lowest-mass systems due to their smaller gas reservoirs and potential wells.This significantly influences L X due to its densitysquare dependence.For example, Balogh et al. (1999) model groups with isentropic cores show a L ∼ T 5 slope, as observed in low-mass groups.According to Balogh et al. (1999); Voit & Bryan (2001a); Babul et al. (2002); McCarthy et al. (2004), this indicates that heating and/or cooling have substantially altered the distribution of the hot X-ray gas.Hence, by probing scaling relations, one can obtain constraints on the group's feedback by accounting for the line emission mechanism.
Fig. 4 shows the rest frame 0.5 − 2.0 keV X-ray luminosity emitted within R 500 against the mean core-excised spectroscopic temperature for the simulated groups in our four simulations at z = 0: (i) The complete Simba-C simulation (blue solid line), (ii) the original Simba simulation (red dashed line), (iii) the Simba NoAGN simulation (green tight dot-dashed line), and finally (iv) the Simba NoFeedback (magenta loosely dot-dashed line).For our two main comparison simulation (Simba-C/Simba) results, we show a light blue and a light red band to display the 1σ-deviation in each temperature bin for each simulation, respectively.We include at least 10 haloes in each temperature bin (see Table 1).In cases where a temperature bin contains less than 10 haloes, the individual halo values are not binned and are represented with the following markers: (i) Simba-C (blue circles), (ii) Simba (red squares), (iii) Simba NoAGN (green stars), and (iv) Simba NoFeedback (magenta diamonds).This approach provides some insight into the emerging trends for the more massive galaxy groups and clusters, even where the trends cannot necessarily be confirmed as statistically significant because there are fewer than 10 haloes per bin.This approach will be consistently applied throughout Sec. 3.
For comparison, we have included low redshift X-ray observations from Pratt et al. (2009), using data from the Representative XMM-Newton Cluster Structure Survey (REXCESS), Eckmiller et al. (2011) using the Highest X-ray FLUx Galaxy Cluster Sample (HIFLUGCS), Lovisari et al. (2015) using XMM-Newton observations for a complete sample of galaxy groups, and finally O'Sullivan et al. ( 2017) using the Complete Local Volume Groups Sample (CLoGS) from XMM-Newton and Chandra observations.9In addition, we include a black dashed line representing the 'self-similar' model results in the 0.5 − 2.0 keV band, incorporating both line radiation and bremsstrahlung.This self-similar model adopts the slopes from Table 1 in Lovisari et al. (2021).
From Fig. 4, we note that the complete Simba-C simulation and its 1σ-error range generally overlap with the observations, although very slightly low for the coldest groups.Following the trends of individual halo values for the most massive groups and clusters, the Simba-C simulations seem to match the observations the most closely.This suggests that Simba-C broadly succeeds in determining L X through a combination of cooling and feedback.Compared to Simba-C, Simba has reduced luminosities in the lower halo mass regime, to the point where Simba obtained too low luminosities compared to observations, although some observations are still in its 1σ-error.This improvement from Simba-C likely stems from the overall lower metal mass fractions in Simba-C when adopting the new chemical evolution model, as demonstrated in Sec.3.5, which reduces metal-line cooling and consequently yields more hot gas.In the higher halo-mass regime, the opposite trend is observed, with individual haloes tending to be slightly too bright on average.However, when taking into account the 1σ-deviations, the improvements obtained by Simba-C over Simba in the L X,0.5−2.0 − T spec,corr relation are only significant in the colder groups, while in the warm groups these differences are insignificant.
Furthermore, Fig. 4 shows that the 'self-similar' model tends to overestimate the L X,0.5−2.0 − T spec,corr relation compared to observations.This well-known result is commonly attributed to feedback and/or cooling, which selectively removes hot gas more from smaller systems.Among our simulations, the NoFeedback run is the most similar to the self-similar model.Therefore, the impact of cooling can be assessed from this model and it alone is not sufficient to explain observations (see also McCarthy et al. 2004McCarthy et al. , 2008)).The NoAGN model has a lower L X value than the NoFeedback model, but still generally exceeds the observations and is still higher than the two main simulation results, suggesting that stellar feedback plays a role in gas removal (see also Liang et al. 2016), but that AGN feedback remains the primary L X reduction source Robson & Davé (2020).1.The number of haloes in each temperature bin is provided for each simulation across all T spec,corr plots.The T spec,corr value shown in the bin names denotes the lower T spec,corr value of the bin.All bins increase with a log T spec,corr = 0.1 bin size.

Mass-Temperature scaling relation
The relationship between mass and temperature is not anticipated to be very sensitive to baryonic processes, as for bound systems, the temperature should predominantly reflect the gravitational potential primarily driven by dark matter.However, modest departures from self-similarity can still occur due to the interaction between feedback heating and gas removal in lower-mass systems, impacting T spec .Fig. 5 shows the M − T spec,corr relation for the mass of the simulated groups within the central R 500 region for our four simulations at z = 0.For comparison, we included low-redshift X-ray observational results from Sun et al. (2009) using Chandra archival data, Eckmiller et al. (2011) using the HIFLUGC Survey, Kettula et al. (2013) using COSMOS results, and Lovisari et al. (2015) with their XMM-Newton observations.We also present the original O'Sullivan et al. ( 2017) CLoGS results based on the M − T spec,corr relation of Sun et al. (2009), as well as the revised CLoGS results utilizing SIMBA-C's M 500 − T spec,corr relationship for −0.6 ≤ log(T spec ) ≤ 0 to determine the mass of the CLoGS groups. 10Consequently, we cannot draw conclusions about our simulations that match the CLoGS results, but it provides insight into the observational trends of the lower-temperature groups.We use our updated results when examining the IGrM entropy in Fig. 6.
As expected, Fig. 5 reveals no strong differences between simulations for this particular scaling relation.The four simulations produced M 500 − T spec,corr curves, which are within 1σ of each other in all temperature ranges.The scaling is also consistent with that obtained by Liang et al. (2016) using a simulation without AGN feedback.They all demonstrate reduced masses for a given T spec on group scales, indicating departures from self-similarity driven primarily by the interplay of radiative cooling (present in all models) and the measurement of T spec .Furthermore, we can also conclude that T spec,corr can be considered as a proxy for the halo mass.This property allowed us to retroactively determine the CLoGS's M 500 values used here.However, temperature is consistently used in further plots due to its direct comparison with observations, as used in the CLoGS paper (O'Sullivan et al. 2017), and other galaxy group X-ray scaling relation property studies, e.g., Lovisari et al. (2021).
All simulations closely follow the observational trends, even though most of the overlapping observational results are from the CLoGS sample, which is to an extent by construction.This holds true even for the NoFeedback run, since this scaling relation is more determined by gravitational processes than by feedback interactions; i.e., the gas expands until its temperature is consistent with the gravitational potential of the group.Hence, while this scaling relation does not provide discriminatory power between models, it is reassuring that the simulations can reproduce the observations down to 10 Previously, the CLoGS masses were estimated using the scaling relations of Tier 1+2 groups from Sun et al. (2009).the smallest systems and support the viability of Simba-C for X-ray group studies.

Entropy-Temperature scaling
As highlighted in Balogh et al. (1999), Babul et al. (2002), Lewis et al. (2000), and Voit & Bryan (2001a), entropy serves as a valuable quantity to study how the IGrM is influenced by cooling and/or heating processes.This is because a significant portion of the Universe's baryons reside in intergalactic space and experience heating through gravitationally driven shocks (Davé et al. 2001).Once heated, they settle into gravitational potential wells and adopt the characteristic temperature of the enclosing dark matter.However, the mean intensity of the X-ray emissions from the baryons reflects the amount of non-gravitational energy.The emissivity of baryons depends on how severely they are compressed and how this injection affects the baryon distribution (Balogh et al. 1999;Voit & Bryan 2001a;McCarthy et al. 2004McCarthy et al. , 2008)).Stellar and AGN feedback can restrict this compression (Babul et al. 2002;McCarthy et al. 2004), thus reducing the X-ray luminosity.These processes are essential because gravitational-only processes would excessively produce the 0.5 − 2.0 keV X-ray background when establishing the entropy dis-tribution (Pen 1999; Wu et al. 2000).In other words, the lowest entropy (most compressible) gas needs to be eliminated (Voit & Bryan 2001b).Non-gravitational heating and radiative cooling/subsequent condensation are potential physical reasons for this change in the low entropy gas (McCarthy et al. 2004).Furthermore, a lower limit to the entropy of the intracluster increases the L X − T relation, because shallower potential wells of low temperature groups/clusters are less capable of overcoming resistance to compression (Balogh et al. 1999;Voit & Bryan 2001a;Babul et al. 2002).We plot the canonical proxy for entropy11 in Fig. 6, rather than the thermodynamic specific entropy, given by: where k B is the Boltzmann constant and n e (r) is the electron number density within a thin spherical shell at radius r.The gas distribution is organized so that the lowest entropy is in the centre of the halo, while the highest entropy is in the outer limits of the group (Mc-Carthy et al. 2004, 2008).For this scaling relation, we also show S 2500 since the most notable improvement can be seen in this regime.In Fig. 6, we show the gas entropy versus the spectroscopic temperature measured at R 500 and R 2500 (the inner core of the halo) for the simulated groups in our four simulations at z = 0.The entropy was calculated by taking the average on a radial shell between R x and 1.05 × R x .For comparison, we show low-redshift X-ray observational results from Sun et al. (2009), using the Chandra archival data, as well as the CLoGS results from O'Sullivan et al. ( 2017), using XMM-Newton and Chandra observations.Lastly, we show with the dashed black lines, the power-law indices of α = 1 and α = 0.74 corresponding to the best-fit R 500 and R 2500 values as determined by Sun et al. (2009) for the scaling of the group entropy.
From Fig. 6, it is evident that the different physical processes (feedback, dust, enrichment, etc.) impact the entropy at both radii as expected.For Simba/Simba-C, which include both feedback mechanisms that severely inhibit gas compression, a higher entropy is measured.The increased entropy suppresses gas cooling and condensation further until there is no more gas below the cooling threshold to form stars.If feedback is inefficient, condensation will remove this gas from the intracluster medium.If the feedback is highly efficient, it will increase entropy and convect to the outer regions of the cluster, resulting in a decrease in the gas density (Voit & Bryan 2001a).In return, this reduces the luminosity and steepens the L X −T relation (cf.Balogh et al. 1999;McCarthy et al. 2004) also seen in Fig. 4.
When comparing these two simulations with the observations, we note that although the 1σ-error between Simba-C and Simba overlaps with each other at both radii (especially in the inner-core region), the average gas entropy in the haloes tends to be lower in Simba-C than in Simba.This results in the entire 1σ-error range of Simba-C also being lowered.This positively impacts the Simba-C simulation which now closely matches, on average, with the S 2500 gas entropy CLoGS observations (O'Sullivan et al. 2017).The impact on the S 500 results are less pronounced, with the CLoGS observations tending to be, on average, between the Simba and Simba-C results, but still within the 1σ-error.Furthermore, by examining the individual halo entropies, the general trend for Simba-C simulation most closely follows the higher mass groups and the slope of the self-similar model.Although true at both radii, it is more notable in the inner-core region, with the Simba simulation's S 2500 tending to be a bit more flat.From this we can conclude that Simba-C improves the gas entropy in the group inner-cores with respect to Simba, even though their overlapping 1σ-error indicates that this improvement is insignificant.This is an important finding, since Simba's entropy profiles have been shown to flatten in an extended entropy core (Oppenheimer et al. 2021;Altamura et al. 2023).Therefore, with Simba-C obtaining, on average, a lower entropy compared to Simba at both radii, it goes in the direction of relieving some of this tension and requires further investigation of the entropy profiles to provide context for this improvement.This is being done in a follow-up study by Padawer-Blatt et al. (in prep).
Furthermore, Fig. 6 also shows the impact of the inefficient metal cooling that Simba-C obtained.In Hough et al. (2023) §2.6, it was found that with the introduction of the new stellar feedback and metals model (Chem5), the Simba-C simulation was underproducing metals.This led to a shallow mass-metallicity relation (MZR).The Chem5 model inherently produced fewer metals, resulting in less efficient metal cooling.To address this problem, the strength of the AGN feedback (a heating process) in Simba-C was reduced.12This reduction is now reflected in the slightly reduced gas entropy obtained from Simba-C.This decrease in entropy allowed the Simba-C simulation to partially resolve some of the overcorrection/overestimation of the gas entropy observed after the AGN feedback was included into the Simba simulation, i.e., going from Simba NoAGN to Simba, which increased the entropy by approximately 0.3 dex).
Simulations lacking feedback, especially AGN feedback, dominated by gravitational-only shocks, exhibit lower entropy values at R 2500 than the other simulations.This, in turn, leads to an overproduction of the X-ray luminosity (see Fig. 4 and McCarthy et al. 2004).It should also be noted that the NoFeedback model (cooling only) is closely aligned with the self-similar model; however, this alignment does not have physical significance.When gas cooling is allowed, it initially becomes more centrally concentrated, and the system's luminosity increases (McCarthy et al. 2004).However, as the gas continues to cool, stars start to form, reducing the hot gas fraction.Consequently, this reduces the luminosity of the system at a given temperature.The Simba NoAGN result (already studied in Robson & Davé 2020) is consistent with the results of Liang et al. (2016).

Total baryon, stellar and gas fractions within the galaxy groups
As expected, AGN feedback has the greatest impact on the scaling relations.It lowers the gas content of the groups, thereby reducing the luminosity and increasing the entropy, with only a slight impact on the temperature.However, the addition of stellar feedback (NoFeedback to NoAGN) had a similar impact on the scaling relations, albeit less intensive.This process should also be reflected in the baryonic mass fractions in groups.Also, as shown in Fig. 6, Simba-C, with its lower amounts of metals, requiring a 'weaker' AGN feedback model, resulted in a slight but valuable improvement to the entropy, by lowering the observed gas entropy within the inner-core region of the group.In this section, we examine whether the simulated group's mass components are realistic by exploring various baryonic mass components as functions of the halo mass and compare among our different simulations and observations.This process has been extensively studied in previous literature with a specific focus on the impact of AGN on the original Simba simulation by comparing it with the Simba NoAGN simulation (see Robson & Davé 2020), where they found that the addition of AGN feedback not only greatly impacted the mass fractions, but also resulted in the Simba simulation matching the observations quite well.Therefore, Simba-C is not expected to result in greatly differing mass fractions, but is expected to show, at most, a mild departure from the Simba simulation due to the alteration of the AGN feedback strength.Our focus will therefore be largely on the validation of the Simba-C's mass fractions, and show the three other simulations only for comparison.Fig. 7 shows the mass fractions for four different properties: (i) The baryonic mass fraction (top left panel), (ii) the hot T > 5×10 5 K diffuse IGrM gas mass fraction (top right panel), (iii) the stellar mass fraction (bottom left panel), and finally (iv) the cold IGrM gas fraction (bottom right panel) -within R 500 for the simulated groups in our four simulations at z = 0 versus their halo mass.For comparison, we include low redshift X-ray observational results from Eckmiller et al. (2011) using the HIFLUGC Survey, Laganá et al. (2013) using XMM-Newton, Chandra and the Sloan Digital Sky Survey (SDSS) observations, Gonzalez et al. (2013); Lovisari et al. (2015) both using XMM-Newton data, and Loubser et al. (2018) using Brightest Group Galaxies (BGGs) from the CLoGS sample (O'Sullivan et al. 2012;Kolokythas et al. 2022).It must be noted that the BGG results only provide a lower limit on the CLoGS mass scales, although it is expected that the BGG's stellar mass would dominate the CLoGS group's stellar mass.Furthermore, the CLoGS sample is chosen to have the early-type galaxies as the BGGs, which is not necessarily the case for the simulations.Therefore, the CLoGS BGGs are not a direct comparison but provide insight into the lower-mass groups.The cosmological baryon fraction assumed in these simulations of Ω b /Ω m = 0.16 is indicated by the black dashed line.
From Fig. 7, we see distinct groupings of models.The two simulations that include the AGN feedback produce low baryonic mass fractions in all components (also seen in figure 6 by Robson & Davé 2020), while those lacking AGN feedback produce high baryonic contents.This split persists through every mass fraction for all mass ranges.Simba-C shows only a small increase in the different mass fractions compared to Simba, validating our hypothesis that Simba-C should be at most a slight departure from Simba.Examining the small differences between Simba-C and Simba in more detail, we note that Simba simulation predicts slightly lower baryonic and hot diffuse gas mass fractions than Simba-C at lower masses, but converges for halo masses log M 500 > 13.5 M ⊙ .This is responsible for driving the higher entropies seen in Simba (Oppenheimer et al. 2021).
In summary, the primary physics module necessary to obtain realistic mass fractions is AGN feedback for Simba, as determined by Robson & Davé (2020).This was also found by Henden et al. (2018) using the FABLE simulation when AGN feedback is included, indicating that baryon fractions are a strong constraint in this process (see e.g.Oppenheimer et al. 2021).Furthermore, as shown in Cui et al. (2022), a consistent jet velocity implemented in Simba for the AGN feedback is more efficient at reducing the gas fractions in galaxy groups than clusters.Therefore, agreement of Simba-C with observations in the group regime is non-trivial and has been difficult to obtain in other models (e.g.Barnes et al. 2017;McCarthy et al. 2018).Simba-C has an increase, albeit slight, in the mass fractions compared to Simba, owing to the lowering of the AGN feedback strength.This motivates us to examine the hot-gas metallicities in more detail.

Metal enrichment of the IGrM
Given the vital role played by the chemical enrichment model in star formation through metal cooling and establishing the abundance ratios of various elements in these groups, it is interesting to delve into the chemical enrichment of these groups.
The metal content in the intergalactic medium originates from the transport out of the ISM primarily through large-scale galactic outflows (e.g.Aguirre et al. 2001;Davé et al. 2008;Oppenheimer et al. 2012;Veilleux et al. 2013).These outflows simultaneously establish the mass-metallicity relation in galaxies (Finlator & Davé 2008;Davé et al. 2011;Hirschmann et al. 2013;Somerville & Davé 2015;Liang et al. 2016).However, in the group and cluster environment, it remains less clear whether these enriching outflows are driven by stellar or AGN feedback and to what extent gas stripping processes contribute.In this context, we focus on the observed abundances and abundance ratios of silicon, iron, and oxygen.This allows us to constrain the underlying driver of IGrM enrichment and evaluate how well our models match with observations.In Fig. 8, we show the global Fe abundance (top row) and the global Si abundance (bottom row) for the simulated groups within R 500 in our four simulations at z = 0 versus the halo mass.We make a distinction between the mass-weighted (left column) and emissionweighted (right column) abundances.In addition, these abundance calculations only involve hot diffuse gas.For comparison, massweighted X-ray [Fe/H] abundances from Yates et al. (2021) using the Mapping Nearby Galaxies at APO (MaNGA) and the Multi-Unit Spectroscopic Explorer (MUSE) observations are included.Emission-weighted X-ray [Si/H] and [Fe/H] abundances from Hels- From Fig. 8, it is evident that the emission-weighted abundance ratios are, on average, higher than the mass-weighted X-ray abundances.This is expected because much of the gas mass resides on the outskirts of groups where the metallicities are lower, while a significant portion of the emission originates from the central region where the metallicities are higher.This underscores the importance of creating realistic mock X-ray images, for example, using methods such as Mock Observations of X-ray Halos and Analysis (MOXHA) (Jennings & Davé 2023), to ensure accurate comparisons with metallicity measurements.While we defer a detailed investigation of this aspect to future work, for now, we consider the emission-weighted measures as a reasonable proxy for what would be observed.
Mass-weighted abundances establish a more direct link to the underlying physical processes of group-wide enrichment.In this context, it is evident that Simba-C exhibits lower abundances than Simba in both iron and silicon, by approximately ∼ 0.2 dex.Recall that the overall metal production rate is lower in Simba-C, but was re-tuned to match the galaxy mass-metallicity relation (Hough et al. 2023).With a higher fraction of metals retained within galaxies, this worsens the differences in IGrM enrichment.
The Simba NoAGN model exhibits substantially lower metallicity than Simba.This is despite the fact that the stellar mass and hence metal production in NoAGN are substantially higher.Therefore, AGN feedback plays an major role in the ejection of metals into the IGrM.This is due to Simba's AGN kinetic feedback being hydrodynamically decoupled for some time upon ejection; hence, they explicitly cannot retain ISM metals.Therefore, the impact of AGN feedback arises from quenching galaxies via the heating of ISM gas, allowing this gas, along with its associated metals, to join the hot IGrM.
In the case of NoFeedback, the IGrM metal enrichment can re- sult only from tidal stripping.This establishes a baseline for other models, although it is essential to consider that NoFeedback produces significantly more metals overall because of the formation of a larger number of stars.
We now focus on the emission-weighted X-ray Fe and Si abundances in the right panels.Generally, most models show an abundance of [Fe/H] of about ∼ 0.2−0.4,with the NoAGN model slightly below this.However, considering that both the models and the observations show a large scatter in the [Fe/H] abundance, we cannot draw conclusive conclusions.[Si/H] predictions in all models appear to be broadly similar to the observations, but the data do not extend into the group regime to allow direct comparisons.Interestingly, there is not much difference between the original Simba simulation and the Simba-C simulation, or even the simulations without feedback.This is in contrast with the mass-weighted abundances, illustrating the high bias when measuring metallicities only from the X-ray-emitting gas and cautioning against over-interpretation of such data in terms of metal formation mechanisms and timescales.The results can be reconciled by noting that much of the X-ray emission comes from the central region, so potentially the central metallicities could be similar even if the overall mass-weighted ones (dominated by mass in the outskirts) are different.In future work, we will examine group metal profiles.
Fig. 9 is similar to Fig. 8, but shows the global abundance ratios [Si/O] (top row) and [Si/Fe] (bottom row) instead of the chemical abundances. 13These abundance ratios are scaled to the solar level to Anders & Grevesse (1989).
From Fig. 9, we observe notable differences between the four simulations.First, unlike in Fig. 8, we do not see significant differences between the mass-weighted and emission-weighted calculations.Second, three clear trends emerge: (i) The Simba-C simulation have a significantly higher [Si/O] abundance ratios with respect From these differences, it is clear that with the introduction of a stellar feedback system and its chemical enrichment process (NoFeedback → NoAGN), a specific trend/ lope is obtained as a function of temperature.The slope height was only significantly altered again with the introduction of the updated chemical enrichment model and the stellar feedback system.Combining this with the fact that the abundance ratios seem to be insensitive to AGN feedback, it shows that it is crucial that the simulation contains an accurate stellar feedback and chemical enrichment model for the creation of realistic groups and clusters.

REDSHIFT EVOLUTION OF VARIOUS PHYSICAL PROPERTIES AS A FUNCTION OF GROUP TEMPERATURE
In this section, we briefly study the evolution of the 1 keVtemperature groups at all redshifts.We track a specific type of halo as a function of redshift (z = 2 to 0) for all four of our simulations.
We discuss the same plots as in the previous section (Sec.3).However, since this is an evolution study of very particular simulated galaxies groups, we will not show any observations, which are only for z ≈ 0. Furthermore, this section only concentrates on available trends to see if interesting topics can emerge for future studies, and our main goal of determining how the Simba-C simulation with its updated physics modules compares with observations and to Simba remains.

Scaling relations
Here, we present the three scaling relations in specifically chosen temperature bins with equally spaced bin sizes of log T spec = 0.075 keV (log M 500 = 0.208 M ⊙ ), at redshifts z = 2, 1, 0.5 and z = 0. Table 2. Amount of haloes in the 1 keV group temperature bin used in Fig. 10.
group luminosity L X,0.5−2.0 , (ii) middle panel -the group mass M 500 , and (iii) bottom panel -the inner-group region gas entropy S 2500 . 14able 2 shows the number of haloes in the 1 keV group temperature bin used in Fig. 10.Due to the smaller box volume of the Simba NoAGN and Simba NoFeedback the 10 halo limit had to be lowered to 8 haloes to produce complete evolutionary tracks only for these two simulations.This reduces the statistical reliability for these two simulations, but it is a necessary modification.
We present these results for each of our four simulations: (i) The full Simba-C simulation (blue circle solid line), (ii) the original Simba simulation (red square dashed line), (iii) the Simba NoAGN simulation (green star tight dot-dashed line), and finally (iv) the Simba NoFeedback simulation (purple diamond loosely dotdashed line).In addition, we show the 1σ-deviation for each calculated quantity at every redshift for all four of our simulations, indicated by the error bars.We applied a small offset at each redshift to separate the error bars, making it easier to distinguish them in the plot.We use this approach throughout the entirety of Sec. 4. We consider the following redshifts for each plot: z = 2 (Cosmic Noon -highest star formation rate), z = 1, z = 0.5 (start of the accelerated cosmic expansion), and z = 0 (present-day).
From Fig. 10 (upper panel), we observe that without feedback, the evolution of L X,0.5−2.0 for the 1 keV groups consistently exhibits increasing luminosity, resulting, on average, in the brightest haloes.However, when feedback (first stellar and then AGN) is included, the luminosity is reduced by almost an order of magnitude by z = 1, and then followed by a continued lowered luminosity up to z = 0 for Simba-C and Simba.Interestingly, the increasing luminosity appears to start only from z = 1 to 0, while without stellar feedback, it increases throughout the evolution between z = 2 and z = 0.However, when the large errors are taken into account, only the general increasing luminosity as a function of redshift trend and the fact that feedback lowers the luminosity in the galaxy groups' early evolution are significant.The differences between each individual simulation provide only potential nonsignificant patterns, e.g.stellar feedback impacting the period when the group's luminosity begins to evolve.
From Fig. 10 (middle panel), we observe that the evolution of the M 500 for the 1 keV groups in our four simulations is similar to the trends found in the L X,0.5−2.0 − T spec,corr plot (upper panel).The masses of the haloes increase with time, similar to the increase of the luminosity with time.However, the four simulations yielded more closely aligned trends, particularly between Simba-C and the original Simba simulation, which resulted in nearly identical evolutionary tracks for M 500 − T spec,corr .Therefore, the Chem5 model plays no role in the evolution of this scaling relation.
From Fig. 10 (bottom panel), we observe that the gas entropy in the haloes' inner core region -R 2500 -for the 1 keV groups seems to have minimal evolution.Furthermore, AGN feedback, as expected, appears to have the largest impact on the outcome of the gas entropy, resulting in the only significant difference between the four simulation's entropy.The only other interesting effect is that both Simba and Simba-C appear to experience a decrease in gas entropy between z = 0.5 and z = 0 for the former and between z = 1 and z = 0 for the latter in the inner core region.This change is very small, and when the large errors are taken into account, this effect is negligible in this study but deserves a more detailed investigation to determine the origin of this effect.

Evolution of physical properties within galaxy groups
Similarly to the previous section, we discuss the evolution of the various galaxy group properties; however, here we focus on the evolution of the physical properties that govern the group structure.We show the following plots in the specifically chosen T spec,corr = 1 keV group temperature bin: (i) mass fractions are shown in Fig. 11, ( From Fig. 11, we note that the trends appear to be similar to the group mass fractions in Fig. 7.However, these trends show only the evolution of a specific type of galaxy group with a halo spectroscopic temperature of T spec,corr = 1 keV.Interestingly, between Simba-C and Simba there are minor notable differences that arise in the evolution of the mass fractions.Specifically, it seems that Simba-C has more hot diffuse gas at z = 0, while Simba has more cold gas at z = 2.This results in a slight difference in the baryonic mass over the evolution period, but with virtually no difference in the stellar mass.
Since AGN feedback plays an important role in mass fraction calculations, these minor differences can be a direct result of the re-calibrated AGN feedback strength in Simba-C.The differences between Simba NoFeedback and Simba NoAGN simulations gives further indication that stellar feedback can affect the fraction of the mass component and to some extent the evolution and should therefore be taken into account, but AGN feedback remains the primary contributor to the observed differences seen in the mass fractions, as expected.
From Fig. 12, we note that the evolution of the global Si and Fe abundances within these galaxy groups is remarkably similar for three of the simulations.The only major exception to these trends appears to be the Simba NoFeedback simulation, which experiences a substantial increase in both these elements between z = 2 and z = 0; however, this could be the result of the low number of haloes in 50 Mpc h −1 .The other three simulations also show an increase in Fe and Si, although only slightly.This increase is expected to occur as a result of the creation of new metals through stellar feedback.
If the increase experienced by the Simba NoFeedback is not related to the low amount of haloes, it is interesting to note that even without stellar feedback this simulation is still able to obtain similar global Si and Fe abundances at z = 0 solely through stellar evolution and an abundance of star formation from the cooling gas.Therefore, the inclusion of metals with star formation through cooling will still evolve these abundances, regardless of the accompanying stellar feedback model.Hence, the yields should be as accurate as possible to allow the accompanying stellar feedback model to obtain the resulting abundance trends at the correct stage in the galaxy group's evolution.This again motivates the need for the newly updated chemical enrichment model when studying the chemical evolution of the systems.Fig. 13 shows that the global abundance ratios evolve significantly and that they are very sensitive to the stellar feedback physics models included in the simulations.This is one of our most significant differences between the four simulations on group scales.We note three interesting phenomena within the ratios and weighted temperatures for the galaxy groups with T spec,corr = 1 keV.First, the trends/slopes remain relatively consistent regardless of the simulation.Only the values of the abundance ratios differ (either all of them systematically increase or decrease).Second, in both abundance ratios, a ≳ 0.2 dex decrease/increase is obtained as a result of the addition of the instantaneous recycling of metals approximation and its stellar feedback model.This is then reversed by ∼ 0.1 dex with the updated stellar feedback and chemical enrichment model in Simba-C.Third, it seems that AGN feedback plays a minimal role in the evolution of the abundance ratios, especially for [Si/O], confirming Fig. 9's conclusion.
If we focus only on the average, since the 1σ-error range indicates a non-significant difference between Simba-C and Simba, we find an interesting phenomenon between these two simulations.On average, Simba-C obtained a higher abundance of Si compared to the abundance of O, as seen in the increase in the abundance ratio of [Si/O] (Fig. 9) compared to Simba.On the other hand, they both obtained a similar [Si/Fe] abundance ratio at z = 0 (also in Fig. 9).This indicates that Si and Fe scales similarly between the two simulations.Therefore, we have an increase in the abundance of Si and Fe relative to the abundance of O in Simba-C.This corresponds to the frequently used [α/Fe] trend (Wallerstein 1962;Kobayashi et al. 2020b;Kobayashi & Taylor 2023).This trend has a plateau for high [α/Fe] abundance ratios at [Fe/H]< −1 values and then decreases to [α/Fe]∼ 0 after [Fe/H]∼ −1 (owing to an increase in Fe from SNe Ia).Simba-C can successfully reproduce this pattern unlike Simba, as shown in Hough et al. (2023).However, the interplay between the Si, O, and Fe abundances is complex and may not be fully understood from these graphs.A possible explanation for the lowering of O relative to Si and Fe could be the result of the introduction of 'failed' SNe, as shown in Kobayashi et al. (2020b), where all three of these elements were affected.
From the above discussion, we can conclude that to obtain the most realistic abundance ratios, an accurate chemical enrichment and its corresponding stellar feedback model are crucial.

SUMMARY
In this paper, we examine the halo and galaxy group X-ray properties in detail for the cosmological simulation known as Simba-C, the most recent and up-to-date version of Simba.A significant fraction of the baryons within these galaxy groups exist in the form of hot diffuse gas, enabling the study of galaxy groups through X-ray observations.
To identify the haloes containing the galaxy groups within our simulation, we employed the Amiga Halo Finder to generate a catalogue containing information about gas, dark matter, and star particles, along with their corresponding host galaxies, located within each halo.These catalogues were then utilized to analyze the halo Xray properties, including radius, temperature (both T X and T spec,corr ), luminosity, mass, entropy, and metallicity, using the X-ray property of the IntraGroup Medium Python package (XIGrM).
Utilizing these catalogues, we presented and discussed general halo properties such as the halo mass function, the galaxy stellar mass function, and the virial mass and X-ray luminosity as a function of T X .
Furthermore, due to the complexities surrounding the various feedback mechanisms (i.e.stellar and AGN feedback) and metal content (i.e., metal yields) and their influence on simulations, we used different versions of the Simba simulation in an attempt to understand how each physics module's implementation influenced the X-ray properties, with a specific focus on the updated chemical enrichment and stellar feedback physics.These simulations include two 100 Mpc h −1 box simulations: (i) The original published Simba simulation, with an instantaneous recycling of the metals model approximation, an AGN feedback model, and a dust model.(ii) Our main result, the complete Simba-C simulation with its updated chemical enrichment, a self-consistent stellar feedback model, a recalibrated AGN feedback strength, and the reintegrated dust model from Simba.These two simulations are accompanied by two 50 Mpc h −1 box simulations for comparison: (iii) A no feedback Simba simulation, and a (iv) Simba simulation with only the simplistic instantaneous recycling of the metals model.
In Sec. 3, we presented the present-day (z = 0) X-ray scaling relations, namely L X,0.5−2.0 − T spec,corr , M 500 − T spec,corr , and S 500/2500 − T spec,corr (Figs. 4,5,and 6,respectively).We compared the scaling relations with various low-redshift X-ray observations.The first notable result is that the complete Simba-C simulation (our main simulation) appears to be the most consistent in matching the observations of the three scaling relations.This demonstrates that the new chemical enrichment, with its accompanying stellar feedback model, as well as the re-calibration process regarding the AGN feedback strength, is a necessary addition to the Simba simulation.
Second, as expected, the AGN feedback is the most important mechanism for obtaining realistic scaling relations in these simulations, matching the findings from Robson & Davé (2020), while the dust had minimal effects on the outcome of these relations.Simba-C reduces some of the overcorrected halo X-ray properties that came about with the introduction of AGN feedback (see Fig. 6).This could be due to the Chem5 model, which originally produced fewer metals, leading to a weak metal cooling function, which resulted in the recalibration of the AGN feedback strength.
Physical properties, namely, mass fractions and abundance ratios, also showed that AGN feedback played an important role in determining physical properties (Figs. 7,8,and 9,respectively).In fact, for the mass fractions, AGN feedback is the only necessary physics module to obtain simulations that can match the observations.From the abundance ratios, we observed that the Simba-C simulation resulted in a [Si/O] trend that differs from Simba.Interestingly, only the introduction of stellar feedback (NoFeedback to NoAGN) or the update to the chemical enrichment model (Simba to Simba-C) changed the abundance ratios, with the latest change partially reverting some of the changes obtained in the first update.
In Sec.4.1, we presented the same properties; however, we showed the evolution of the 1 keV temperature groups at redshifts z = 2, z = 1, z = 0.5, and z = 0 (Figs. 10 to 13).In most cases, we reach the same conclusions; for instance, AGN feedback being the largest contributor to the differences shown in the simulations.However, the global [Si/O] and [Si/Fe] results (Fig. 13) showed that the abundance ratios are sensitive to stellar feedback.Furthermore, the metal yields of the Chem5 model produced an increase in the abundance of Si and Fe, relative to O (an α-element).This is expected due to the pattern emerging from the frequently-used [α/Fe] ratios, which Simba-C can successfully reproduce, unlike Simba.Therefore, even though AGN feedback is crucial for the simulation to obtain realistic galaxy/galaxy groups, an accurate stellar feedback and its chemical enrichment model are needed to produce realistic abundance ratios.
Lastly, two minor interesting patterns were noted: (i) Simba-C did not change the already correctly modelled M 500 − T spec,corr scaling relation and the mass fractions.Both of which Simba already managed to match the observations.(ii) Simba-C also has more hot diffuse gas at z = 0, while Simba has more cold gas at z = 2, slightly impacting the baryonic mass during the evolution period, but with virtually no difference in the stellar mass.
Future work will include a follow-up study based on the X-ray profiles for each galaxy group property with comparisons to the recent findings of Altamura et al. (2023), where they found that EAGLE-like simulation models do not solve the entropy core problem by studying the S /S 500 profiles, which were revealed to be to flat (Padawer-Blatt et al. in prep.).This would contextualize our im- proved gas entropy S 500 results, since the original Simba simulation also showed a flat gas entropy profile.We are also investigating a follow-up study based on the impact that the Chem5 module has on simulated clusters and their scaling relationships and global properties.This will allow us to determine the Chem5 model's impact on various scales, starting with individual stellar populations to galaxies, groups of galaxies, and clusters.
ST/V000594/1, the Atracción de Talento Contract no.2020-T1/TIC-19882 was granted by the Comunidad de Madrid in Spain, and the science research grants were from the China Manned Space Project.He also thanks the Ministerio de Ciencia e Innovación (Spain) for financial support under Project grant PID2021-122603NB-C21 and HORIZON EUROPE Marie Sklodowska-Curie Actions for supporting the LACEGAL-III project with grant number 101086388.CK acknowledges funding from the UK Science and Technology Facility Council through grant ST/R000905/1 and ST/V000632/1.The Flatiron Institute is supported by the Simons Foundation.Any opinion, finding, and conclusion or recommendation expressed in this material is that of the author(s), and the NRF does not accept any liability in this regard.

Figure 1 .
Figure 1.Comparison of the dust mass function between the Simba-C simulation and the published version of Simba at z = 0, compared to observations from Dunne et al. (2011) and Clemens et al. (2013).The Simba-C simulation's median results are shown by the blue line with its spread in the light blue band, while the red line displays the median Simba results for comparison.
Collaboration et al. (2018) ΛCDM cosmology of Ω m = 0.3, Ω Λ = 0.7, Ω b = 0.048, and H 0 = 68 km s −1 Mpc −1 .The different versions of the simulations are as follows: (i) The Simba simulation without any included feedback mechanisms, but metal injection from stellar evolution is still present -Simba NoFeedback.(ii) The Simba simulation utilizing the instantaneous recycling of the metals stellar feedback model, but excluding AGN feedback.This configuration resembles the simulation used in Liang et al. (2016) -Simba NoAGN.

(
iii) The original Simba simulation, as described inDavé et al. (2019), which contains the instantaneous metal recycling model, the updated AGN feedback model, and the described dust model -Simba.(iv) The complete Simba-C simulation incorporating the newly updated chemical enrichment and stellar feedback model, the recalibrated AGN feedback model from the original Simba simulation, and the re-integrated dust model.We consider this to be our main simulation/result of the study -Simba-C.

Figure 2 .
Figure2.The halo mass function for haloes with at least three (red), two (blue), or one (magenta) luminous galaxies, as well as the complete halo population in the simulation (black).The dashed vertical line shows the halo mass (10 11 M ⊙ ) cut-off introduced in the XIGrM halo analysis.All Simbabased simulations have a halo mass resolution of 6.8×10 9 M ⊙ , corresponding to 64 dark matter particles.

Figure 6 .
Figure 6.Gas entropy at R 500 and R 2500 for the various simulations with their lines as described in Fig. 4. Observations of the following low redshift group data are included for comparison: Sun et al. (2009) (diamonds) and O'Sullivan et al. (2017) (stars).The dashed lines are the best-fit power-law indices of α = 1 and α = 0.74 for the S-T relation at R 500 and R 2500 , respectively, for the full group+cluster sample from Sun et al. (2009).

Figure 7 .
Figure 7. Stellar and gas mass fractions within R 500 for the various different simulations with their lines as described in Fig. 4. For comparison, observations of the following low-redshift group data are included: Eckmiller et al. (2011) (squares), Laganá et al. (2013) (open stars), Gonzalez et al. (2013) (pentagons), Lovisari et al. (2015) (crosses), and Loubser et al. (2018) (filled stars).The upper left panel shows the total baryonic fraction.The black line indicates the cosmological value of the simulation, Ω b /Ωm = 0.16.The upper right panel shows the hot gas fraction.The bottom left panel shows the stellar mass fraction.The bottom right panel shows the cold gas fraction (i.e.diffuse gas with T < 5 × 10 5 K and the galactic ISM).The simulation results include stars in the galaxies and those comprising the diffuse intragroup stars (IGS) component.OnlyGonzalez et al. (2013) accounts for the IGS.

Figure 8 .
Figure 8. Global Fe (top) and Si (bottom) abundances within R 500 for the various different simulations with their lines as described in Fig. 4. The left column shows the mass-weighted abundances, while the right column shows the emission-weighted abundance in the IGrM.For comparison, observations of the following low redshift group data are included: Helsdon & Ponman (2000) (grey diamonds), Peterson et al. (2003) (black triangles) and Yates et al. (2021) (squares).

Figure 9 .
Figure 9. Global [Si/O] (top) and [Si/Fe] (bottom) abundance ratios within R 500 for the various different simulations with their lines as described in Fig. 4. The left column shows the mass-weighted abundances, while the right column shows the emission-weighted abundance in the IGrM.

Figure 10 .
Figure 10.The average: a) L X,0.5−2.0 -values (upper panel) and b) M 500value (middle panel) both within R 500 , as well as the c) gas entropy S 2500 within R 2500 (bottom panel), for the simulated group haloes with T spec,corr = 1 keV over redshift.The 1σ-error bars are shown.The simulations included for comparison are: Simba-C (blue circle/solid line), Simba (red square dashed line), Simba NoAGN (green star/tight dot-dashed line), and Simba NoFeedback (purple diamond/loosely dot-dashed line).We include a small offset at each redshift to separate the error bars for visibility.
ii) global Fe and Si abundances are shown in Fig. 12 , and finally (iii) global abundance ratios [Si/O] and [Si/Fe] are shown in Fig. 13.Table 2 shows the number of haloes in the 1 keV group temperature bin used in Figs.11 -13.

Figure 11 .
Figure 11.The average stellar and gas mass fractions within R 500 for the simulated group haloes with log T spec,corr = 0 keV over redshift.The 1σ-error bars are shown.The simulations included for comparisons are shown with the same lines as described in Fig. 10.The upper left panel shows the total baryonic fraction.The black line indicates the simulation's cosmological value, Ω b /Ω m = 0.16.The upper right panel shows the hot gas fraction.The bottom left panel shows the stellar mass fraction.The bottom right panel shows the cold gas fraction (i.e.diffuse gas with T < 5 × 10 5 K and the galactic ISM).

Figure 12 .
Figure 12.Global Fe (top) and Si (bottom) abundances within R 500 for the simulated group haloes with T spec,corr = 1 keV over redshift.The 1σ-error bars are shown.The simulations included for comparisons are shown with the same lines as described in Fig. 10.The left column shows the mass-weighted abundances, while the right column shows the emission-weighted abundances in the IGrM.

Figure 13 .
Figure 13.Global Si-O (top) and Si-Fe (bottom) abundance ratios within R 500 for the simulated group haloes with T spec,corr = 1 keV over redshift.The 1σ-error bars are shown.The simulations included for comparisons are shown with the same lines as described in Fig. 10.The left column shows the mass-weighted abundance ratios, while the right column shows the emission-weighted abundance ratios in the IGrM.