text:
# Direct comparison of phase-space distributions of K⁻ and K⁺ mesons in heavy-ion collisions at SIS energies - evidence for in-medium modifications of kaons ?

## Abstract

The ratio of K$`^{}`$ to K$`^+`$ meson yields has been measured in the systems 96Ru+96Ru at 1.69 A GeV, 96Ru+96Zr at 1.69 A GeV, and 58Ni+58Ni at 1.93 A GeV incident beam kinetic energy. The yield ratio is observed to vary across the measured phase space. Relativistic transport-model calculations indicate that the data are best understood if in-medium modifications of the kaons are taken into account.

Recently there has been considerable effort, both experimentally and theoretically, to investigate changes of hadron properties in a hot and dense nuclear medium. In particular, a variety of theoretical approaches consistently predict that the effective mass of kaons increases slightly with increasing baryon density, while the mass of antikaons is expected to drop substantially . This phenomenon could lead to the formation of a kaon condensate in a dense hadronic environment , which in turn would effect the nuclear equation of state, and have consequences for the physics of neutron stars . The modifications of the properties of kaons in a hadronic medium might originate from the partial restoration of the chiral symmetry of QCD .

The question whether the kaon masses are modified in a dense hadronic environment can be addressed experimentally with studies on kaons produced in heavy-ion collisions at bombarding energies around 1-2 GeV per nucleon, which is close to the production threshold in elementary, nucleon-nucleon reactions (1.6 and 2.5 GeV for K+ and K- mesons, respectively). At these beam energies, kaons are most likely produced in the early stage and in the central region of the collision , where densities of up to 3 times the normal nuclear matter density and temperatures in the order of 100 MeV can be reached . The production rate can be influenced not only by the surrounding medium (e.g., its density) or the properties of the nucleons and their resonances, but also by possible changes of the kaon properties themselves.

The observation of the enhanced K- meson yield at midrapidity in heavy-ion collisions with respect to the elementary, nucleon-nucleon reactions  is a very interesting signature, possibly related to a substantial in-medium drop of the effective mass of K- mesons. A change of the effective mass of a particle can be understood as an effect of a density dependent potential. Gradients of this potential cause forces that act on particles. While antikaons are attracted into regions of high baryon density, kaons are repelled from these regions. This effect offers an explanation  for the characteristic directed side-flow patterns of K+ mesons that are observed experimentally . It also predicts characteristic changes of the final-state phase-space distributions of kaons and antikaons .

In this Note we report on measurements of K+ and K- mesons produced in 96Ru+96Ru/96Zr collisions at 1.69 A GeV and in 58Ni+58Ni collisions at 1.93 A GeV incident beam kinetic energy. The experiments were performed at SIS/GSI, using the FOPI experimental setup  which allows for a simultaneous measurement of all charged reaction products. Thus final-state distributions of the particles can be directly compared within the same event sample and with the same acceptance. Results on the phase-space population of $`\pi ^\pm `$, $`p`$, $`d`$, and K+ in the Ni+Ni experiment were reported elsewhere . Here, we show, for the first time in this energy regime, the ratio of K- to K+ meson yields in the backward hemisphere. We extract the ratio across a relatively wide region of phase space, which provides high sensitivity to the dynamics of the propagation of kaons through the medium. We observe that the phase-space distributions of K- and K+ mesons differ and discuss the origin of this effect.

An ensemble of events biased to small impact parameters has been selected by requiring high charged-particle multiplicity in the polar-angle range $`7^{}<\mathrm{\Theta }_{Lab}<30^{}`$ on the trigger level. For the Ni+Ni experiment, 4.7$``$106 events were selected, corresponding to the centralmost 11% of the total geometrical cross section. Since no difference in strangeness production was found between the Ru+Ru and Ru+Zr systems , the accumulated statistics was combined for a total of 7.7$``$106 events, corresponding to the centralmost 14% of the geometrical cross section.

The identification of K+ and K- mesons with the FOPI detector relies upon the information on specific energy loss and track curvature in the Central Drift Chamber (CDC), and on the measurements of Time-of-Flight (ToF) in the surrounding Plastic Scintillator Barrel detector . The acceptance is thus restricted to the polar-angle range in the laboratory reference frame 39$`{}_{}{}^{}<\mathrm{\Theta }_{lab}<`$130 and to transverse momenta $`p_t>`$ 0.1 GeV/c. The finite detector resolution limits the identification of K+ mesons to laboratory momenta $`p_{lab}<0.5`$ GeV/c. Due to the much lower K- yield compared to that of K+ meson, the background contamination in the former increases more rapidly with momentum. In order to eliminate possible distortions due to a momentum dependent background contribution, the considered momentum range for K- meson identification is restricted to $`p_{lab}`$ below 0.32 and 0.34 GeV/c in case of the Ru+Ru/Zr and Ni+Ni experiments, respectively. At $`p_{lab}`$ $``$ 0.32 (0.34) GeV/c the measured K- meson yields have less than 20% background contamination and the K+ meson yields are practically background free (the contamination is less than 5%).

With the different upper momentum limits for positive and negative kaons, mentioned above, around 26000 K+ and 240 K- mesons have been identified in the Ni+Ni experiment. The combined statistics of the Ru+Ru/Zr experiment is about 40000 K+ and 220 K- mesons. The mass spectrum of particles with charge $`\pm `$1 measured in the Ru+Ru/Zr experiment is shown in Fig.1.a. Peaks from K+ and K- mesons are clearly visible. The portion of the phase space populated by K+ mesons measured in the Ru+Ru/Zr experiment is shown in Fig.1.b in terms of normalized rapidity ($`y^{(0)}=y^{lab}/y^{CM}1`$, where $`y^{CM}`$ is half of the beam rapidity) and transverse momentum. In this representation -1 and 0 on the rapidity axis correspond to the target and the midrapidity, respectively. The yield of K+ mesons is depicted by the contour lines on a linear scale. The geometrical limit at $`\mathrm{\Theta }_{lab}`$=39 and the upper $`p_{lab}`$ limits for K+ and K- meson identification are depicted by dashed lines. The solid lines show the polar-angle range in the center-of-mass (c.m.) reference frame 150$`{}_{}{}^{}<\mathrm{\Theta }_{cm}<`$165, which will be referred to later.

To quantitatively examine and compare the phase-space distributions of K+ and K- mesons, where for the latter case low statistics does not allow to extrapolate the measured yields to experimentally not accessible regions of the phase space, we study the ratio of K- to K+ meson yields in the limited phase-space region defined by the K- meson identification. This offers two advantages with respect to analysing the single particle distributions. (i) Experimental difficulties, like detection efficiencies and acceptance deficiencies, cancel to a large extent . (ii) In-medium effects act in opposite ways on K- and K+ mesons, hence the ratio should reveal these more clearly.

In FOPI, the efficiency for particle detection is given by the tracking efficiency in the CDC and the matching efficiency with the ToF Barrel. Possible systematic bias on the measured K-/K+ ratio was estimated using a Monte-Carlo simulation in which the full detector response was modelled with the GEANT package . The simulated data were analysed in the same way as the experimental data. Comparing the output of the simulation to its input, the final K-/K+ ratio was found to be overestimated by 15%, independently of transverse momentum and rapidity. This effect is attributed to different efficiencies of the track finding for positively and negatively charged particles due to the geometry of the CDC. A similar asymmetry of the efficiency was reported in  for $`\pi ^+`$ and $`\pi ^{}`$ mesons. All data points shown in the following figures are corrected for this systematic bias, i.e., all ratios are reduced to 87% of the directly measured value. Furthermore, the systematic uncertainties due to the identification criteria ($`<`$30%) and the background contamination ($`<`$20%) were estimated by varying the conditions that were imposed on track parameters in order to select K+ and K- mesons. These systematic errors are depicted by light-grey shaded areas in Fig.2 and 3. The systematic distortion due to the nuclear scattering of K- mesons in the target material is neglected, since for the targets used in the experients (1% interaction length in the beam direction and a transverse diameter of about 12 mm), this interaction probability is estimated to be below 3%.

