text:
# Small Numerators Canceling Small Denominators: Is Dyson’s Hierarchical Model Solvable?

## Abstract

We present an analytical method to solve Dyson’s hierarchical model, involving the scaling variables near the high-temperature fixed point. The procedure seems plagued by the presence of small denominators as in perturbative expansions near integrable systems in Hamiltonian mechanics. However, in all cases considered, a zero denominator always comes with a zero numerator. We conjecture that these cancellations occur in general, suggesting that the model has remarkable features reminiscent of the integrable systems.

PACS: 05.50.+q, 11.10.Hi, 64.60.Ak, 75.40.Cx

In many physical problems involving nonlinear flows, a common strategy consists in constructing a system of coordinates where the flow becomes linear. The action-angle variables in Hamiltonian mechanics provide well-know examples of such a procedure. Whenever the angle variables can be constructed, they evolve linearly with time and expressing the original variables in terms of the new ones solves completely the original problem. The problems for which well-defined angle-action variables can be constructed (e.g. Kepler’s problem or the free rigid body) are very distinguished and called integrable systems. A large number of numerical experiment has lead us to believe that in a generic way, small perturbations destroy integrability. This point of view was was first inferred by H. Poincaré who pointed out the existence of small denominators in the canonical transformation designed to eliminate the angle dependence of a perturbed Hamiltonian.

In this letter, we discuss the question of small denominators for renormalization group (RG) flows. The variables which play the role of angle variables are the scaling variables introduced by Wegner . Near a fixed point, the RG flows can be linearized. The problem of expressing the physical quantities in terms of variables which transform as in the linear approximation when the non-linear terms are taken into account, is analogous to removing the angle dependence of a perturbed Hamiltonian. If the task can be carried through, one obtains analytical expressions for the RG flows. However, as we will show, small denominators appear. Does this mean that, as in Hamiltonian mechanics, in generic situations the construction will fail?

Surprisingly, we found in a numerical calculation performed with Dyson’s hierarchical model , that zero denominators were systematically canceled by zero numerators. These remarkable cancellations suggest that either the model considered is as distinguished as the integrable systems of classical mechanics or that there exists a general mechanism that allows us to circumvent the small denominator problem for RG flows. The example of the two-dimensional Ising model shows the importance of having a non-trivial model which can solved in closed form. The results presented below indicate that Dyson’s hierarchical can be solved analytically.

During the last decades, the RG method has been successfully applied to many important problems in field theory and statistical mechanics: the critical behavior of ferromagnets and superconductors, the confinement of quarks or the generation of mass for the W and Z bosons. However, its practical implementation is still a formidable enterprise. Typically, expansions are often available near fixed points, but not to orders large enough to allow one to extrapolate between fixed points. Unfortunately, the calculation of physical quantities (e. g. the magnetic susceptibility) beyond an order of magnitude estimation, requires a calculation of the flows in crossover regions. One has then to rely to Monte Carlo simulations to achieve this goal. One also needs to select a small set of interactions which closes reasonably well under RG transformation near both fixed points. Interesting examples of such lattice Monte Carlo calculations are given in Refs.  for scalar theory and  for gauge theories.

In order to get analytical results, further approximations are needed. One possibility consists in using hierarchical approximations such as the one derived by Wilson in Ref.  and resulted into the “approximate recursion formula”. In this approximation, only the local interactions get renormalized and the flow can be calculated from a simple integral formula. Retracing Wilson’s construction from the beginning, one can restore the other renormalizations perturbatively. In order to perform this task, one would like to have a closed form solution in the hierarchical approximation. In order to achieve this goal, we have proposed  to use the Fourier representation of the integral formula and to construct the scaling variables in this basis. In this process, we identified the existence of small denominators possibly ruining the whole approach. This problem can be avoided in special circumstances, for instance for flows starting exactly along the unstable direction of a non-trivial fixed point. With this restriction, we found  expansions with overlapping domains of convergence in the crossover region.

In the following, we discuss the problem of small denominators in the construction of the scaling variables near the high-temperature (HT) fixed point of Dyson’s hierarchical model. The treatment of the this model is mathematically similar to the one of the approximate recursion formula. However there exists a large literature on Dyson’s model , the non-trivial fixed point is known very precisely  and the HT expansion studied to a very large order .

For a description of this model as a spin model, we refer to Ref. , while the details of the derivation of the RG flows as expressed below can be found in . To make a long story short, all the information regarding the local interactions after $`n`$ RG transformations is encoded in a function $`R_n(k)=1+a_{n,1}k^2+a_{n,2}k^4+\mathrm{}`$. The logarithm of this function generates the connected zero-momentum Green’s functions at finite volume. The recursion formula reads

$$R_{n+1}(k)=C_{n+1}\mathrm{exp}\left[\frac{1}{2}\frac{^2}{k^2}\right]\left[R_n\left(\frac{\sqrt{c}k}{2}\right)\right]^2.$$
(1)
We fix the normalization constant $`C_n`$ so that $`R_n(0)=1`$. We use the parametrization $`c=2^{12/D}`$ which implies that a free massless field scales in the same way as in a usual $`D`$-dimensional theory. In this formulation, the temperature dependence has been absorbed in the initial $`R_0(k)`$. For an Ising measure, $`R_0(k)=\mathrm{cos}(\sqrt{\beta }k)`$, while in general, we have to numerically integrate the Fourier transform of the local measure to determine the coefficients of $`R_0(k)`$ expanded in terms of $`k`$. In general, $`a_{n,l}`$ is of order $`\beta ^l`$ in the HT expansion.

In the HT phase, polynomial truncations of order $`l_{max}`$ in $`k^2`$ provide rapidly converging approximations . The RG flows can be expressed in terms of the quadratic map

$$a_{n+1,l}=\frac{u_{n,l}}{u_{n,0}},$$
(2)
with

$$u_{n,\sigma }=\mathrm{\Gamma }_\sigma ^{\mu \nu }a_{n,\mu }a_{n,\nu },$$
(3)
and

$$\mathrm{\Gamma }_\sigma ^{\mu \nu }=(c/4)^{\mu +\nu }\frac{(1/2)^{\mu +\nu \sigma }(2(\mu +\nu ))!}{(\mu +\nu \sigma )!(2\sigma )!},$$
(4)
for $`\mu +\nu \sigma `$ and zero otherwise. We use “relativistic” notations. Repeated indices mean summation. The greek indices $`\mu `$ and $`\nu `$ go from $`0`$ to $`l_{max}`$, while latin indices $`i`$, $`j`$ go from 1 to $`l_{max}`$.

The diagonalization of the linear RG transformation near the HT fixed point is quite simple because it is of the upper triangular form. From Eq. (4), one finds the spectrum

$$\lambda _{(r)}=2(c/4)^r.$$
(5)
in agreement with Ref. . Using the matrix of right eigenvectors, $`_l^i\psi _i^r=\lambda _{(r)}\psi _l^r`$, we introduce new coordinates such that $`a_{n,l}=\psi _l^rh_{n,r}`$. This diagonalizes the linear RG transformation. Note that the form of the eigenvectors guarantees that $`h_{n,l}`$ is also of order $`\beta ^l`$. More details are given Ref. , where all the quantities related to the HT fixed point are dressed with a “tilde” omitted here . In summary, the RG flows in the new coordinates can be written as

$$h_{n+1,l}=\frac{\lambda _{(l)}h_{n,l}+\mathrm{\Delta }_l^{pq}h_{n,p}h_{n,q}}{1+\mathrm{\Lambda }^ph_{n,p}+\mathrm{\Delta }_0^{pq}h_{n,p}h_{n,q}},$$
(6)
with coefficients calculable from Eq. (4).

We now express the $`h_{n,l}`$ in terms of the scaling variables $`y_{n,1},\mathrm{},y_{n,l_{max}}`$ which transform under a RG transformation as $`y_{n+1,i}=\lambda _iy_{n,i}`$. If we can construct functions $`h_l`$ and $`y_l`$ such that $`h_{n,l}=h_l(𝐲_n)`$ and $`y_{n,l}=y_l(𝐡_n)`$, then we get a complete analytical expression of $`h_{n,l}`$ (which contains all the thermodynamical quantities) in terms of $`h_{0,l}`$ (which depends on the initial energy density):

$$h_{n,l}=h_l(\lambda _1^ny_1(𝐡_0),\lambda _2^ny_2(𝐡_0),\mathrm{}).$$
(7)
The feasibility of this approach is demonstrated for a one-dimensional example in Ref. .

