Photoelectric heating effects on the evolution of luminous disk galaxies
MMon. Not. R. Astron. Soc. , 1– ?? (2005) Printed 18 September 2020 (MN L A TEX style file v2.2)
Photoelectric heating effects on the evolution of luminousdisk galaxies
Omima Osman , (cid:63) , Kenji Bekki , and Luca Cortese ICRAR M468 The University of Western Australia 35 Stirling Hwy, Crawley Western Australia 6009, Australia University of Khartoum - Department of Physics. Al-Gamaa Ave, Khartoum 11115, Sudan. P.O.Box 321
Accepted, Received 2005 February 20; in original form
ABSTRACT
Photoelectric heating (PEH) influences the temperature and density of the interstellarmedium (ISM), and potentially also affecting star formation. PEH is expected to havea stronger effect on massive galaxies, as they host larger dust reservoirs compared todwarf systems. Accordingly, in this paper, we study PEH effects in Milky Way–likegalaxies using smoothed particle hydrodynamics (SPH) code which self–consistentlyimplements the evolution of the gas, dust, and interstellar radiation field (ISRF).Dust evolution includes dust formation by stars, destruction by SNe, and growth indense media. We find that PEH suppresses star formation due to the excess heatingthat reduces the ISM density. This suppression is seen across the entire range of gasfractions, star formation recipes, dust models, and PEH efficiencies investigated byour code. The suppression ranges from negligible values to approximately a factorof five depending on the specific implementation. Galaxy models having higher gasfraction experience higher star formation suppression. The adopted dust model alsoalters the extent of star formation suppression. Moreover, when PEH is switched on,galaxy models show higher gas outflow rates and have higher loading factors indicativeof enhanced SNe feedback. In gas–rich models (i.e. a gas fraction of 0.5), we also findthat PEH suppresses the formation of disk clumps via violent disk instabilities, andthus suppresses bulge formation via clumps migration to the central regions.
Key words:
ISM: dust – galaxies: ISM – galaxies: evolution – stars: formation
For several decades researchers have been studying the im-portance of the photoelectric heating (PEH) of the gas bydust in the thermodynamical balance of the ISM (Watson1972; Draine 1978; Bakes & Tielens 1994; Wolfire et al. 2003;Weingartner & Draine 2001; Hill et al. 2018). Electronsejected from dust grain surfaces by far–ultraviolet (FUV)photons are primarily responsible for the heating of the coldneutral medium (CNM) and diffuse atomic hydrogen (HI)regions (Wolfire et al. 1995; Ingalls, Reach, & Bania 2002).Those electrons store their kinetic energy in the gas throughcollision with the different chemical species therein. PEH isvery efficient in coupling FUV radiation (6 − . eV ) to thegas, and it is mainly caused by the smallest size dust grains,i.e. Polycyclic Aromatic Hydrocarbons (PAH, Bakes & Tie-lens 1994; Okada et al. 2013). The rate and efficiency ( F e )of the PEH depend on the detailed ISM and dust proper-ties such as the interstellar radiation field (ISRF) strength,ISM density and temperature, electron density, and grain (cid:63) E-mail: [email protected] sizes (see models proposed by, e.g. Bakes & Tielens 1994and Weingartner & Draine 2001).Furthermore, PEH influences the dynamical evolutionof galaxies through suppressing star formation as a result ofthe drop in gas density. However, the importance of PEH inregulating star formation is still arguable compared to theother feedback mechanisms such as SNe. For instance, us-ing numerical simulation, Forbes et al. (2016) investigatedSNe feedback and PEH in dwarf galaxies and concluded thatPEH is the dominant process in regulating star formation.On the other hand, Hu et al. (2017) found the contrary andreported that PEH is unable to create outflows and pro-duce the multiphase structure of the ISM, unlike SNe. Huet al. (2017) rigorously justified the discrepancy betweentheir results and those of Forbes et al. (2016) by attributingit mainly to inconsistent cooling implementation by Forbeset al. (2016). For Milky Way–like galaxies, Tasker (2011)who studied the influence of star formation and PEH on thelifetime of giant molecular clouds (GMCs) found that PEHchanges the ISM structure and suppresses the initial star for-mation (in the period between 50 and 200 Myrs). After 200Myrs, the GMCs where PEH is switched on have higher star c (cid:13) a r X i v : . [ a s t r o - ph . GA ] S e p O. Osman, K. Bekki, and L. Cortese formation rate (SFR) since they return enough gas to con-tinue forming stars. It should be noted that Tasker’s (2011)study did not include SNe feedback. Bekki’s (2015) studyon disk galaxies, which involved a sophisticated implemen-tation of the dust hydrodynamical interaction with the restof the galaxy components, also predicted SFR suppressionwhen PEH is included in simulations.PEH is commonly studied in the context of feedbackin dwarf galaxies (e.g. Forbes et al. 2016; Hu et al. 2017;Emerick et al. 2019). This is justified by their small sizeand mass, which allow for high–resolution simulations, andthe fact that they are chemically young objects (i.e. simpleISM to model). However, dwarf galaxies, have lower dustcontent compared to massive galaxies (Remy–Ruyer et al.2014; Grossi et al. 2015) as the bulk of their stars (initialsource of dust) formed at later epoch compared to massivegalaxies (Cowie et al. 1996; Gavazzi et al. 1996). Moreover,dust growth in the ISM of dwarf galaxies requires a moreextended time period to overcome dust production by stars(the critical metallicity to set on dust growth is low, but thegalaxy takes a long time to reach it, Inoue 2011; Asano etal. 2013). Zhukovska (2014) found that dust growth in dwarfgalaxies becomes important in a timescale of 0.1 to 1 Gyr.We expect PEH to have a significant influence on theevolution of more massive, luminous disk galaxies since theyefficiently form dust through star formation and dust growth(see, e.g. Dwek 1998; Hirashita 2013; Aoyama et al. 2017 fordiscussions on the dust processing in the ISM). In particular,the PEH rate is directly proportional to the dust–to–gas ra-tio (Eq. 5 section 2.2.2 of this paper, Bakes & Tielens 1994;Wolfire et al. 2003; Bergin et al. 2004) which has highervalues in massive galaxies. Dust and gas are closely linked,and they interact through several processes that influenceboth (e.g. PEH, dust growth and destruction, molecular hy-drogen formation, Hirashita 2000; Nozawa, Kozasa & Habe2006; Wakelam et al. 2017). Thus, we expect PEH also to in-fluence the spatial distribution of the ISM components andregulate dust evolution. These two last points can not beresolved assuming a constant dust–to–gas ratio since dustevolution (varying dust–to–metal ratio) is proven to be rel-evant (e.g. Dwek & Scalo 1980; Hirashita 2000; Asano et al.2013; Bekki 2015; Chiang et al. 2018).Thus the purpose of this paper is to investigate therole of PEH in suppressing star formation in Milky way–likegalaxies. We also study the dependence of this suppressionon the gas fraction, the PEH efficiency, and the galaxy mass.Moreover, the paper investigates how PEH self–regulatesdust evolution, gas and H consumption, and their spatialdistribution (the ISM structure). To accurately model PEH,one needs to simultaneously model the evolution of the gas,dust, and ISRF (Table 1 shows the dust physics includedin the present models). Although several theoretical stud-ies succeeded in modelling the time and space varying ISRF(e.g. Forbes et al. 2016; Hu et al. 2017; Emerick et al. 2019)using adaptive mesh refinement hydrodynamical codes andsmoothed particle hydrodynamical (SPH) simulations, noneof them modelled (i) dust formation and evolution, (ii) for-mation of H on dust grains, and (iii) self–shielding of H effects, all of which could be important for the ISM of lumi-nous disk galaxies (Osman et al. 2020).Here we present a SPH simulation code that self– Table 1.
Dust physics included in the models.Dust physics This work Reference a Dust formation (cid:51)
1, 2Dust destruction (cid:51)
1, 2, 3Dust growth (cid:51)
1, 2, 3Dust shattering and coagulation (cid:55) (cid:51)
3, 4, 5Radiation Pressure on dust (cid:55)
3, 6Electron density of the ISM (cid:55) formation on dust grains (cid:51)
1, 3Grain size distribution (cid:55)
2, 8Grain composition (cid:55) (cid:51)
10, 3Varying dust sticking coefficient (cid:55)
11, 12Varying destruction efficiency (cid:55) (cid:55) (cid:51)
1, 3Grain charge distribution (cid:55) a consistently models the evolution of the gas, dust, and ISRF.However, our code is still limited in its ability to model thePEH process adequately since it does not resolve the elec-tron density in the ISM. The electron density is one of thequantities that determine F e along with the gas temperatureand radiation field. For this reason, we are forced to treat F e as a parameter that ranges from 0.003 to 0.05 (Bakes &Tielens 1994). From an observational point of view, the ratiobetween [CII] 158 µ m line and the total infrared emissionor between [CII] 158 µ m + [OI] 63 µ m and PAH emissionis used as an estimate of the PEH efficiency (Croxall et al.2012; Beirao et al. 2012; Herrera–Camus et al. 2018). Obser-vations of NGC 1097 and NGC 4559 galaxies show a PEHefficiency of up to 0.06 (Croxall et al. 2012). Rollig et al.(2006) derived an analytical formula for electron density inphoton–dominated regions using KOSMA − τ model (Storzeret al. 1996).The rest of this paper is organized as follows: sections 2and 3 describe the simulations set up and the main resultsof the study, respectively. Discussion on the results is givenin section 4, while section 5 contains the conclusions andsummary of the paper. The SPH chemodynamical model used in the present studyis a modified version of the model presented in Bekki (2013)and Osman et al. (2020) (hereafter B13 and O20). In thisupdated version, we advantage from the detailed modellingof the dust evolution, and the time and space varying ISRFpresented in O20 to introduce PEH into the model self–consistently. Hence, we briefly describe the model and referthe reader to B13 and O20 for further details. This study c (cid:13) , 1– ?? ust processes models an isolated Milky Way–like (MW–like) disk withNFW (Navarro, Frenk & White 1996) dark matter halo den-sity profile. Initially, the gaseous halo is in hydrostatic equi-librium and has a NFW density profile. Tables 2 and 3 showsome of the underlying parameters and physics included inthe models. The chemical abundance of eleven elements such as He, C,N, O, Mg and Ca is followed in time after their ejectionby SNIa, SNII, and AGB stars (B13; Bekki 2015). Metalsejected are equally distributed among the neighbouring gasparticles and non–instantaneous recycling is assumed (e.g.,Bekki & Shioya 1998). A fraction of the ejected metals con-denses to form dust grains or accretes onto dust grain sur-faces in the ISM (B13). We adopt the dust model proposedby Dwek (1998) here. The model accounts for: (a) dust for-mation by AGB stars, SNIa, and SNII, (b) dust growth byaccretion of the ISM gas–phase metals, (c) dust destructionby SNe blasts, and (d) formation of polycyclic aromatic hy-drocarbons (PAHs). The dependence of the dust growth anddestruction processes on the ISM properties (temperature,density, and metallicity) is accounted for, i.e. the timescaleof those processes is not constant throughout the galaxy(O20). H formation on dust grains is also self–consistentlyimplemented with the dust evolution (B13). For star formation, we mainly apply a H–dependent starformation recipe (SFR depends on the total density of thegas), however, a H –dependent recipe (SFR depends onthe H density, see B13) is used for comparison with theH–dependent recipe results in one model. In H–dependentrecipe, a gas particle is transformed into a stellar particlewith Salpeter IMF if: the local density exceeds a thresh-old density ( ρ th = 1 cm − in the present study), the localvelocity field is consistent with gravitational collapse (i.e.div v < E sn = 10 erg per SNa) is divided into 90% thermal feedbackand 10% kinematic feedback consistent with the numericalsimulations by Thornton et al. (1998). This energy is dis-tributed equally among the neighbouring particles (Bekkiet al. 2013).The radiative cooling is implemented using thecooling curve by Rosen & Bregman (1995) for 100 (cid:54) T < K and the MAPPING III code for T (cid:62) K (Suther-land & Dopita 1993). The gas also cools down beyond 100 Kvia H cooling to a floor of 10 K. The implementation of thePEH associated with the dust is described in the followingsubsection. We implement the same PEH model implemented in Bekki(2015) where the analytic formula for the photoelectric heat-ing rate (nΓ pe ) proposed by Bakes & Tielens (1994) isadopted: Table 2.
Description of the basic parameters valuesfor the MW–like models.Physical properties Parameter valuesTotal Mass a M h = 10 M (cid:12) Structure b r vir = 245 kpc, c = 10Initial H fraction 0.01Initial metallicity [Fe / H] = 0 .
30 dexMetallicity gradient − .
04 dex/kpcInitial dust/metal ratio 0.4SF c ISRF, ρ th = 1 cm − IMF Salpeter ( α = 2 . d (cid:15) dm = 935 pc, (cid:15) g = 94 pcGas mass resolution m g ∼ × M (cid:12) a M h = M dm + M g , where M dm and M g are the to-tal masses of dark matter halo and gas in a galaxy,respectively. b For the structure of the dark matter halo NFWprofile with a virial radius ( r vir ) and a c –parameter isadopted. c ρ th is the threshold gas density for star formation.The interstellar radiation field (ISRF) is included inthe estimation of H mass fraction in this model. d (cid:15) dm and (cid:15) g are for the dark matter and gas, respec-tively. Table 3.
The grid of the models analyzed for this study. f g , f b , F e , and SF recipe are the gas fraction, baryonic fraction, PEHefficiency, and the star formation recipe dependent on either thetotal gas (H) or molecular gas (H ) densities. Dust evolution (dueto formation by stars, growth in dense media, and destruction bySNe) is switched on in the model if indicated (cid:51) . M14 and M15 aremodels of interacting pair of galaxies. M19 and M20 are modelswith baryon fraction 50% less than the rest of the models (i.e., f b = 0.03).Model ID f g f b F e SF recipe Dust evolutionM1 0.1 0.06 0.05 H (cid:51)
M2 0.1 0.06 0.0 H (cid:51)
M3 0.5 0.06 0.05 H (cid:51)
M4 0.5 0.06 0.0 H (cid:51)
M5 0.3 0.06 0.05 H (cid:51)
M6 0.3 0.06 0.0 H (cid:51)
M7 0.03 0.06 0.05 H (cid:51)
M8 0.03 0.06 0.0 H (cid:51)
M9 0.1 0.06 0.03 H (cid:51)
M10 0.1 0.06 0.01 H (cid:51)
M11 0.1 0.06 0.003 H (cid:51)
M12 0.1 0.06 0.05 H (cid:51) M13 0.1 0.06 0.05 H (cid:55)
M14 0.1 0.06 0.05 H (cid:51)
M15 0.1 0.06 0.0 H (cid:51)
M16 0.5 0.06 0.03 H (cid:51)
M17 0.5 0.06 0.01 H (cid:51)
M18 0.5 0.06 0.003 H (cid:51)
M19 0.5 0.03 0.05 H (cid:51)
M20 0.5 0.03 0.0 H (cid:51) c (cid:13) , 1– ?? O. Osman, K. Bekki, and L. Cortese n Γ pe = βF e nG ergcm − s − (1)where β = 1 × − , n, F e , and G are the number densityof the ISM, PEH efficiency, and the intensity of the FUVfield in units of Habing (1968), respectively.The PEH rate and efficiency ( F e ) depend on the under-lying ISM and dust properties such as the ISRF strength,ISM density and temperature, electron density, and grainsizes (Bakes & Tielens 1994; Weingartner & Draine 2001).To estimate the time and space varying G in the presentmodel, we follow the following steps:1- estimate the FUV radiation strength . For each stellarparticle, the stellar population synthesis codes by Bruzual& Charlot (2003) are used to estimate its SED accordingto its age and metallicity, and hence estimate the strengthof the FUV part of the ISRF. Not all the FUV radiation isused in PEH, part of which is exhausted by dust extinction.Thus the flux at a wavelength λ in the FUV part from the i th stellar particle around the j th gas particle is given bythe following equation according to the screen model: f λ,i = f λ, ,i e − τ λ,j r i,j /h j , (2)where f λ, ,i , τ λ,j , r i,j , and h j are the original flux of thestellar particle, FUV optical depth, distance between the gasand the stellar particles, and SPH smoothing length of thegas particle.2- estimate the fraction of light absorbed by dust , thisfraction is calculated using the following equation: F ex,i,j = (cid:90) e − τ λ,j r i,j dr, (3)where for the adopted SPH kernel (equals 0 at h j ), the in-tegration in this equation ranges between 0 and 1. To esti-mate τ λ,j for each gas particle, first the optical dust extinc-tion ( A V,j ) is estimated using the gas column density andthe dust–to–gas ratio. Then, A FUV,j is estimated using theCalzetti’s extinction law which relates A V,j and A FUV,j (i.e.Eq. 4 in Calzetti et al. 2000), and thus τ λ,j ( τ λ ≈ . A λ ).The gas column density is estimated using the 3D hydrogendensity and h j , i.e. ρ j ( H ) h j . Hence G ,i,j is given by thefollowing equation: G ,i,j = F ex,i,j g ,i,j , (4)where g ,i,j is G ,i,j assuming no dust extinction.Estimating the ISM electron density involves followinga network of chemical reactions of the abundant elementsin the ISM (Rollig et al. 2006) as well as the grain size andcharge distributions (Bakes & Tielens 1994; Weingartner &Draine 2001; Wolfire et al. 2003) which are outside the scopeof this study. Accordingly, we are forced to treat the F e as aparameter that ranges from 0.003 to 0.05 (Bakes & Tielens1994). Several recent hydrodynamical simulations implementedPEH in the study of the ISM cloud evolution, feedback in dwarf galaxies, and evolution of disk galaxies. In the follow-ing, we review the PEH implementation in a few of thesestudies. Most of these studies used the analytic formula pro-posed by Bakes & Tielens (1994) in one form or another.In dwarf galaxies, Hu et al. (2017) and Emerick et al.(2019) adopted the formula in Eq. 2, following Bakes & Tie-lens (1994), Wolfire et al. (2003), and Bergin et al. (2004).Both studies were able to track free electrons in their mod-els. However, Emerick et al. (2019) applied F e that varieswith the gas density in the form F e = 0.0148n . (becausetheir chemical network does not track electrons from all es-sential elements), and they used Remy–Ruyer et al. (2014)fit to estimate the dust–to–gas ratio. Hu et al. (2017) foundPEH to be a subdominant feedback mechanism in dwarfgalaxies, while Emerick et al. (2019) reported that the mul-tichannel feedback (including PEH) in their models resultedin a realistic evolution of a dwarf galaxy that is consistentwith observations in terms of star formation and metallicityof outflows. Forbes et al. (2016) used a recipe in which thePEH rate is directly proportional to the flux of the FUVphotons and the density of metals ( Zn H ) with the approxi-mation that temperature and electron density have negligi-ble effects in cold/dense media.Γ pe = 1 . × βF e nG eff D ergcm − s − (5)where G eff = G e ( − . × − DN H,tot ) is the attenuated ra-diation field strength in units of Habing and D is the dust–to–gas ratio. Although Eqs 1 and 2 appear similar, the un-derlying physics of the dust and ISM is different.In disk galaxies, Butler et al. (2017) used the same for-mula used by Forbes et al. (2016) for the PEH rate in theirkpc study on how stellar feedback regulates SFR. The au-thors concluded that including only SNe feedback results insimilar SFRs to the SFRs in their models where PEH andionization due to the FUV and extreme ultraviolet (EUV)radiation, respectively, are included as well. However, theISM has very different temperature and chemical states, andthe young stars have a different distribution. Bekki (2015)implementation is the same as the implementation we arepresenting here with F e = 0.003. Tasker (2011) studied theeffect of the PEH and star formation on the GMCs forma-tion and evolution using Eq. 3 for the PEH rate followingWolfire et al. (2003).Γ pe = βF e G (cid:40) exp − ( R − R ) /H R ergs − r (cid:62) kpc exp − (4 − R ) /H R ergs − r < kpc (6)where F e = 0.05, H R , and R are the scale length (= 4.1kpc, Wolfire et al. 2003) and the radial scale length at 8kpc, respectively. Other implementations can also be foundin Choi et al. (2017) and Hill et al. (2018). In all of thosestudies, dust evolution was not explicitly followed. Fig. 1 shows the time evolution of the SFR (upper row) andthe total gas mass (bottom row) in models M3, M4 (f g (gasfraction) = 0.5), M5, M6 (f g = 0.3), M1, M2 (f g = 0.1),and M7, M8 (f g = 0.03). The red dashed, and black solid c (cid:13) , 1– ?? ust processes −1.5−1.0−0.50.00.51.0 l o g ( S F R ) ( M ⊙ ⊙ y r ) f g = 0.5 f g = 0.3 f g = 0.1 f g = 0.03 l o g ( M g ) ( M ⊙ ) W⊙ PEHW⊙O PEH
Ti e (Gyr)
Figure 1.
The time evolution of the SFR (upper row) and the total gas mass (bottom row) in models M3, M4 (f g (gas fraction) = 0.5),M5, M6 (f g = 0.3), M1, M2 (f g = 0.1), and M7, M8 (f g = 0.03). The red dashed and black solid lines represent models with F e = 0.05(W/ PEH) and F e = 0.0 (W/O PEH), respectively. Time (Gyr) −0.7−0.6−0.5−0.4−0.3−0.2−0.1 l o g ( S F R ) ( M ⊙ ⊙ y r ) F e = 0.05F e = 0.03F e = 0.01F e = 0.003F e = 0.0 Figure 2.
The time evolution of the SFR in models M1, M9,M10, M11, and M2 with F e = 0.05, 0.03, 0.01, 0.003, and 0.0 andf g = 0.1, respectively. lines represent models with F e = 0.05 and 0.0, respectively.When F e is greater than zero, the PEH is switched on (W/PEH) in the model with efficiency F e , and it is switched off(W/O PEH) when F e = 0.0. Table 3 gives F e values for thedifferent models.The diffuse heating caused by the PEH results in an ISMwith average temperature one order of magnitude higher inM1, M3, and M5 compared to M2, M4, and M6, and twiceas hot in M7 compared to M8. Thus, the density of theISM drops in models with PEH compared to models with-out PEH. The rise in temperature and decrease in densityincrease the stability of the gas particles in M1, M3, M5,and M7 against gravitational collapse to form stars (lower density and higher Jeans’ mass) compared to the gas parti-cles in M2, M4, M6, and M8 models . Accordingly, the SFR(upper row) and gas consumption (bottom row) in modelswith PEH are lower than in models without PEH. In the ab-sence of late episodic or continuous gas accretion (isolateddisk galaxy), the SFR gradually declines by 0.6 to 0.2 dex asthe gas fraction decreases from 0.5 to 0.03, respectively. Thisdecline is due to the gas consumption by star formation.In addition to the SFRs suppression in models withPEH compared to models without PEH, the magnitudeof suppression (how much SFR is reduced when PEH isswitched on) increases with the gas fraction, which can alsobe seen in the gas evolution (bottom row). For instance,when the gas fraction is increased from 0.1 in M1 and M2(middle–right column) to 0.5 in M3 and M4 (left column),the SFRs rose considerably, especially in the case where PEHis switched off. The average SFRs increased by 0.65 dex inM3 compared to M1 and by 1.14 dex in M4 compared to M2.Moreover, PEH suppressed star formation by 0.74 dex in M3compared to M4, while a suppression of 0.26 dex occurredin M1 compared to M2.One of the models’ limitations is the unresolved elec-tron density which forced us to treat F e as a parameter thatranges from 0.05 to 0.003 (Bakes & Tielens 1994). Fig. 2shows the time evolution of the SFR in models M1, M9,M10, M11, and M2 with F e = 0.05, 0.03, 0.01, 0.003, and0.0 and f g = 0.1, respectively. The suppression of the SFR,as in the case of gas fraction, is present in all models withPEH regardless of F e value, however, its magnitude dependson the actual value of F e adopted (increases with F e ). F e is not constant and varies according to the ISM conditions,hence adopting one value for F e throughout the galaxy re-sults in error in estimating the mean SFR corresponding to c (cid:13) , 1–, 1–
The time evolution of the SFR in models M1, M9,M10, M11, and M2 with F e = 0.05, 0.03, 0.01, 0.003, and 0.0 andf g = 0.1, respectively. lines represent models with F e = 0.05 and 0.0, respectively.When F e is greater than zero, the PEH is switched on (W/PEH) in the model with efficiency F e , and it is switched off(W/O PEH) when F e = 0.0. Table 3 gives F e values for thedifferent models.The diffuse heating caused by the PEH results in an ISMwith average temperature one order of magnitude higher inM1, M3, and M5 compared to M2, M4, and M6, and twiceas hot in M7 compared to M8. Thus, the density of theISM drops in models with PEH compared to models with-out PEH. The rise in temperature and decrease in densityincrease the stability of the gas particles in M1, M3, M5,and M7 against gravitational collapse to form stars (lower density and higher Jeans’ mass) compared to the gas parti-cles in M2, M4, M6, and M8 models . Accordingly, the SFR(upper row) and gas consumption (bottom row) in modelswith PEH are lower than in models without PEH. In the ab-sence of late episodic or continuous gas accretion (isolateddisk galaxy), the SFR gradually declines by 0.6 to 0.2 dex asthe gas fraction decreases from 0.5 to 0.03, respectively. Thisdecline is due to the gas consumption by star formation.In addition to the SFRs suppression in models withPEH compared to models without PEH, the magnitudeof suppression (how much SFR is reduced when PEH isswitched on) increases with the gas fraction, which can alsobe seen in the gas evolution (bottom row). For instance,when the gas fraction is increased from 0.1 in M1 and M2(middle–right column) to 0.5 in M3 and M4 (left column),the SFRs rose considerably, especially in the case where PEHis switched off. The average SFRs increased by 0.65 dex inM3 compared to M1 and by 1.14 dex in M4 compared to M2.Moreover, PEH suppressed star formation by 0.74 dex in M3compared to M4, while a suppression of 0.26 dex occurredin M1 compared to M2.One of the models’ limitations is the unresolved elec-tron density which forced us to treat F e as a parameter thatranges from 0.05 to 0.003 (Bakes & Tielens 1994). Fig. 2shows the time evolution of the SFR in models M1, M9,M10, M11, and M2 with F e = 0.05, 0.03, 0.01, 0.003, and0.0 and f g = 0.1, respectively. The suppression of the SFR,as in the case of gas fraction, is present in all models withPEH regardless of F e value, however, its magnitude dependson the actual value of F e adopted (increases with F e ). F e is not constant and varies according to the ISM conditions,hence adopting one value for F e throughout the galaxy re-sults in error in estimating the mean SFR corresponding to c (cid:13) , 1–, 1– ?? O. Osman, K. Bekki, and L. Cortese −0.8−0.7−0.6−0.5−0.4
H-dependentH -dependent W/ DEW/O DE
Time (Gyr) l o g ( S F R ) ( M ⊙ / y r ) Figure 3.
The time evolution of the SFR in models M1 (reddashed, F e = 0.05) and M12 (black solid, F e = 0.05) with H andH dependent SF recipes, respectively, top panel. The bottompanel shows the same in models M1 (red dashed, W/ DE) andM13 (black solid, W/O DE) with and without dust evolution,respectively. a standard deviation of up to 0.1 dex. The same argumentholds for models M16, M17, and M18 with gas fraction 0.5,where the estimation of the SFR is associated with an errorcorresponding to 0.2 dex standard deviations.The physics included in the models influences the re-sultant SFRs as well. In Fig. 3 we examine the influence ofthe SF recipe and the dust model on the previous results.The top panel shows the time evolution of the SFR in mod-els M1 (red dashed) and M12 (black solid) with H and H dependent SF recipes, respectively. In both models, PEH isswitched on ( F e = 0.05). Adopting SFR that scales with thetotal gas density rather than with the H density results in aslight overestimation of the SFR (0.02 dex on average). Thedensity of the total gas is more often around the thresholddensity for the star formation, unlike H density, which re-sults in a reduction of the star formation. The two modelsshare similar SFRs only when the gas is mostly molecular.The slightly lower SFR in M12 beyond 500 Myrs resulted indust build–up which, in turn, suppressed the star formationfurther through PEH. With a slight delay, the build–up ofdust enhanced H abundance. In this case, where the ISMhas significantly high molecular hydrogen fraction, adoptingeither of the star formation recipes results in a small dif-ference. On the contrary, in the case of the model withoutPEH, the difference is more significant since it has a lowerhydrogen fraction throughout the course of the evolution.The bottom panel of Fig. 3 shows the time evolutionof the SFR in models M1 (red dashed, W/ DE) and M13(black solid, W/O DE). Dust evolution due to formation bystars, destruction by SNe, and growth in dense media is im-plemented in M1. In M13, constant dust–to–gas ratio (D) isused instead. PEH is implemented in both models with F e =0.05. Implementing constant D results in an overestimationof the SFR (0.04 dex on average) since D decreases steadily as the galaxy evolves and consumes its gas which suppressesthe PEH effect. In the model with dust evolution, dust abun-dance increases gradually with time until it reaches a peakbefore declining at later epochs. The adopted dust modelwould not influence the SFR in models without PEH andwith H–dependent SF recipe since these models are hardlysensitive to how much dust is present in the ISM. This im-plies that models that implement constant D underestimateor overestimate PEH effect in suppressing SFRs dependingon the adopted D. Fig. 4 shows the xy projections of the mass surface density ofthe total gas (upper row), H (middle–upper row), metals(middle–bottom row), and dust (bottom row) (Σ G , Σ H ,Σ M , and Σ D , respectively) in logarithmic scale at T = 1 Gyr.The left column shows models M1 (left subfigures in eachpanel) and M2 (right subfigures in each panel). The rightcolumn shows models M3 (left subfigures in each panel) andM4 (right subfigures in each panel).The extra heat supplied to the ISM through PEH causesthe 3D density of the individual gas particles ( n ) to dropwhich in turn suppresses star formation. Accordingly, M1and M3 models after 1 Gyr of evolution have 6% and 31%higher amount of gas in total compared to the gas in M2 andM4 models, respectively. The process also influences H anddust abundances, as the relatively quiet environment (lessstar formation and SNe going off) allows H to continueforming on dust grains (Cazaux & Tielens 2004; Fukui &Kawamura 2010) and dust grains to continue growing inmolecular clouds (Savage & Sembach 1996; Jones 2000; Sofia2004). The result is 55% and 43% higher amount of H and33% and 46% higher amount of dust in M1 and M3 models,respectively. Metals are affected in such a way that modelswithout PEH have a higher amount of metals. This is a resultof the higher SFRs (more metals produced by stars) and theless/more efficient dust growth/destruction in these modelscompared to models with PEH. Thus, the gas, H , and dustsurface densities are higher in models M1 and M3 comparedto models M2 and M4, while M2 and M4 have higher metalssurface densities. The increase/decrease in the total mass ofthe H , dust, and total gas, in 1 Gyr (the mass growth rate)is given in Table 4.Furthermore, the xy projection of the total gas, H , anddust in M1 model show multiple high–density cores apartfrom the central region, while metals projection does notshow such cores since at the position of these cores met-als are used up by dust growth. The area and density con-trast between these cores and their background are prob-ably smaller than the characteristic values of clumps. Onthe contrary, the xy projections in M2 do not show as manyhigh–density regions. Although gas distribution in M1 modelshows a rather small bar or a boxy/peanut–shaped bulge(Berentzen et al. 2007) compared to the bar seen in M2,the stellar distributions show spiral arms in both models, abar in M1, and a boxy–shaped bulge in M2. It is difficult toconclude that one galaxy/galaxy–model is clumpier than theother just by visual inspection. In case we adopt the simpledefinition that clumps are high–density contrasts embeddedin a low–density background covering an area of at leastone square kpc (consistent with clump size in high–redshift c (cid:13) , 1– ?? ust processes W/ PEH
W/O PEH log(Σ G ) (M ⊙ /kpc ) W/ PEH W/O PEH log(Σ G ) (M ⊙ /kpc ) W/ PEH
W/O PEH log(Σ H ) (M ⊙ /kpc ) W/ PEH W/O PEH log(Σ H ) (M ⊙ /kpc ) W/ PEH
W/O PEH log(Σ M ) (M ⊙ /kpc ) W/ PEH
W/O PEH log(Σ M ) (M ⊙ /kpc ) W/ PEH
W/O PEH log(Σ D ) (M ⊙ /kpc ) W/ PEH W/O PEH log(Σ D ) (M ⊙ /kpc ) Figure 4.
The xy projection of the mass surface density of the total gas (upper row), H (middle–upper row), metals (middle–bottomrow), and dust (bottom row) (Σ G , Σ H , Σ M , and Σ D , respectively) in logarithmic scale at T = 1 Gyr. The left column shows M1 (leftsubfigures in each panel) and M2 (right subfigures in each panel) models. The right column shows M3 (left subfigures in each panel) andM4 (right subfigures in each panel) models.c (cid:13) , 1– ?? O. Osman, K. Bekki, and L. Cortese −10 −8 −6 −4 −2 0 2 4246
W/ PEH
Par icles coun −10 −8 −6 −4 −2 0 2 4
W/O PEH
Ho diffuse Ho denseWarmCold10 l o g ( T ) ( K ) log(n) (cm −3 ) −10 −8 −6 −4 −2 0 2 4246 W/ PEH
Par icles coun −10 −8 −6 −4 −2 0 2 4
W/O PEH
Ho diffuse Ho denseWarmCold10 l o g ( T ) ( K ) log(n) (cm −3 ) Figure 5.
2D histograms of the temperature versus density at T = 1 Gyr. Left panel shows models M1 (left subfigure) and M2 (rightsubfigure). Right panel shows models M3 (left subfigure) and M4 (right subfigure). The vertical and horizontal dashed lines roughlyillustrate regions of the phase diagram where the different ISM phases are contained according to their temperature. clumpy galaxies), we find that only the H projection in M2shows distinct clumps.All the different components of the ISM (HI, H , met-als, and dust) show clumpy structure in M3 and M4 modelswith M4 model showing more prominent clumps comparedto their background. Young stellar populations associatedwith these gas clumps indicate the efficient star formationtherein. M3 model also shows a disrupted bar with a bulgecompared to the bulge seen in M4 model. The bulge in M4could, possibly, be a remnant of a bar that formed earlierand was destroyed due to the accretion of clumpy gas (barself–destruction, Pfenniger & Norman 1990; Norman et al.1996; Kormendy & Kennicutt 2004; Kormendy 2013) whichis indeed present. The bulge could also be formed by migra-tion of the gas clumps to the central regions (e.g., Elmegreenet al. 2008) without passing through the bar formation stage.It is worth noting that the disk in M4 is a disrupted diskwhich means the inclusion of PEH (with efficiency as lowas 0.003) helps in stabilizing gas–rich (f g (cid:62) > K and n > cm − ) seen in models M1, M3, andM4 is SNe shocked–heated gas located in the disk, whilethe low–density gas (n < − cm − ) resides in the halo,mostly, ( | z | > < | z | < K and the warm/coldgas with a temperature of around 10 K is smoother in thepresence of PEH.ISM in other simulations without PEH (e.g. Forbes etal. 2016; Hu et al. 2017) has hot gas fraction that is likelyhigher compared to our models without PEH. Moreover, inother studies with and without PEH there is no to littledense hot gas (
T > K and n > cm − ) formed inthe models (e.g. Hu et al. 2016; Forbes et al. 2016; Emericket al. 2019) apart from the model in Capelo et al. (2018)study. Capelo et al. (2018) noted that this gas is kept arti-ficially hot because of the cooling scheme they implement.One interesting diagram is shown in Aoyama et al. (2017;Fig. 2), which is similar to our models with PEH althoughtheir model is without PEH. It is also interesting to see thesuppression of hot gas formation ( T > K) when PEH isincluded in Forbes et al. (2016) (however they have a cool-ing problem in their models; see Hu et al. 2017) models. InForbes et al. (2016) and Tasker (2011) studies, PEH seemsto suppress the scatter in the phase diagram.Fig. 6 illustrates normalized 1D histograms of the totalgas and H surface densities (top row), metals and dust sur-face densities (bottom row) at T = 1 Gyr. The left columndepicts histograms of models M1(red–empty histograms)and M2 (blue filled–hatched histograms). The right columndepicts models M3 (same as M1) and M4 (same as M2).Histograms are normalized in such a way that the total areaunder the histogram is unity. The total gas and metals sur-face densities have similar distributions in models with andwithout PEH and f g = 0.1 (left column), while H and dustsurface densities distributions are slightly shifted towardshigher values in M1 (with PEH). The differences in the ISMproperties between models with and without PEH becomemore considerable in models with f g = 0.5 (right column).In this case, distributions of the total gas and metals surfacedensities are also different. Models with f g = 0.03 and 0.3behave similar to the models in the left and right columns, c (cid:13) , 1– ?? ust processes Figure 6.
1D normalized histograms of the total gas and H surface densities (top row), metals and dust surface densities (bottom row)at T = 1 Gyr. Left column shows models M1 (red–empty histograms) and M2 (blue filled–hatched histograms). Right column showsmodels M3 (same as M1) and M4 (same as M2). respectively. Generally, the differences in the ISM proper-ties between each pair of models (with and without PEH)diminish as the gas fraction decreases. Clumpy disk galaxies are observed not only in the high red-shift universe but also in the local universe with less fre-quency though (Shibuya et al. 2016; Buck et al. 2017). Ac-cordingly, we ran a pair of models that share the same char-acteristics with M3 and M4 but with baryon fraction 50%less than what is present in M3 and M4 (i.e. f b = 0.03) to im-itate present–day gas–rich galaxies, namely, M19 and M20models. Fig. 7 shows the xy projection, xz projection, andthe radial profiles of the total gas (Σ G , left column from topto bottom, respectively) and the young stars (Σ NS , rightcolumn from top to bottom, respectively). M19 is shown inthe left subfigure of each panel in the top and middle rowsand represented by the red dashed lines in the bottom row( F e = 0.05). M20 is shown in the right subfigure of eachpanel in the top and middle rows and represented by theblack solid lines in the bottom row ( F e = 0.0).M19 shows smoother disk structure compared to theclumpy disk in M20 which is a result of violent disk insta-bility. The massive gas clumps in M20 are associated with massive formation of young stars (right column, a commonfeature of clumps, Genzel et al. 2008; Forster Schreiber et al.2009), unlike M19 where stellar clumps are less pronounced.While M20 shows clumps that do not follow a specific pat-tern and a slightly asymmetric disk (a feature of clumpygalaxies, Conselice et al. 2004), M19 shows a symmetricdisk with clumps that are fragments of spiral arms (Inoue& Yoshida 2018; 2019). Viewed edge–on, M19 disk resem-bles present–day galaxies with a triangular bulge, and M20resembles chain galaxies (Cowie, Hu, & Songaila 1995). Theradial profiles show that M19 returns more gas in most ofthe disk compared to M20 since it forms stars at a lowerrate. The clumpy structure is evident in both the gas andyoung stars profiles.Migration of disk clumps to the central regions providesa viable mechanism for bulge formation (Noguchi 1999; Im-meli et al. 2004; Dekel et al. 2009). After evolving M19 andM20 for a period of 3 Gyrs, M20 formed a bulge one orderof magnitude more massive than the bulge in M19. Fig. 8shows the cumulative histograms of the initial location onthe disk (at T = 0) of the bulge gas particles (at T = 3Gyr) excluding particles in the bulge at T = 0, i.e. wherethe bulge particles located in the disk initially. Red dashedand black solid lines represent models M19 and M20, respec-tively. In the case of M20, gas particles from farther parts of c (cid:13) , 1–, 1–
1D normalized histograms of the total gas and H surface densities (top row), metals and dust surface densities (bottom row)at T = 1 Gyr. Left column shows models M1 (red–empty histograms) and M2 (blue filled–hatched histograms). Right column showsmodels M3 (same as M1) and M4 (same as M2). respectively. Generally, the differences in the ISM proper-ties between each pair of models (with and without PEH)diminish as the gas fraction decreases. Clumpy disk galaxies are observed not only in the high red-shift universe but also in the local universe with less fre-quency though (Shibuya et al. 2016; Buck et al. 2017). Ac-cordingly, we ran a pair of models that share the same char-acteristics with M3 and M4 but with baryon fraction 50%less than what is present in M3 and M4 (i.e. f b = 0.03) to im-itate present–day gas–rich galaxies, namely, M19 and M20models. Fig. 7 shows the xy projection, xz projection, andthe radial profiles of the total gas (Σ G , left column from topto bottom, respectively) and the young stars (Σ NS , rightcolumn from top to bottom, respectively). M19 is shown inthe left subfigure of each panel in the top and middle rowsand represented by the red dashed lines in the bottom row( F e = 0.05). M20 is shown in the right subfigure of eachpanel in the top and middle rows and represented by theblack solid lines in the bottom row ( F e = 0.0).M19 shows smoother disk structure compared to theclumpy disk in M20 which is a result of violent disk insta-bility. The massive gas clumps in M20 are associated with massive formation of young stars (right column, a commonfeature of clumps, Genzel et al. 2008; Forster Schreiber et al.2009), unlike M19 where stellar clumps are less pronounced.While M20 shows clumps that do not follow a specific pat-tern and a slightly asymmetric disk (a feature of clumpygalaxies, Conselice et al. 2004), M19 shows a symmetricdisk with clumps that are fragments of spiral arms (Inoue& Yoshida 2018; 2019). Viewed edge–on, M19 disk resem-bles present–day galaxies with a triangular bulge, and M20resembles chain galaxies (Cowie, Hu, & Songaila 1995). Theradial profiles show that M19 returns more gas in most ofthe disk compared to M20 since it forms stars at a lowerrate. The clumpy structure is evident in both the gas andyoung stars profiles.Migration of disk clumps to the central regions providesa viable mechanism for bulge formation (Noguchi 1999; Im-meli et al. 2004; Dekel et al. 2009). After evolving M19 andM20 for a period of 3 Gyrs, M20 formed a bulge one orderof magnitude more massive than the bulge in M19. Fig. 8shows the cumulative histograms of the initial location onthe disk (at T = 0) of the bulge gas particles (at T = 3Gyr) excluding particles in the bulge at T = 0, i.e. wherethe bulge particles located in the disk initially. Red dashedand black solid lines represent models M19 and M20, respec-tively. In the case of M20, gas particles from farther parts of c (cid:13) , 1–, 1– ?? O. Osman, K. Bekki, and L. Cortese
W/ PEH
W/O PEH log(Σ G ) (M ⊙ /kpc ⊙ ) W/ PEH W/O PEH log(Σ NS ) (M ⊙ /kpc ⊙ )W/ PEH W/O PEH log(Σ G ) (M ⊙ /kpc ⊙ ) W/ PEH W/O PEH log(Σ NS ) (M ⊙ /kpc ⊙ ) R (kpc) l o g ( Σ G ) ( M ⊙ ⊙ p c ) W⊙PEHW⊙O PEH
R(kpc) −2.0−1.5−1.0−0.50.00.51.01.5 l o g ( Σ N S ) ( M ⊙ ⊙ p c ) Figure 7. xy projection, xz projection, and the radial profiles of the total gas (Σ G , left column from top to bottom, respectively) andthe young stars (Σ NS , right column from top to bottom, respectively). M19 is shown in the left subfigure of each panel in the top andmiddle rows and represented by the red dashed lines in the bottom row ( F e = 0.05). M20 is shown in the right subfigure of each panelin the top and middle rows and represented by the black solid lines in the bottom row ( F e = 0.0). the disk contribute to the bulge formation, unlike gas par-ticles in model M19. The fraction of the disk particles thatend up in the bulge at 3 Gyr is also higher in the case ofM20. Gas driven out to the halo by star formation activity and itssubsequent SNe events (White & Rees 1978; Dekel & Silk1986; White & Frenk 1991; see Veilleux et al. 2005; Naab &Ostriker 2017 for overviews) is shown in Fig. 9. Fig. 9 showsthe xz projection of the gas surface density at T = 1 Gyr for models M1 (left) and M2 (right). The disk where PEH inswitched on (left) is slightly thicker and shorter comparedto the disk where PEH is switched off, however, what is in-teresting is the amount of gas above the disk in each model.The higher amount of gas in the halo of M1 indicates ahigher gas loading factor; accordingly, the time evolutionof the gas loading factor ( µ ) is investigated and is shownin Fig. 10 for models M1 and M2 (left panel), and M3 andM4 (right panel). Red dashed and black solid lines representmodels M1 and M3, and M2 and M4, respectively. µ is theratio between the gas outflow rate and SFR. The outflowrate ( ˙ M out ) is defined as the mass flux crossing a surface c (cid:13) , 1– ?? ust processes R (kpc) N W/ PEHW/O PEH
Figure 8.
The cumulative histograms of the initial location onthe disk (at T = 0) of the bulge gas particles (at T = 3 Gyr)excluding particles in the bulge at T = 0. The red dashed andblack solid lines represent models M19 and M20, respectively.
Table 4.
A subset of the models presented in Table 3 where ∆ M g ,∆ M H , and ∆ M D are the mass growth rate ( M final − M initial M initial )of the total gas, H , and dust, respectively, in 1 Gyr. R m is thechange in the gas amount in the halo ( | z | > M g ∆ M H ∆ M D R m ( × M (cid:12) )M1 − .
070 31.8 0.34 1.40 ± − .
130 13.8 − .
12 0.02 ± − .
065 47.0 1.80 135.00 ± − .
350 26.4 0.49 7.30 ± − .
072 54.5 1.90 34.00 ± − .
240 37.0 1.90 0.57 ± − .
037 31.9 0.25 0.014 ± − .
045 22.7 0.16 0.015 ± − .
086 38.6 0.60 0.43 ± − .
112 20.7 0.0006 0.05 ± − .
121 13.3 − .
11 0.003 ± located at | z | = 2 kpc per unit time and is calculated usingthe sum: ˙ M out = (cid:80) i m g,i −→ v i . ˆ n/ ∆ x where ˆ n is the normaldirection to the disk, m g,i and −→ v i are the mass and velocityof the i th gas particle. Only particles with −→ v i . ˆ n > x is taken to be 0.5 kpc. M1 and M3have persistently higher gas outflow rates over time com-pared with M2 and M4, and µ shows that these outflowrates are quite significant compared to the SFR and couldinfluence it. Although the lower star formation activity inM1 and M3 compared with M2 and M4 implies that M2and M4 have higher SNe rate which drives gas outflows, theISM in M1 and M3 is hotter and less dense compared to theISM in M2 and M4. For SNe to drive gas outflows efficiently,the ambient ISM needs to have low density (SNe efficiencydecreases with increasing ambient density, see equation 10in Naab & Ostriker 2017; Hu et al. 2017; Hu 2019). Accord-ingly, we argue that PEH could influence the efficiency ofSNe feedback since it influences the density of the ISM. Ad-ditionally, the gas driven out to the halo takes longer to cooldown and fall back to the disk. Furthermore, the efficiencyof SNe events in farther, less enriched, regions of the disk isalso altered, resulting in less enriched gas outflows.To validate those arguments we investigated the ISM W/ PEH
W/O PEH log(Σ G ) (M ⊙ /kpc ⊙ ) Figure 9. xz projection of the gas surface density at T = 1 Gyrfor models M1 (left) and M2 (right). density ambient to SNe, metals and dust distributions inthe halo, and tracked halo gas particles back to the disk.In tracking halo gas particles, we followed all the particlesthat are in the halo ( | z | > R m ( M (cid:12) )) inthe halo ( | z | > R m values account for the extent of the haloadopted (calculated for | z | >
1, 2, and 3 kpc). The extentof the disk is taken to be 24, 35, 40, and 50 kpc for modelswith f g = 0.03, 0.1, 0.3, and 0.5, respectively. The amountof gas in the halo increases with the gas fraction and PEHefficiency because of both influence SFRs in the models. R m is defined as follows: R m = M g,h,f − M g,h,i M (cid:12) (7)where M g,h,f and M g,h,i are the total gas in the halo ( | z | > c (cid:13) , 1– ?? O. Osman, K. Bekki, and L. Cortese
Time (Gyr) μ W/ PEHW/O PEH
Time (Gyr) μ Figure 10.
Gas loading factor (ratio of the outflow rate to the SFR) as a function of time for models M1 and M2 (left panel), and M3and M4 (right panel). The red dashed lines represent models M1 and M3, while the black solid lines represent models M2 and M4.
Figure 11.
Distributions of the density of the ISM ambient to SNe in two pairs of models M1 and M2 (left: as in Fig. 6) and M3 andM4 (right: as in Fig. 6).
Several processes lead to star formation suppression ingalaxies such as stellar feedback that results in galacticwinds (e.g. White & Rees 1978; Dekel & Silk 1986; White &Frenk 1991; Ceverino & Klypin 2009; Muratov et al. 2015;El–Badry et al. 2016), AGN feedback (Nulsen et al. 2005;McNamara & Nulsen 2007; Fabian 2012; Cicone et al. 2014),and environmental effects (e.g. tidal or ram–pressure strip-ping, Gunn & Gott 1972; Moore et al. 1996; Abadi et al.1999; Bekki 2014; Poggianti et al. 2017). Here, we haveshown that photoelectric heating (PEH) of the gas by elec-trons ejected from dust grains suppresses star formation byreducing/increasing the density/temperature of the gas andprobably creating an environment for SNe feedback to bemore efficient. We tried to demonstrate the second point byinvestigating the gas amount in the halo, gas loading factors,time spent in the halo by gas particles, dust and metals dis-tributions in the halo, and the ISM ambient to SNe in mod-els with gas fractions 0.1 and 0.5. All these investigations,except for the ISM ambient to SNe in models with gas frac-tion 0.1 (see section 3.4 however), show that SNe efficiencyis altered somehow in models with PEH. Models with PEH have a higher quantity of gas in their halos ( | z | > F e throughout the galaxy is an oversim- c (cid:13) , 1– ?? ust processes plification which results in errors in estimating SFRs thatare corresponding to standard deviation of 0.1 to 0.2 dex formodels with gas fractions 0.1 and 0.5, respectively. Thoseerrors become less significant compare to the suppression asthe gas fraction increases. PEH and the subsequent reduction of the star formationactivity influence the ISM in a non–linear manner whereabundances, radial profiles, and 2D distributions of the dif-ferent components of the ISM are altered. H abundanceis enhanced because of the lower consumption by star for-mation and dissociation by UV radiation from young stars.Moreover, the lower SNe rate enhances dust abundance byeffectively reducing dust destruction. Thus, H continues toform on dust grains, and dust grains continue to grow inmolecular clouds. Most of the H in those models, however, isnot star–forming molecular gas (i.e. warm molecular hydro-gen). On the contrary, metals are found in lower quantitiesin models with PEH, not only because of the lower amountof metals supplied by stars (lower SFRs, Hu et al. 2016) butalso because of efficient dust growth which exhausted a highfraction of metals. Additionally, the inclusion of PEH leadsto a flattening in the metallicity gradients.PEH effects on the structure and properties of the ISMare found to be more significant as the gas fraction increases.In particular, models without PEH and with gas fraction of0.5 developed pronounced gaseous and stellar clumps (seeFigs 4 and 7). These clumps form gravitating systems to-wards which gas and disk stars gravitate, and they last forseveral dynamical times. Models with PEH formed less pro-nounced, short–lived ( <
100 Myrs) gaseous clumps along thespiral arms. These models can be linked to the clumpy diskgalaxies at high and low redshifts (Cowie, Hu, & Songaila1995; van den Bergh et al. 1996; Conselice et al. 2004; Genzelet al. 2008; Tadaki et al. 2014; Murata et al. 2014; Guo et al.2015; Shibuya et al. 2016; Buck et al. 2017; Guo et al. 2018)where clumps in rotationally supported systems form viadisk instabilities such as Toomre (1964) instability (Noguchi1998; 1999; Shapiro et al. 2008; Dekel et al. 2009; Bournaud& Elmegreen 2009) and spiral arms fragmentation (Inoue& Yoshida 2018; 2019). The fate of these clumpy galaxiesis still under debate, however, PEH could have had an im-portant role in suppressing the formation, shortening the lifetime of these clumps, and further stabilize the disks that arethought to be the progenitors of present–day galaxies (e.g.Shlosman & Noguchi 1993; Noguchi 1998; 1999; Hopkins,Quataert & Murray 2012).Disk clumps are, arguably, able to migrate to the centralregions of the galaxy and contribute to the bulge formation(Noguchi 1999; Immeli et al. 2004; Dekel et al. 2009). Afterevolving models M19 and M20 (Fig. 7) for 3 Gyrs, we reporttentative evidence of bulge formation suppression via PEH.The model without PEH (M20) has a bulge one order ofmagnitude more massive than the bulge in the model withPEH (M19). Our forthcoming studies will investigate PEHeffects on clumps formation and evolution.
Time (Gy ) −0.65−0.60−0.55−0.50−0.45−0.40−0.35−0.30 l o g ( S F R ) ( M ⊙ / y ) Inte -W/ PEHInter-W/O PEH
Figure 12.
The time evolution of the SFR in models: M14 (reddashed, F e = 0.05) and M15 (black solid, F e = 0.0) of interactingpair of galaxies. When PEH is implemented in the models, H density be-comes sufficiently high for a H –dependent star formationrecipe to maintain SFR similar to the SFR when a H–dependent recipe is implemented. This is not the case formodels without PEH since their H density is considerablylow compared to the total gas density which results in sup-pression of the SFR. Hence, with PEH, implementing eitherof the star formation recipes makes little difference. Thecontrary occurs when a constant dust–to–gas ratio is imple-mented (e.g. Hu et al. 2017; Emerick et al. 2019) instead ofimplementing explicit dust evolution. SFRs are influencedin models with PEH and not in models without PEH. SFRsare overestimated in the case of models with PEH. Thus,the magnitude of SFR suppression is underestimated in thiscase, however, it could be also overestimated depending onthe adopted dust–to–gas ratio. Another important property that influences feedback pro-cesses in general and photoelectric heating in particular isthe galaxy mass. Several studies showed that PEH sup-presses star formation in dwarf galaxies (Forbes et al. 2016;Hu et al. 2017) and massive galaxies (Tasker 2011; Bekki2015). Hence, we have run an extra set of models that in-cludes Milky Way, M33, and SMC like (M h = 10 , 10 ,and 10 M (cid:12) ) with a gas resolution of about 2 × M (cid:12) .The models share the same parameters with the MW–likemodels presented here except for the mass, size, and reso-lution. We confirm the previous results and report that themagnitude of suppression depends on the gas fraction aswell as the galaxy mass. The magnitude decreases as themass decreases, however, it is not entirely clear when thegas fraction is small (f g ∼ In Fig. 12, we explore whether or not our findings hold inthe case of interacting galaxies. Models M14 (red dashed, F e c (cid:13) , 1– ?? O. Osman, K. Bekki, and L. Cortese = 0.05) and M15 (black solid, F e = 0.0) contain interactingpair of galaxies. At the beginning (in the first ∼
500 Myrs)when the two galaxies are far apart, they act as isolatedgalaxies in terms of the PEH suppression of SFR, i.e. themodel with PEH has SFR that is about 0.14 dex lower thanthe SFR in the model without PEH. At the time of theirclosest encounter ( ∼
600 Myrs), the model with PEH re-tained enough gas to starburst with higher magnitude/SFRand longer duration compared to the model without PEH.Afterwards, the SFR in the model with PEH drops slightlybelow the SFR in the model without PEH.
In this paper, we have used our original SPH code to studythe influence of photoelectric heating (PEH) of the gas byelectrons ejected from dust grains on the evolution of lu-minous disk galaxies. The evolution of the gas, dust, andISRF are self–consistently implemented in the code. Dustevolution includes dust formation by stars, destruction bySNe, and growth in dense media. However, the code doesnot reach the necessary resolution to resolve the electrondensity in the ISM, hence, we treat the PEH efficiency ( F e )as a parameter that ranges between 0.05 and 0.003 (Bakes& Tielens 1994). Our main results are summarized in thefollowing.(i) The diffuse heating caused by the PEH results in anISM with an average temperature that is a factor of two toone order of magnitude higher compared to the ISM wherePEH is switched off. This causes the density of the ISMto drop by a few and in turn, SFRs are suppressed. Thesuppression occurs in models with different gas fractions,star formation recipes, dust models, and PEH efficiencies.However, the magnitude of suppression (how much SFR isreduced when PEH is switched on) depends on the specificparameters adopted.(ii) On the other hand, PEH enhances SNe feedback byaltering the gas desnity. This enhancement increases the gasfraction in the halo. The gas driven out to the halo takesa long time (longer than the simulation timescale) to cooldown and fall back onto the disk which makes it unavailablefor star formation. The study of the gas mass loading factorsin a few models indicates the efficiency of the outflows inaffecting star formation.(iii) The moderate consumption of the gas by star for-mation in models with PEH results in a higher abundanceof all the different ISM components except for the gas–phasemetals. Metals are less abundant in models with PEH com-pared to models without PEH because of their lower produc-tion rate by stars and their consumption by dust growth.Accordingly, PEH (through dust evolution) influences theglobal properties, radial profiles, and spatial distributions ofthe ISM. Perhaps one of the most important radial profilesare the metallicity profiles, we find that PEH flattens thoseprofiles.(v) Gas–rich disk galaxies in high redshift universe ( z ∼ ACKNOWLEDGEMENT
OO is a recipient of an Australian Government ResearchTraining Program (RTP) Scholarship. LC is a recipi-ent of an Australian Research Council Future Fellowship(FT180100066) funded by the Australian Government.
DATA AVAILABILITY
The data underlying this article will be shared on a justifiedrequest to the corresponding author.
REFERENCES
Abadi M. G., Moore B., Bower R. G., 1999, MNRAS, 308,947Anderson M. E., Bregman J. N., 2010, ApJ, 714, 320Aoyama S. et al., 2017, MNRAS, 466, 105Asano R. S. et al., 2013, MNRAS, 432, 637Bakes E. L. O., Tielens A. G. G. M., 1994, ApJ, 427, 822Behroozi P. S., Conroy C., Wechsler R. H., 2010, ApJ, 717,379Beirao P. et al., 2012, ApJ, 751, 144Bekki k., Shioya Y., 1998, ApJ, 497, 108Bekki K., 2013, MNRAS, 432, 2298Bekki K., Shigeyama T., Tsujimoto T., 2013, MNRAS, 428,L31Bekki K., 2014, MNRAS, 438, 444Bekki K., 2015, MNRAS, 449, 1625Berentzen I., Shlosman I., Martinez-Valpuesta I., Heller C.H., 2007, ApJ, 666, 189Bergin E. A. et al., 2004, ApJ, 612, 921Bournaud F., Elmegreen B. G., 2009, ApJ, 694, L158Bruzual G., Charlot S., 2003, MNRAS, 344, 1000Buck T. et al., 2017, MNRAS, 468, 3628Butler M. J. et al., 2017, 841, 82Cacciato M., Dekel A., Genel S., 2012, MNRAS, 421, 818Calzetti D. et al., 2000, ApJ, 533, 682Capelo P. R. et al., 2018, MNRAS, 475, 3283Cazaux S., Tielens A. G. G. M., 2004, ApJ, 604, 222Ceverino D., Klypin A. 2009. ApJ, 695, 292Chiang I-Da et al., 2018, ApJ, 865, 23Choi E. et al., 2017, ApJ, 844, 31Cicone C. et al., 2014, A&A, 562, 25Conselice C. J. et al., 2004, ApJ, 600, L139Cowie L., Hu E., Songaila, A. 1995, AJ, 110, 1576Cowie L. L., Songaila A., Hu E. M., Cohen J. G., 1996, AJ,112, 839Croxall K. V. et al., 2012, ApJ, 747, 81Daddi E. et al., 2008, ApJ, 673, L21Daddi E. et al., 2010, ApJ, 713, 686Dekel A., Silk J., 1986, ApJ, 303, 39Dekel A. et al., 2009, Nature, 457, 451Draine B. T., 1978, ApJSS, 36, 595Dwek E., Scalo J. M., 1980, ApJ, 239, 193Dwek E., 1998, ApJ, 501, 643 c (cid:13) , 1– ?? ust processes El-Badry K. et al., 2016, ApJ, 820, 17Elmegreen B. G., Bournaud F., Elmegreen D. M., 2008, ApJ,688, 67Emerick A., Bryan G. L., Mac Low Mordecai–Mark, 2019,MNRAS, 482, 1304Fabian A. C., 2012, ARA&A, 50, 455Fisher D. B. et al., 2014, ApJ, 790, L30Forbes J. C. et al., 2016, Nature, 535, 523Forster Schreiber N. M. et al., 2009, ApJ, 706, 1364Fukugita M., Hogan C. J., Peebles P. J. E., 1998, ApJ, 503,518Fukui Y., Kawamura A., 2010, ARA&A, 48, 547Gavazzi G., Pierini D., Boselli A., 1996, A&A, 312, 397Genzel R. et al., 2008, ApJ, 687, 59Grossi M. et al., 2015, A&A, 574, A126Gunn J. E., Gott J. Richard III, 1972, ApJ, 176, 1Guo Y. et al., 2015, ApJ, 800, 39Guo Y. et al., 2018, ApJ, 853, 108Habing H. J., 1968, Bull. Astron. Inst. Neth., 19, 421Herrera–Camus R. et al., 2018, ApJ, 861, 95Hill A. S., Mac Low Mordecai–Mark, Gatto A., Ibanez-MejiaJ. C., 2018, ApJ, 862, 55Hirashita H., 2000, PASJ, 52, 585Hirashita H., 2013, POS (LCDU 2013) 027Hitschfeld M., Kramer C., Schuster K. F., Garcia-Burillo S.,Stutzki J., 2009, A&A, 495, 795Hopkins P. F., Quataert E., Murray N., 2012, MNRAS, 421,3488Hu Chia-Yu et al., 2016, MNRAS, 458, 3528Hu Chia-Yu et al., 2017, MNRAS, 471, 2151Hu Chia–Yu, 2019, MNRAS, 483, 3363Immeli A., Samland M., Westera P., Gerhard O., 2004, ApJ,611, 20Ingalls J. G., Reach W. T., Bania T. M., 2002, 579, 289Inoue A. K., 2011, Earth Planets and Space, 63, 1027Inoue S., Yoshida N., 2018, MNRAS, 474, 3466Inoue S., Yoshida N., 2019, MNRAS, 485, 3024Jones A. P., 2000, JGR (Journal of Geophysical Research),105, 10257Komatsu E. et al., 2011, ApJS, 192, 18Kormendy J., 2013, In XXIII Canary Islands Winter Schoolof Astrophysics, Secular Evolution of Galaxies eds. J.Falcon-Barroso & J. H. Knapen (Cambridge: CambridgeUniversity Press), p. 1Kormendy J., Kennicutt R. C., 2004, ARAA, 42, 603McNamara B. R., Nulsen P. E. J., 2007, ARA&A, 45, 117Mo H. J. et al., 2005, MNRAS, 363, 1155Moore B. et al., 1996, Nature, 379, 613Moster B. P., Naab T., White S. D. M., 2013, MNRAS, 428,3121Murata K. L. et al., 2014, ApJ, 786, 15Muratov A. L. et al., 2015, MNRAS, 454, 2691Naab T., Ostriker J. P., 2017, ARA&A, 55, 59Navarro J. F., Frenk C. S., White S. D. M., 1996, ApJ, 462,563Noguchi M., 1998, Nature, 392, 253Noguchi M., 1999, ApJ, 514, 77Norman C. A., Sellwood J. A., Hasan H., 1996, ApJ, 462,114Nozawa T., Kozasa T., Habe A., 2006, ApJ, 648, 435Nulsen, P. E. J. et al., 2005, ApJ, 628, 629Okada Y. et al., 2013, A&A, 553, A2 Osman O., Bekki K., Cortese L., 2020, MNRAS, submittedOstriker J. P., Peebles P. J. E., 1973, ApJ, 186, 467Pfenniger D., Norman C., 1990, ApJ, 363, 391Poggianti B. M. et al., 2017, ApJ, 844, 21Puech M., 2010, MNRAS, 406, 535Remy–Ruyer A. et al., 2014, A&A, 563, A31Ribeiro B. et al., 2017, A&A, 608, A16Robertson B. E., Bullock J. S., 2008, ApJ, 685, L27Rollig M. et al., 2006, A&A, 451, 917Rosen A., Bregman J. N., 1995, ApJ, 440, 634Savage B. D., Sembach K. R., 1996, AR A&A, 34, 279, 329Schneider R., Ferrara A., Natarajan P., Omukai K., 2002,ApJ, 571, 30Schneider R., Ferrara A., Salvaterra R., Omukai K., BrommV., 2003, Nature, 422, 869Shapiro K. L. et al., 2008, ApJ, 682, 231Shibuya T., Ouchi M., Kubo M., Harikane Y., 2016, ApJ,821, 72Shlosman I., Noguchi M., 1993, ApJ, 414, 474Sofia U. J., 2004, ASP, 309, 393Storzer H., Stutzki J., Sternberg A., 1996, A&A, 310, 592Sutherland R. S., Dopita M. A., 1993, ApJS, 88, 253Tacconi L. J. et al., 2008, ApJ, 680, 246Tadaki Ken-ichi et al., 2014, ApJ, 780, 77Tasker E. J., 2011, ApJ, 730, 11Thornton K., Gaudlitz M., Janka H.-Th., Steinmetz M.,1998, ApJ, 500, 95Toomre A., 1964, ApJ, 139, 1217van den Bergh S. et al., 1996, AJ, 112, 2Veilleux S., Cecil G., Bland–Hawthorn J., 2005, ARA&A,43, 769Wakelam V. et al., 2017, Molecular Astrophysics, 9, 1Watson W. D., 1972, ApJ, 176, 103Weiner B. J. et al., 2006, ApJ, 653, 1027Weingartner J. C., Draine B. T., 2001, ApJSS, 134, 263White S. D. M., Rees M. J., 1978, MNRAS, 183, 341White S. D. M., Frenk C. S., 1991, ApJ, 379, 52Wolfire M. G. et al., 1995, ApJ, 443, 152Wolfire M. G., Mckee C. F., Tielens A. G. G. M., 2003, ApJ,587, 278Zhukovska S., 2014, A&A, 562, A76 c (cid:13) , 1–, 1–