Since produced particles, and especially kaons, are found to be emitted almost isotropically in the c.m. system , we plot in Fig.2 the K-/K+ ratio as a function of the kinetic energy in the c.m. reference frame ($`E_{cm}^{kin}`$) for both the Ru+Ru/Zr (a) and the Ni+Ni (b) experiments. The polar-angle range 150$`{}_{}{}^{}<\mathrm{\Theta }_{cm}<`$165 has been chosen since there the kinetic-energy acceptance is largest (see Fig.1.b). We observe that the ratio rises towards small $`E_{cm}^{kin}`$.

In order to test whether the effect can be caused by the electric field of the hadronic fireball formed in the collision, we performed numerical calculations following the studies on the influence of the Coulomb potential on charged-particle spectra reported in . We followed the propagation of K- and K+ mesons in an electric field of a net charge $`Z`$, which was originally homogeneously distributed in a sphere of a $`R_f`$ radius, and expanded radially with a $`\beta _{rad}`$ flow velocity. Initially, the K- and K+ mesons were evenly distributed within the expanding volume, and had identical energy spectra that corresponded to an isotropic emission from a midrapidity source of a $`T_f`$ temperature. The yields of K- and K+ mesons were arbitrarily normalized in order to allow a direct comparison with the experimental data. The parameters of the simulation (the total charge $`Z`$, the expansion velocity $`\beta _{rad}`$, the temperature $`T_f`$, and the radius $`R_f`$) were varied in a reasonable range in order to model possible freeze-out conditions. The dark-grey shaded areas in Fig.2 correspond to the results obtained with different sets of parameters. Comparing this to the data, we conclude that the influence of the Coulomb potential of the net positive charge of colliding ions on the K-/K+ ratio is too small to account for the observed relative narrowing (widening) of the K- (K+) meson final-state phase-space distributions.

In Fig.3 we plot the K-/K+ ratio as a function of $`y^{(0)}`$, i.e., in the direction parallel to the beam axis. The result is biased by the detector acceptance, but has the best statistical significance, which is given by the error bars attached to the data points. We observe that the K-/K+ ratio rises towards midrapidity. How representative the trend in the data is for the unbiased rapidity density distribution depends on the variation of the ratio as function of transverse momentum for fixed rapidity. For a meaningful comparison with the model predictions (see below), an acceptance filter was applied to the results of the model, taking into account the angular boundaries and the maximum laboratory momentum. The distortion due to acceptance turned out to be smaller than 20% in the case of the rapidity dependence of the ratio shown in Fig.3. The kinetic-energy dependence of the ratio presented in Fig.2 is not biased by the acceptance in the considered range of $`E_{cm}^{kin}`$.

Since the elementary reaction threshold is higher for K- than it is for K+ mesons (2.5 and 1.6 GeV respectively), one might expect that at the time of production, the average kinetic energies are lower for K- than for K+ mesons. However, when the incident beam energy is far below threshold, the production of K- is dominated by channels involving intermediate baryon resonances and/or pions. According to Relativistic Boltzmann-Uehling-Uhlenbeck (RBUU) transport-model calculations  owing to this mechanism the initial momentum-space distribution of K- is even wider than that of K+ mesons. In addition, particles are further rescattered, which tends to equalize their phase-space distributions. The solid lines in Fig.2 and 3 show the values of the K-/K+ ratio predicted by the RBUU transport model  when describing the kaon scattering in a fashion consistent with the free particle properties. No significant dependence of the ratio on $`y^{(0)}`$ and $`E_{cm}^{kin}`$ is found, which results in a manifest contradiction to the data. The same conclusion has been drawn with another realization of a RBUU-type model .

Finally, we try to explain the variation of the K-/K+ ratio in the phase space by in-medium modifications of kaon properties. Dashed and dotted lines in Fig.2 and 3 show the values of the ratio predicted by the RBUU model when in-medium effects are taken into account by a linear dependence of the in-medium potential on density . Two scenarios with different strength $`\mathrm{U}(\rho =\rho _{})`$ of the in-medium kaon potentials at normal nuclear matter density are presented.