We now discuss the construction of the $`h_l`$. We use the expansion

$$h_l=y_l+\underset{i_1,i_2,\mathrm{}}{}s_{l,i_1i_2\mathrm{}}y_1^{i_1}y_2^{i_2}\mathrm{},$$
(8)
where the sums over the $`i`$’s run from $`0`$ to infinity in each variable with at least two non-zero indices. In the following, we use the notation $`𝐢`$ for $`(i_1,i_2,\mathrm{})`$. Plugging the expansion into Eq. (6), and requiring that the one step advance is obtained by rescaling the scaling variables by their associated eigenvalue, we obtain

$$s_{l,𝐢}=\frac{N_{l,𝐢}}{\left(_m\lambda _{(m)}^{i_m}\lambda _r\right)}.$$
(9)
with

$`N_{l,𝐢}=`$ $`{\displaystyle \underset{𝐣+𝐤=𝐢}{}}(\mathrm{\Delta }_l^{pq}s_{p,𝐣}s_{q,𝐤}s_{l,𝐣}{\displaystyle \underset{m}{}}\lambda _{(m)}^{j_m}\mathrm{\Lambda }^ps_{p,𝐤})`$ (10)
$`{\displaystyle \underset{𝐣+𝐤+𝐫=𝐢}{}}s_{l,𝐣}{\displaystyle \underset{m}{}}\lambda _{(m)}^{j_m}\mathrm{\Delta }_0^{pq}s_{p,𝐤}s_{q,𝐫}.`$ (11)
For a given set of indices $`𝐢`$, we introduce the notation

$$_q(𝐢)=\underset{m}{}i_mm^q.$$
(12)
One sees that $`_0`$ is the degree of the associated monomial and $`_1`$ its order in the HT expansion (since $`y_l`$ is also of order $`\beta ^l`$). Given that all the indices are positive and that at least one index is not zero, one can see that if $`𝐣+𝐤=𝐢`$ then $`_q(𝐣)<_q(𝐢)`$ and $`_q(𝐤)<_q(𝐢)`$. Consequently, Eq. (10) provides a solution order by order in $`_0`$ or in $`_1`$ (since the r.h.s is always of lower order) assuming that the small denominator problem can be avoided.

Using the parametrization $`c=2^{12/D}`$, the zero denominators in Eq. (9) appear when

$$Dl(D+2)=D_0(D+2)_1.$$
(13)
If $`D`$ and $`D+2`$ have no common factors, this equation has non-trivial solutions for indices such that $`_0=(2+D)q+1`$, for $`q`$ a strictly positive integer and $`l`$ such that $`_1=Dq+l`$. If $`D`$ and $`D+2`$ have common factors, we can proceed in the same way but with these common factors removed from both $`D`$ and $`D+2`$.

We have investigated numerically the well-studied case $`D=3`$ with a polynomial truncation $`l_{max}=25`$. For practical reasons, we have limited our study the case where $`h_l`$ depends only on the two leading variables $`y_1`$ and $`y_2`$. This restriction is self-consistent since if we start for instance with $`y_3=0`$, the multiplicative renormalization of the scaling variables implies that this condition stays valid after $`n`$ iterations. Since $`i_3=i_4=\mathrm{}=0`$, we will use notations such as $`s_{l,i_1,i_2}`$ or $`h_l(y_1,y_2)`$. In addition, we have concentrated our attention to the large order behavior in the leading variable $`y_1`$ and calculated the first order corrections in the subleading variable $`y_2`$ If $`i_2=0`$, we have the following sets of small denominators: $`i_1=5q+1`$ and $`l=2q+1`$. If $`i_2=1`$, we have $`i_1=5q`$ and $`l=2q+2`$. We have calculated the numerators corresponding to these zero denominators for $`q=1`$, 2, 3 and 4.

We have checked our calculation of the $`h_l`$ with two different methods. First we have considered a configuration $`h_{n,l}=h_l(y_1,y_2)`$ for particular values of $`y_1`$ and $`y_2`$ and then calculated the one step backward configuration $`h_{n1,l}=h_l(\lambda _1^1y_1,\lambda _2^1y_2)`$. We then calculated the one step forward using these $`h_{n1,l}`$ and the exact formula Eq. (6). Comparison between the two $`h_{n,l}`$ for values of $`y_1`$ and $`y_2`$ varying between 0.1 and 0.001 and for various truncations in the power of $`y_1`$ and $`y_2`$ considered, showed errors scaling like the powers neglected with coefficients compatible with the order of magnitude of the coefficients involved in the expansion. Second, we used the expression $`y_1(h_l)`$ calculated in Ref.  up to order 11 in $`\beta `$ and checked $`y_1(h_l(y_1,y_2))=y_1`$ with errors ranging between $`10^{14}`$ and $`10^{16}`$ for the various coefficients of the higher order terms.

In Fig. 1 we show the absolute value of $`N_{3,j,0}`$. The calculations have been performed with different arithmetic precisions. We have used $`Mathematica`$ 4.0 with $`\$MaxPrecision`$ and $`\$MinPrecision`$ both set to the same value $`pr`$. We have considered the cases $`pr`$=20, 30 and 40. The graphs shows $`N_{3,6,0}`$ is more than twenty orders of magnitudes smaller than its peers when $`pr=20`$ and then drops by ten orders of magnitude each time the precision is increased by 10, while all the other values stay stable. This is a strong evidence for $`N_{3,6,0}=0`$

In Figs. 2 and 3, we show $`N_{l,i_1,0}`$ and $`N_{l,i_1,1}`$ calculated with $`pr`$=20. One sees that $`N_{2q+1,5q+1,0}`$ and $`N_{2q+2,5q,1}`$ are more than 20 orders of magnitude smaller than the naive interpolation. We have checked in each of these cases that the small value drops by ten order of magnitudes each time the precision is increased just as in Fig. 1.

In summary, for each zero denominator considered, we found a zero numerator. We thus conjecture that all numerators corresponding to zero denominators are zero. If this conjecture is correct, the $`s_{l,𝐢}`$ corresponding to non-zero denominators are calculable from Eq. (10) and those corresponding to zero denominators are undetermined. The calculations performed above have been done with these coefficients treated as undetermined. The figures have been drawn with these coefficients set to zero. This choice is not essential, we have considered other choices where the undetermined coefficients have been set to values of the same order as the coefficients above and below (in $`i_1`$) and reached identical conclusions. This is a non-trivial statement. For instance, the undetermined coefficient $`s_{3,6,0}`$ appears explicitly in the numerator of the equation for $`s_{5,11,0}`$, however it is multiplied by a very small number which drops when $`pr`$ is increased.

We have checked that the individual terms in the zero numerators were not zero. The fact that we were able to obtain very precise cancellations without any fine-tuning suggests the existence of closed form formulas or of a symmetry forbidding these terms. More generally, the equally spaced spectrum (on a log scale), the number of terms at a given level $`m`$ (the order in HT expansion) equal to the number of partitions of $`m`$  and the existence of “nested” ambiguities are somehow reminiscent of string theory. A possible starting point to discover these hypothetical symmetries would be to exploit the fact that since for instance $`y_3`$ and $`y_1^6`$ transform the same way under a RG transformation, there exist ambiguities in the construction of $`h_l`$ in terms of the scaling variables.

We should say a few words about the small denominators near other fixed points. The spectrum of the linearized RG transformation near the non-trivial fixed point  can be calculated numerically . Investigating the small denominators of the form $`\lambda _l^k\lambda _m^r`$ for the first twenty eigenvalues and for powers not larger than 100. The best solution we found was $`\lambda _2^9\lambda _4`$ with two parts in a thousand. A more detailed study is necessary to decide if the rate of decay of the coefficients is sufficient to take care of the smaller denominator which will appear at larger orders. On the other hand, the spectrum at the Gaussian fixed point  is $`\lambda _j=2c^j`$. There are many zero denominators in integer dimensions, e.g., $`\lambda _1=\lambda _2^2`$ for $`D=3`$. If the numerators are not zero, one can “repair”  the situation by considering $`n`$-dependent coefficients. This excludes a solution of the form of Eq. (7) but it generates logarithmic corrections which are necessary. These can even be observed in the large order of the HT temperature expansion in Ref. .

In conclusion, we have found remarkable cancellations of small denominators by small numerators. If these cancellations occur in general, it is possible to calculate analytically the thermodynamical quantities of Dyson’s hierarchical model in the HT phase using Eq. (10). Our results suggest that this model has features analogous to the integrable systems in Hamiltonian mechanics. As such, the model would stand out as a first approximation to be used in situations where the conventional perturbative expansions are not reliable.

