{'fraction_non_alphanumeric': 0.03485144018906674, 'fraction_numerical': 0.019371084760106093, 'mean_word_length': 4.188925779158188, 'pattern_counts': {'":': 0, '<': 30, '<?xml version=': 0, '>': 21, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Clump Survival and Migration in VDI Galaxies: an Analytic Model versus Simulations and Observations

We address the nature of the giant clumps in high-z galaxies that undergo Violent Disc Instability, attempting to distinguish between long-lived clumps that migrate inward and short-lived clumps that disrupt by feedback. We study the evolution of clumps as they migrate through the disc using an analytic model tested by simulations and confront theory with CANDELS observations. The clump ``bathtub"model, which considers gas and stellar gain and loss, is characterized by four parameters: the accretion efficiency, the star-formation-rate (SFR) efficiency, and the outflow mass-loading factors for gas and stars. The relevant timescales are all comparable to the migration time, two-three orbital times. A clump differs from a galaxy by the internal dependence of the accretion rate on the varying clump mass. The analytic solution, involving exponential growing and decaying modes, reveals a main evolution phase during the migration, where the SFR and gas mass are constant and the stellar mass is rising linearly with time. This makes the inverse of the specific SFR an observable proxy for clump age. Later, the masses and SFR approach an exponential growth with a constant specific SFR, but this phase is hypothetical as the clump disappears in the galaxy center. The model matches simulations with different, moderate feedback, both in isolated and cosmological settings. The observed clumps agree with our theoretical predictions, indicating that the massive clumps are long-lived and migrating. A non-trivial challenge is to model feedback that is non-disruptive in massive clumps but suppresses SFR to match the galactic stellar-to-halo mass ratio.

The common wisdom is that most of these clumps are formed in-situ by gravitational disc instability (Toomre 1964). This has been studied in simulations of isolated discs (Noguchi 1999;Immeli et al. 2004b,a;Bournaud, Elmegreen & Elmegreen 2007;Genzel et al. 2008) and in cosmological simulations (Dekel, Sari & Ceverino 2009;Agertz, Teyssier & Moore 2009;Ceverino, Dekel & Bournaud 2010;Ceverino et al. 2012), in particular Mandelker et al. (2014, M14) and Mandelker et al. (2017, M17). With the high gas fraction at high redshift (e.g., Daddi et al. 2010;Tacconi et al. 2018), the characteristic Toomre clump mass is as large as a few percent of the disc mass, thus associated with "violent" disc instability (VDI), where the clumps play a significant dynamical role in the galaxy evolution (e.g., Dekel, Sari & Ceverino 2009;. It may be argued that the standard analysis of linear Toomre instability, assuming small perturbations in an otherwise uniform disc, may be invalid in its pure form in high-z cosmological discs, as they are continuously perturbed in a nonlinear manner and are constantly stimulated by external perturbations, such as mergers, tidal interactions and intense instreaming. Indeed, we find in cosmological simulations that Lagrangian proto-clump regions, in their linear regime, sometimes have high Toomre-Q values in the range 2 − 5, as opposed to the Q ∼ 1 expected in the standard Toomre analysis (Inoue et al. 2016, Mandelker et al., in prep.). This calls for a modified instability theory for VDI. Other works have called into question the validity of linear Toomre analysis in the context of high-z disc galaxies for other reasons (e.g., Behrendt, Burkert & Schartmann 2015;Tamburello et al. 2015;Mayer et al. 2016). Nevertheless, the successes of the standard Toomre analysis in explaining and predicting many features of high-z clumpy discs (see also Fisher et al. 2017) indicates that it can serve as a useful tool for obtaining a qualitative understanding. According to VDI, bound clumps that survive disruptive feedback are predicted to form with a Toomre mass and a high gas fraction preferably in the outer disc ( §6 below) and to migrate due to VDI-driven torques into the disc centre in a few disc orbital times, which is a few hundred Megayears at z ∼ 2 (Dekel, Sari & Ceverino 2009;Ceverino, Dekel & Bournaud 2010;Cacciato, Dekel & Genel 2012;Forbes, Krumholz & Burkert 2012;Forbes et al. 2014). We note, however, that this migration may become slower after a wet compaction event (Zolotov et al. 2015), when the disc becomes an extended ring whose migration rate is suppressed by a massive central bulge (Dekel et al. 2020b). This paper is partly motivated by an ongoing controversy concerning whether the massive clumps disrupt on a disc dynamical timescale due to stellar feedback or survive longer and possibly complete their migration to the central bulge in ten dynamical times or more. Dekel, Sari & Ceverino (2009), based on Dekel & Silk (1986), estimated that supernova feedback by itself may not have enough power to disrupt the massive clumps. Murray, Quataert & Thompson (2010) argued that momentum-driven radiative stellar feedback, enhanced by infra-red photon trapping, could disrupt the clumps on a dynamical timescale, as it does in the local giant molecular clouds. However, Krumholz & Dekel (2010) pointed out that such an explosive disruption would be possible in the high-redshift giant clumps only if the star-formation rate (SFR) efficiency in a free-fall time is ff ∼ 0.1, significantly larger than what is implied by the observed Kennicutt-Schmidt relation in different types of galaxies at different redshifts, namely ff of about one to a few percent (Tacconi et al. 2010;Krumholz, Dekel & McKee 2012;Freundlich et al. 2013). Dekel & Krumholz (2013) proposed instead that outflows from high-redshift clumps, with mass-loading factors η of order unity to a few, are driven by steady momentum-driven outflows from stars over many tens of free-fall times. Their analysis is based on the finding from high-resolution 2D simulations that radiation trapping is negligible because it destabilizes the wind (Krumholz & Thompson 2012. Each photon can therefore contribute to the wind momentum only once, so the radiative force is limited to ∼ L/c, where L is the luminosity of the clump. Combining radiation, protostellar plus main-sequence winds, and supernovae, Dekel & Krumholz (2013) estimated the total direct injection rate of momentum into the outflow to be ∼ 2.5L/c. The adiabatic phase of clustered supernovae and main-sequence winds may double this rate (Gentry et al. 2017). The predicted values are thus η ∼ 1 − 5, shown in Dekel & Krumholz (2013) to be consistent with the values deduced from the pioneering observations of outflows from giant clumps (Genzel et al. 2011;Newman et al. 2012). They concluded that most massive clumps are expected to complete their migration prior to gas depletion. With the additional gas accretion onto the clumps, they argued that the clumps are actually expected to grow in mass as they migrate inward.
Simulations that put in very strong winds by hand at the sub-grid level , or that include enhanced radiative feedback of 30 − 50L/c due to strong trapping (Hopkins, Quataert & Murray 2012;Oklopčić et al. 2017), indeed produced very short-lived clumps (SL clumps) that are typically unbound and disrupt on timescales in the ball park of the disc dynamical timescales, ten to a few tens of million years at z ∼ 2. On the other hand, simulations with more modest momentum driving of 2 − 3L/c Bournaud et al. 2014;Mandelker et al. 2017) reveal the formation of clumps of two types. In particular, the VELA-3 cosmological simulations that were analyzed in M14 and M17 reveal that a significant fraction of the giant clumps above ∼ 10 8.5 M are long-lived clumps (LL clumps) that are compact, round and bound and they survive feedback during their inward migration that lasts ∼ 300 Myr or more. This is while the less massive clumps tend to be SL clumps that are diffuse, elongated and unbound and they disrupt on one to a few disc dynamical timescale. These SL clumps may in certain ways be similar to the clumps detected in the simulations with stronger feedback, though some of the low-mass VELA-3 clumps do survive for a few dynamical times before they more gradually disrupt by the accumulated supernovae and possible tidal stripping. Cos-mological simulations of a similar nature to VELA-3, of the same suite of galaxies but with a stronger supernova feedback, turn out to produce more SL clumps at the expense of LL clumps even in the massive end (VELA-6, Ceverino et al. 2021, in preparation). In turn, isolated galaxy simulations have emphasized the dependence of clump formation, mass and longevity on the gas fraction in the disc where they form Renaud, Romeo & Agertz 2021). Thus, one may consider two extreme general types of giant clumps, one where the clumps are long-lived and manage to migrate through a fair fraction of the disc radius, and the other where the clumps form and die young not too far from where they formed. When they are very short lived, they may poetically resemble popcorn or the lights of a Xmass tree.
The main difference between the observable predictions of these extreme types of clumps is the distribution of stellar ages among the massive in-situ clumps, which are typically expected to be less than ∼ 50 Myr if the clumps disrupt on a disc dynamical timescale, while they are expected to last for up to a few hundred Megayears if the clumps survive till the end of migration. The distribution of radial distances may be another distinguishing feature, where it is expected to be broader for migrating climps. The time evolution of LL clumps is expected to be observable as correlations between clump properties and clump age. Due to the migration, this evolution could be translated to gradients of clump properties as a function of radial distance from the galaxy center. In contrast, the SL clumps are not expected to show a wide span of ages and migration-driven radial gradients. For both types of clumps, certain radial gradients might have been implemented at clump formation, reflecting possible radial gradients in disc conditions (see §7 below). The radial gradients in clump properties for LL and SL clumps, as well as for clumps that were formed ex-situ to the disc and for the background disc itself, were measured in cosmological simulations (M14,M17). They confirmed and partly quantified the above expectations for the LL clumps. These theoretical expectations for the two scenarios should be sharpened and physically understood better in order to allow a meaningful comparison with observations.
We focus here on a theoretical understanding of the evolution of clump properties during migration through the disc, as a function of the physical ingredients involved. For this we first follow the evolution of the giant clumps using a simple analytic model. This idealized model is in the general spirit of the bathtub model used extensively for the whole galaxy or the disc (Bouché et al. 2010;Davé, Finlator & Oppenheimer 2012;Cacciato, Dekel & Genel 2012;Lilly et al. 2013;Krumholz et al. 2018), except that it applies here to a single giant clump. The feature unique to clumps is the proportionality of the gas accretion rate onto a clump to its total clump mass, which can vary on a short timescale. For a whole galaxy, the gas (baryonic) accretion rate onto the galaxy can be crudely related to the gas (baryonic) mass, as the cosmological specific accretion rate of gas onto the central galaxy is expected to be roughly comparable to the total specific accretion rate, (e.g., . The accretion onto a clump, unlike a whole galaxy, is not directly constrained by the evolution of the cosmological accretion rate . During clump migration, clump gas turns into stars and the clump exchanges gas and stars with the disc through accretion, feedback-driven outflows and tidal stripping (Bournaud, Elmegreen & Elmegreen 2007;Bournaud et al. 2014). The evolution of clump properties is governed by the balance between these processes under conservation of gas and stellar mass. The evolution of clump properties is obtained here by analytically solving the two corresponding continuity equations, and this evolution in time can be translated to radial gradients within the disc. The analytic model allows a study of how the clump properties and evolution depend on the assumed efficiency of accretion, star formation and outflows. It also points at connections between the different clump properties, and suggests alternative observable proxies to quantities that suffer from large observational uncertainties, such as stellar age.
The model predictions are compared to the results of hydrodynamical AMR simulations, both of isolated gas-rich galaxies mimicking high-redshift discs at high resolution and of zoomed-in high-redshift galaxies in their cosmological setting, spanning a range of moderate feedback strengths.
Preliminary observational estimates of clump ages and radial gradients using a relatively small number of clumps Guo et al. 2012Guo et al. , 2015Zanella et al. 2019) seem to be consistent with the scenario of LL migrating clumps. The same is true for low-redshift gas-rich analogs of high-redshift discs (Fisher et al. 2019;Lenkić et al. 2021). However, the individual estimates of clump properties carry large uncertainties, especially the age estimates but also the mass and SFR estimates (e.g., Huertas-Company et al. 2020;Ginzburg et al. 2021). It is therefore important to use large samples of clumps, and to consider systematic errors in their properties (e.g., as done in Huertas-Company et al. 2020). A major step in the direction of a large sample has been made by Guo et al. (2018), who assembled and studied a very large sample consisting of thousands of clumps in massive galaxies from the CANDELS survey at z = 0.5 − 3. In a first test of our model against observations, attempting to distinguish between the competing scenarios for clump survival, we perform below a comparison of our theoretical predictions to the clump properties and their radial gradients in a large subsample of the observed CANDELS clumps, focusing on star-forming clumps at z = 1.5−3. In this comparison, we allow for systematic errors in the clump masses. Comparisons of our model to the larger sample of CANDELS clumps, with improved mass estimates via machine learning (Huertas-Company et al. 2020), are deferred to a future study. Using a machinelearning method for comparing theory to observations, The two types of clumps as identified in the VELA-3 sim-ulations, LL and SL clumps, were found in the CAN-DELS sample, with consitent relative distributions of clump masses, radial positions and host-galaxy masses (Ginzburg et al. 2021).
A potential practical difficulty in the comparison to observations may arise due to contamination of the observed clumps by stars from the underlying disc. This may smear out the differences in clump properties, such as their age, and possibly generate false radial gradients even in the SL scenario. If the density contrast in clumps is significant, the contamination effects are weak (M17). If, however, the contrast is low, these effects could be significant (Oklopčić et al. 2017). One should therefore seek distinguishing observable clump properties that are less sensitive to contamination by the disc, or that can be corrected for it, and attempt to estimate the possible effects of contamination both in the simulations and in the observations. As pointed out, the distribution of clump ages is the key distinguishing feature between the competing scenarios for clump survival. However, the observational age estimates are very crude and potentially biased, e.g., by an assumed SFR that exponentially decays in time in the SED fitting. We use our model and simulations to evaluate the actual star-formation history (SFH) in clumps, and attempt to derive an observable proxy for clump age based on the actual SFH. Motivated by the results of our toy model, we test the hypothesis that the SFR in a clump is rather constant with time, such that the inverse of the sSFR is a proxy for time since clump formation, capable of distinguishing between the competing scenarios. We also evaluate, using isolated-galaxy simulations with low gas fraction and strong feedback as well as short-lived clumps in the cosmological simulations, to what extent the inverse of the sSFR is a proxy for age also for SL clumps.
The paper is organized as follows. The following two sections present the bathtub toy model for clumps. In §2 we summarize and parameterize the processes of clump mass gain and loss, and in §3 we put the ingredients together in the continuity equations for gas and stars and solve them analytically for the evolution of clump properties. In passing, we propose a method for observationally estimating clump ages based on their sSFR. The subsequent two sections present the results from the simulations in comparison to the model predictions. In §4 we describe the isolated-galaxy and cosmological simulations used and the way we identify clumps in them and evaluate their properties. In §5 we present the simulation results for clump evolution and compare them to the toy-model predictions in order to test the validity of the model and interpret the simulation results. In §6 we summarize the observational results reported in Guo et al. (2018) from CANDELS and compare them to the theory predictions. In §7 we use the observations and theoretical predictions to constrain the nature of the giant clumps. In §8 we summarize our conclusions and discuss them.

ANALYTIC MODEL: MASS GAIN & LOSS
In order to spell out a bathtub model for clump evolution, we first address here each of the main processes of mass gain and loss in the gas and stellar components of the clumps, and derive the corresponding rates of mass exchange and the associated timescales. These ingredients will be inserted in the conservation equations for gas and stars in §3.

Clump migration
The giant clumps in VDI discs migrate toward the disc centre following angular-momentum and energy loss by several processes such as torques from the perturbed disc, clump-clump interactions and dynamical friction. The migration time can be estimated in several different ways (e.g., Dekel, Sari & Ceverino 2009;Krumholz & Burkert 2010;, all leading to estimates in the same ball park of ten to twenty disc dynamical times. Assuming that the clump maintains its initial Toomre mass, we adopt here the estimate based on clump encounters by (Dekel, Sari & Ceverino 2009), where t d = R d /V d is the disc crossing time, R d is the characteristic disc radius and V d is the characteristic disc circular velocity. The quantity δ is the mass fraction in the cold disc within the disc radius, which at the cosmological steady state has been argued to be δ ∼ 0.3 (Dekel, Sari & Ceverino 2009). The Toomre parameter is Q ∼ 1 1 . The migration time is thus comparable to 2−3 orbital times at the outer disc, t orb = 2πt d . With a higher value of Q (e.g., Inoue et al. 2016), the migration time would be somewhat longer. On the other hand, if the clump changes its mass during migration, the migration time scales inversely proportional to the clump mass M c (t) [as in dynamical friction, or in clump encounters, equations 16-18 of Dekel, Sari & Ceverino (2009) where M T is the initial mass of the clump, assumed to be the Toomre mass. A corrected estimate for t mig can be obtained once M c (t) is known. The inward radial migration velocity when the clump mass is M c can be crudely estimated from eq. (1) and eq. (2) as Assuming that the clump forms at r ∼ R d (leading to an upper limit for the migration time), one can integrate V mig in time to estimate the clump position r(t), and evaluate t mig when r(t mig ) = 0. We find that the change 1 In a thick disc the Toomre Q is expected to be slightly smaller than unity (Goldreich & Lynden-Bell 1965), while in cosmological simulations we sometimes find Q values that are somewhat larger than unity (Inoue et al. 2016). The migration time in our simulations therefore ranges from ten to a few tens of dynamical times.
in total mass during migration is typically not large, so the initial estimate of t mig , which is crude anyway, is a fair approximation for our purpose here. However, eq.
(3) turns out not to be a good approximation for translating the time dependence of the clump mass to radial gradients within the disc.

Gas accretion onto clumps
As the clump spirals in toward the disc centre, it accretes matter from the surrounding disc. An estimate of the accretion rate is provided by the entry rate into the tidal (Hill) sphere of the clump in the galaxy, R T , Here ρ d is the density in the cold disc (gas or young stars), πR 2 T is the cross section for entry into the tidal sphere, and σ d is the velocity dispersion in the disc representing here the relative velocity of the clump with respect to the rest of the ring within which the clump is orbiting. We crudely assume here that ρ d and σ d as well as R T are cosntant within the disc throughout the clump migration. The parameter α represents the fraction of the mass entering the tidal radius that is actually bound to the clump, and is crudely expected to be on the order of 0.2-0.4 (Dekel & Krumholz 2013;Mandelker et al. 2014). If only gas that is inflowing with respect to the clump center becomes bound to the clump, the value of α would be on the low side. We assume here that α is constant throughout the clump lifetime, and test this assumption using simulations in §5.
The tidal or Hill radius R T about the clump is where the self-gravity force by the clump balances the tidal force exerted by the total mass distribution in the galaxy along the galactic radial direction. If the disc is in marginal Toomre instability with Q ∼ 1, this is the same as the Toomre radius of the proto-clump patch that eventually contracts to form the clump (Dekel, Sari & Ceverino 2009), where the clump mass is given by with M d referring to the mass of the cold disc. Also when Q ∼ 1, the disc half-thickness H d is comparable to and with the numerical factor corresponding to a flat rotation curve. Note that eq. (4) is valid only when R T is not larger than H d , such that the clump accretes from a 3D distribution around it.
We can now evaluate the rate of clump growth by accretion using eq. (4). We insert R T from eq. (5), write , and use eq. (7) for H d to obtain the desired expression for the accretion rate, The timescale for accretion is thus, This turns out to be about half the migration time of eq. (1), already indicating that the clump mass may grow by a factor of a few during the migration.

Star Formation Rate
For a clump of gas mass M g the SFR can be modeled as where t ff is the clump free-fall time and ff the corresponding SFR efficiency. The free-fall time is where M c and R c are the clump mass and radius. The SFR efficiency per free-fall time has to be of order ff ∼ 0.02 in order to match the observed Kennicutt-Schmidt law at different environments and redshifts (e.g., Krumholz, Dekel & McKee 2012). We thus assume ff to be a constant during the clump lifetime. The effective value of ff , as defined for the clump mass and free-fall time, may be larger if star formation actually occurs in dense sub-regions within the clump where the free-fall timescale is shorter. For our purpose here, it is convenient to replace ff by d , where t ff in the clump is replaced by the disc t d , as in eq. (11), and where we assume a constant density contrast between the clump interior and the disc over all. The constancy of ff , and therefore of d , is expected in the Toomre regime of star formation, that is relevant for VDI highz discs (Krumholz, Dekel & McKee 2012), so we adopt it here. The ratio t d /t ff is given by the square root of the density contrast between the star-forming region in the clump and the disc, namely t d /t ff ∼ 3 − 10 if the density contrast is 10−100 (e.g., Ceverino et al. 2012). We thus expect d ∼ 0.1 or higher. We assume that a fraction µ of the mass in forming stars remains in stars, the rest assumed to be instantaneously lost from the stars due to supernovae and stellar winds and deposited back in the inter-stellar medium within the clump, namely in M g . The value of µ is estimated to be between 0.5 and 0.8, depending on the duration over which it is computed . For the relevant timescales for clump evolution, on the order of ∼ 100 Myr, we adopt hereafter µ = 0.8.
The two observable timescales of interest are the depletion time and star-formation time. The depletion time is where the actual timescale for gas consumption by star formation is only slightly different, µ −1 t dep . The starformation time is defined as the inverse of the specific SFR (sSFR), where f sg is the stellar-to-gas mass ratio, related to the conventional gas fraction f g by

Gas Outflow
The star formation in the clump produces gas outflow from the clump via various stellar feedback processes, including supernovae feedback and radiative pressure feedback (e.g., Dekel & Krumholz 2013). One can approximate the outflow rate aṡ where η is the mass-loading factor. We assume here that when averaged over a dynamical time (or more) η is a constant throughout the clump lifetime. It is estimated to be of order unity, or a few, both from observations (e.g., Genzel et al. 2011;Newman et al. 2012) and based on theoretical considerations (e.g., Dekel & Krumholz 2013). The timescale for gas outflow is thus For η ∼ 1, this is comparable to the depletion timescale, which for d ∼ 0.1 is comparable to the migration timescale.

Stellar mass exchange
Accretion of stars from the disc into clumps is expected to be less efficient than the gas accretion because even when they enter the Hill sphere they are less likely to become bound to the clump, as they tend to have a larger velocity dispersion and they do not dissipate. On the other hand, tidal stripping of older stars from the clump is expected to be more efficient than gas stripping because these stars tend to be on less bound orbits. We may therefore expect a net effect of stellar mass loss, that may partly compensate for the stellar mass gain by star formation. The rate of tidal stripping after clump formation is not expected to be large, and is likely to be smaller than both the rate of gas accretion and the rate of gas outflow driven by feedback. For a net stellar mass loss from the clumpṀ s,loss , stripping minus accretion, we define Negative values would refer to a net accretion rather than stripping. Since η s is expected to be relatively unimportant, we allow it to be constant during the clump migration, and use its value only for estimating the qualitative effect of stellar exchange between the clump and the disc. 2

ANALYTIC CLUMP BATHTUB MODEL
We note that the timescales for changes in clump properties as evaluated in §2 (t ac , t sf , t dep , t out ), with the fiducial values of the model parameters (α ∼ 0.2, d ∼ 0.1, and η ∼ 1, with f sg ∼ 1), are all in the ball park of ∼ 10t d , which is on the order of the migration time t mig . This implies that the clump properties are expected to evolve during the migration, significantly but not by orders of magnitude. The ingredients discussed in §2 are now put together in mass conservation equations for each clump, to be solved analytically and to be compared to simulations ( §5) and observations ( §6).

Gas mass conservation
The gas mass in the clump varies due to the net effect of growth by accretion and decrease by gas consumption into stars and gas loss by outflows, Using the expressions of the previous section, this be-comesṀ Eq. (20) could be re-written in short aṡ (21) However, the timescale t g = g −1 is not necessarily constant, as f g may evolve in time even if α, η and d are roughly constant. In fact, g(t) may change sign as f g evolves, so the associated timescale g −1 may diverge at this moment. Thus, unfortunately, M g (t) cannot be obtained by a simple integration of eq. (21) with a constant g, and one has to work harder on solving eq. (20).
Eq. (20) is reminiscent of the analogous continuity equation for the gas mass in the whole galaxy, with source and draining terms on the right hand side (eq. 1 of . For a whole galaxy, while the draining term is similar, the source term is external, driven by the cosmological gas accretion rate and its penetration into the central galaxy. In the case of a constant source term, the draining term being proportional to negative M g drives the system to a quasi-steady state, where M g is constant on a timescale shorter than the long timescales for variations in the accretion time and in the disc dynamical time, and therefore the SFR follows the accretion rate, as it varies slowly with cosmological time. For a clump, however, the source term gas outflow mass-loading factorṀout = η SFR eq. 16 1 ηs stellar exchange rate factorṀ s,loss = ηs SFR eq. 18 1 µ fraction of formed stellar mass left in stars Model solution t 1 , t 2 timescales of growing and decaying modes depends on the total clump mass, which involves both M g and M s and could vary on a short timescale. This leads to a different solution. In order to solve eq. (20), it is convenient to separate M c into M s +M g and rearrange eq. (20) into terms that are proportional separately to M s and M g , The first term is proportional to M s , which in general, like M g , varies in time. This implies that in order to follow the evolution of M g one has to also take into account the variation in stellar mass. This reflects the main difference between the bathtub model for a clump and for a whole galaxy.

Stellar mass conservation
The stellar mass evolves at a ratė A slight overestimate of the clump stellar mass can be obtained for the simple case of η s = 0, where the stellar mass grows primarily by star formation, ignoring stellar mass exchange with the disc. An underestimate of the clump mass would be obtained with η s = µ, where the stellar mass loss exactly compensates for the star formation, such that the stellar mass remains constant in time.

Analytic solution for clump evolution
Eq. (22) and eq. (23) can be solved analytically for the evolution of M g (t), M s (t) and the associated quantities of interest. We assume that the model parameters are constant during the migration of the clump, including α, d , η and η s (µ = 0.8 is naturally assumed throughout), as well as t d of the disc. We set the mass units such that the initial total clump mass at t = 0, assumed to be the Toomre mass, is M ci = 1. The initial condition is characterized by the initial gas fraction f gi , giving rise to the initial gas mass and stellar mass M gi and M si . Taking a time derivative of eq. (22) we obtain We try a solution of the sort Inserting each term of eq. (28) in eq. (25), we obtain a quadratic equation for each β, independent of the corresponding m, The two roots (which always exist as c < 0) are These represent the two characteristic timescales The large root, β 1 = β 2 + d, is commonly positive, representing an exponentially growing mode that dominates at t > t 1 . The small root, β 2 , could be positive as well, such that t 1 < t 2 . Alternatively it could be negative, thus representing an exponentially decaying mode, with t 2 = |β 2 | −1 either smaller or larger than t 1 , which may in principle dominate at t < t 1 and t < t 2 .
The coefficients m 1 and m 2 are determined by the initial condition M gi , via The first equality in eq. (33) is based on the definition of g in eq. (21), and the second equality is based on a time derivative of M g (t) from eq. (28). Solving eq. (32) and eq. (33) yields The dependence of these coefficients on the initial condition is both through the explicit M gi and through g i . The stellar mass is obtained by a time integration of the SFR, as given in eq. (23), namely integrating M g (t) from eq. (28). We obtain Eq. (28) and eq. (35) provide the evolution of the stellarto-gas mass ratio, (36) The conventional gas fraction is given by f g (t) = [1 + f sg (t)] −1 , and the timescale for star formation is t sf (t) = d −1 f sg (t)t d from eq. (14). Figure 1 shows the time evolution of the different clump quantities of interest for a given choice of the model parameters. Time is plotted with respect to 10t d , assumed to represent t mig . Mass is in units of the initial clump mass M ci . Shown are M g (blue), M s (red), f g (green), t sf (black) and t dep (magenta). The parameters chosen for this figure are α = 0.2, d = 0.1, η = 1 and η s = 0. The two panels show the solutions for two different initial conditions, f gi = 1 (left) and f gi = 0.7 with an initial stellar component (right). The exponents in this case, independent of f gi , are (β 1 , β 2 ) = (0.058, −0.138) in units of t −1 d , thus representing a growing and a decaying mode. The corresponding times t = |β| −1 , in units of t d , are (t 1 , t 2 ) = (17.3, 7.25), namely both on the order of t mig and with t 1 somewhat larger. We identify three stages of evolution, as follows.

Exponential growth at late times
At t > t 1 (and t > t 2 ) the system approaches a quasi steady state. In this stage the growing mode dominates, so both M g (and the SFR) and M s grow exponentially, but their ratio converges to an asymptotic value, This makes t sf saturate at a maximum level of For the choice of parameters in Fig. 1, the asymptotic values are f sg = 1.38, corresponding to f g = 0.42, and t sf = 13.8t d . The exponential growth of M g (and therefore M s ) in this regime is because the accretion of gas from the disc into the clump is proportional to the total clump mass itself, which is driven by the growing stellar mass that has become comparable to and slightly larger than the gas mass. This results in a quasi steady state with a coherent exponential growth of both M g and M s at a fixed f g .
Another way to understand this asymptotic exponential growth regime is by noting that with a constant f g the time derivative of g(t) in eq. (21) vanishes, so g is a constant andṀ g ∝ M g , yielding an exponential growth M g ∝ exp(gt).
We note, however, that if the migration is as rapid as estimated in eq. (1), this asymptotic exponential phase is largely hypothetical as it is typically supposed to occur after the clump has completed its migration and merged into the bulge.

The main stage of migration: a constant SFR
Prior to t 1 , at t < t 1 (and t < t 2 ), the migrating clump is at it's main stage of evolution. As t → 0, the gas mass approaches a constant M gi , and the stellar mass approaches M si accordingly, In this whole regime prior to t 1 , M g is roughly constant, therefore the SFR is roughly constant, so the growing term of M s is proportional to t. From eq. (14) and eq. (36), we obtain during this period For a small enough M si (f gi ∼ 1 or f sgi 1), this expression simplifies to and 0.7 (right). Shown are clump gas mass Mg (blue), stellar mass Ms (red), gas fraction fg (green), depletion time t dep (magenta), star-formation time t sf =sSFR −1 (black) and its corrected version t sfc for f gi < 1, eq. (45) (dashed black), which differs from t sf only when f gi < 1. We see the phase during migration where the SFR, proportional to Mg, is roughly constant, such that t sfc is a proxy for the clump lifetime (dotted black). At late times the clump converges to a quasi-steady state with an exponential mass growth and a constant t sf , but this is supposed to 'occur' only after the migration is complete. For f gi < 1 there is an early phase dominated by the initial gas fraction, where t sf is an overestimate of the clump time, but this is properly corrected for by t sfc , which is a close lower bound for the clump time.
This linear growth with t is seen in Fig. 1 in the case f gi = 1 throughout the whole stage t < t 1 .
For the above choice of parameters and f gi = 1, the coefficients for M g (t) are (m 1 , m 2 ) (0.3, 0.7) in units of M gi . The fact that in this case m 2 > m 1 > 0 means that the decaying mode is dominant at early times, t < t 2 , such that M g slowly declines with time until the growing mode takes over after t 1 . However, the gas mass is roughly constant until t ∼ t 1 , as the growing and decaying modes roughly balance each other. For f gi ∼ 1, the minimum of M g (t), obtained between t 2 and t 1 , is at a mass value not much smaller than the initial M gi ∼ 1. It turns out that for all values of f gi the curves for M g (t) go through a fixed point, near t ∼ 10t d (t 9t d and M g 0.7 for the above choice of parameter values).
We note that for this choice of parameters, the total clump mass M c is roughly constant until the end of the migration at t ∼ 10t d . This is because the clump mass is dominated by the gas mass that is roughly constant. The constant M c makes the accretion rate into the clumps roughly constant, so the solution tends toward a constant M g , with the input balanced by the output, similar to the quasi-steady-state solution for a whole galaxy under a constant cosmological accretion rate ).

Effect of initial stellar mass at early times
With an initial stellar component such that f gi < 1, there is a third, early stage, where the star-formation time approaches an asymptotic minimal value valid at t < t i , defined by This follows from eq. (36), with |β 1 t| 1 and |β 2 t| 1 and demanding that the second term in the numerator is smaller than the first. As long as t i is significantly smaller than t mig , t 1 and t 2 , there is a period between these two times where the evolution is of the sort discussed in §3.4.2, namely a roughly constant SFR and t sf t. For example, if min{t 1 , t mig } ∼ 10t d , this requires f sgi 10 d (µ − η s ) ∼ 1, namely an initial gas fraction close to unity.
However, once M si is not negligible, the early evolution of t sf deviates from the pure linear t sf ∝ t. As can be seen in Fig. 1, already for f gi = 0.7, the growth of t sf (t) throughout the period (t d , t mig ) is significantly shallower than a linear growth, with t sf a severe overestimate of t at early times.
If one wishes to use t sf as an observable proxy for the clump age t, and make sure that it is not a severe overestimate of the age at early times, one could define a corrected estimator for the star-formation time by As seen in the right panel of Fig. 1, for f gi = 0.7, t sfc is indeed growing roughly linearly with time prior to ∼ 10t d , where it deviates from t by only ∼ 20%. The initial stellar mass is M si = (1 − f gi )M ci , so it can be estimated for a given initial (Toomre) clump mass and an initial gas fraction that can be estimated from the gas fraction in the cold component of the clump-forming disc.

Dependence on model parameters
One can show that for the parameters within a factor of ∼ 2 about their fiducial values η ∼ 1, α ∼ 0.2 and d ∼ 0.1, the value of t 1 is robustly of order t 1 > ∼ 10t d , quite insensitive to the actual values of the parameters. If the estimated rapid migration of eq. (1) is valid, this implies that t 1 ∼ t mig , so the asymptotic exponential growth phase occurs only after the end of migration. Thus, the main stage where t sf ∼ t extends all the way till the end of migration. This is provided that the initial clump is gas rich such that the onset of this regime, at t i , is long before t 1 . The corrected t sfc is a good proxy for t even prior to t i .
Short-lived clumps may be modeled by a strong feedback. When η is larger than unity, say η ∼ 3, while β 1 remains similar to the case of η ∼ 1 discussed above, the value of |β 2 | becomes significantly larger than β 1 , so the gas mass declines significantly during t < t 1 . Still, t sf ∝ t in this regime, though it is an overestimate of t if f gi < 1. However, the corrected t sfc is an excellent approximation to the clump time for any value of f gi , especially during the first few dynamical times after which the short-lived clump is expected to be disrupted.
We note in passing that in some local Giant Molecular Clouds (GMC) there are certain indications for an early exponential growth while the clump is still gas dominated (Palla & Stahler 1999). Our toy model would allow an early exponential growth when t 1 is sufficiently small compared to 10t d , or when the migration is slower. This requires significant deviations from the fiducial values of the model parameters assumed above for the high-z giant clumps in their main phase of evolution, in particular very efficient accretion with α closer to unity, and/or negligible outflows with η 1, as well as a low SFR efficiency d < 0.1. This may be valid in the early, collapsing phase of the GMC. The exponential phase may alternatively be reached if the migration is significantly slower than assumed here. This would be the case if the disc becomes a ring about a large central mass, which suppresses the ring's inward mass transport rate (Dekel et al. 2020b). Such conditions are expected after a generic wet-compaction event that typically occurs when the galaxy becomes larger than the golden mass of M s,gal ∼ 10 10 M . We use here simulations similar to those described in Bournaud et al. (2014, e.g., model G2) and in Perret et al. (2014). The RAMSES AMR code (Teyssier 2002) is utilized with a maximum resolution of 3.5 pc and an SFR efficiency of ff = 0.02 per free-fall time on the local cell size. Six gas-rich disc galaxies were simulated in isolation, mimicking the evolution of z ∼ 2 discs for a duration of ∼ 1 Gyr. The galaxy baryonic masses are in the range 4×10 10 −4×10 11 M and the initial gas fraction is 0.6. The typical disc crossing time at the radius where most clumps form is t d ∼ 24 ± 3 Myr.
These simulations include subgrid models for supernova feedback, and four of them also incorporate models for stellar feedback by photo-ionization and by radiation pressure (RP), as described in Renaud et al. (2013). The radiative pressure is assuming a momentum supply to the outflow at a rate of ∼ 2.5L/c for a luminosity L, namely the photons are weakly trapped with each photon contributing 2.5 times. This radiative contribution is comparable to the contribution by supernova feedback and is comparable to the momentum driving by all stellar sources as predicted by Dekel & Krumholz (2013). A few clumps from these simulations were analyzed in Bournaud et al. (2014), and the relevant results are summarized in figures 6 and 8 there. The simulations with and without radiationpressure feedback, both including supernova feedback, are termed hereafter RP and NoRP respectively. Figure 3 shows face on views of the surface cold-gas density in NoRP and RP simulated disks during their main clumpy phase. They both show giant star-forming clumps with masses M c > 10 8 M , tending toward larger masses when the RP is not implemented.

Isolated clump analysis
We briefly describe here our method for detecting clumps in the isolated galaxy simulations, and refer to Bournaud et al. (2014) for further details. The clumps are first detected visually in 2D face-on projections of the baryonic mass. The center of each clump is then determined in 3D through an iterative procedure as follows. An initial guess for the clump center is obtained visually, assuming that it lies in the disk mid-plane. The center of mass of the baryons within a sphere of radius 200 pc about this initial center is determined, and the procedure is repeated about this new center until the correction is smaller than the resolution of the simulation. The clump radius and mass about this centre is determined as follows. Using a radius R 1 , initially set to 100 pc, we compute the mass in a sphere of radius R 1 and in a sphere with a volume twice as large (i.e. a radius 2 1/3 R 1 ). If the mass increase is larger than 30%, we iterate the process with a larger radius R 1 +100 pc until the mass increase is by less than 30%. In practice, the clump radii range from 200 pc to 500 pc. These positions and radii are used to measure gas and stellar masses, mass flows, and SFR in the clumps. Given the method used, clump masses are determined with an uncertainty better than 30%. Background subtraction was performed by subtracting the average density near a spherical contour around the clump at twice the clump radius. Another technique, based on the average density at the same galacto-centric radius as the clump center but outside of the clump itself, gave very similar results. This background subtraction may be an over correction, because some of the subtracted background may belong to the clump (see §A).
Here we analyze all the giant clumps with total masses M c 10 8 M , which amount to 17 clumps from the two NoRP galaxies and 31 clumps form the four galaxies with RP.
The SFR in each clump is directly measured from the simulations. The gas outflow rate from a clump, being mostly in the polar direction perpendicular to the disc, is measured through parallel disc surfaces of radius 1 kpc at a height ±1 kpc above and below the clump center and the galactic disc plane. This turned out to be a less noisy measurement than the outflow through the spherical boundary of the clump, which suffers from diffuse disk gas entering and leaving the clump ). The former method was found to lead on average to only a small, 23 percent underestimate of the outflow. The gas accretion rate onto the clump, allowed to come from all directions, is then evaluated through mass conservation, eq. (19) (with µ = 0.8), where the net rate of change of M g is measured from M g in two successive output times separated by ∼ 10 Myr. The net stellar mass-loss rate is measured by identifying stellar particles moving in and out from the clump.

The simulations
We use here two zoom-in cosmological simulations from the VELA suite, out of 34 galaxies with halo masses 10 11 − 10 12 M at z ∼ 2, which have been used to explore many aspects of galaxy formation at high redshift (e.g., Ceverino et al. 2014;Moody et al. 2014;Zolotov et al. 2015;Ceverino, Primack & Dekel 2015;Inoue et al. 2016;Tacchella et al. 2016a,b;Tomassetti et al. 2016;Ceverino et al. 2016b,a;Mandelker et al. 2017;Dekel et al. 2020a,b). The simulations utilize the ART code (Kravtsov, Klypin & Khokhlov 1997;Kravtsov 2003;Ceverino & Klypin 2009), which follows the evolution of a gravitating N -body system and the Eulerian gas dynamics with an AMR maximum resolution of 17.5 − 35 pc in physical units at all times. The darkmatter particle mass is 8.3 × 10 4 M and the minimum mass of stellar particles is 10 3 M . The code incorporates gas and metal cooling, UV-background photoionization and self-shielding in dense gas, stochastic star formation, stellar winds and metal enrichment, thermal feedback from supernovae (Ceverino, Dekel & Bournaud 2010), and feedback from radiation pressure. In the current implementation of radiation pressure feedback, referred to as model "RadPre" in Ceverino et al. (2014), "RP" in M17, and VELA-3 in other places, the momentum driving efficiency is ∼ 3L/c. This is comparable to the theoretical estimates and to the radiative pressure contribution in the RP isolated simulations (which is comparable to the supernova contribution), but is on the low side compared to certain other strongfeedback simulations (e.g., Genel et al. 2012;Oklopčić et al. 2017). Further details regarding the feedback and other sub-grid models can be found in Ceverino et al. (2014) and M17.
M17 have analyzed the clumps in the VELA-3 galaxies as well as in counterparts that lacked the radiation pressure feedback (referred to as RP and NoRP, or VELA-3 and VELA-2 , respectively). They found that with RP a significant fraction of the clumps of M c < 2 × 10 8 M are disrupted on timescales in the ball park of one to a few disc dynamical times, while most of the more massive clumps survive and migrate to the disc centre. The migration inward of these longlived clumps produced galacto-centric radial gradients in their properties, such as mass, stellar age, gas fraction and sSFR, which distinguished them from the short-lived clumps. However, in that study the output timesteps were ∼ 100 Myr, which prevented a detailed study of the individual clump evolution. In this study, we have re-simulated two of the VELA galaxies with RP (V07 and V19, see table 1 in M17), saving about 10 snapshots per disc crossing time. Both V07 and V19 have extended clumpy phases of VDI, in the redshift ranges z = 2.5 − 1 and 5 − 3 respectively, where massive, long-lived, star-forming clumps are continuously formed.
As in M17, the disc is defined as a cylinder with radius R d and half height H d , containing 85% of the cold gas (T < 1.5 × 10 4 K) and the young stars (ages < 100 Myr) within 0.15R v . Thus, if the disc surfacedensity profile is exponential, R d corresponds to twice the exponential scale radius. The disc stars are also required to obey a kinematic criterion, that their specific angular momentum parallel to that of the disc, j z , is at least 70% of the maximal possible value it could have given its galacto-centric distance, r, and orbital velocity, v, namely j max = vr.
The disc crossing time in the main body of the disc is computed as t d = R rot /V rot , where R rot and V rot are the mass-weighted radius and azimuthal velocity in the cylindrical ring of height ±0.5 kpc and radii (0.5−1)R d . For an exponential disc R rot 0.78R d . The evolution of t d in the two simulations in the relevant redshift ranges is shown Fig. 2. In V07, the main VDI phase occurs at z ∼ 2.5−1, following a dramatic phase of wet compaction into a star-forming "blue nugget" (Zolotov et al. 2015, figure  2). During this time, the disc radius steadily increases from R d ∼ 9 to 20 kpc while H d /R d ∼ 0.1−0.2. The disc crossing time increases from t d ∼ 32 to 57 Myr, the baryonic disc mass increases from M d ∼ 2.7 to 6.5×10 10 M , and the gas fraction in the disc decreases from f g ∼ 0.28 to 0.13. The "cold" mass fraction in the disc, referring to cold gas and young stars, declines from f cold ∼ 0.4 to 0.14.
In V19, the VDI phase occurs at higher redshifts, z ∼ 5 − 3. During this phase the disc radius is R d ∼ 3−4 kpc, H d /R d ∼ 0.1−0.2 and t d ∼ 10 Myr. The disc mass increases from M d ∼ 1.0 to 1.5×10 10 M , the gas fraction in the disc decreases from f g ∼ 0.23 to 0.07, and the cold mass fraction is decreases from f cold ∼ 0.45 to 0.18. V19 thus represents a galaxy that is quite different from V07; being at higher redshift it is more compact and with a shorter dynamical time, where the specific accretion rate, sSFR and gas fraction are significantly higher. These two galaxies allow us to explore cosmological galaxies at extremely different conditions, and hopefully identify robust features in the evolution of their clumps that are captured by our toy model.
The gas fractions in the discs of these two simulations are somewhat lower than estimated in typical observed galaxies at similar redshifts. This has been discussed in detail in several papers which used the VELA suite (e.g., Zolotov et al. 2015;Tacchella et al. 2016b;Mandelker et al. 2017). As highlighted in those papers, while they are state-of-the-art in terms of high-resolution AMR hydrodynamics and the treatment of key physical processes at the subgrid level, the VELA simulations are not perfect in terms of their treatment of star-formation and feedback, much like other simulations. Star-formation tends to occur too early, leading to lower gas fractions later on. The stellar masses at z ∼ 2 are a factor of ∼ 1.5−2 higher than inferred for similar mass halos from abundance matching (e.g., Rodríguez-Puebla et al. 2017;Moster, Naab & White 2018;Behroozi et al. 2019). However, for the purposes of the present study, the relatively low gas fractions during the peak VDI phase would only underestimate the actual accretion of fresh gas onto clumps during their migration, providing a lower limit on clump survival. The effect of gas fraction on clump properties and survival in simulated isolated discs is further discussed in . Regardless, the toy model can still be tested with the relevant parameters. Figure 3 shows face on views of V07 at z = 2.33 and V19 at z = 4.56, shortly after the start of their VDI phase. The figure shows the surface density of the cold mass, integrated over ±R d perpendicular to the disc, where R d 12 and 3.2 kpc for V07 and V19 respectively. While these two discs are at very different redshifts, with different masses and sizes, both show giant starforming clumps with masses M c > 10 8 M .

Cosmological clump analysis
Clumps are identified in 3D and followed through time following the method detailed in M17. Here we briefly summarize the main features.
We search for clumps within a cube of side L = 3R d centered on the main galaxy centre. Using a cloud-incell interpolation, we deposit the mass in a uniform grid with a cell size of ∆ = 70 pc, 2 − 4 times the maximal AMR resolution. We then smooth the density in the grid cells, ρ, into a smoothed density ρ W , using a spherical Gaussian filter of FWHM of min(2.5 kpc, 0.5R d ).
The density fluctuation at each grid point is δ ρ = (ρ−ρ W )/ρ W . This is performed separately for the cold mass and the stellar mass, and at each point we adopt the maximum of the two residual values. We define clumps as connected regions containing at least 8 grid cells with a density fluctuation above δ min ρ = 10. No attempt has been made to remove unbound mass from the clump. We define the clump centre as the baryonic density peak and the clump radius, R c , as the radius of a sphere with the same volume as the clump. Ex situ clumps, which joined the disc as minor mergers, are identified by their dark matter content and the birth place of their stellar particles and are not considered further here.
The SFR in the clumps is calculated using the mass in stars younger than 30 Myr, which is long enough for good statistics and short enough for stellar mass loss to be negligible. Outflow rates from the clumps are measured using each of three spherical shells of width ∆r = 140 pc centered at radii r = R c , 1.5R c and 2R c . The gas outflow rate isṀ = ∆ −1 i V r ρ, where the sum is over cells within the shell with V r > 0 and a 3D velocity larger than the escape velocity from the shell, V 2 > 2GM c /r. We adopt the average of the outflow rates in the three shells, which may vary by ∼ 0.1 dex. This is different from how we measure outflows in the isolated galaxy simulations, through parallel discs. This is because the cosmological disc planes are thicker and highly perturbed and the outflows through slabs are noisier due to ongoing accretion and mergers. We verified that the two methods do not lead to systematically different results. The gas accretion rates are then computed as in the isolated galaxies, namely from the gas continuity equation, eq. (19) (using µ = 0.8), after obtaining the net rate of change from the clump gas mass  ). Both show giant clumps, indicating less massive clumps when RP is incorporated. Bottom: Cosmological simulations, V07 at z = 2.33 (left) and V19 at z = 4.56 (right) from the VELA-3 suite. Shown here is the surface density of the "cold" mass, comprised of cold gas (T < 1.5×10 4 K) and young stars (< 100 Myr). The disc radii are R d = 11.8 and 3.2 kpc respectively. In both cases the disc half-thickness is H d ∼ R d /6. The baryonic masses of the discs are M d /(10 10 M ) ∼ 3.50 and 1.35, and the fractions of cold-mass are 28% and 33% respectively. Both discs show giant clumps with masses Mc > 10 8 M . The gas density in and around the ∼ 1 kpc vicinity of the clumps is similar in the isolated simulations and the cosmological simulations at z ∼ 2, while the background density is higher in the cosmological simulations where gas is continuously accreted.
in successive output times. Using instead baryonic-mass conservation yields nearly identical results. Stellar mass loss or gain is measured by identifying stellar particles moving in and out from the clump, as for the isolated galaxy simulations.
Individual clumps, that contain at least 10 stellar particles, are traced through time based on their stellar particles. For each such clump at a given snapshot, we search for all "progenitor clumps" in the preceding snapshot, defined as clumps that contributed at least 25% of their stellar particles to the current clump. If a given clump has more than one progenitor, we consider the most massive one as the main progenitor and the others as having merged, thus creating a clump merger tree. If a clump in snapshot i has no progenitors in snapshot i−1, we search the previous snapshots back to two disc crossing times before snapshot i. If no progenitor is found in this period, snapshot i is declared the initial, formation time of the clump, and for that clump t is set to zero at that time. 3 We select for analysis here clumps that formed at z = 2.5−1 and 5−3 for V07 and V19 respectively. Based on the disc crossing times, the expected migration time for clumps, which is about 10−20 dynamical times, is in the range 400−1000 Myr for V07 and only ∼ 150 Myr for V19. In both cases, we find this to be a good match to the actual migration times of massive clumps. We only consider for analysis here clumps that survive for longer than 50 Myr. Similar to our analysis of the isolated simulations, we also limit our sample of clumps based on their masses at 50 Myr (averaged from 40 to 60 Myr), to be more massive than 10 8 M . In order to relate to the observed clumps that are detected in UV (Guo et al. 2018), we also apply a lower SFR threshold of 0.1M yr −1 at t = 200 Myr. This turns out to have no effect in V19, but it removes about half the clumps that meet the mass and lifetime criteria in V07, all in the mass range 10 8 −10 8.5 M (see also M17, Figure 10). We note that these passive clumps that were removed, most of which have formed ex-situ, tend to move on perturbed orbits that deviate from the main disc, such that the gas density along their orbits fluctuates, invalidating the toy model assumption of a constant density. The orbits of the clumps selected by the criteria of mass and SFR above are all confined to the disc such that the accretion model is valid. Our final clump sample contains 37 clumps from V07 and 12 clumps from V19. We note that we have not applied here an explicit selection based on the M17 characterization of the clump as LL or SL. While most of the clumps selected here are naturally LL in-situ clumps, some of them may be SL clumps, though the 50 Myr age threshold implies that they live more than one-two disc dynamical times (SL clumps a la M17 can live for up to 20 clump free-fall times, which could reach a few disc dynamical times).
In summary, our sample consists of massive, starforming, bound and mostly long-lived clumps from four types of simulations. First, two RAMSES isolated-disc simulations, NoRP and RP, with 17 and 31 clumps from two and four galaxies respectively. Second, two ART-VELA cosmological simulations (with RP), V07 in the redshift range z = 2.5 − 1.0 and V19 at z = 5 − 3, with 37 and 12 clumps respectively. In all these simulations the feedback is not overwhelmingly strong, such that the massive clumps above 10 8.5 M tend to survive feedback during their migration.

Model validity and parameters
Figures 4 to 7 display the evolution of clump properties in the four types of simulations during the clump lifetimes, starting from their formation time and evolving for tens of disc dynamical times. The typical clumps form in the outer disc, migrate inwards (see, e.g., Fig. 9), and merge into the bulge at the galaxy center after t mig ∼ (10 − 20)t d , as predicted in eq. (1). A fraction of the clumps live longer, e.g., because they have low masses compared to the Toomre mass, or because they form at very large radii, or for other reasons. Shown in the figure are the medians over the clumps of the given simulation type, and the 68% scatter about them at time t since clump formation, measured with respect to the disc dynamical time, t d . The values of t d (z) for the different galaxies are shown in Fig. 2. Shown in comparison are the results from the toy model with the parameter values as derived separately from each of the four simulation types. The four top panels refer to the model parameters α, d , η and η s . The third row from top shows the clump gas and stellar masses. Before stacking, the total mass of each clump, M c , is normalized to unity at t n = 3t d , and M g and M s are normalized accordingly, thus they are given with respect to M c (t n ). The actual median clump masses at t n are M c (t n ) = 10 8.43 , 10 8.36 , 10 8.35 , 10 8.15 M for NoRP, RP, V07 and V19 respectively. The SFR in the bottom-left panel is normalized to M c /t d at t n . The bottom-right panel shows t sf , the inverse of sSFR. Recall that t dep /t d is given by d −1 .
These figures allow an evaluation of the validity of the toy model by verifying the degree of constancy of the assumed model parameters and then estimating their values in comparison to their a priori estimates. They then serve for exploring the validity of the robust model predictions in the main stage of clump evolution.
We note that in certain ways the clumps go through an initial adjustment phase of a duration on the order of a disc dynamical time after which they settle to a more stable behavior. This is expected in the isolated simulations, where the clumps form out of an otherwise unperturbed disc, but it is evident also in the cosmological simulations. For example, the specific accretion rate tends to start high and settle into a rather constant lower value. Similarly, the SFR tends to climb to a peak value and then decline and settle into a rather constant value. We crudely identify the time when the clump settles in as t n , which we choose to be t n = 3t d . This t n is used as the anchor of the normalization of the clump masses before stacking the clumps and determining the median for the plot. We focus hereafter on the more stable phase after t n , and evaluate the model parameters in this phase.
We next inspect the evolution of the model parameters one by one in the four simulation types. The median efficiency of gas accretion, α, as defined in eq. (4), is roughly constant after the first disc dynamical time, except in the isolated simulations NoRP, where it is high till t n and lower but rather flat in the range t/t d 4−15. The fluctuations are on the order of ±0.1dex. This makes the model assumption of a constant α crudely acceptable in the main stage of evolution. In the isolated simulations we estimate α 0.2 and 0.25 in NoRP and RP respectively. In the cosmological simulations we estimate at late times α 0.25 and 0.4 for V07 and V19 respectively. These are in the ball park of the values assumed a priori.
The median SFR efficiency per disc dynamical time, d , as defined in eq. (11), is rather constant, with fluctuations on the order of ±0.1dex, validating the model assumption of a constant d . This is partly by construction, as the simulations assume a local SFR law which is roughly consistent with a constant local ff in eq. (11). However, the fact that the global properties of M g anḋ M sf are related to t d via a roughly constant d does not follow in a trivial way. In the isolated simulations we deduce d 0.25 and 0.2 for NoRP and RP respectively. In the cosmological simulations we estimate d ∼ 0.5 and 0.15 for V07 and V19 respectively. We note that with a  . Evolution of simulated clump properties versus toy model, here for the isolated galaxies with RP. The details are the same as in Fig. 4, here with 31 clumps. The bottom-right panel shows t sf also for short-lived clumps from simulations of isolated galaxies with a low gas fraction and stronger feedback.
higher density threshold for star formation within the clumps, corresponding to lower values of t ff , one can obtain higher values of d while the values of ff are the same universally, as d / ff = t d /t ff . The higher value of d for V07 compared to V19 indeed reflects a higher value of t d /t ff ∼ 10 compared to ∼ 5 in V19, stemming from the different density contrasts of the clumps with respect to the disc in the two simulated galaxies. Similarly, the slightly higher value for NoRP versus RP re- Figure 6. Evolution of simulated clump properties versus toy model, here for the cosmological galaxy V07 at z = 2.5−1. The details are the same as in Fig. 4, here with 37 clumps.
flects the tendency of the star-forming region to be of higher density when the feedback is weaker.
The median gas outflow mass-loading factor, η, is roughly constant in the cosmological simulations after ∼ t d , but in the isolated simulations it starts lower and flattens off only after ∼ t n = 3t d . The fluctuations are on the order of ±0.2dex or smaller. We note that the instantaneous mass-loading factor can fluctuate by a factor of a few over a disc dynamical time (e.g., Fig. 13 of M17), for example following a sudden boost of the SFR due to a clump merger, but these peaks become significantly lower when smoothed over 0.1dex in t/t d , validating the assumption of a constant η after sufficient smoothing. In the isolated simulations we estimate η 0.1 and 1 for NoRP and RP respectively. The low η is a unique feature of NoRP, set by construction by the weak feedback assumed. In the cosmological simulations we estimate η 1 and 1.7 for V07 and V19 respectively. These val-ues of order unity are similar to the results from Dekel & Krumholz (2013), where the outflows from simulated clumps are compared to observed clumps. Similar massloading values are indicated by observations (Schroetter et al. 2019;Förster Schreiber et al. 2019) The median stellar exchange parameter, η s , as de-fined in eq. (18), carries large uncertainties, but it is low in all cases except V19. In the isolated simulations we estimate η s ∼ 0.1 and 0.25 for NoRP and RP respectively, with fluctuations on the order of ±0.2dex or smaller. The weaker stellar stripping in NoRP is likely caused by the clumps being more compact when the feedback is weaker. In V07 we deduce η s ∼ 0.1 in the range t/t d ∼ 2 − 15, and small negative values outside this range. In V19 there are actually negative values of η s during the main stage of evolution, namely a net gain of stellar mass. This is mostly accretion of young stars due to their high density in the disc of V19 at the high redshifts analyzed. A word of caution is that the value of η s may be partly contaminated by unbound disc stars in the clump identification, given that the density contrast of relatively old stars in the clumps is typically lower than that of the gas and SFR. In summary, the stellar mass exchange in our simulations tends to be small, and in most cases less important than the other ingredients of the bathtub model.
We conclude that the three important model parameters, α, d and η, typically fluctuate by only ±0.2dex or less during the main stage of evolution, and can therefore be approximated as constants during the main phase of the clump lifetime. This makes the model a valid approximation in this regime, with η s negligible except in V19. The values of these parameters as estimated from the simulations are in the ball park of the a priori expectations.

Simulated clump evolution versus model
As a result of the approximate constancy of the basic parameters, the main phase of clump evolution, after one or a few disc dynamical times, is characterized by a fairly constant SFR and M g and an associated linear growth of M s . The gas fraction is declining accordingly. The total clump mass is initially roughly constant, as long as the clump is gas dominated, and it is growing with an increasing pace as the stellar mass becomes important, approaching a linear growth.
Inspecting the evolution of masses in the different simulation types, we see that in the isolated runs, as the clumps grow from a no-clump uniform disc, the clump gas mass is growing from zero until it converges to roughly the model-predicted level. In the cosmological runs, as expected, the gas mass roughly matches the model predictions starting at an earlier time. The stellar mass matches the prediction pretty well at all times. The constancy of the SFR after t n is apparent in all the simulation types, at the level of ±0.1dex.
The star-formation time t sf , the inverse of the sSFR, approaches a linear growth with time t since clump formation, reflecting the linear growth of M s and the constancy of the SFR. Being an overestimate of t at early times, t sf becomes a good proxy for t during the main stage of clump evolution. In most cases, the corrected t sfc is a better estimate and sometimes a close lower limit for t. This is true for all the simulation types.
The depletion time t dep in units of t d , which according to eq. (13) is given by the inverse of d , is rather constant. During the main stage of evolution (and more so in the initial phase), t dep serves as an upper limit for t, while near ∼ 10t d it becomes comparable to t sf and a fair estimate of t.
At late times, toward the expected asymptotic regime, some of the few clumps that survive for more than 20t d indeed seem to show marginal indications for an accelerated growth of stellar mass at a constant gas fraction and t sf , as expected by the analytic solution. This, however, is not apparent in the figures that show the medians of stacked samples of clumps.
After evaluating the validity of the model through the constancy of its parameters, and crudely estimating the model parameters for each of the simulation types, we also show in Figs. 4 to 7 the specific model predictions for the evolution of the medians of the simulated quantities M g , M s , SFR, t sf and t sfc relevant to each simulation. The specific values of the model parameters for each simulation type, as marked in the figures, were estimated from each simulation ( §5.1, the four upper panels of the figures). They tend to be in the ball park of the a priori expectations, though they are not necessarily identical in the best fits of the different simulations.
The model seems to provide fair qualitative fits to the simulations of all four types. This is especially true during the main phase of clump evolution, at (2−10)t d , after the initial settling phase and before the end of migration, where both the model and simulations show the phase of roughly constant SFR and M g and linear growth of M s , and demonstrate that t sf is a good proxy for time since clump formation. The NoRP simulation seems to show marginal hints for the late-time asymptotic phase, which is supposed to occur after the end of migration, where the masses grow exponentially and f g and t sf saturate to a constant value.

Time indicators for short-lived clumps
While the analysis so far focused on the clumps that survive for more than 50 Myr, we wish to verify to what extent t sf , t sfc and t dep can serve as observable indicators of time since formation also for the short-lived clumps. This will be useful in the attempt to distinguish between the scenarios of SL and LL clumps.
The SL clumps typically disrupt before or during the early stages of the main phase of evolution, by ∼ 3t d . However, their evolution during their short lifetimes may be similar to that of the LL clumps. In particular, similarly to the LL clumps in their initial phase, we expect t sf for SL clumps to be smaller than t dep , as t sf = f sg t dep , and the gas fraction in SL clumps is expected to be above 50% (as in the simulations of Oklopčić et al. 2017). After a short initial phase, the t sf of SL clumps is expected to grow linearly with time, much like the LL clumps, but the SL clumps are disrupted before their t sf reaches the levels of t dep .
In order to crudely test the above predictions, we inspect simulated SL clumps. For this we appeal to isolated simulations similar to the ones used here (Fensch Figure 8. Age as a proxy for time. Shown are the medians and 68% scatter of the average stellar age in clumps versus time since clump formation, for the clumps in V07. The age is a proxy for time after the first disc dynamical time. & Bournaud 2021), but where most of the clumps are short-lived. This was obtained by starting with discs of a low gas fraction, f g ∼ 0.2, and adopting a strong feedback including thermal supernova feedback, radiation pressure and photoionization. These simulations typically show 3-6 clumps of baryonic mass ∼ 10 8 M and larger for at least half the time. The key difference from the other simulations used here is the resulting shorter clump lifetimes. Figure 5, bottom-right panel, also shows the median and 68% scatter for t sf in these SL clumps as a function of time from clump formation. The median for the SL clumps is indeed similar to the median for the LL clumps from the gas-rich simulations with RP. After the short initial phase of duration ∼ t d , the star-formation time t sf becomes a good proxy for time t also for SL clumps at (1−3)t d .

Clumps in CANDELS
In Guo et al. (2018, hereafter G18), we have analyzed 3187 clumps identified in 1269 galaxies from the CANDELS/GOODS-S field in the redshift range z = 0.5 − 3, based on the sample described in Guo et al. (2015). The sample, clump detection and estimates of clump properties are described in these papers, while we bring only a very brief summary here.
For each galaxy, the stellar mass and SFR were determined by SED fitting, while the SFR was also determined directly from the UV flux. Star-forming galaxies were selected at 0.5 < z < 3 to have global M s > 10 9 M and sSFR> 10 −1 Gyr −1 . An apparent magnitude cut of HF160W < 24.5 AB has been applied to ensure reliable morphology and size measurements. In order to be able to identify and measure clumps given the resolution, the sample has been limited to galaxies whose effective radii along the galaxy semi-major axis (SMA) is larger than 0.2". To minimize the effect of dust extinction and clump blending, only galaxies with axial ratio q > 0.5 were included. This sample consists of 1655 galaxies.
The clumps were detected in rest-frame Near-UV. A blob of Near-UV excess is identified as a clump if its UV luminosity is at least 3% of the total UV luminosity of the host galaxy. A total of 3187 clumps have been identified in 1269 galaxies. In the following first comparison we use a subsample of these galaxies, with log M s,gal = 9.8 − 11.4 in the redshift-range z = 1.5 − 3 where UV fluxes are available, a total of 607 clumps in 284 galaxies.
The flux in each of nine HST bands was measured for each clump using an aperture of 0.18", and an aperture correction has been applied. Background disc subtraction has been applied in each band, using different levels of correction, ranging from no subtraction to very aggressive subtraction (G18, Table 2). The clump properties, including M s , age and dust reddening, were derived by fitting their HST SEDs to stellar population synthesis models (Bruzual & Charlot 2003) with Chabrier IMF (Chabrier 2003), as described in Guo et al. (2012). The SFR used here, for galaxies at z > 1.5, was estimated directly from the UV flux. The SFR estimates from UV are independent of the mass, age and SFR estimates from SED fitting, so we prefer to use the UV SFR here.
We use below the clump properties as derived after a moderate correction for disc contamination, the fiducial background subtraction denoted aperi v3 in Table 2 of G18. The radial gradients turn out to show the same qualitative features for the different levels of correction, from a very defensive subtraction (aper3) to a very aggressive subtraction (aperi v1), except for the case of no background subtraction at all (aperi v4).
As a note of caution, we learn from Huertas-Company et al. (2020) that using SED fitting on aperture as in G18 seems to overestimate the stellar mass of clumps by a factor of a few to several. This factor, however, is rather constant for all clump masses and for clumps in all galactic radii (and therefore likely for clumps of all ages). In our comparisons to the G18 observed clump masses (and SFR, which could be subject to a similar overestimate), this overestimate may thus affect the overall normalization, but it has little effect on the trends. To accommodate this uncertainty in mass normalization, we will normalize our model relations to match the masses at young ages.

Simulated time vs age and distance
The available observed quantity for time is the mean stellar age in a clump, while the main quantity of interest for dynamical clump evolution is the time since its formation. The age may be larger due to older stars that formed outside the clump and became part of it at or after clump formation. These could be a real bound part of the clump, or unbound stars that are temporarily within the clump boundaries, or foreground/background disc stars that are associated with the clump by the 2D observational identification method. In order to find out to what extent the stellar age can serve as a proxy for time since clump formation, we inspect in Fig. 8 the median age as a function of t, both in units of t d , for the clumps in the cosmological simulation V07. We learn that after the first disc dynamical time the age is a good proxy for time since formation, indicating that the contamination by background stars is small. In V19 the contamination is more severe. This could be because at this high redshift the accreted disc stars are < ∼ 100 Myr old (see the time range of negative η s in Fig. 7) while the migration time is only ∼ 150 Myr, making the clump stellar ages dominated by the young background stars.
Another observed quantity is the projected distance of the clump from the galaxy center relative to the semimajor axis of the galaxy, d/a. Since the clumps tend to form in the outer disc and migrate inward, d/a may also serve as a proxy for time since formation, and gradients of clump properties as a function of d/a are expected to be related to clump evolution in time. In order to obtain an impression of this, we show in Fig. 9 the median of d/R d as a function of t/t d for the clumps in the cosmological simulations. We see that after ∼ 3t d , the clump distance is declining with time, as expected, reflecting the inward migration. The effective trend in this regime seems to be roughly d ∝ t −1/2 for the two cosmological simulations. A fair fit is provided by d/a e −t/tmig . This indicates that the translation from t to d cannot be modeled in detail by assuming a constant radial velocity during migration. The migration time is (10−20)t d , in the ball park of the prediction in eq. (1). Figure 10 puts together the observed CANDELS clumps from G18, in the galaxy mass range log M sgal = 9.8−11.4 and in the redshift range z = 1.5 − 3 (where UV-based SFR estimates are available), a total of 607 clumps in 284 galaxies. Shown are the clump stellar mass and UVbased SFR with respect to the whole galaxy, as well as the star-formation time within the clump t sf = sSFR −1 . These properties are plotted against clump age (left panels) and against distance from the galaxy center d with respect to the semi-major axis a (right panels). The medians and 68% scatter are shown in log-spaced bins.

Observed clump properties vs age & distance
The left panels of Fig. 10 show the clump properties as a function of its age. While the t sf values derived from the SED-based SFR are correlated with the clump ages by construction, as they both arise from the same SED fitting, the t sf derived from the UV-based SFR is not expected to be degenerate with the age. We notice in the top-left panel that t sf and clump age are still significantly correlated. If the age is a proxy for time since clump formation, as in V07 (Fig. 8), we learn from the observed correlation that the star-formation time t sf , much like the age, can indeed serve as a proxy for t.
The median clump SFR (bottom-left panel) is rather constant as a function of clump age. This is consistent with the constant SFR as a function of time predicted by the toy model during the main clump evolution phase. The median stellar mass is growing roughly linearly with age (middle-left panel). This is as expected from the toy model once the age of the clump stars reflects the time since clump formation.
Our bathtub toy-model results for t sf , M s and SF R are shown in comparison assuming that stellar age and time since clump formation are the same. The parameter choice for a good fit is α = 0.3, d = 0.1, η = 1, η s = 0.3 and f gi = 0.97. The stellar mass and SFR are normalized to match the observational results at young ages. The model reveals remarkable fits to the rates of evolution of the medians of the observed clump properties as a function of age.
The right panels of Fig. 10 show the gradients of the same clump properties as a function of projected distance from the center d/a. The top-right panel shows t sf vs d/a. The anti-correlation between the two is consistent with the expected migration inwards and assertion that t sf is a proxy for time since clump formation. As expected, the observed median SFR (bottom-right panel) hints to only a weak correlation with distance, consistent with the roughly constant SFR with age, and the expected constant SFR with time during the main stage of evolution. The median of the observed M s (middleright panel) shows a strong anti-correlation with distance, reflecting the strong correlation of M s with age and the anti-correlation between age and d/a. This is consistent with the predicted growth of M s in the main stage of clump evolution combined with the inward migration.
In principle, one could extend the toy-model predictions to the gradients as a function of d/a by attempt- Figure 10. Observed clump properties against clump stellar age (left) and against galacto-centric distance with respect to the semimajor axis (right). Top: star-formation time t sf =sSFR −1 . Middle-bottom: clump stellar mass and SFR with respect to the whole galaxy. The points represent 607 clumps in 284 CANDELS galaxies with log Ms gal = 9.8−11.4 at z = 1.5−3, from the sample analyzed in G18. The medians and 68% scatter are shown in log-spaced bins. A model fit is shown as a function of time, with the parameter choice α = 0.3, d = 0.1, η = 1, ηs = 0.3 and fg0 = 0.97, and with a normalization that matches the properties at young ages. The median radial gradients from the simulations (M17) are shown in comparison to the observed gradients as a function of d/a, normalized to match the observations when necessary.
ing to model the migration process and thus translate t/t d to galacto-centric distance and to its observable counterpart d/a. This requires accurate modeling of the migration process which we do not attempt here. Instead, we display in comparison to the observed radial gradients the medians from the radial gradients of the simulated clumps (following M17). We learn that the slopes of the simulated and observed radial gradients are similar, consistent with clump migration. The normalizations of M s and SFR in the simulations and observations are different, and were adjusted in the figures as marked in the labels. These offsets largely result from the limited observed resolution of more than 1 kpc, which causes an overestimate in the masses and SFR assigned to clumps (Huertas-Company et al. 2020;Ginzburg et al. 2021).

Other observations
While we have learned that the large sample of CAN-DELS clumps from G18 provides convincing evidence in favor of long-lived clumps, smaller samples provide supporting evidence. For example, several other studies have attempted to estimate clump ages using stellar population synthesis modeling (Elmegreen & Elmegreen 2005;Elmegreen et al. 2009;Genzel et al. 2011;Förster Schreiber et al. 2011;Wuyts et al. 2012;Guo et al. 2012;Zanella et al. 2015;Soto et al. 2017). These estimates carry large systematic and random uncertainties, partly due measurement errors and partly due to uncertainties in the SED modeling including the star-formation histories, dust extinction and background subtraction. Given these uncertainties, the estimated ages span a wide range, from ∼ 10 Myr to < ∼ 1 Gyr, with most estimates indicating average ages of a few hundred Myr, consistent with the results of G18 used above and with the predicted LL clumps. Cibinel et al. (2017) used ALMA to observe six clumps in a main-sequence clumpy galaxy at z ∼ 1.5 in the HUDF. From the stellar masses and SFRs they estimate t sf ∼ 70−900 Myr. Using ALMA CO(5-4) observations, they also constrain the gas mass and derive crude upper limits for the depletion times, t dep < (0.1−2) Gyr. When stacking the six clumps together, they estimate t dep ∼ 200 Myr, comparable to and slightly smaller than the typical value of t sf for these clumps, with f g ∼ 20 − 50%. This argues for long-lived clumps, that had time to consume gas into stars. Such gas fractions are consistent with the typical gas fractions in simulated LL clumps as seen in Figs. 4 to 7 and lower than expected for SL clumps whose gas fractions are predicted to be > 50% (Oklopčić et al. 2017).
A few other studies have examined radial gradients of clump properties with galacto-centric distance in small samples of clumps and galaxies Guo et al. 2012;Soto et al. 2017;Zanella et al. 2019). They all find evidence for older, redder, more massive clumps closer to the disc centre, in qualitative agreement with the finding of G18, and the predicted migration-induced radial gradients of LL clumps.

DISTINGUISHING FEATURES OF LL VS SL
With the theory predictions in hand one can attempt to identify the observable properties that may allow a distinction between survival versus disruption of massive in-situ clumps. We refer to long-lived (LL) versus shortlived (SL) in-situ clumps, ex-situ clumps (Ex), and the background disc (D). We base the predictions on the model of the current paper and the simulation results in M14 and M17. The distinguishing properties, listed by observables, are as follows: • Position in the disc: The migrating LL are expected at all radii while the SL are expected where they form, preferentially in the outer disc. Ex tend to populate the outer disc and they typically migrate less because they tend to oscillate aboou the disc.
• Star-formation time: The SFR time t sf , the inverse of the sSFR, is a proxy for age for in-situ clumps (LL and SL), but not necessarily for Ex and D. The LL were crudely expected in M17 to show a declining radial gradient, crudely t sf ∝ d −1 , but no such radial gradient is expected for SL. A similar radial gradient is expected for Ex and D, with the t sf of D higher by ∼ 0.6dex. • Stellar mass: For LL typically M s > 10 8 M , while for SL typically M s < 10 8.5 M , on average smaller than LL by ∼ 0.5dex. For Ex M s > 0.07M d , on average larger than LL by 0.5dex. A declining radial gradient is expected based on M17, for LL crudely M s ∝ d −1 and for SL M s ∝ d −0.5 , with the slope for Ex in between.
• SFR: The SFR is expected to be roughly constant with age and distance for the in-situ clumps (LL and SL).
• Gas fraction: LL are expected to show a range of f g values, possibly above and below ∼ 0.5, while SL would tend to have only high values, f g > 0.5, and Ex would tend to show very low values, f g < 0.1. Based on M17 we expect a radial gradient for LL, crudely f g ∝ t −1 , namely decreasing with decreasing d, but roughly constant values for SL and Ex.
Based on the above, the distinguishing properties listed by clump type are as follows: • The LL clumps and the SL clumps are different from one another in age, t sf , M s , and f g , but not in SFR. In the VELA simulations, they are also distinct by the distribution of distances from the galaxy center and by the value of M s . • The LL and Ex clumps are distinct by age, M s , f g , and more weakly by M s and SFR.
• However, as discussed in appendix §A, contamination of the clumps by disc may give false values and generate false radial gradients, especially for the intensive quantities such as age and f g .
• The disc contamination of the SFR is expected to be relatively small, given the high density contrast between clump and disc, so the SFR could be corrected reliably. Hopefully, the correction to M s can also be manageable, estimated in our simulations to be at the level of a few tens of percent. If so, then t sf can be estimated, and t sfc can serve as a proxy for clump age.
• Most of the ex-situ clumps can be removed based on age (> 500 Myr), f g (< 0.1), and maybe also mass (> 0.07M d ).
Our practical conclusions for distinguishing between LL and SL clumps are thus as follows: • The radial gradient of the age proxy t sf is the most promising distinguishing feature among the quantities that can be corrected for disc contamination.
• The radial gradient of M s and its amplitude are potential distinguishing features that may be corrected for disc contamination, under the assumption that any gradient that might have been established at formation is weaker.
• Age and f g are distinguishing features in principle, but they may be useful only if the disc contamination is small.
The observed radial gradients of CANDELS clumps (G18), as summarized in Fig. 10, indicate consistency with the scenario of long-lived, migrating clumps, and inconsistency with the clump population being dominated by short-lived clumps. In particular, the distribution of clump ages, either as estimated from SED fitting or from t sf , indicates that most clumps are older than 100 Myr. The clumps span a broad range of distances from the galaxy center, and both age and t sf show significant gradients with d/a. Furthermore, the clump stellar mass shows a strong gradient with d/a.

CONCLUSION
Motivated by the desire to understand the nature of the so-prominent giant clumps in high-z star-forming disc galaxies, and in particular attempting to distinguish between long-lived and short-lived clumps, we studied the evolution of migrating clumps, using an analytic toy model and different types of simulations, and confronted theory with CANDELS observations. Our theoretical results for the migrating long-lived clump evolution can be summarized as follows: • The unique feature of the bathtub toy model of clumps, distinguishing it from the bathtub model for a whole galaxy, is that the gas accretion rate into the clump scales with the total clump mass, which can vary on the relevant timescale and is not directly constrained by the cosmological accretion rate. This allows an analytic solution specific for clumps, a sum of exponential growing and decaying modes.
• The analytic model reveals a main phase of clump evolution, in which the SFR and gas mass are constant, and therefore the stellar mass is growing linearly in time, such that the inverse of the sSFR is a very useful proxy for clump age. This phase is valid after an initial period of ∼ t d , the disc dynamical time, and until after ∼ 10 t d , namely until the clump has likely completed its inward migration. This constant SFR is very different from the commonly assumed "tau-model" exponential decay behavior in a closed box, and from the cosmological evolution of a whole galaxy where the sSFR is systematically decreasing in time in proportion to (1 + z) 5/2 following the cosmological specific accretion rate. The tau-model, which is commonly assumed in clump mass estimates, is likely to lead to wrong results.
• Later on, the solution is intended to approach an asymptotic exponential growth of mass and SFR with constant specific SFR and gas fraction, but this phase is unlikely to materialize as the typical migrating clumps are expected to disappear at the galaxy center by ∼ 10−20t d . • The main model parameters are the accretion efficiency per dynamical time α ∼ 0.2 − 0.4, the SFR efficiency in a disc dynamical time d ∼ 0.1 − 0.3 (corresponding to ff ∼ 0.02), and the gas outflow massloading factor η ∼ 1 − 3. The stellar exchange between the clump and the disc has a less-important role.
• The evolution of clumps in simulations, both of gasrich isolated galaxies and of zoomed-in high-z galaxies in a cosmological setting, with different feedback strengths (though all at a moderate level), is found to be in qualitative agreement with the analytic predictions. The simulations confirm the validity of the toy model as a sensible approximation, they yield model parameter values in the expected range, and they match the model-predicted evolution patterns. The model simplifying assumptions such as the constant disc properties at different radii and constant model parameters as the clump migrates and evolves, turn out to be sensible in view of the qualitative match of the model to the simulations. The toy model thus provides a very useful tool for further studies of giant clumps in high-z galaxies.
In view of the above theoretical expectations, the observational results for massive clumps in a sample of more than 600 clumps at z = 1.5 − 3 from CANDLES can be summarized as follows: • The observed clumps show a broad distribution of clump ages, extending from a few to ∼ 500 Myr and beyond. This is true for ages as estimated from SED fitting or UV as well as from the star-formation time t sf , the inverse of sSFR, which is confirmed to be a useful age proxy, as predicted. They also show a wide spread of distances from the galaxy center, which are anti-correlated with the ages, as expected from inward migration.
• The observed clumps show evolution patterns and radial gradients of clump stellar mass and SFR that resemble the theoretical predictions for migrating clumps. In particular, the constant SFR over time and radial distance d, and the corresponding linear rise of stellar mass with age and inversely with d, resemble the predicted main phase of clump evolution.
• The observations thus favor the scenario where the massive clumps, roughly above one per cent of the galaxy stellar mass, largely survive feedback, remain long-lived and bound, and migrate inward on a timescale of a few disc orbital times. The feedback strength as incorporated in the range of simulations used here, with a mass-loading factor of order unity, is successful in reproducing the observations, ruling out much stronger, disruptive feedback. • It is likely, though, that lower-mass short-lived clumps co-exist, as well as ex-situ clumps that have merged into the disc. We pointed out features that could help distinguishing them from the in-situ clumps.
The radial gradient of the age proxy t sf , the inverse of sSFR, is the most promising distinguishing feature between LL and SL clumps, among the quantities that can be corrected for disc contamination. The radial gradient of the clump stellar mass M s and its amplitude are potential distinguishing features that may be corrected for disc contamination. The other estimates of age and gas fraction f g are distinguishing features in principle, but they may be useful only in cases where the disc contamination is small.
A potential caveat is contamination of the signal of migration by gradients in clump properties that were implanted at clump formation by radial variations in disc properties. Our simulations indicate that the initial radial gradients of clump properties are weaker than the gradients induced by clump migration and the strong gradients observed, but this requires further investigation.
Preliminary kinematical analysis of observed z = 1−2 star-forming discs (Genzel et al., private communication) indicates detections of radial velocities inward of gas emitting Hα and CO, as residual velocities after subtracting a rotating disc and a dispersive bulge. The indicated inward velocities seem to be in the ball park of the predicted migration rate, yet to be properly observed and analyzed. If these residual velocities are largely associated with clumps, they may provide direct evidence for clump migration, thus supporting their long-lived nature.
The need for survival of massive clumps does not permit a strong ejective and disruptive feedback in these clumps. On the other hand, the indicated low overall stellar-to-halo mass ratio, based on abundance matching of observed galaxies with simulated ΛCDM haloes (e.g., Rodríguez-Puebla et al. 2017;Moster, Naab & White 2018;Behroozi et al. 2019), requires strong preventive feedback that would suppress the total stellar mass. The subgrid models for supernova and stellar feedback, as incorporated in current simulations, fail to simultaneously fulfill these two conflicting requirements. Thus, understanding the complex nature of feedback presents a non-trivial challenge.