The attractive K- potential influences the results shown in Fig.2 and 3 in a systematic fashion: with increasing depth of the K- potential, the K-/K+ ratio increases on average, and in addition the slope of the ratio with respect to rapidity and kinetic energy rises. While the first effect is caused by the corresponding drop of the effective mass, the rapidity and kinetic-energy dependencies are generated by the forces originating from the gradients of the potentials. Similar observations have been made on the $`p_t`$ dependence of the K-/K+ ratio . Varying the repulsive potential for K+ mesons does not modify the ratios significantly, as also shown in . It has to be noted, however, that the calculations presented above underestimate the production yield of K+ in the Ni+Ni collisions by a factor of two, while with very similar parameters ($`\mathrm{U}(\rho =\rho _{})=`$+20 MeV for K+ the production rate is reproduced in . Despite this discrepancy, in-medium modifications of kaon masses are presently the only mechanism explaining the trends found in the data.

Recently very similar observations of the K-/K+ ratios in heavy-ion collisions have been obtained at AGS energies (11.6 A GeV) . Discrepancies to transport-model calculations were interpreted as signal for multi-body collisions owing to the very high densities that are reached in this energy regime. While the sensitivity of kaon-production ratios to in-medium properties is enhanced in the vicinity of the elementary reaction thresholds, it remains to be seen by careful model analysis whether both mechanisms can be clearly separated.

In summary, we measured the ratio of K- to K+ meson yields as a function of different kinematic variables in the experiments 96Ru+96Ru96/Zr at 1.69 A GeV and 58Ni+58Ni at 1.93 A GeV incident beam kinetic energy. We found the K- final-state phase-space distribution narrower than that of K+ mesons. It is unlikely that this effect is due to the different kinematical conditions of meson production or to the influence of the Coulomb potential of the net positive charge of colliding ions on meson propagation. However, it can be explained when assuming modifications of kaon properties in a dense nuclear medium. Effective theories that implement and exploit chiral symmetry breaking patterns of QCD argue that these modifications are a consequence of the restoration of the symmetry in a hot and dense nuclear matter . However, the drop of the K- effective mass can be also explained as an effect of the in-medium modifications of the $`\mathrm{\Lambda }(1405)`$ spectral function due to the Pauli blocking of the proton . In addition, it was suggested in  that in heavy ion experiments the chaotic as well as the coherent movement of the baryonic matter may additionally mask the influence of the in-medium potential on the measured antikaon yields. In order to address the outstanding problems, more systematic studies of various systems and energies and with a better acceptance is needed. A further theoretical clarifications of the origin of the in-medium effects on kaons is also necessary.

We are grateful to E.L.Bratkovskaya, W.Cassing and G.Q.Li for providing us with the RBUU events and for useful communications. This work has been supported by the German BMBF under contracts POL-119-95, RUM-005-95, UNG-021-96 and RUS-676-98 and by the Deutsche Forschungsgemeinschaft (DFG) under projects 436 RUS-113/143/2 and 446 KOR-113/76/0. Support has also been received from the Polish State Committee of Scientific Research (KBN) under grants 2P302-011-04 and 2P03B 05716, from the Korean Science and Engineering Foundation (KOSEF) under grant 985-0200-004-2, from the Hungarian OMFB under contract D-86/96 and from the Hungarian OTKA under grant T029379.


format: text

text:
# CHEMICAL EVOLUTION OF THE GALAXY BASED ON THE OSCILLATORY STAR FORMATION HISTORY

## 1 INTRODUCTION

Revealing the star formation histories (SFHs) of galaxies is essential in understanding the galaxy formation and evolution. The SFH of the Galaxy (Milky Way) is worth studying, since a large number of stars are observed individually and the SFH is inferred directly from the age distribution of the stars. The SFH is closely related to the chemical evolution of the Galaxy. For example, the age–metallicity relation (e.g., Pagel 1997) of Galactic stars is generally believed to originate from the chemical enrichment of the Galaxy as a result of star formation.

Eggen, Lynden-Bell, & Sandage (1962) have pioneered the modeling of the Galactic SFH and its chemical evolution. From the correlation between the ultraviolet excess and orbital eccentricity of stars, they have concluded that the Galaxy formed by collapse on a free-fall timescale from a single protogalactic cloud. An alternative picture of halo formation has been proposed by Searle & Zinn (1978). They have argued that the Galactic system formed from the capture of fragments such as dwarf galaxies over a longer timescale than that proposed by Eggen et al. In any case, determining the timescale of the infall of matter and the chemical enrichment is an important issue to resolve the formation history of the Galaxy.

A number of papers have investigated the formation (e.g., Burkert, Truran, & Hensler 1992) and chemical evolution (e.g., Matteucci & François 1989) of the Galaxy and other spiral galaxies (e.g., Lynden-Bell 1975; Sommer-Larsen 1996). Many models of the SFH of the Galaxy have treated the formation of the Galactic disk through gas infall from the halo. This scenario (the so-called infall model) can be consistent with the age–metallicity relation of the disk stars (e.g., Twarog 1980), if a reasonable SFH is used. Moreover, the infall model provides a physically reasonable way of solving the G-dwarf problem (e.g., Pagel 1997, p.236), contrary to the closed-box model which tends to overpredict the number of the low-metallicity stars.

Since stars are formed from interstellar medium (ISM), one of the factors that determine the star formation rate (SFR) is the gas content of galaxies. Indeed, the SFR and the gas density is closely related (Kennicutt, 1998). The most commonly used relation is called the Schmidt law (Schmidt, 1959). It assumes that $`\mathrm{SFR}\rho ^n`$, where $`\rho `$ is the gas density and $`n=1`$–2. As long as such a law is assumed and the infall of gas occurs continuously as a smooth function of time, the predicted SFH is also a smooth function of time.

Though the “classical” (i.e., smooth) infall model is widely accepted, there are observational data that suggest intermittent or oscillatory star formation activities in spiral galaxies. This means that the SFH is not a smooth function of time. Kennicutt, Tamblyn, & Congdon (1994) have shown that the ratio of present-to-past SFR in spiral sample has a significant scatter. More recently, Tomita, Tomita & Saitō (1996) have analyzed the far-infrared to $`B`$-band flux ratio $`f_{\mathrm{FIR}}/f_B`$ of 1681 spiral galaxies (see also Devereux & Hameed 1997). The indicator $`f_{\mathrm{FIR}}/f_B`$ represents the ratio between the present SFR and the averaged SFR over the recent Gyr. They have shown order-of-magnitude spread of $`f_{\mathrm{FIR}}/f_B`$ and suggested a violent temporal variation of the SFR.

The intermittence of the SFH in the Galactic disk is recently suggested by Rocha-Pinto et al. (2000a, hereafter R00). They have provided the SFH of the Galaxy inferred from the stellar age of the solar neighborhood, using 552 late-type stars. The age of each star has been estimated from the chromospheric emission in the Caii H and K lines (Soderblom, Duncan, & Johnson, 1991). After metallicity-dependent age correction, completeness correction, and scale-height correction111The scale height is dependent on the age of stars., they have derived the age distribution of the stars. Then, after correcting for evolved stars, they have derived the SFH. They have also asserted that their SFH derived from the stars in the solar neighborhood is representative of the SFH in the whole disk, since the diffusion timescale of stars is much shorter than the Galactic age. The discussion in this paper is based on Figure 2 of R00. Based on their data, R00 have suggested that the star formation activity of the disk is intermittent or variates violently. Their suggestion has been statistically confirmed by Takeuchi & Hirashita (2000, hereafter TH00), who have also shown that the typical timescale of the variation is 2 Gyr. We note that Hernandez, Valls-Gabaud, & Gilmore (2000) have also found an oscillatory component of the SFH in the solar neighborhood.

Theoretically, the intermittent or oscillatory SFH is easily reproduced if we treat the ISM as a nonlinear open system (Ikeuchi, 1988). Ikeuchi & Tomita (1983, hereafter IT83) have considered the ISM composed of three phases (cold, warm, and hot) as suggested by McKee & Ostriker (1977) and modeled the time evolution of the fractional masses of the three components (see also Habe, Ikeuchi, & Tanaka 1981). Since the mass exchange among the three components is a nonlinear process, the limit-cycle evolution of the fractional masses can emerge (see also Scalo & Struck-Marcell 1986; Korchagin, Ryabstev, & Vorobyov 1994). The limit-cycle behavior is supported by Kamaya & Takeuchi (1997), who have interpreted the various levels of the star formation activities in spiral galaxies shown observationally by Tomita et al. (1996) in the framework of the nonlinear open system model. Their interpretation is based on the galaxy-wide limit-cycle behavior of the ISM.

Another interesting topic is the chemical evolution in such a limit-cycle ISM. If the oscillatory SFH is considered, we may find an oscillation in a chemical enrichment process. For example, the age–metallicity relation of the stars in the Galactic disk may scatter because of the oscillation. We will examine quantitatively such a scatter caused by the limit-cycle evolution.

In this paper, we model the oscillatory SFH in the Galactic disk proposed by R00 by combining the infall model and the limit-cycle model. The chemical evolution in the oscillatory SFH is also investigated. This paper is organized as follows. First, in § 2 we model the chemical evolution of the Galactic disk by using the infall model. Then, in § 3 we review the limit-cycle model of ISM. Some results derived from the equations are described in § 4. Finally, we discuss the SFH and the chemical evolution in the limit-cycle ISM in § 5.

## 2 CHEMICAL EVOLUTION MODEL

The chemical evolution model of the Galactic disk is constructed here. The model is based on the infall model, which is characterized by the gradual gas infall from the halo. We adopt a one-zone model for simplicity. In other words, the phenomena of the ISM are averaged in space. This simple treatment is advantageous because the response of the chemical evolution on the parameters is easy to examine. When we compare the result with the observational data, however, we should be careful whether the data are averaged or not. We comment on the one-zone approximation in § 3.4.

### 2.1 Gas and Metal

The changing rates of the gas mass ($`M_\mathrm{g}`$) and metal mass ($`M_i`$, where $`i`$ denotes the species of the metal; $`i=\mathrm{Fe}`$, O, etc.) in the Galactic disk are described by a set of differential equations

$`{\displaystyle \frac{dM_\mathrm{g}}{dt}}`$ $`=`$ $`\psi +E+F,`$ (1)
$`{\displaystyle \frac{dM_i}{dt}}`$ $`=`$ $`X_i\psi +E_i+FX_i^\mathrm{f},`$ (2)
where $`\psi `$ is the SFR, $`E`$ is the total injection rate of gas from stars, $`F`$ is the rate of gas infall from halo, $`X_i`$ is the abundance of $`i`$ (i.e., $`X_iM_i/M_\mathrm{g}`$), $`E_i`$ is the injection rate of element $`i`$ from stars, and $`X_i^\mathrm{f}`$ is the abundance of $`i`$ in the infall material (see e.g., Tinsley 1980 for the basic treatment of chemical evolution of galaxies). Introducing $`X_i^\mathrm{f}`$ enables us to treat the infall of pre-enriched gas. An early enrichment in the halo may be important for the initial phase of the disk formation (e.g., Ikuta & Arimoto 1999).

In this paper, we choose two tracers for the metallicity, O and Fe. Almost all the oxygen is produced by high-mass stars, while the iron is produced mainly by Type Ia supernovae (SNe) as well as by high-mass stars. Thus, we include the contribution from Type Ia SNe in our formulation. We adopt the combination of the instantaneous recycling approximation and the delayed production approximation as formulated by Pagel & Tautvaišienė (1995). The evolution of Fe abundance based on a model of Type Ia SN has been considered in Kobayashi et al. (1998). With these approximations, $`E`$ and $`E_i`$ at $`t`$ are expressed by using the SFR as a function of time, $`\psi (t)`$:

$`E`$ $`=`$ $`R_{\mathrm{ins}}\psi (t)+R_{\mathrm{del}}\psi (t\tau ),`$ (3)
$`E_i`$ $`=`$ $`[R_{\mathrm{ins}}X_i(t)+Y_{i,\mathrm{ins}}]\psi (t)+[R_{\mathrm{del}}X_i(t\tau )+Y_{i,\mathrm{del}}]\psi (t\tau ),`$ (4)
where $`R`$ and $`Y_i`$ are the returned fraction of gas from stars and the fractional mass of the newly formed element $`i`$, respectively, and the subscripts “ins” and “del” denote the instantaneous recycling and delayed production parts, respectively. The delay time $`\tau `$ is set to 1.3 Gyr, according to § 3 of Pagel & Tautvaišienė (1995).

We assume that the disk begins to form at $`t=0`$. Its age is assumed to be 15 Gyr. When $`t\tau <0`$, all the functions whose arguments depend on $`(t\tau )`$ are set to zero; for example, $`\psi (t\tau )=0`$ when $`t<\tau `$. The initial condition is summarized in § 3.5.

With the above expressions (eqs.  and ), equation (1) becomes

$`{\displaystyle \frac{dM_\mathrm{g}}{dt}}=(1R_{\mathrm{ins}})\psi (t)+R_{\mathrm{del}}\psi (t\tau )+F(t),`$ (5)
and the combination of equations (2) and (5) leads to

$`M_\mathrm{g}{\displaystyle \frac{dX_i}{dt}}`$ $`=`$ $`Y_{i,\mathrm{ins}}\psi (t)+Y_{i,\mathrm{del}}\psi (t\tau )+R_{\mathrm{del}}\psi (t\tau )[X_i(t\tau )X_i(t)]`$ (6)
$``$ $`F(t)[X_i(t)X_i^\mathrm{f}].`$
### 2.2 Infall Rate of Gas

For the infall rate $`F`$, we follow TH00. They assumed an exponential form for the infall rate

$`F(t)={\displaystyle \frac{M_0}{t_{\mathrm{in}}}}\mathrm{exp}(t/t_{\mathrm{in}}),`$ (7)
where $`M_0`$ indicates the total mass that can fall into the galaxy. In other words,

$`{\displaystyle _0^{\mathrm{}}}F(t)𝑑t=M_0.`$ (8)
Normalizing equation (5) by $`M_0`$ leads to

$`{\displaystyle \frac{df_\mathrm{g}}{dt}}=(1R_{\mathrm{ins}})\stackrel{~}{\psi }(t)+R_{\mathrm{del}}\stackrel{~}{\psi }(t\tau )+{\displaystyle \frac{1}{t_{\mathrm{in}}}}\mathrm{exp}(t/t_{\mathrm{in}}),`$ (9)
where

$`f_\mathrm{g}{\displaystyle \frac{M_\mathrm{g}}{M_0}}\text{and}\stackrel{~}{\psi }{\displaystyle \frac{\psi }{M_0}}.`$ (10)
Equation (6) is also normalized by $`M_0`$ as

$`f_\mathrm{g}{\displaystyle \frac{dX_i}{dt}}`$ $`=`$ $`Y_{i,\mathrm{ins}}\stackrel{~}{\psi }(t)+Y_{i,\mathrm{del}}\stackrel{~}{\psi }(t\tau )+R_{\mathrm{del}}\stackrel{~}{\psi }(t\tau )[X_i(t\tau )X_i(t)]`$ (11)
$``$ $`{\displaystyle \frac{X_i(t)X_i^\mathrm{f}}{t_{\mathrm{in}}}}\mathrm{exp}(t/t_{\mathrm{in}}).`$
### 2.3 Star Formation Law

In order to include the three-phase model of the ISM composed of cold, warm and hot components, we modify the Schmidt law with the index $`n=1`$ (Schmidt 1959) as

$`\psi =M_\mathrm{g}X_{\mathrm{cold}}/t_{},`$ (12)
where $`X_{\mathrm{cold}}`$ is the mass ratio of the cold component to the total gas mass ($`M_\mathrm{g}`$) and $`t_{}`$ is the timescale of the cold-gas consumption to form stars. In other words, we consider that stars are formed from cold clouds on a gas consumption timescale of $`t_{}`$. The time evolution of $`X_{\mathrm{cold}}`$ will be modeled based on IT83, which treated the ISM as a nonlinear open system, in § 3. Equation (12) is equivalent to

$`\stackrel{~}{\psi }=f_\mathrm{g}X_{\mathrm{cold}}/t_{}.`$ (13)
### 2.4 Choice of the Parameters of the Chemical Evolution

According to Pagel & Tautvaišienė (1995), we choose the fractional masses of newly formed elements (see eq. ) as follows: $`Y_{\mathrm{O},\mathrm{ins}}/X_\mathrm{O}=0.70`$, $`Y_{\mathrm{O},\mathrm{del}}/X_\mathrm{O}=0.0`$, $`Y_{\mathrm{Fe},\mathrm{ins}}/X_{\mathrm{Fe}}=0.28`$, and $`Y_{\mathrm{Fe},\mathrm{del}}/X_{\mathrm{Fe}}=0.42`$, where the subscript $``$ indicates the solar value. They explained the observed metallicity of Galactic stars along with an infall model. For the abundances in the inflow gas, we examine two cases: one is the primordial case, $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0,\mathrm{\hspace{0.17em}0})`$; the other is the pre-enriched case, $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0.1,\mathrm{\hspace{0.17em}0.25})`$222This means that $`[\mathrm{Fe}/\mathrm{H}]=1`$ and $`[\mathrm{O}/\mathrm{Fe}]=0.4`$ in the inflow material. We note $`[\mathrm{Fe}/\mathrm{H}]=\mathrm{log}(X_{\mathrm{Fe}}/X_{\mathrm{Fe}})`$, $`[\mathrm{O}/\mathrm{H}]=\mathrm{log}(X_\mathrm{O}/X_\mathrm{O})`$, and $`[\mathrm{O}/\mathrm{Fe}]=\mathrm{log}(X_\mathrm{O}/X_\mathrm{O})\mathrm{log}(X_{\mathrm{Fe}}/X_{\mathrm{Fe}})`$.. The parameter in the delayed production approximation, $`\tau `$, is set as $`\tau =1.3`$ Gyr. The returned fractions of gas, $`R_{\mathrm{ins}}`$ and $`R_{\mathrm{del}}`$, are determined as follows: $`R_{\mathrm{ins}}=0.16`$, and $`R_{\mathrm{del}}=0.13`$. The details about these values are described in Appendix A. We will determine $`t_{}`$ and $`t_{\mathrm{in}}`$ in § 3.2.