This research was supported in part by the Department of Energy under Contract No. FG02-91ER40664. Y. M. thanks the Aspen Center for Physics for its hospitality in Summer 2000 while this work was in progress and for a conversation there with L. Kadanoff.


format: text

text:
# Silicate Emission in the TW Hydrae Association

## 1 INTRODUCTION

Since the discovery 15 years ago that a large fraction of all main sequence stars possessed excess infrared emission (Aumann, 1988), the study of these “Vega-type stars” and their precursors has grown enormously. With the detection of the edge-on debris disk in $`\beta `$ Pictoris (Smith and Terille, 1984), and the equally important discovery that the debris exhibited a silicate emission band whose structure was similar to solar system comets (Knacke et al., 1993) it was apparent that at least some of these systems might be true planetary systems, perhaps similar to what our own solar system was like in its youth. Since that time, disks and tori around other stars have been successfully imaged.

Recently, Kastner et al. (1997) identified a previously unrecognized nearby young stellar association that consisted of 5 T Tauri (TT) or T Tauri-like stars, including the classical TT star TW Hya. Designated as the TW Hya Association (TWA), the list of objects in this stellar group is now over a dozen objects and growing (Webb et al., 2000). Located at a distance of about 50 pc, TWA represents the nearest region of recent ($`<`$20 Myr) star formation. It is no longer embedded in a molecular cloud, yet one of its members, TW Hya itself, possesses both a circumstellar molecular cloud (CO, HCN, and CN have been detected; Kastner et al. (1997)) and an infrared excess. HR 4796A, which is also believed to be a member of TWA (Webb et al., 2000) has a dusty torus that has been recently resolved (Jayawardhana et al., 1998; Koerner et al., 1998; Schneider et al., 1999). Another member of TWA is HD 98800, a multiple star system with dust concentrated around one unresolved pair (Gehrz et al., 1999; Koerner et al., 2000). Combining a variety of criteria, Weintraub et al. (2000) derive an age of 5-15 Myr for TWA.

Because TWA is the nearest group of pre-main sequence stars known, it is providing imaging possibilities that are in some ways superior to that of the Taurus-Auriga complex, and at an epoch of great importance to the study of planet formation. This association provides a laboratory for studying a number of stars of the same age but different masses, multiplicities, etc. Although the luminosities of low-mass stars in TWA are substantially less than those of younger stars of the same mass (such as in Taurus-Auriga), the relative proximity of the TW Hya association and the implementation of sensitive IR array detectors on large telescopes opens up a realm of study of PMS disk evolution not possible before.

In this paper we report observations of the 3-13 $`\mu `$m spectra of TW Hya, HR 4796A, and HD 98800. The mid-IR spectrum of HD 98800 has already been reported by Sylvester and Skinner (1996) using ground-based data, and by Walker and Heinrichsen (2000) using ISOPHOT. Our data extends to shorter wavelengths than the Sylvester & Skinner data, while having a better signal-to-noise ratio than the ISOPHOT observations. Both are useful for realistically determining the photospheric contribution to the observed emission.

## 2 OBSERVATIONS

Spectra of the TWA targets were obtained with the Aerospace Corporation’s Broadband Array Spectrograph System (BASS) configured to f/35 at the 3m NASA Infrared Telescope Facility (IRTF). The BASS consists of a pair of cooled prisms which disperse the spectrum onto two 58-element blocked impurity band (BIB) linear arrays that simultaneously span the 3-13.5 $`\mu `$m region. The spectral dispersion ranges from about 30 to 125 over each of the 3-6 and 6.5-13.5 $`\mu `$m regions. At the IRTF, the BASS entrance aperture subtends 3.4 arcsec. The observations of TW Hya were obtained on 21 March 1998 (UT). HD 98800 was observed on 8 February 1998 (UT). HR 4796A was observed on 12 May 1998 (UT).

Details of the observations are listed in Table 1, and the resultant spectra are shown in Figures 1 (TW Hya), 2 (HD 98800), and 3 (HR 4796A), where we have plotted the spectral flux ($`\lambda `$Fλ in W m-2). All BASS spectral data points such as those in the atmospheric CO2 and H2O bands with errors exceeding the observed flux (i.e., signal/noise less than unity) have been rejected. The spectrum of TW Hya was flux-calibrated using observations of $`\alpha `$ CMa (Sirius) at approximately the same airmass, obtained just prior to the TW Hya data. The spectrum of HD 98800 was calibrated using observations of $`\alpha `$ Boo (Arcturus) at the same airmass, and converted to the same flux system using the Arcturus/Sirius flux ratio determined by Russell and Mazuk (1998) for the BASS instrument. For HR 4796A, no calibrator at the same airmass were observed. To establish the spectrum, it was calibrated against one star, $`\alpha `$ Lyr (Vega) at a lower airmass, and another star, $`\beta `$ Gem (Pollux) at higher airmass. For the latter, we also used the Pollux/Sirius ratio of Russell and Mazuk (1998). Outside of those wavelengths that are severely affected by telluric lines, and which are usually rejected from the figures due to their inferior signal-to-noise ratio (SNR), the resultant fluxes are nearly identical. We illustrate this fact by plotting both resultant spectra in Figure 3. Because of the similarity of these spectra, we have taken the mean values and plotted them against the appropriate model atmosphere spectrum in Figure 4. The detail near the 10 $`\mu `$m band is shown in Figure 5.

## 3 Discussion of the Spectra of Individual Objects

### 3.1 TW Hya

Figure 6 shows the spectral energy distribution of TW Hya from 0.44 to 100 $`\mu `$m. The photometric data in the 12, 25, 60, and 100$`\mu `$m bands were obtained by IRAS. In addition, measurements in the B, V, R, I, J, H, K, and L photometric bands by Ruchinski and Krautter (1983) are included. Because TW Hya is a variable star, we show the mean value of the Rucinski and Krautter measurements, while the error bars indicate the range in those values. We also show the Kurucz model photosphere spectrum of a K7 V star, normalized to the K and L-band flux.

As can be seen in Figure 6, TW Hya possesses a significant infrared excess. In fact, the emission longward of 8 $`\mu `$m, which is dominated by the thermal emission of the circumstellar dust, accounts for 23% of the total observed luminosity of the star. Originally, Ruchinski and Krautter (1983) fit the data at 1.25 $`\mu `$m and shortward with a blackbody spectrum, which left a flux excess between 1.25 and 3.5 $`\mu `$m. However, model atmospheres of stars in this temperature range have line opacities that naturally reproduce the shape of the spectrum observed in this region. In Figure 6 we illustrate this by fitting a Kurucz model atmosphere of a K7 star, and normalize it to the 2-4 $`\mu `$m BASS data. The same model fit is included in Figure 1. The fit is not perfect down to the shortest wavelength, but TW Hya is a variable star, the BASS data are not simultaneous with the observations of Rucinski and Krautter, and the star possesses a short-wavelength excess attributed to the disk-star boundary layer (Muzerolle et al., 2000) which we have not included here. In any case, there does not appear to be any significant excess flux between 2 and 8 $`\mu `$m.

In TW Hya (see Figure 1) the emission between 8 and 13 $`\mu `$m is dominated by the silicate band at 10 $`\mu `$m. The three points between 9.5 and 9.8 $`\mu `$m with the large error bars lie in the telluric ozone band and may not be reliable, even though the data on TW Hya were obtained at the same airmass as that of Sirius. The point at 9.4 $`\mu `$m with the smaller error bar sits right at the edge of this band.

### 3.2 HD 98800

The spectrum of HD 98800 between 3 and 13.5 $`\mu `$m is shown in Figure 2. Also shown in the figure are the photometric data of Zuckerman and Becklin (1993) and Koerner et al. (2000). The agreement between absolute fluxes derived from the BASS data and the other two sets is good, although the BASS fluxes are slightly lower. The shape and total flux level is also the same as the spectral data obtained by Sylvester and Skinner (1996), whose observations cover the 8-24 $`\mu `$m region, and who used the photometry from another instrument in their modeling. The value of the BASS observations are that they cover the wavelength regions dominated by both photospheric and dust emission simultaneously and on the same photometric system, allowing the photospheric flux to be extrapolated out into the 10 $`\mu `$m region using a single set of simultaneous data. The ISOPHOT observations of Walker and Heinrichsen (2000) cover the 5.8-11.6 $`\mu `$m region, but suffer from such large systematic fluctuations in the signal that the silicate feature cannot even be seen.

The BASS spectrum shortward of 6 $`\mu `$m is well-fit by the photospheric model normalized near 3 $`\mu `$m. This is consistent with the results of Koerner et al. (2000), who fit similar models to the individual A and B components of the system (each of which is a binary), and find that the flux is dominated by photospheric emission out to 7 $`\mu `$m in both A and B components.

### 3.3 HR 4796A

In Figures 4 and 5 we show the spectrum of HR 4796A and its appropriate model atmosphere. The 8-13.5 $`\mu `$m flux is weak, and half is photospheric. Koerner et al. (1998), using thermal imaging with the MERLIN instrument on Keck, derive a total flux at 12.5 $`\mu `$m of 5.4x10-14 Wm-2, of which it was estimated that half was photospheric. For the BASS observations, the 12-13 $`\mu `$m flux is higher than those derived from MERLIN, but the fraction due to photospheric emission is comparable. The lack of any detectable excess flux shortward of 8 $`\mu `$m indicates that there is little dust hotter than 450 K or closer than about 2 AU (depending somewhat on the optical albedo and infrared emissivity of the grains). Augereau et al. (1999) have modeled the spectral energy distribution of HR 4796A, including the hot inner region, with a 2-component dust model. In their model, the excess emission near 10 $`\mu `$m comes primarily from huge (larger than 100 $`\mu `$m in size) grains located near 9 AU from the star. Unfortunately, the quality of the spectra here are inadequate the detailed spectral features expected from such material.

## 4 COMPARISON WITH SOLAR SYSTEM COMETS AND OTHER DUSTY PMS STARS

### 4.1 General Considerations

Some main sequence and pre-main sequence stars with dusty debris disks possess silicate grains whose spectral features resemble those of long-period comets. Knacke et al. (1993) were the first to demonstrate this by comparing the 10 $`\mu `$m emission feature in $`\beta `$ Pic with 1P/Halley and Levy 1990 (C/1991 L3). All three objects exhibit an emission band with maxima or shoulders near 9.5 $`\mu `$m and 11.2 $`\mu `$m, the latter indicative of crystalline olivine. Subsequently, HAEBEs embedded in star-forming regions have been examined (i.e. Hanner et al. (1995)) and some, such as HD 150193 in the $`\rho `$ Oph cloud, are found to possess features of a similar nature. The two isolated HAEBEs (which might also be called post-HAEBEs since they no longer embedded in nebulosity) HD 31648 and HD 163296, were shown to have spectra similar to that of Hale-Bopp and Levy (Sitko et al., 1999). However, the stellar disk dust in these systems emitted relatively more strongly near 9.4 $`\mu `$m than the comets.

It should be remembered that in the case of a comet, we are dealing with material with a single well-defined distance from the sun, while this is not the case for a disk. Nevertheless, Wooden et al. (2000) have shown that over a range of a factor of three in heliocentric distance, the spectral shape of Hale-Bopp was, with the exception of one subfeature, invariant. Furthermore, a comparison of the ISO spectra of the HAEBE star HD 100546 and Hale-Bopp indicate that, while the overall spectral energy distributions are significantly different in the way expected (more cool dust in the stellar disk), the shapes of the narrow features were still quite similar (Malfait et al., 1998). However, the shape of the 10 $`\mu `$m band is sensitive to the abundance and form of pyroxenes, as evidenced in the variety of spectral shapes of pyroxene-rich interplanetary dust particles (IDPs) that have been analyzed in the laboratory (Sanford and Walker, 1985; Bradley et al., 1992). It is interesting that the one subfeature that was observed to vary with heliocentric distance in Hale-Bopp has been attributed to crystalline pyroxenes (Wooden et al., 2000). This suggests that perhaps the main reason for the excess 9.4 $`\mu `$m emission of the disks compared to comets is due to the presence of dust at temperatures much higher than those generally observed in solar system comets, which allows the crystalline pyroxenes to be enhanced, but on a grander scale than was seen in Hale-Bopp.

### 4.2 TW Hya

In Figure 7 we show the 10 $`\mu `$m band of TW Hya, normalized to the continuum by dividing by a gray body that matches the observed flux at 8 and 13 $`\mu `$m. Also shown are the data for HD 31648, HD 163296, comet Hale-Bopp, and comet Levy from Sitko et al. (1999). Although the data for TW Hya has a modest SNR, it is apparent that TW Hya does not have the same shape as the other objects. In particular, the well-defined crystalline olivine feature seen at 11.2 $`\mu `$m in both the comets and HD 31648 and HD 163296 seems to considerably less-defined in TW Hya. At best, a modest change of spectral slope is present at that wavelength. In Figure 8 we compare the spectrum of TW Hya with three HAEBEs (Sitko et al., 2000) whose emission features are representative of the wide variety that are seen in these PMS objects (excluding objects dominated by organic emission features). That of RY Tau exhibits a silicate feature which is strong and relatively featureless, except for a single maximum near 9.7 $`\mu `$m. In HD 163296 (again) the olivine feature is present, but not as prominent as in comets such as Hale-Bopp and Levy. Finally, the spectrum of HD 35187 has a weak silicate band, but more comet-like (flat-topped) in shape. In Figure 9 we also compare the spectrum of TW Hya to that of RW Aur, a PMS G5e star with circumstellar CO emission in the first vibrational overtone band (Sitko and Hanson, 2000). Because the strength of the silicate band with respect to the underlying continuum in RW Aur is only one-quarter that of TW Hya, for comparison purposes the band strength with respect to the continuum in RW Aur has been scaled up by a factor of 4.

Of the four “comparison” star, the silicate band in TW Hya is most similar in shape and strength to that of RY Tau, a F8Ve PMS star in the Taurus-Auriga cloud. The spectral shape in the silicate band of TW Hya is also nearly identical to that of RW Aur (G5e), but its strength is considerably stronger. This should not be construed to imply that TW Hya and the other late-type PMS stars are systematically different from the earlier-type HAEBEs. Most HAEBEs are not strong 11.2 $`\mu `$m emitters. Furthermore, high-quality mid-IR spectra on a large sample of TTs are lacking.

In Figure 10 we show the flux/continuum ratio of HD 98800 with that of TW Hya. While the two sets of data overlap reasonably well between 11 and 13 $`\mu `$m, much of this may simply be due to the normalization procedure used to get the flux/continuum ratio. Overall, the band strength seems to be weaker in HD 98800, and there may be structural differences as well. However, when the band strength is intrinsically weak, the normalization procedure tends to magnify the systematic effects in the data (such as small errors in cancellation of telluric lines between target and calibration star) so that spectral “features” may not be real. In such cases, a direct comparison of the observed fluxes may be more useful. In Figure 11 we compare the observed fluxes of these two stars directly (after suitable re-scaling). Here, the differences are more readily apparent.

Due to the weakness of the dust emission in HR 4796A compared to the extrapolated photospheric flux, the flux/continuum ratio was not determined.

## 5 DISCUSSION

In Figure 11 we compare the spectra of all three stars, normalized to the same flux between 3 and 5 $`\mu `$m, where all three are dominated by photospheric emission. It is apparent that the net strength of the silicate emission feature of the dust surrounding the three TWA stars is different, and it cannot simply be the result of differing ages, since the stars are believed to be approximately coeval. Thus the degree to which small silicate grains survive and/or are modified in these systems must be dependent on other parameters as well, such as stellar mass, presence and location of companions, and the like. With only three TWA stars observed so far, it is premature to draw any specific conclusions in this regard. Furthermore, the quality of the spectra of the dust in HR 4796A and HD 98800 preclude any sweeping generalizations concerning their detailed mineralogies. Instead, we examine some of the inferences that can be made concerning the general nature of the material surrounding TW Hya specifically.

The 11.2 $`\mu `$m crystalline olivine feature has not been detected with the ISO Short Wavelength Spectrometer (SWS) in the emission spectra of the youngest (i.e. embedded) PMS stars, but is seen observed in a number of more evolved objects (Waelkens et al., 2000). Because the majority of the objects observable with ISO SWS must be bright, most of the targets have been HAEBEs. The TTs were too faint for SWS (van den Ancker, 2000). Only a few were observed at more modest spectral resolution with ISOPHOT (see Natta et al. (2000)), but much of these data have inadequate signal-to-noise to make any conclusions regarding the presence of crystalline material.