## 3 LIMIT-CYCLE MODEL OF THE ISM

We model the oscillatory behavior of the Galactic SFH proposed by R00. IT83 have shown that an oscillatory behavior of the fractional masses of three components (cold, warm, and hot) emerges if one considers the ISM to be a nonlinear open system. An introduction and details concerning nonlinear open systems are found in Nicolis & Prigogine (1977). We adopt the model by IT83 to explain the oscillatory SFH by R00.

As long as the infall timescale ($`t_{\mathrm{in}}`$) is much longer than the oscillatory timescale, the effect of the infall on the limit-cycle evolution is not significant. Indeed, as shown in Table 1, for the present case the case where $`t_{\mathrm{in}}9`$ Gyr, which is much longer than the period of oscillation proposed by R00 ($`1`$ Gyr; see also Takeuchi & Hirashita 2000). Thus, we can apply the original model by IT83, which did not include the effect of infall.

### 3.1 Model Equations

We review the model by IT83. This model is used to calculate the time evolution of the mass fraction of the three ISM phases (see also § 1). The result is used to calculate the SFR through equation (13).

The ISM is assumed to consist of three components (McKee & Ostriker 1977); the hot rarefied gas ($`T10^6`$ K, $`n10^3`$ cm-3), the warm gas ($`T10^4`$ K, $`n10^1`$ cm-3), and cold clouds ($`T10^2`$ K, $`n10`$ cm-3). The fractional masses of the three components are $`X_{\mathrm{hot}}`$, $`X_{\mathrm{warm}}`$, and $`X_{\mathrm{cold}}`$, respectively. A trivial relation holds:

$`X_{\mathrm{hot}}+X_{\mathrm{warm}}+X_{\mathrm{cold}}=1.`$ (14)
The following three processes are considered (see IT83 and Ikeuchi 1988 for the details):  the sweeping of a warm gas into a cold component at the rate of $`aX_{\mathrm{warm}}`$ ($`a5\times 10^8`$ yr-1);  the evaporation of cold clouds embedded in a hot gas at the rate of $`bX_{\mathrm{cold}}X_{\mathrm{hot}}^2`$ ($`b10^7`$$`10^8`$ yr-1);  the radiative cooling of a hot gas through collisions with a warm gas at the rate of $`cX_{\mathrm{warm}}X_{\mathrm{hot}}`$ ($`c10^6`$$`10^7`$ yr-1). Writing down the rate equations and using equation (14), we obtain

$`{\displaystyle \frac{dX_{\mathrm{cold}}}{d\tau }}`$ $`=`$ $`BX_{\mathrm{cold}}X_{\mathrm{hot}}^2+A(1X_{\mathrm{cold}}X_{\mathrm{hot}}),`$ (15)
$`{\displaystyle \frac{dX_{\mathrm{hot}}}{d\tau }}`$ $`=`$ $`X_{\mathrm{hot}}(1X_{\mathrm{cold}}X_{\mathrm{hot}})+BX_{\mathrm{cold}}X_{\mathrm{hot}}^2,`$ (16)
where $`\tau ct`$, $`Aa/c`$, and $`Bb/c`$.

The solutions of equations (15) and (16) are classified into the following three types, according to the values of $`A`$ and $`B`$ (IT83):

1. $`A>1`$; all the orbits in the $`(X_{\mathrm{cold}},X_{\mathrm{hot}})`$-plane reduce to the node (0, 1),

2. $`A<1`$ and $`B>B_{\mathrm{cr}}`$; all the orbits reduce to a stable focus $`[(1A)/(AB+1),A]`$,

3. $`A<1`$ and $`B<B_{\mathrm{cr}}`$; all the orbits converge on a limit-cycle orbit,

where $`B_{\mathrm{cr}}(12A)/A^2`$. Apparently, case 3 is important in the interpretation of the oscillatory SFH shown in R00. Thus, we choose the parameters that satisfy case 3 as will be described in the next subsection.

### 3.2 Choice of Parameters of the Limit-Cycle Model

Since the timescale of the variation of SFR derived by R00 is $`1`$ Gyr (see also TH00), we first investigate whether the period of the limit-cycle can be the order of 1 Gyr. Indeed, a Gyr-timescale cycle is possible in the natural parameter range. According to Figure 3 of IT83, the period can be $`10^2/c`$ when we choose $`A=0.3`$ and $`A=0.5`$. Since $`1/c`$ is of the order of $`10^6`$$`10^7`$ yr, $`10^2/c1`$ Gyr is possible. Thus, we choose $`A=0.3`$ or $`A=0.5`$ and $`1/c=10^7`$ yr to demonstrate the Gyr-scale oscillation of the Galactic SFH. The subsequent discussions are unchanged if we adopt another set of the parameters that satisfies the oscillation period of $`1`$ Gyr.

Here we confirm that the adopted parameters are within the reasonable range of the physical properties of the ISM in the Galactic disk. First, Ikeuchi & Tomita (1983) estimated $`a`$ from the SN rate and the maximum radius of an SN remnant (SNR) and obtained $`a5\times 10^8`$ yr-1. Next, we estimate $`b`$ as the reciprocal of the evaporation timescale of a cold cloud. The evaporation timescale estimated in Hirashita (2000a) may be applicable in the present case and we obtain $`b10^7`$ yr-1. Finally, $`c`$ is estimated from the collision rate of a cold cloud with SNRs. The collision rate $`t_{\mathrm{col}}`$ is estimated as $`t_{\mathrm{col}}(\pi R_{\mathrm{SNR}}^2n_{\mathrm{SNR}}v)^1`$, where $`R_{\mathrm{SNR}}`$ is the typical size of a SNR, $`n_{\mathrm{SNR}}`$ is the number density of SNRs in the interstellar space, and $`v`$ is the typical relative velocity between a SNR and a cloud. If we put $`R_{\mathrm{SNR}}=50`$ pc, $`n_{\mathrm{SNR}}=10^6`$ pc-3,333The values of $`R_{\mathrm{SNR}}`$ and $`n_{\mathrm{SNR}}`$ are estimated according to Ikeuchi & Tomita (1983). We assume that the typical lifetime of a SNR is $`10^7`$ yr. and $`v=100`$ km s-1, we obtain $`t_{\mathrm{col}}10^7`$ yr. This means that $`ct_{\mathrm{col}}^110^7`$ yr-1. $`A=0.3`$ and $`B=0.5`$ are easily satisfied if we assume for example $`a=6\times 10^8`$ yr, $`b=10^7`$ yr, and $`c=2\times 10^7`$ yr, all of which are consistent with the above order-of-magnitude estimates.

TH00 adopted a star formation law $`\stackrel{~}{\psi }=f_\mathrm{g}/t_{\mathrm{sf}}`$ and did not consider the effect of a multi-component medium with phase changes. In order to use their choices of the parameter values $`A`$ and $`1/c`$ we must relate their $`t_{\mathrm{sf}}`$ (the timescale in the classical smooth infall model) to our $`t_{}`$ (the one in the oscillatory infall model). This is achieved by averaging the gas consumption rate defined as $`\stackrel{~}{\psi }/f_\mathrm{g}=X_{\mathrm{cold}}/t_{}`$ (eq. ) over the whole galactic age where the oscillatory part of the efficiency is smoothed out and becomes $`1/t_{\mathrm{sf}}`$. In other words,