The dust in TW Hya itself does not appear to exhibit the degree of crystallinity, as judged by the strength of the 11.2 $`\mu `$m olivine feature, that some other objects with resolved debris disks have, and which is seen in many, but not all, solar system comets (Hanner et al., 1994). For example, in the HAEBE star HD 163296, a star whose age is believed to be between 4 and 6 Myr, the feature is definitely present (Sitko et al., 1999; van den Ancker et al., 2000). About half (or more) of long-period comets studied so far also possess the 11.2 $`\mu `$m feature. These comets have their recent origin in the Oort cloud, and were presumably originally formed in the Uranus-Neptune region of the early solar system. By contrast, short-period comets, whose origin is probably in the Kuiper belt, tend to have weaker silicate features (Lynch et al., 1992, 1995, 2000; Hanner et al., 1996). The mean particle sizes must be greater (or the feature otherwise hidden in some fashion) in these latter objects, and the degree of crystallinity cannot be assessed.

When does the crystalline material form in these disk systems, if at all? Bradley et al. (1999) have demonstrated that the most primitive unequilibrated interplanetary dust particles, of presumed cometary origin, exhibit a smooth silicate band essentially devoid of crystalline features. This is consistent with their glassy internal structure. The spectra of these particles resemble those of disk systems that do not exhibit the 11.2 $`\mu `$m feature. However, the large porous grains in which they are embedded, along with significant amounts of crystalline material, have spectra which resemble those older disk systems and solar system comets that do exhibit this crystalline feature. It may be that much of the crystalline material is produced early in the development of the protosolar nebula, and is later agglomerated with the more pristine material (Bradley et al., 2000). One could even envision a scenario whereby disk or comet dust shed and heated near the star is recycled over the entire lifetime of the formation of the bulk of the cometesimals in the nebula. In such a scenario, “pristine” objects might be highly “contaminated” by dust processed at a variety of temperatures.

The lack of the feature in TW Hya, whose age is greater than that of some other stars exhibiting the feature, could be due to a number of factors. Perhaps it is simply too young, and hasn’t entered the stage where this material is produced. We simply do not have enough data on lower-mass stars to know for sure if they are slower than the higher-mass stars in entering this phase (and not all of them do, to the best of our knowledge). Or perhaps the material is currently not being delivered close enough to the star yet for the material to emit strongly. In addition, the sheer volume of space surrounding the higher-mass stars that is capable of thermally processing grains is larger than in lower-mass stars.

## 6 CONCLUSION

We have presented 3-13 $`\mu `$m spectra of three stars in the TW Hya Association: TW Hya, HR 4796A, and HD 98800. TW Hya possesses a strong silicate emission band that is generally similar in strength and shape to many other pre-main sequence stars, such as the Herbig Ae/Be stars (HAEBEs) and T Tauri stars (TTs). However, the feature does not show the unambiguous presence of crystalline material, as evidenced by the 11.2 $`\mu `$m feature, as is observes in some HAEBEs and in some long-period solar system comets. HR 4796A has only a very weak excess emission above photospheric levels, and nothing can be said about the presence of a significant silicate band from these data. HD 98800 is intermediate between these two in terms of the strength of its silicate band, but the data are inadequate to say anything about the mineralogy or degree of crystallinity.

For the object where we have the best data, TW Hya, the lack of obvious crystallinity is not entirely surprising, as it is absent in all of the youngest pre-main sequence stars, and in most but not all of the older ones. It is present in many of the long-period comets of our own solar system, however, and our failure to detect it so far in a younger star of roughly solar mass is intriguing, though not unexpected. We are far from understanding how the crystalline material is produced and incorporated into the grains of our own solar system and in the disk systems of some of the higher-mass stars.

One of the main obstacles preventing more definitive conclusions is the lack of high-quality mid-IR spectra of TTs. Many can be reached using ground-based instruments, but require long (in excess of an hour or two) integration times to achieve the necessary signal-to-noise. A truly comprehensive study will require an instrument like the Space Infrared Telescope Facility (SIRTF).

Support for this work was provided for MLS through NASA’s Origins of Solar Systems Program grant NAG5-9475, and the University Research Council and Physics Department of the University of Cincinnati. Support for DKL and RWR was provided by the The Aerospace Corporation’s Independent Research and Development program. We would like to thank Ann Mazuk, and Ted Tessensohn for technical support and Bill Golisch and Dave Griep of the IRTF for expert telescope operation. The authors would also like to thank Robert Joseph for being flexible in the scheduling of the IRTF time for this project with our other programs. We also thank Carol Grady for her comments on the manuscript. The IRAS broadband fluxes were obtained via SIMBAD.


text:
# Scaling regimes for second layer nucleation

## I Introduction

Atomistic nucleation theory was developed in the sixties and seventies in response to the realization that critical nuclei of atomic dimensions are common in far from equilibrium thin film growth . During the past decade the subject has experienced a vigorous renaissance, driven by the invention of atomic scale imaging techniques as well as by the availability of realistic large scale simulations . While this has resulted in a wide range of new theoretical developments, the classical theory has stood the test of time remarkably well. This is surprising because classical nucleation theory relies, in the jargon of statistical mechanics, on an approximation of mean field type, in the sense that it considers only the evolution of spatially averaged concentrations of the different surface species – adatoms and clusters – and thereby neglects effects of spatial correlations and fluctuations.

It was recently pointed out that the classical theory manifestly fails in describing nucleation events on top of islands bounded by strong step edge barriers . The onset of second layer nucleation on top of the first layer islands is an important factor determining the growth mode and morphology of multilayer films. Tersoff and coworkers  developed a theory of second layer nucleation based on the classical mean field estimate

$$\omega _{\mathrm{mf}}\nu An^{i^{}+1}$$
(1)
for the rate $`\omega `$ of nucleation events on the island, in terms of the in-layer hopping rate $`\nu `$, the island area $`A`$, the adatom density $`n`$ and the size $`i^{}`$ of the critical (= largest unstable) cluster. In Refs. it was shown that for $`i^{}=1`$, 2 the mean field estimate disagrees with microscopic considerations  and computer simulations. The failure of mean field theory was qualitatively attributed to the small number of atoms typically present on the island, and the associated large fluctuations in the occupancy of the island.

The question then arises whether mean field behavior may be recovered at larger values of the critical nucleus size, where nucleation necessarily involves a greater number of atoms, and thus fluctuation effects should be reduced. A partial answer to this question was provided in recent work by Heinrichs, Rottler and Maass (HRM) . They distinguish between two different nucleation modes: In the fluctuation-dominated mode the mean number of adatoms on the island (in the quasi-stationary state prior to nucleation) is much smaller than the number $`i^{}+1`$ required for nucleation, and therefore nucleation is a rare, large fluctuation, while in the mean field regime, which corresponds to the classical approach of , the mean number of adatoms exceeds $`i^{}+1`$. Consistency considerations then show that second layer nucleation on top of first layer islands is fluctuation-dominated for $`i^{}2`$ and mean-field-like for $`i^{}3`$.

In this paper I elaborate on the observations of HRM, with the goal of further clarifying the role of fluctuations in second layer nucleation. I begin by rederiving the main results of HRM within the approach of , which relies on expressing the nucleation rate in terms of the basic time scales of the process. While less quantitative than the treatment of , this provides, perhaps, a somewhat more transparent and unified view. A simple, universal criterion for fluctuation-dominated second layer nucleation is found which involves the first layer island density exponent and the critical nucleus size (Eq.(28)). I then extend the analysis to $`d`$-dimensional nucleation, finding the same general structure as in the two-dimensional case, but with a decreasing influence of fluctuations with increasing $`d`$; as a byproduct, I recover the expression for the one-dimensional first layer island density exponent recently derived in Ref.. Finally, I generalize the theory of nucleation on top of multilayer mounds developed in  to $`i^{}>1`$. In contrast to the case of nucleation on first layer islands, here fluctuation-dominated and mean field like nucleation regimes are found to coexist for arbitrary values of $`i^{}`$. This illustrates that the importance of fluctuations depends crucially on the growth situation.

## II Stationary rates of fluctuation-dominated nucleation

We consider a two-dimensional island of fixed linear extension $`L`$. Atoms are deposited onto the island at rate $`F`$, they migrate on the island with an in-layer hopping rate $`\nu `$ and descend from the island at the inter-layer hopping rate $`\nu ^{}`$. The mean time interval between subsequent arrivals of atoms on the island is then

$$\mathrm{\Delta }t\frac{1}{FL^2}$$
(2)
and the mean residence time of an atom on the island in steady state is of the form