$`\stackrel{~}{\psi }/f_\mathrm{g}=X_{\mathrm{cold}}/t_{}=1/t_{\mathrm{sf}},`$ (17)
where $``$ indicates the time average of the quantity over the Galactic history ($`0<t<15`$ Gyr). Since $`X_{\mathrm{cold}}=0.65`$ for $`(A,B)=(0.3,\mathrm{\hspace{0.17em}0.5})`$ (Appendix B), $`t_{}=0.65t_{\mathrm{sf}}`$. Following TH00, we examine the three sets of $`(t_{\mathrm{sf}},t_{\mathrm{in}})`$ as summarized in Table 1. They determined the parameters by fitting their infall model to the observed SFH proposed by R00. All the three models provide an almost identical SFH and reproduce the smoothed trend of the SFH (T00). Thus, it is meaningful to examine all the three cases. However, Model A gives an infall timescale much larger than that in e.g., Matteucci & François (1989), although they assumed the same time dependence of infall as that in this paper. This long timescale indicates that the infall rate is almost constant over the history of the Galactic disk. We note that Model A predicts the highest metallicity of the three models (§ 4).

### 3.3 Treatment of the Delayed Production

Here we should comment on the delayed production approximation. It assumes that all the delayed production at $`t`$ is determined by the SFR at $`t\tau `$. In reality, $`\tau `$ differs among Type Ia SNe. Thus, the delayed production is determined by the averaged SFR around $`t\tau `$. Expecting that the lifetimes of Type Ia progenitors are comparable to, or longer than, 1 Gyr ($`\text{the period of the SFR oscillation in the Galaxy}`$) the averaged SFR around $`t\tau `$ is described as $`M_\mathrm{g}(t\tau )X_{\mathrm{cold}}/t_{}=M_\mathrm{g}(t\tau )/t_{\mathrm{sf}}`$, where $`M_\mathrm{g}(t\tau )`$ is the gas mass at $`t\tau `$ and equation (17) is used. Thus, we hereafter assume that

$`\stackrel{~}{\psi }(t\tau )={\displaystyle \frac{f_\mathrm{g}(t\tau )}{t_{\mathrm{sf}}}}.`$ (18)
### 3.4 Comment on One-Zone Treatment

Following IT83, the structure of a model galaxy is approximated by one zone. The simplicity of the one-zone approximation gives the advantage that the background physical processes are easy to see. In this subsection, we discuss the one-zone treatment.

Habe et al. (1981) stated in their § 7 that for the one-zone assumption to be acceptable it is necessary that the mean distance between supernova remnants (SNRs) be less than 100 pc if a characteristic lifetime of SNRs of $`\tau _{\mathrm{life}}10^7`$ yr and a mean expansion velocity of 10 km s-1 are adopted. This is because the SNRs should affect the whole disk for the one-zone treatment. A distance of less than 100 pc means that there are $`N10^4`$ SNRs in a galaxy disk, if the disk size of 10 kpc is assumed. This number is reasonable if SNe occur every $`10^3`$ yr ($`\tau _{\mathrm{life}}/N10^7`$ \[yr\]/$`10^4`$). Considering that the SN rate in a spiral galaxy is typically 1/100–1/50 yr-1 (Cappellaro et al. 1993), the mean distance between SNRs is less than 100 pc even if 10–20 massive stars are clustered in a region.

The information of a place can travel over a distance of $`c_{\mathrm{s},\mathrm{eff}}t_{\mathrm{cyc}}`$ in one period of the limit cycle, where $`c_{\mathrm{s},\mathrm{eff}}`$ is the effective sound speed in the ISM and $`t_{\mathrm{cyc}}`$ is the oscillatory period. Estimating these quantities as $`c_{\mathrm{s},\mathrm{eff}}10`$ km s-1 (a typical velocity dispersion in the interstellar space) and $`t_{\mathrm{cyc}}1`$ Gyr, we obtain $`c_{\mathrm{s},\mathrm{eff}}t_{\mathrm{cyc}}=10`$ kpc. Thus, the information of the limit-cycle behavior can propagate over the whole disk. Thus, the assumption of the limit-cycle behavior over the whole disk may be good. In this paper, as a first step, we treat the model galaxy as being a one-zone object. Here we should note that the gas transport in the radial direction is difficult if we consider the angular momentum conservation. This difficulty may be a cause of the radial gradient of the metallicity and the gas-to-star fraction in spiral galaxies.

The above discussions are not a satisfactory “proof” for the limit-cycle oscillation on the scale of the whole Galactic disk. (But it is important that it is not rejected.) At present, thus, the cyclic star formation over the whole disk is an assumption that easily explains the variation of the star formation activity observed in the solar neighborhood by R00. We note that it also explains the variety of the star formation activity of spiral galaxies (Kamaya & Takeuchi 1997). Hence, in this paper, we base our discussion on the limit-cycle behavior on a whole-disk scale.

When we compare our result with the observational data, we should carefully consider to what extent the data is averaged. The range of Galactocentric radii that enter in the average depends on the age. Considering that the diffusion of stellar orbits on a scale of 1 kpc occurs in 0.2 Gyr (Wielen, 1977), it is reasonable to assume that the observational quantities are averaged on a scale of more than 1 kpc. Thus, as a first step, we adopt the one-zone treatment for the Galactic disk to see the chemical evolution in an oscillatory SFH. We should extend our model to multi-zone treatment as Chiappini, Matteucci, & Gratton (1997) (see also Romano et al. 2000 for a recent work) in the future. Observationally, the formalism in Meusinger (1991) may be useful in order to link the global observed SFH with the SFH in different Galactocentric annuli.

### 3.5 Initial Conditions

The initial condition is set as follows: $`X_\mathrm{O}(t=0)=0`$, $`X_{\mathrm{Fe}}(t=0)=0`$, $`f_\mathrm{g}(t=0)=0`$, $`X_{\mathrm{cold}}(t=0)=0.1`$, and $`X_{\mathrm{hot}}(t=0)=0.7`$. The convergence to the limit-cycle occurs on the timescale of a few periods. The results are however not strongly dependent on the choice of the initial conditions for $`X_{\mathrm{cold}}`$ and $`X_{\mathrm{hot}}`$.

## 4 RESULTS

In this section, the results calculated from the equations above are presented. They are compared with the observational data. Before displaying the results, we review the solving processes of the equations. First, the mass fraction of the cold gas, $`X_{\mathrm{cold}}`$, is calculated by equations (15) and (16). $`X_{\mathrm{cold}}`$ is used to determine the SFR at $`t`$ through equation (13). As for the delayed contribution expressed in $`\stackrel{~}{\psi }(t\tau )`$, we take into account the scatter of the lifetimes of progenitors of Type Ia SNe (§ 3.3) and assume equation (18). The SFH and the chemical evolution are modeled by the infall model. The evolution of the gas mass normalized by the total available gas mass $`M_0`$ is calculated by equations (9). For $`t_{\mathrm{sf}}`$ and $`t_{\mathrm{in}}`$, we examine three cases listed in Table 1. These three cases are also examined in TH00. The chemical evolution is calculated by equation (11).

### 4.1 Star Formation History

The SFH calculated by our model is presented in Figure 1a. Since the three models predict almost the same SFH as indicated by TH00, we present only the result of Model A in Table 1. We also show the SFH observationally determined by R00 in Figure 1b in order to demonstrate the qualitative similarity between the model prediction and the observation.

### 4.2 Metallicity Evolution

We test the model with the metallicity data of Galactic stars. Age–metallicity relation, G-dwarf metallicity distribution, and \[Fe/O\]–\[Fe/H\] relation are examined. First of all, we should note that the yields may be uncertain because of the treatment of convection, nuclear reaction rates, mass loss in the asymptotic giant branch phase, etc. If the yields are systematically larger/smaller than assumed in this paper, the metallicities predicted by our model should be systematically larger/smaller. Thus, quantitative agreement by fine tuning of the parameters might be meaningless. However, the qualitative behavior of the metallicity evolution in the limit-cycle ISM is not altered even if the yield changes.