$$\tau =\frac{aL^2}{\nu }+\frac{bL}{\nu ^{}},$$
(3)
where $`a`$ and $`b`$ are geometry-dependent constants. The first term on the right hand side of (3) is of the order of the traversal time $`\tau _{\mathrm{tr}}L^2/\nu `$ required for the atom to migrate across the island, while the second term describes the increase in the residence time due to the step edge barriers. The second term dominates if $`Ll_{\mathrm{ES}}`$, where

$$l_{\mathrm{ES}}\nu /\nu ^{}$$
(4)
is the Ehrlich-Schwoebel length .

The stationary mean number of atoms on the island is given by $`\tau /\mathrm{\Delta }t`$ , and hence the condition for the fluctuation-dominated mode of HRM reads

$$\frac{\tau }{\mathrm{\Delta }t}1.$$
(5)
It will turn out to be useful to write this in the form $`LL_F`$, where the flux scale $`L_F`$ is given by

$$L_F\{\begin{array}{cc}(\nu ^{}/F)^{1/3}:\hfill & Ll_{\mathrm{ES}}\hfill \\ (\nu /F)^{1/4}:\hfill & Ll_{\mathrm{ES}}.\hfill \end{array}$$
(6)
For given $`L`$, $`\nu `$ and $`\nu ^{}`$, Eq.(5) can always be satisfied by making the flux sufficiently small.

Nucleation with a critical nucleus size $`i^{}`$ is treated within the “noninteracting particle model” , in which the lifetimes of all unstable clusters are assumed comparable to the inverse hopping rate $`1/\nu `$. Then the third relevant time scale in the problem, in addition to (2) and (3), is given by the encounter time $`\tau _{\mathrm{enc}}`$ required for the $`i^{}+1`$ atoms to find each other on the island. In order of magnitude we have

$$\tau _{\mathrm{enc}}L^{2i^{}}/\nu .$$
(7)
A simple derivation is given below in Section IV.

Under the condition (5), each nucleation event is uniquely associated with the deposition of the $`i^{}+1`$’th atom. The number of nucleation events per unit time can then be written as

$$\omega =p_{\mathrm{nuc}}/\mathrm{\Delta }t,$$
(8)
where $`p_{\mathrm{nuc}}`$ is the probability for a freshly deposited atom to participate in a nucleation event. The latter is generally of the form

$$p_{\mathrm{nuc}}=p_{i^{}+1}p_{\mathrm{enc}},$$
(9)
where $`p_{i^{}+1}`$ denotes the probability that $`i^{}`$ atoms were present on the island prior to the arrival of the $`i^{}+1`$’th, and $`p_{\mathrm{enc}}`$ is the probability that the $`i^{}+1`$ atoms encounter each other before one of them leaves the island again. The quantity $`p_{i^{}+1}`$ is given by

$$p_{i^{}+1}=\frac{1}{i^{}!}\left(\frac{\tau }{\mathrm{\Delta }t}\right)^i^{}.$$
(10)
The encounter probability can be estimated as

$$p_{\mathrm{enc}}1e^{\tau /\tau _{\mathrm{enc}}}$$
(11)
where for simplicity we have equated the residence time of the $`i^{}+1`$ atoms with that of a single atom (in fact the two differ by a factor $`i^{}+1`$). Thus $`p_{\mathrm{enc}}`$ is of order unity for $`\tau /\tau _{\mathrm{enc}}1`$ and of the order $`p_{\mathrm{enc}}\tau /\tau _{\mathrm{enc}}`$ when $`\tau /\tau _{\mathrm{enc}}1`$.

In the case $`i^{}=1`$ the encounter time (7) is of the same order as the traversal time $`L^2/\nu `$, whereas the residence time (3) is of the order of $`\tau _{\mathrm{tr}}`$ or larger. Thus one always has $`p_{\mathrm{enc}}=𝒪(1)`$ and $`p_{\mathrm{nuc}}p_2`$ . In contrast, for $`i^{}>1`$ it is possible to have a transition from $`\tau \tau _{\mathrm{enc}}`$ to $`\tau \tau _{\mathrm{enc}}`$ with increasing island size. This obviously requires that $`\nu /\nu ^{}1`$, so that $`\tau L/\nu ^{}\tau _{\mathrm{enc}}`$ for small islands, and occurs at the characteristic length scale

$$L_\times (\nu /\nu ^{})^{1/(2i^{}1)}.$$
(12)
We note that this is small compared to the Ehrlich-Schwoebel length, and therefore the transition to encounter-limited nucleation occurs before the step edge barriers become irrelevant, in the sense that $`\tau `$ becomes of the same order as $`\tau _{\mathrm{tr}}`$.

In the encounter-limited regime the expression (8) agrees in order of magnitude with the mean field estimate (1). To see this, note that the adatom density on the island is of order $`\tau /(\mathrm{\Delta }tL^2)=F\tau `$, and hence

$$\omega _{\mathrm{mf}}\nu L^2\left(\frac{\tau }{L^2\mathrm{\Delta }t}\right)^{i^{}+1}\frac{1}{\mathrm{\Delta }t}\left(\frac{\tau }{\mathrm{\Delta }t}\right)^i^{}\frac{\tau }{\tau _{\mathrm{enc}}}\frac{p_{i^{}+1}p_{\mathrm{enc}}}{\mathrm{\Delta }t}.$$
(13)
In this sense the mean field approach remains applicable also outside of its strict range of validity .

We summarize the preceding considerations in the following three expressions for the fluctuation-dominated nucleation rate, which apply with increasing island size $`L`$:

$$\omega \frac{1}{\mathrm{\Delta }t}\left(\frac{\tau }{\mathrm{\Delta }t}\right)^i^{}F(F/\nu ^{})^i^{}L^{3i^{}+2}LL_\times (\mathrm{regime}\mathrm{II})$$
(14)
$$\omega \frac{1}{\tau _{\mathrm{enc}}}\left(\frac{\tau }{\mathrm{\Delta }t}\right)^{i^{}+1}F(F/\nu ^{})^i^{}(\nu /\nu ^{})L^{i^{}+3}L_\times Ll_{\mathrm{ES}}(\mathrm{regime}\mathrm{III})$$
(15)
$$\omega \frac{1}{\tau _{\mathrm{enc}}}\left(\frac{\tau }{\mathrm{\Delta }t}\right)^{i^{}+1}F(F/\nu )^i^{}L^{2i^{}+4}Ll_{\mathrm{ES}}(\mathrm{regime}\mathrm{IV}).$$
(16)
Regime III does not exist for $`i^{}=1`$. For $`i^{}>1`$, the encounter-limited regimes III and IV extend also into the mean field realm, while regime II occurs only in the fluctuation-dominated nucleation mode.

Note that Eqs.(14) and (15) correspond to the regimes II and I identified in  for $`i^{}=2`$. Here we prefer to enumerate the regimes in accordance with HRM . Their regime I will be dealt with below in Section III B.

## III Consistent scaling regimes in two dimensions

While Eqs.(14,15,16) provide, in the specified ranges of parameters, valid expressions for the nucleation rate on top of an island of given size, in their application to a specific growth situation additional consistency requirements arise : When a given expression for $`\omega `$ is used to compute the island size at which nucleation occurs, it has to be checked that this island size does indeed lie in the appropriate range.

In this section we derive the resulting scaling regimes for second layer nucleation in two dimensions. Following earlier work  we assume a population of equal sized first layer islands of density $`N`$. The island size increases with time $`t`$ or coverage $`\theta =Ft`$ according to the growth law

$$L(t)\sqrt{\theta /N}.$$
(17)
For a nucleation rate of the form

$$\omega =F\mathrm{\Omega }L^k$$
(18)
the critical island size at which a second layer nucleus forms is then given by

$$L_c(N\mathrm{\Omega })^{1/(k+2)}.$$
(19)
We further introduce the notation

$$\mathrm{\Gamma }\nu /F,\alpha \nu ^{}/\nu =l_{\mathrm{ES}}^1,$$
(20)
where usually $`\mathrm{\Gamma }1`$, $`\alpha 1`$. The first layer island density is then of the order of

$$N\mathrm{\Gamma }^\chi ,$$
(21)
which defines the island density exponent $`\chi `$. In standard nucleation theory  it is given by the expression

$$\chi =\frac{i^{}}{i^{}+2}.$$
(22)
At the moment we will however leave $`\chi `$ unspecified to allow for the possibility of different first layer nucleation mechanisms. Generalizations of (22) will be derived in Section IV.

### A Stationary nucleation regimes

Inserting Eqs.(14,15,16) into (19) one obtains expressions for the critical island size of the general form

$$L_c\mathrm{\Gamma }^\gamma \alpha ^\mu ,$$
(23)
where the exponents $`\mu `$ and $`\gamma `$ are given by

$$\gamma =\frac{\chi +i^{}}{3i^{}+4},\mu =\frac{i^{}}{3i^{}+4}(\mathrm{regime}\mathrm{II})$$
(24)
$$\gamma =\frac{\chi +i^{}}{i^{}+5},\mu =\frac{i^{}+1}{i^{}+5}(\mathrm{regime}\mathrm{III})$$
(25)
$$\gamma =\frac{\chi +i^{}}{2i^{}+6},\mu =0(\mathrm{regime}\mathrm{IV}).$$
(26)
Consider first regime II. Consistency requires the inequalities $`L_c\alpha ^\mu \mathrm{\Gamma }^\gamma L_\times \alpha ^{1/(2i^{}1)}`$ and $`\alpha ^\mu \mathrm{\Gamma }^\gamma L_F(\alpha \mathrm{\Gamma })^{1/3}`$ to be satisfied simultaneously, which implies

$$\mathrm{\Gamma }^{(13\gamma )/(13\mu )}\alpha \mathrm{\Gamma }^{\gamma /(\mu +1/(2i^{}1))}.$$
(27)
Since $`\alpha 1`$ and $`\mathrm{\Gamma }1`$, these inequalities can be satisfied only if the exponent of $`\mathrm{\Gamma }`$ on the left hand side is smaller than the one on the right hand side. Inserting the expressions (24) this is found to lead to the condition

$$\chi <\frac{2}{i^{}+1}.$$
(28)
The analysis of regime III, based on the requirement $`L_\times L_cL_F`$, yields the same condition (28), while in regime IV the relevant inequalities $`l_{\mathrm{ES}}L_cL_F`$ translate into

$$\chi <\frac{3i^{}}{2}.$$
(29)
Using the expression (22) for $`\chi `$, we see that the fluctuation-dominated regimes II and III become inconsistent at $`i^{}=i_{c1}^{}=(1+\sqrt{17})/22.56`$. For larger $`i^{}`$ regime III survives in its mean field version, while regime II disappears. In regime IV the transition from fluctuation-dominated to mean-field behavior occurs at $`i^{}=i_{c2}^{}=2`$; in this regime the critical island size is of the order of the first layer island distance, $`\gamma =\chi /2`$, and hence the problem really reduces to that of first layer nucleation.

The general form of the conditions (28,29) illustrates that the upper critical nucleus size, at which the mean-field nucleation mode sets it, is intimately linked to the expression for the first layer island density exponent. Scenarios different from the standard case (22) can arise if, due to e.g. metastable clusters  or surface reconstruction , the critical nucleus on the substrate is smaller or larger than that on the island. Another example of nonstandard behavior is kinetically limited first layer nucleation, where Eq.(22) is replaced by  $`\chi =2i^{}/(i^{}+3)`$ and the upper critical nucleus sizes shift to $`i_{c1}^{}=\sqrt{3}1.732`$ and $`i_{c2}^{}=\sqrt{13}21.6055`$. In addition the growth law (17) plays an important role, as will be discussed below.

### B Nonstationary nucleation

For very strong step edge barriers the residence time $`\tau `$ becomes the largest time scale in the problem, and may effectively be set to infinity . In this situation all atoms deposited on the island remain there, and a stationary state where deposition and loss of adatoms compensate each other is never reached. A distinction between fluctuation-dominated and mean-field nucleation scenarios is nevertheless possible : In the first case nucleation occurs as soon as the $`i^{}+1`$-th atom arrives on the island, while in the second case many more atoms have to be deposited.

Using the growth law (17) it is easy to show that the time required to deposit a few (i.e., $`i^{}+1`$) atoms on the island is

$$\tau _{\mathrm{dep}}N^{1/2}/F.$$
(30)
At this point the island size is of the order $`N^{1/4}`$, corresponding to exponents

$$\gamma =\chi /4,\mu =0,(\mathrm{regime}\mathrm{I}/\mathrm{fl})$$
(31)
for the fluctuation-dominated part of regime I. To see whether nucleation is actually fluctuation dominated, the time scale (30) is to be compared to the encounter time (7), evaluated at the island size $`N^{1/4}`$; fluctuation-dominated nucleation occurs if $`\tau _{\mathrm{enc}}\tau _{\mathrm{dep}}`$. Once more this yields the condition (28), which can therefore be regarded as a universal criterion for fluctuation-dominated second layer nucleation. For $`\chi >2/(i^{}+1)`$ mean field theory can be used to compute the nucleation rate, and one arrives at the scaling exponents

$$\gamma =\frac{\chi (i^{}+2)1}{2i^{}+6},\mu =0(\mathrm{regime}\mathrm{I}/\mathrm{mf}).$$
(32)
A graphical representation of the sequence of scaling regimes as a function of $`i^{}`$ is provided in Figure 1.

While in the stationary nucleation regimes of Section III A the criterion (28) appears in a rather indirect way, in the nonstationary case it has a straightforward interpretation in terms of a comparison of the time scales $`\tau _{\mathrm{dep}}`$ and $`\tau _{\mathrm{enc}}`$. This is easily generalized to first layer island growth laws of the form

$$L(t)N^{1/2}\theta ^\beta $$
(33)
which lead to

$$\chi <\frac{2\beta +1}{i^{}+1}$$
(34)
instead of (28).

## IV Dimensionality dependence

In statistical physics fluctuation effects are typically more prominent in low spatial dimensionalities. It is therefore of some interest to repeat the above considerations for growth of $`d`$-dimensional islands on a $`d`$-dimensional substrate. The residence time is then still given by an expression of the form (3), while the interarrival time becomes

$$\mathrm{\Delta }t\frac{1}{FL^d}.$$
(35)
To estimate the encounter time, it is useful to visualize the trajectories of the $`i^{}+1`$ diffusing adatoms in a configuration space of dimensionality $`d(i^{}+1)`$, where the motion is restricted to a region of linear extent $`L`$. Nucleation occurs on a $`d`$-dimensional submanifold, the codimension of which is $`di^{}`$. We are thus dealing with a first passage problem in $`di^{}`$ dimensions, where the encounter time plays the role of the mean first passage time. This implies that

$$\tau _{\mathrm{enc}}\{\begin{array}{cc}L^2/\nu :\hfill & di^{}<2\hfill \\ L^{di^{}}/\nu :\hfill & di^{}2.\hfill \end{array}$$
(36)
Since the residence time is always at least equal to the traversal time $`L^2/\nu `$, it follows that encounter-limited nucleation in the sense of $`\tau _{\mathrm{enc}}\tau `$ is possible only for dimensionalities $`d>2/i^{}`$. This recovers the upper critical dimension $`d_c=2/i^{}`$ for diffusion-limited reactions of $`i^{}+1`$ particles .

### A Island density exponent in $`d`$ dimensions

As a warm-up excercise we derive here a formula for the exponent $`\chi `$ of the first layer island density valid for general $`d`$. The appropriate expression for the nucleation rate is (16), which is now to be interpreted as the rate of nucleation in a region of size $`L`$. The distance $`l_DN^{1/d}`$ between first layer nuclei is determined by the condition that there should be one nucleation event per area $`l_D^d`$ during the monolayer deposition time $`1/F`$, i.e.

$$\omega (l_D)/F=𝒪(1).$$
(37)
For $`di^{}2`$ this yields

$$\chi =\frac{di^{}}{d+2+2i^{}}(di^{}2),$$
(38)
which agrees with (22) for $`d=2`$ and with the expression $`\chi =i^{}/(3+2i^{})`$ derived by Kallabis, Krapivsky and Wolf for $`d=1`$, $`i^{}2`$ . For $`di^{}<2`$ we obtain instead

$$\chi =\frac{di^{}}{(d+2)i^{}+d}(di^{}<2).$$
(39)
This case is physically realized only for $`d=i^{}=1`$, where (39) reduces to the well-known value $`\chi =1/4`$ .

These results are easily generalized to the growth of islands of (fractal or integer) dimensionality $`d^{}`$ on a $`d`$-dimensional substrate. In that case the time of interest is not the monolayer growth time $`1/F`$, but rather the time required to cover a fraction of the substrate which is of order unity. While the true coverage is related to island density and island size through $`\theta NL^d^{}`$, the covered fraction of the $`d`$-dimensional substrate is $`\theta _dNL^dN^{1d/d^{}}\theta ^{d/d^{}}`$. Setting the inverse nucleation rate equal to the time in which $`\theta _d=𝒪(1)`$ is reached then yields the condition

$$\omega (l_D)/Fl_D^{dd^{}},$$
(40)
which for $`di^{}2`$ implies the island density exponent

$$\chi =\frac{di^{}}{2(i^{}+1)+d^{}}(di^{}2).$$
(41)
For $`d=1`$, 2 this agrees with the general expression given in .

### B Second layer nucleation in $`d`$ dimensions

For compact $`d`$-dimensional first layer islands the growth law (17) generalizes to

$$L(t)(\theta /N)^{1/d},$$
(42)
and the critical island size for a nucleation rate of the form (18) becomes

$$L_c(\mathrm{\Omega }N)^{1/(d+k)}.$$
(43)
Using the expressions (14,15,16) for the nucleation rates together with (35,36) and (42) it is then straightforward to repeat the considerations of Section III. Rather than developing in detail the general case, we focus here on the transition between fluctuation-dominated and mean field nucleation regimes.

The analysis is most transparent in the nonstationary case. One finds that the deposition time is given by (30) for all $`d`$, and the critical island size for second layer nucleation is

$$L_cN^{1/2d}\sqrt{l_D}\mathrm{\Gamma }^{\chi /2d}(\mathrm{regime}\mathrm{I}/\mathrm{fl}).$$
(44)
Using this to evaluate the encounter time (7), it turns out that nucleation is always fluctuation-dominated for $`di^{}<2`$, while for $`di^{}2`$ the condition for the fluctuation-dominated regime is again given by the universal relation (28). The same relation also governs the disappearance of the fluctuation-dominated regime II for general $`d`$.

Setting the island density exponent (38) equal to $`2/(i^{}+1)`$ then yields the expression

$$i_{c1}^{}=\frac{1}{2d}(4d+\sqrt{8d^2+(d+4)^2})$$
(45)
for the upper critical nucleus size at which the fluctuation-dominated regimes I/fl, II and III/fl disappear. This is a decreasing function of $`d`$ which approaches $`i_{c1}^{}=1`$ for $`d\mathrm{}`$. As expected, the fluctuation-dominated regime becomes smaller in higher dimensions.

To analyse the transition to mean field nucleation in regime IV we use the fact that the critical island size for second layer nucleation is of the order of the first layer island spacing, i.e. $`\gamma =\chi /d`$ and $`L_cl_D`$. The condition for fluctuation-dominated nucleation therefore simply reads $`\tau /\mathrm{\Delta }tFl_D^{2+d}/\nu 1`$, or $`\chi <d/(2+d)`$. This clearly always holds for the expression (39), while for (38) it is true provided $`i^{}<i_{c2}^{}`$ with

$$i_{c2}^{}=1+2/d.$$
(46)
Again $`i_{c2}^{}1`$ for $`d\mathrm{}`$, and $`i_{c2}^{}<i_{c1}^{}`$ for all $`d`$.

The expressions for the scaling exponents $`\gamma `$ and $`\mu `$ in the remaining regimes read

$$\gamma =\frac{\chi (i^{}+2)1}{d(i^{}+3)},\mu =0(\mathrm{regime}\mathrm{I}/\mathrm{mf})$$
(47)
$$\gamma =\frac{\chi +i^{}}{(d+1)i^{}+2d},\mu =\frac{i^{}}{(d+1)i^{}+2d}(\mathrm{regime}\mathrm{II})$$
(48)
$$\gamma =\frac{\chi +i^{}}{i^{}+2d+1},\mu =\frac{i^{}+1}{i^{}+2d+1}(\mathrm{regime}\mathrm{III}).$$
(49)
The scaling regimes for the one-dimensional case are shown in Figure 2. The upper critical nucleus sizes are $`i_{c1}^{}=(3+\sqrt{33})/24.372\mathrm{}`$ and $`i_{c2}^{}=3`$.

## V Nucleation on the top terrace of a mound

In multilayer growth, the suppression of interlayer transport leads to the formation of pyramidal mounds . For strong step edge barriers, in the sense that $`l_{\mathrm{ES}}`$ exceeds the distance $`l_D`$ between first layer islands, the mound separation is set by $`l_D`$ and remains constant during growth . In this growth regime the mounds are wedding-cake-like stacks of islands upon islands, with a characteristic up-down-asymmetry: While the valleys between mounds are deep crevices, at the hilltops flat terraces of a characteristic size $`L_{\mathrm{top}}<l_D`$ are found . $`L_{\mathrm{top}}`$ is determined by the nucleation rate on the top terrace through the requirement that on average one nucleation event should occur during the growth of a monolayer, so that

$$\omega (L_{\mathrm{top}})/F=𝒪(1).$$
(50)
For a nucleation rate of the general form (18) this implies

$$L_{\mathrm{top}}\mathrm{\Omega }^{1/k}.$$
(51)
A more precise calculation the factor of proportionality can be found in .

As in the case of second layer nucleation, application of the expressions (14, 15, 16) for the nucleation rate yields estimates for the top terrace size of the form

$$L_{\mathrm{top}}\mathrm{\Gamma }^\gamma ^{}\alpha ^\mu ^{}$$
(52)
where the scaling exponents $`\gamma ^{}`$ and $`\mu ^{}`$ depend on $`i^{}`$ and on the scaling regime in question. Inspection shows that the equivalent of the nonstationary regime I discussed above in Section III B is trivial in the case of mound growth: When $`\alpha \mathrm{\Gamma }^1`$, the top terrace size becomes comparable to the lattice constant, $`L_{\mathrm{top}}=𝒪(1)`$, corresponding to pure statistical Poisson growth .

The fluctuation-dominated regime II with scaling exponents

$$\gamma ^{}=\mu ^{}=\frac{i^{}}{3i^{}+2},(\mathrm{regime}\mathrm{II},\mathrm{mounds})$$
(53)
covers the parameter range $`\mathrm{\Gamma }^1\alpha \mathrm{\Gamma }^{\delta _1}`$, where $`\delta _1=i^{}(2i^{}1)/(2i^{}(i^{}+1)+2)<1`$ for all $`i^{}`$. In marked contrast to second layer nucleation, here regime II exists for all $`i^{}`$.

For $`\mathrm{\Gamma }^{\delta _1}\alpha \mathrm{\Gamma }^{\chi /2}`$ nucleation is governed by the exponents of regime III,

$$\gamma ^{}=\frac{i^{}}{i^{}+3},\mu ^{}=\frac{i^{}+1}{i^{}+3}(\mathrm{regime}\mathrm{III},\mathrm{mounds}).$$
(54)
This regime is fluctuation-dominated for $`i^{}<2`$, while for $`i^{}2`$ it splits up into fluctuation-dominated and mean field-like subregimes with identical scaling exponents. Both subregimes exist for all $`i^{}2`$, though the fluctuation-dominated part becomes very small for large $`i^{}`$. For $`\alpha >\mathrm{\Gamma }^{\chi /2}`$ one enters the regime of weak step edge barriers, where $`l_{\mathrm{ES}}<l_D`$ and the shape and evolution of the mound morphology are qualitatively different from the wedding-cake regime discussed here . A pictorial representation of the scaling regimes is given in Figure 3.

## VI Summary

The goal of this paper has been to extend the analysis of  to arbitrary critical nucleus sizes, as well as to clarify the relationship between the recent stochastic approaches to second layer nucleation  and the earlier mean field theory of Ref.. The essential new ingredient that enters the problem for $`i^{}>1`$ is the encounter time $`\tau _{\mathrm{enc}}`$, which may become larger than the residence time $`\tau `$. The mean field approach is valid, at least as far as the values of scaling exponents are concerned, when $`\tau \tau _{\mathrm{enc}}`$, but it fails when $`\tau \tau _{\mathrm{enc}}`$. Apart from recovering the results of , the analysis has revealed a remarkably simple and universal criterion for fluctuation-dominated second layer nucleation, Eq.(28), which may be useful in situations where nucleation is governed by different mechanisms in the first and second layer. Furthermore I have verified the expectation that fluctuation-dominated nucleation should be more prominent in low dimensionalities, and I have shown that no upper critical nucleus size (beyond which fluctuations can be neglected) exists for nucleation on the top terrace of a mound.

## Acknowledgements

I am grateful to S. Heinrichs, P. Maass and T. Michely for useful discussions. This work was supported by DFG within SFB237.