#### 4.2.1 Age–metallicity relation

The age–metallicity relation of stars in the Galactic disk provides us with information on its chemical enrichment history. Thus, our model is worth testing by using the age–metallicity relation of the stars in the solar neighborhood. The sample is provided by Rocha-Pinto et al. (2000b), which used the same sample as R00.

First, we examine the case where the infall gas is of primordial abundance \[i.e., ($`X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},`$ $`X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0,\mathrm{\hspace{0.17em}0})`$\]. The age–metallicity relation predicted by our model is shown in Figure 2a. The solid, dotted, and dashed lines represent the results in Models A, B, and C, respectively. We also present the observational data of the age–metallicity relation by Rocha-Pinto et al. (2000b) (see their Table 3). Model A predicts the highest present metallicity, since its short gas consumption timescale leads to the most efficient chemical enrichment. Even in Model A, however, the discrepancy between the model prediction and data is significant in the low-metallicity range. Rocha-Pinto et al. (2000b) noted the high initial metallicity of the disk and attributed it to the pre-enrichment of the gas before the formation of the first stars in the disk. They also have shown that the age–metallicity relation determined from the chromospheric age can deviate upward (i.e., metallicity is overestimated) for larger ages from the real relation because of the uncertainty in the age estimation. Thus, we also examine the case where the infalling gas in enriched with metal \[$`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0.1,\mathrm{\hspace{0.17em}0.25})`$\]. The result is shown in Figure 2b. We see that the discrepancy between the model prediction and the observational data is reduced. From the viewpoint of modeling, the yield is also uncertain. Thus, we do not try any fine-tuning the age–metallicity relation.

In spite of such uncertainty, we can discuss the qualitative behavior of the age–metallicity relation. The amplitude of the oscillation of metallicity is smaller than the typical scatter of the observed data points ($`0.3`$ dex). The small amplitude is natural, because metallicity is determined by all the past history of star formation, and thus the present oscillatory star formation does not significantly contribute to the metallicity.

Rocha-Pinto et al.’s data shown in Figure 2 seem to show a recent increase in \[Fe/H\]. The pre-enriched infall cannot solve this increase, since the infalling gas has too low a metallicity. However, we should carefully examine whether the most recent data point in Figure 2 represents the SFH of the whole Galactic disk, because if the recent chemical enrichment rate in the solar neighborhood is significantly higher than that in the whole Galaxy, the most recent data point would naturally show a systematically higher metallicity. Moreover, data sets shown by other authors do not necessarily show such an increase in the recent metallicity (e.g., Twarog 1980). It is necessary to analyze the observational data further before we construct a theoretical model for the recent \[$`15t1`$ (Gyr)\] increase in metallicity.

We can expect that the oscillation behavior is more prominent for oxygen than for iron, because all the oxygen is produced from the “instantaneous” part (i.e., $`Y_{\mathrm{O},\mathrm{ins}}Y_{\mathrm{O},\mathrm{del}}`$). We present the result for the case of pre-enriched infall \[i.e., $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0.1,\mathrm{\hspace{0.17em}0.25})`$\] in Figure 3. Indeed, the amplitude of the oscillation is larger than Figure 2b. However, the oscillation would not explain the scatter of the oxygen abundance of the stars in the Galactic disk.

#### 4.2.2 Metallicity distribution

The G-dwarf metallicity distribution is also tested along with our model, since the primary motivation for the infall model is to solve the G-dwarf problem (e.g., Pagel 1997). The probability distribution function $`P(\mathrm{log}X_i)`$ of the metallicity is calculated from our model as

$`P(\mathrm{log}X_i)d\mathrm{log}X_i=C\stackrel{~}{\psi }dt,`$ (19)
where the constant $`C`$ is the normalization so that

$`C{\displaystyle _0^{15\mathrm{Gyr}}}\stackrel{~}{\psi }𝑑t=1.`$ (20)
From equation (19), we obtain the following analytical expression for $`P`$:

$`P(\mathrm{log}X_i)=C(\mathrm{ln}10)\stackrel{~}{\psi }(t)X_i(t)\left({\displaystyle \frac{dX_i}{dt}}\right)^1,`$ (21)
where $`dX_i/dt`$ is calculated from equation (11). In comparing the distribution function with the observational data, we should take into account the scatter of the data. Here, we simply convolve $`P`$ with a Gaussian kernel as

$`P_{\mathrm{conv}}(\mathrm{log}X_\mathrm{O}){\displaystyle _{\mathrm{}}^{\mathrm{}}}P(u){\displaystyle \frac{1}{\sqrt{2\pi }\sigma }}\mathrm{exp}\left[{\displaystyle \frac{(\mathrm{log}X_\mathrm{O}u)^2}{2\sigma ^2}}\right]𝑑u,`$ (22)
where we adopt $`\sigma =0.1`$ to compare with Rocha-Pinto & Maciel (1996). We adopt these data because we would like to use a sample of G-dwarfs, whose lifetimes are comparable to the age of the universe. In Figure 4, we show $`P_{\mathrm{conv}}`$ as a function of \[Fe/H\]. The solid, dotted, and dashed lines represent the result in Models A, B, and C, respectively. The histogram shows the data by Rocha-Pinto & Maciel (1996). The two figures (a) and (b) correspond to $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0,\mathrm{\hspace{0.17em}0})`$ and $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0.1,\mathrm{\hspace{0.17em}0.25})`$, respectively (same as Fig. 2a and b, respectively). We see that Model A in Figure 4b seems to be the best of all the models. However, considering the uncertainty in the yields, we do not try any fine tuning. The excess of the observed number of stars around $`[\mathrm{Fe}/\mathrm{H}]0.0`$ is consistent with the data by Rocha-Pinto et al. (2000b). As stated in § 4.2.1, this may be due to the recent significant enrichment in the solar neighborhood.

#### 4.2.3 Evolution of \[Fe/O\]

In order to test the effect of the limit-cycle behavior on the \[Fe/O\] ratio, we examine the relation between \[Fe/O\] and \[Fe/H\]. Since the oxygen is mainly produced by stars with short lifetimes, the effect of the limit-cycle ISM is reflected by the time evolution of the oxygen abundance. On the other hand, the iron is also produced by stars with long lifetimes and the information of the oscillation of ISM phase is lost in the iron abundance. Thus, we expect that \[Fe/O\] oscillates as the limit-cycle evolution of ISM.

In Figure 5, we show \[Fe/O\]–\[Fe/H\] relation. The solid, dotted, and dashed lines represent the results in Models A, B, and C, respectively. The two figures (a) and (b) correspond to $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0,\mathrm{\hspace{0.17em}0})`$ and $`(X_{\mathrm{Fe}}^\mathrm{f}/X_{\mathrm{Fe}},X_\mathrm{O}^\mathrm{f}/X_\mathrm{O})=(0.1,\mathrm{\hspace{0.17em}0.25})`$, respectively (same as Fig. 2a and b, respectively). Indeed, we see that \[Fe/O\] oscillates. However, the amplitude of the oscillation is not large. This is consistent with the age–metallicity relation (Fig. 2). As mentioned in § 4.2.1, the amount of metallicity reflects the past-integrated SFR; thus, as long as the mass of newly formed stars is not dominated in the total stellar mass, the present oscillation has little influence on the metallicity evolution.

## 5 DISCUSSIONS

In this paper, we have modeled the oscillatory SFH proposed observationally by R00. Our model is a combination of an infall model developed in the field of chemical evolution and the limit-cycle model proposed by Ikeuchi (1988) and his collaborators. We discuss our result in the following two subsections.

### 5.1 Limit-Cycle Star Formation History

The oscillatory behavior of the Galactic SFH proposed by R00 is modeled by using the limit-cycle model of SFH. The limit-cycle behavior of the three-phase ISM is suggested by Ikeuchi (1988) and his collaborators. Since the period of a limit-cycle orbit can be $`1`$ Gyr within the framework of Ikeuchi (1988), the Galactic SFH is explained by the limit-cycle model of the ISM.

Recently, Hirashita & Kamaya (2000) have explained the observed scatter of star formation activity of a sample of spiral galaxies by using the limit-cycle model. They provided a consistent modeling that explains the variation in the scatter of star formation activity among the morphological types (Sa–Sc) as shown in Kennicutt et al. (1994). Since the Galaxy is a spiral galaxy, an oscillatory SFH is consistent with the variation of star formation activity seen in other spirals.

Kennicutt et al. (1994) have also presented the ratio of the present SFR to the past averaged SFR (indicated as $`b`$ there). From Figure 2 of R00, we see that $`b`$ (denoted as $`\mathrm{SFR}/\mathrm{SFR}`$ in R00) can be as large as 2–3 for the Galactic SFH. Since the morphological type of the Galaxy is Sbc (e.g., Binney & Merrifield 1998, p. 171), $`b=2`$–3 is within the range of the Sbc/Sc sample in Kennicutt et al. (1994, their Fig. 6). This consistency implies that the oscillatory SFH may be a common nature for all the spiral galaxies.

### 5.2 Chemical Evolution in Limit-Cycle ISM

The chemical evolution of the Galaxy has been investigated in the framework of the limit-cycle ISM model. We have found that the amplitude of the oscillatory behavior of the metallicity is smaller than the observed scatter (Figs. 2 and 5). This indicates that the observed scatter is not attributed to the limit-cycle behavior. The scatter might be explained by chemical inhomogeneity in the Galactic disk.

The oscillatory behavior of the metallicity is not prominent because the metallicity reflects all the past SFH. The integrated contribution from all the past SFH smoothes out the oscillatory behavior of SFR. Thus, from the viewpoint of chemical evolution we conclude that we cannot distinguish between the “smooth” infall model without an oscillatory behavior and the oscillatory infall model proposed in this paper.

Contrary to the metallicity, the dust-to-gas ratio can show a oscillatory behavior (Hirashita, 2000b). This is because the efficiency of the dust formation changes according to phase changes in the gas. Moreover, dust is efficiently destroyed when the mass fraction of the cold gas is small. This oscillation of dust amount may be important for the evolution of infrared luminosity of galaxies.

### 5.3 Another Possible Mechanism for the Variation of SFR

Rocha-Pinto et al. (2000c) gave some indication that the Magellanic Clouds could play a role in the SFH of the Galaxy. It is meaningful to explore their idea. Since our discussion is based on the limit-cycle behavior inherent in the ISM as stated in § 1, we do not include the external force into our formulation. Fortunately, Ikeuchi & Tomita (1983) have investigated the behavior of the ISM in the presence of such an external force. Thus, we discuss the influence of the Magellanic Clouds based on the discussion in Ikeuchi & Tomita (1983).

If a perturbation of the external force to the limit-cycle ISM exists, both the period and the amplitude of the oscillation are affected. Thus, it is not necessary that the period proper to the ISM is about 1 Gyr. Even in the stable-focus case (§ 3.1), an oscillation emerges. Furthermore, the qualitative behavior can be changed: the system can show a chaotic orbit. Since the observational data do not reject such a chaotic orbit affected by the perturbation of the Magellanic Clouds, the strongly variable SFH of the Galaxy might be due to the interaction with the Magellanic Clouds. However, we note that the large scatter of the star formation activity of the spiral sample (e.g., Tomita et al. 1996) implies a general oscillatory behavior of ISM in spiral galaxies. This is naturally explained if such an oscillatory behavior is caused by the limit-cycle evolution inherent in the ISM.

Finally, we would like to note that Chiappini et al. (1997) have also considered SFR variation due to a different mechanism. Their rapidly variable SFR is caused by a density threshold for star formation: They assumed that star formation occurs only if the surface density of the gas exceeds a critical value, while we have assumed no threshold. However, the timescale of the SFR variation in Chiappini et al. is much shorter than 1 Gyr. Although Chiappini et al.’s mechanism may indeed present on a timescale much shorter than 1 Gyr, it is necessary to introduce a mechanism different from Chiappini et al. in order to explain the SFR variability presented by R00. Thus, we have proposed a limit-cycle scenario for the SFR in the Galactic disk.

We thank H. Rocha-Pinto, the referee, for invaluable comments and suggestions that improved this paper very much. We are grateful to S. Mineshige and H. Shibai for continuous encouragement, H. Kamaya and A. Ibukiyama for useful comments, and K. Yoshikawa for an excellent computational environment. H. H. was supported by the Research Fellowship of the Japan Society for the Promotion of Science for Young Scientists. A. B. acknowledges the hospitality of the department of astronomy at Kyoto University where this project was started and financial support from the Japan Society for the Promotion of Science. We fully utilized the NASA’s Astrophysics Data System Abstract Service (ADS).

## Appendix A A. DETERMINATION OF $`R_{\mathrm{ins}}`$ AND $`R_{\mathrm{del}}`$

$`R_{\mathrm{ins}}`$ and $`R_{\mathrm{del}}`$ are described as

$`R_{\mathrm{ins}}`$ $`=`$ $`{\displaystyle _{m_{\mathrm{l},\mathrm{ins}}}^{m_\mathrm{u}}}(mw_m)\varphi (m)𝑑m,`$ (A1)
$`R_{\mathrm{del}}`$ $`=`$ $`{\displaystyle _{m_{\mathrm{l},\mathrm{del}}}^{m_{\mathrm{l},\mathrm{ins}}}}(mw_m)\varphi (m)𝑑m,`$ (A2)
where $`\varphi (m)`$ is the initial mass function (IMF), which is normalized so that the integral of $`m\varphi (m)`$ in the full range of the stellar mass (0.1–100$`M_{}`$ in this paper) becomes unity; $`m_\mathrm{u}`$ is the upper mass cut-off of the stellar mass, and we here adopt $`m_\mathrm{u}=100M_{}`$; $`m_{\mathrm{l},\mathrm{ins}}`$ and $`m_{\mathrm{l},\mathrm{del}}`$ are set as $`5M_{}`$ and $`1M_{}`$, corresponding to the stellar lifetime of $`\tau `$ (the parameter for the delayed production) and the age of galaxies. If we adopt the Salpeter’s IMF ($`\varphi (m)m^{2.35}`$), we obtain $`R_{\mathrm{ins}}=0.16`$ and $`R_{\mathrm{del}}=0.13`$.

## Appendix B B. DETERMINATION OF $`X_{\mathrm{cold}}`$

The cold gas mass fraction averaged over the Galactic lifetime, $`X_{\mathrm{cold}}`$, is used in § 3.2. It is estimated as follows. First, the time evolution of the cold gas is calculated based on the limit-cycle model in § 3 by adopting the parameters as $`A=0.3`$, $`B=0.5`$, and $`1/c=10^7`$ yr. The time evolution of $`X_{\mathrm{cold}}`$ is shown in Figure B1. Next, we estimate $`X_{\mathrm{cold}}`$ by averaging $`X_{\mathrm{cold}}(t)`$ over the Galactic age ($`T_\mathrm{G}`$) as

$`X_{\mathrm{cold}}={\displaystyle \frac{1}{T_\mathrm{G}}}{\displaystyle _0^{T_\mathrm{G}}}X_{\mathrm{cold}}(t)𝑑t.`$ (B1)