text:
# Random strategies with historical memory for the Robin Hood game

## 1. The Robin Hood game

The *Robin Hood* game $`RH(r,s,A)`$ is defined for functions $`r,s:`$ such that $`1r(i)<s(i)`$ for each $`i`$, and for a set $`A`$, as follows:

1. On day $`i`$, the *Sheriff* (of Nottingham) puts $`s(i)`$ bags of gold in the cave, each labelled by an element of $`A`$. No label is used twice (over the course of the entire game).
2. On night $`i`$, *Robin* (Hood) removes $`r(i)`$ bags from the cave.

The game is played for each $`i`$. Robin wins if each bag which was put in the cave is eventually removed from it; otherwise the Sheriff wins.

It is easy to see that if Robin has an unlimited historical memory (knowing at each night $`i`$ which of the bags in the cave appeared first), then he has a winning strategy: On night $`i`$ pick $`r(i)`$ bags out of those which arrived first.

Deterministic strategies for this game were studied, from the set-theoretic point of view, in . Gasarch, Golub and Srinivasan  consider the case where Robin has *no historical memory*, that is, he cannot distinguish between the days where the bags were put in the cave. They suggest the following *probabilistic* strategy for Robin: On night $`i`$, remove random $`r(i)`$ bags out of the cave (with uniform distribution). They say that Robin wins *almost surely* if for each bag put in the cave, its probability of being eventually removed is $`1`$. The probability is taken over Robin’s coin tosses. More precisely, the probability that a bag is not removed is the product of all probabilities $`p_i`$ that the bag is not removed at night $`i`$, and the Sheriff’s winning (or Robin’s loosing) probability is the supremum of all probabilities $`p_x`$ that a bag $`x`$ put in the cave is never removed.111Our definition of the Sheriff’s winning probability is simpler than the one given in , but both our and the proofs of  work for both definitions and actually imply that the definitions are equivalent. Let $`L(i)=_{j=1}^is(j)r(j)`$ denote the number of bags in the cave after night $`i`$. The main result in  is that Robin wins almost surely if, and only if, the series

$$\underset{i=1}{\overset{\mathrm{}}{}}\frac{r(i)}{L(i)+r(i)}$$

diverges; otherwise Robin loses almost surely.

## 2. Strategies with bounded historical memory

We generalize the above mentioned result to the case of bounded historical memory. The typical case is that Robin can, on each night, identify the bags put in the cave on the last $`k`$ days, where $`k`$ is constant. It will turn out that already the natural strategy for historical memory $`k=1`$ is strictly stronger than the natural strategy in the memoryless case. Moreover $`k=2`$ yields a strictly stronger strategy than $`k=1`$, etc. (Theorem 3.1). In fact, the game can be analyzed in a much broader family of cases, as will be shown in the sequel.

The most general case is that Robin can, on night $`i`$, identify the bags put in the cave on the last $`b(i)`$ days, where $`b:\{0\}`$ is a function with $`b(i)i`$ for all $`i`$.222$`b(i)`$ stands for the *bound* on Robin’s historical memory on night $`i`$. (So that $`b(i)0`$ is the memoryless case studied in ). It is natural to denote this game by $`RH(r,s,b,A)`$, but our analysis below is independent of the set $`A`$, so we will simply write $`RH(r,s,b)`$. A key observation is that the following natural restriction leads to a substantial simplification of the analysis of the generalized games.

###### Restriction 2.1.

We pose the restriction that Robin cannot remember anything that he forgot earlier, that is, $`b(i+1)b(i)+1`$ for each $`i`$; equivalently, the function $`ib(i)`$ is nondecreasing.

We suggest the following deterministic and random strategies for Robin, motivated by the strategy given in  for the memoryless case: Call a bag *very old* if Robin cannot tell the day it was put in the cave. An important observation is that Robin can identify the very old bags since he can identify the bags which are *not* very old.

* On night $`i`$ Robin chooses any $`r(i)`$ many bags out of the very old bags. If there are less than $`r(i)`$ many very old bags, then Robin also chooses some of the oldest bags among the ones he remembers, so as to choose $`r(i)`$ bags in total.333To put this more precisely, on night $`i`$ Robin has a partition of the bags in the cave into disjoint (possibly empty) sets $`S_0,S_1,\mathrm{},S_{r(i)}`$ such that $`S_0`$ is the set of very old bags, and for $`k=1,\mathrm{},r(i)`$, $`S_k`$ is the set of bags put in the cave $`r(i)k`$ days ago. If $`|S_0|r(i)`$, then Robin chooses any $`r(i)`$ many bags out of the bags in $`S_0`$. Otherwise, let $`m`$ be the minimal such that $`r(i)<_{k=0}^m|S_k|`$. Then Robin takes all bags in the sets $`S_0,\mathrm{},S_{m1}`$, as well as $`r(i)_{k=0}^{m1}|S_k|`$ many bags from $`S_m`$.
* Same as $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝐷𝐸𝑇}}`$, but the $`r(i)`$ bags are chosen at random, with uniform probability, out of the older bags.444Using the notation of the previous footnote, If $`|S_0|r(i)`$, then Robin chooses at random (with uniform probablity) $`r(i)`$ many bags out of the bags in $`S_0`$. Otherwise Robin takes all bags in the sets $`S_0,\mathrm{},S_{m1}`$, as well as $`r(i)_{k=0}^{m1}|S_k|`$ many random bags from $`S_m`$.

Observe that if $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝐷𝐸𝑇}}`$ is a winning strategy for Robin, then so is $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$.

Write

$$\stackrel{~}{L}(i)=\mathrm{max}\{0,\underset{j=1}{\overset{ib(i)}{}}s(j)\underset{j=1}{\overset{i1}{}}r(j)\}$$

Then $`\stackrel{~}{L}(i)`$ is the number of bags put in the cave on days $`1,2,\mathrm{},ib(i)`$ and not removed until day $`i`$.

###### Proposition 2.2.

1. If $`ib(i)`$ is bounded, then $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝐷𝐸𝑇}}`$ is a winning strategy for Robin in $`RH(r,s,b)`$ (for each $`r`$, and $`s`$).
2. If there exist infinitely many $`i`$ such that $`\stackrel{~}{L}(i)r(i)`$, then $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝐷𝐸𝑇}}`$ is a winning strategy in $`RH(r,s,b)`$.

###### Proof.

(1) is easy. To prove (2), assume that a bag was put in the cave on day $`d`$. By (1) we may assume that $`ib(i)`$ is unbounded. Let $`i`$ be such that $`d<ib(i)`$. Restriction 2.1 ensures that this will also hold for all larger $`i`$’s, so we may assume further that $`\stackrel{~}{L}(i)r(i)`$. This means that on *night* $`i`$, all bags put on days $`ib(i)`$ (in particular, those put on day $`d`$) were removed from the cave. ∎

We may therefore make the following additional restriction.

###### Restriction 2.3.

$`\stackrel{~}{L}(i)>r(i)`$ for all but finitely many $`i`$.

###### Theorem 2.4.

Assume that $`r`$, $`s`$, and $`b`$ satisfy Restrictions 2.1 and 2.3, and Robin uses the strategy $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$.

1. If $`r(i)/\stackrel{~}{L}(i)=\mathrm{}`$, then Robin wins almost surely.
2. If $`r(i)/\stackrel{~}{L}(i)<\mathrm{}`$, then the Sheriff wins almost surely.

###### Proof.

If there is $`i`$ such that $`\stackrel{~}{L}(i)r(i)`$, let $`i^{}`$ be the maximal such $`i`$. Then on night $`i^{}`$, all bags put in the first $`i^{}b(i^{})`$ days are removed, and $`\stackrel{~}{L}(i)>r(i)`$ for all $`i>i^{}`$. Since the convergence of the series in question does not depend on the first few elements, this shows that we may assume that $`\stackrel{~}{L}(i)>r(i)`$ for all $`i`$.

Observe that if $`ib(i)`$ is bounded, then $`\stackrel{~}{L}(i)`$ is eventually equal to $`0`$, contradicting our assumption, thus $`ib(i)`$ is unbounded. Assume that a bag was put in the cave on day $`d`$, then for each large enough $`i`$, $`d<ib(i)`$ so that on night $`i`$, the probability that the bag in question is removed is $`r(i)/\stackrel{~}{L}(i)`$.

Now, the probability that a bag put in the cave on day $`d`$ is *not* eventually removed is

(1)
$$\underset{i=d}{\overset{\mathrm{}}{}}\left(1\frac{r(i)}{\stackrel{~}{L}(i)}\right).$$

The product (1) converges to $`0`$ (i.e., Robin wins almost surely) if $`r(i)/\stackrel{~}{L}(i)=\mathrm{}`$. If $`r(i)/\stackrel{~}{L}(i)<\mathrm{}`$, then the product (1) is positive for $`d=1`$. Thus, its limit when $`d\mathrm{}`$ is $`1`$. Consequently, the Sheriff wins almost surely. ∎

Note that in Theorem 2.4, (1) implies that the other direction of (2) also holds, and (2) implies that the other direction of (1) also holds. Thus, the theorem gives an exact characterization of when Robin wins almost surely and when the Sheriff wins almost surely, and shows that it is always the case that one of them wins almost surely. In the case $`b(i)0`$, Theorem 2.4 reduces to the main result of , described at the end of Section 1.

## 3. Historical memory helps

To make sure that the generalization made in 2.4 is not trivial, we must find instances where additional historical memory changes the strategy’s status from a loosing strategy to a winning one.

Assume that $`b,c:`$. We say that $`c`$ *eventually dominates* $`b`$ if there exists $`m`$ such that for all $`i>m`$, $`b(i)<c(i)`$.

###### Theorem 3.1.

Assume that $`c`$ eventually dominates $`b`$, $`ic(i)`$ is unbounded, and $`b`$ and $`c`$ satisfy Restrictions 2.1 and 2.3. Then there exist functions $`r,s:`$ such that for Robin, $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is a (surely) winning strategy in $`RH(r,s,c)`$ and an almost-surely loosing strategy in $`RH(r,s,b)`$.

###### Proof.

Since convergence of series does not depend on the first few terms, we may assume that for each $`i`$, $`b(i)<c(i)`$. Since $`\stackrel{~}{L}(i)`$ also depends on the bounding function $`b`$, let us denote it here by $`\stackrel{~}{L}_b(i)`$. It follows that for each $`r`$ and $`s`$, $`\stackrel{~}{L}_c(i)<\stackrel{~}{L}_b(i)`$ for all $`i`$. It thus suffices to consider the case where $`c(i)=b(i)+1`$ for each $`i`$, and therefore $`\stackrel{~}{L}_c(i)+s(ib(i))=\stackrel{~}{L}_b(i)`$.

At step $`i`$ of the construction we have the definition of $`s`$ at $`1,\mathrm{},ic(i)`$ and $`r`$ at $`1,\mathrm{},i1`$, and therefore $`\stackrel{~}{L}_c(i)`$ is defined. Define $`r(i)=\mathrm{max}\{i,\stackrel{~}{L}_c(i)\}`$, and define $`s`$ on $`ic(i)+1,\mathrm{},i+1c(i+1)`$ to be $`r(i)^3`$ (so that $`r(i)/\stackrel{~}{L}_b(i)=r(i)/(\stackrel{~}{L}_c(i)+s(ib(i)))=r(i)/(r(i)+r(i)^3)<1/r(i)^21/i^2`$).

Since $`r(i)\stackrel{~}{L}_c(i)`$ for each $`i`$, we have by Proposition 2.2 that $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is a winning strategy in $`RH(r,s,c)`$. Now, $`r(i)/\stackrel{~}{L}_b(i)1/i^2<\mathrm{}`$, so by Theorem 2.4, $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is an almost-surely loosing strategy in $`RH(r,s,b)`$. ∎

In particular, we have the following.

###### Corollary 3.2.

For each $`n=0,1,2,`$, there exist functions $`r,s:`$ such that for Robin, $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is a (surely) winning strategy in $`RH(r,s,n+1)`$ and an almost-surely loosing strategy in $`RH(r,s,n)`$.

## 4. Random Sheriff

The authors of  pose the question of the behavior of the Robin Hood game when the Sheriff’s strategy is random as well. Our analysis in Section 2, being independent of the Sheriff’s moves, shows that the results apply to this case as well. In addition to the first strategy of  which was described in the introduction to the present paper (Section 1), a second random strategy for Robin is sketched in , and is conjectured to be an almost-surely winning strategy in the game $`RH(1,s,[0,1])`$ where $`s`$ is constant.

While we are unable to analyze the second strategy for lack of some details, we can see that the first strategy already works. Here $`b(i)=0`$ (no historical memory) and $`L(i)=_{j=1}^is(j)r(j)=_{j=1}^is1=i(s1)`$, therefore

$$\underset{i=1}{\overset{\mathrm{}}{}}\frac{r(i)}{\stackrel{~}{L}(i)}=\underset{i=1}{\overset{\mathrm{}}{}}\frac{r(i)}{L(i)+r(i)}=\underset{i=1}{\overset{\mathrm{}}{}}\frac{1}{L(i)+1}=\underset{i=1}{\overset{\mathrm{}}{}}\frac{1}{i(s1)}=\mathrm{},$$

thus by Theorem 2.4, $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is an almost-surely winning strategy in this game. Since $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ coincides with the first strategy of , we have that the first strategy of  is also an almost-surely winning strategy against a *random* Sheriff in $`RH(1,s,[0,1])`$ (as well as any game $`RH(r,s,A)`$ with $`r(i)/(L(i)+r(i))`$ diverging). Theorem 2.4 is an extension of this phenomenon to the case of nonzero historical memory.

## 5. Open problems

Among the problems which naturally arise, the following two seem to be the most interesting.

###### Conjecture 5.1.

$`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ is the best random strategy among the strategies which are independent of the index set $`A`$.

Our analysis repeatedly uses Restriction 2.1 posed on the bounding function $`b`$.

###### Problem 5.2.

Analyze the the case where $`b`$ does not satisfy Restriction 2.1.

Finally, our strategy $`\mathrm{𝖮𝗅𝖽𝖾𝗌𝗍}_{\mathrm{𝑅𝑁𝐷}}`$ uses unboundedly many “random bits”. The referee suggests the problem whether good strategies exist, in which the number of random bits used at each specific step is bounded by some constant.


format: text

text:
# A partially quenched analysis of the 𝜂-𝜂' system in 𝑁_𝑓=2 QCD

## 1 INTRODUCTION

Flavour singlet pseudoscalar states bear witness to the topological structure of the QCD vacuum and hence constitute objects of great physical interest. In lattice gauge theory, considerable progress has been achieved recently  in directly accessing such $`\eta ^{}`$ like states from their propagators in $`N_f=2`$ QCD with Wilson fermion action, based on applying spectral approximations to two-loop correlators (hairpin diagrams).

In this note, we extend this work to the more realistic setting of partially quenched QCD with two active sea quark flavours, where valence and sea quarks are allowed to carry different masses. In our comprehensive study we analyzed some 800 independent SESAM configurations on $`16^3\times 32`$ lattices at $`\beta =5.6`$ , with four different sea quark masses, as indicated in Fig.1 111For further details of the analysis, see Refs. ..

## 2 The partially quenched scenario

Partial quenching in the flavour singlet sector of QCD is described in the framework of chiral perturbation theory by a nonsymmetric piece of the effective chiral Lagrangian

$$_0=+\frac{N_f}{2}\left[\mu ^2(\eta ^{})^2+\alpha (_\mu \eta ^{})^2\right].$$
(1)
This expression introduces two effective ‘low energy constants’: the mass gap $`\mu `$ and an interaction parameter accounting for momentum dependence, $`\alpha `$. The Lagrangian translates into two-point correlation functions (propagators) for the neutral pseudoscalar quark bilinears

$$\mathrm{\Phi }_{ii}(x)=\overline{q}_i(x)\gamma _5q_i(x),$$
(2)
where the index $`i`$ runs over the quark and pseudoquarks d.o.f., the latter being geared to act as ‘determinant eaters’ on the valence quarks .

In momentum space the disconnected part of the ps-correlators has a compact form  with the infamous double pole structure

$$D_{ij}(p)=𝒩\frac{(\mu ^2+\alpha p^2)(p^2+M_d^2)}{(p^2+M_{ii}^2)(p^2+M_{jj}^2)(p^2+M_\eta ^{}^2)}$$
(3)
The prefactor, $`𝒩=1/(1+N_f\alpha )`$ and the singlet mass, $`M_\eta ^{}^2=(M_d^2+N_f\mu ^2)/(1+N_f\alpha )`$ carry explicit $`\alpha `$-dependencies. Rememeber that the octet masses ($`M_{ii}`$,$`M_{jj}`$, and $`M_d`$) from valence and sea quarks, respectively, are readily computed from connected diagrams on the lattice. Our present concern is the determination of the mass gap, $`\mu `$, and $`\alpha `$ from our lattice data.

## 3 Results

As it has not been possible to determine the effective constant $`\alpha `$ in the past we explore the potential of our data for achieving progress on this issue by applying our spectral techniques  in two steps.

### 3.1 $`\alpha =0`$ analysis

In a first attempt, we choose $`\alpha =0`$ and study effective mass gaps as determined from the hairpin diagrams, $`D_{ij}(t)`$, by fitting, w.r.t. $`m_{0ij}(t):=\mu `$ at each value of $`t`$, the zero-momentum Fourier transforms of ansatz Eq. 3 to our data. In broken SU(3) with two active quark flavours, we introduce two isosinglet pseudoscalars, designated by the indices $`n`$ (for ‘light’, nonstrange) and $`s`$ (for strange). In this manner, we can study three types of hairpin diagrams, $`D_{nn}`$ , $`D_{ns}`$ , $`D_{ss}`$, with the strange quark mass determined on the lattice from the kaon mass (for $`\kappa `$ locations, see Fig.1). As a result, we find satisfying plateau formations in $`m_{0ij}(t)`$ as illustrated in Fig. 2. In fact, the numerical quality of our signals comfortably allows for the consecutive extrapolations (i) to light valence and then (ii) to light sea quark masses; the chiral sea quark extrapolation is exhibited in Fig. 3.

From the consecutive chiral extrapolations, we obtain for the quadratic mass matrix in the quark-flavour basis

$`:=\left(\begin{array}{cc}m_{nn}^2& m_{0,ns}^2\\ m_{0,sn}^2& m_{ss}^2\end{array}\right),`$ (6)
where $`m_{nn}^2:=M_{nn}^2+m_{0nn}^2`$. With the $`\rho `$-meson mass scale  we find

$`=\left(\begin{array}{cc}(306\pm 91)^2& \sqrt{2}(201\pm 34)^2\\ \sqrt{2}(201\pm 34)^2& (680\pm 15)^2\end{array}\right)\text{MeV}^2.`$ (9)
Diagonalization of $``$ renders

$$M_\eta =292\pm 31\text{MeV},M_\eta ^{}=686\pm 31\text{MeV}.$$
(10)
These numbers from $`N_f=2`$ look promising: the individual mass values – while being low w.r.t. the real three flavour world – yield a mass splitting that compares very well to phenomenology.

### 3.2 First determination of $`\alpha `$

We remark that the outcome of the analysis at this stage is not in accord with the tree approximation to the Lagrangian, Eq. 2, as the latter predicts the mass plateaus to be independent of the valence quark mass, $`m_v`$, contrary to the findings in Fig. 2. Therefore, we repeated the above mass plateau study on a set of nonvanishing $`\alpha `$-parameters, in an attempt to verify the compliance of our data with such independence: inspecting the plot of $`m_{0nn}`$ vs. $`m_v`$ (Fig. 4) it becomes obvious (i) that $`m_{0nn}(m_v)`$ shows a simple, linear behaviour and (ii) that $`\alpha `$ can indeed be adjusted to produce a zero slope of $`m_{0nn}(m_v)`$ 222The situation is less favourable for the gap values extracted from $`D_{ns}\text{[5]}.`$!

In this way it is straightforward to find an optimal value, $`\alpha _{opt}`$, that eliminates (for a given sea quark mass) the $`m_v`$-dependence from $`m_{0nn}`$, resulting in a universal plateau level. This collapse of data into universal plateaus (of heights $`\stackrel{~}{\mu }`$) is exemplified in Fig. 5. The comparision with the situation encountered with $`\alpha =0`$ (Fig. 2) provides clear evidence for our sensitivity in determining $`\alpha _{opt}`$.

After chiral sea quark extrapolation we arrive at the values

$$\alpha =0.028\pm 0.013\mu =203\pm 34\text{MeV}.$$
(11)
Note that the numerical value for the mass gap at the chiral point is robust w.r.t. to the two different approaches presented here. One might say that the allowance for the $`\alpha `$ parameter anticipates the chiral valence quark limits.

In conclusion: the application of spectral methods to the hairpin diagrams allows rather detailed partially quenched studies in the flavour singlet sector, with standard Wilson fermions. This provides a strong motivation to extend this work to the case of overlap fermions.

Acknowledgements We thank G. Rossi and M. Golterman for interesting discussions. The project was supported by NERSC, FZ-Jülich and EU (project HPRN-CT-2000-00145).


text:
# The abundance of massive galaxies 3–6 billion years after the Big Bang.

## Methods

### Stellar Mass Fitting

We adopt a cosmology2 of $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`H_0=70`$ $`\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1`$ and use Vega magnitudes. For the mass-function analysis we use GDDS galaxies with $`K<20.6`$ and $`0.8<z<2`$ (150 galaxies; 89% spectroscopic completeness based on Abraham et al. 17 identification confidence classes $`>=2`$; galaxies without spectroscopic redshifts are assigned to their photometric redshifts). To account for the higher priority given to red galaxies in the slit-mask design we use the ‘sampling weights’ from Abraham et al. 17 which are the selection probability for spectroscopy as a function of $`IK`$ and $`K`$.

For each galaxy we derive the most likely stellar mass and a range of uncertainty. These are derived by evaluating against model optical-IR spectral-energy distributions (SED) to determine the $`M/L_K`$ ratio, and hence the mass. We note that a variety of approaches to accomplish this have been described in the literature which vary in the level of detail in which they treat star formation. For example, Cole et al. 20 used a simple set of monotonic star-formation histories (SFH) with a varying e-folding time-scale and a fixed dust law. A potential problem is posed by the fact that real galaxies have more complex SFHs, for example a recent starburst can make an old galaxy temporarily bluer and lead to an underestimate of $`M/L`$ using this method. One approach to account for this effect is to introduce a second young SED component5 superimposed on the old population; this can be computationally expensive depending on the amount of freedom allowed for the second component.

In a spectroscopic sample there is no degeneracy between SED model fitting and photometric redshift (which is also based on SED fitting). We adopt two-component modelling (using PEGASE.2 21 to calculate spectra) in order to be able to assess the possible biases due to starbursts on the calculated masses and allow a range of dust extinction ($`0A_V2`$ mag) and metallicity ($`0.0004Z0.02`$). The primary component is modelled using a star-formation rate $`\text{SFR}\mathrm{exp}(t/\tau )`$ with $`\tau =`$0.1, 0.2, 0.5, 1, 2, 4, 8 & 500 Gyr (the first approximates an instantaneous starburst and the last a constant SFR). The secondary component is a starburst, modelled with a $`\tau =0.1`$ Gyr exponential, which can occur at any time and have a mass between $`10^4`$ and twice that of the primary component. The model age is constrained to be less than that of the Universe at the appropriate redshift, but any formation epoch is allowed. Like all studies of high redshift star-formation, we must assume a Universal Initial Mass Function (IMF). We primarily use the BG03 IMF27 which fits local cosmic luminosity densities well; it has a similar high-mass slope as the classical Salpeter IMF28 but with a more realistic break at 0.5$`M_{}`$. Stars with masses below the break never contribute significantly to optical/IR light so the overall effect is to re-scale the total stellar masses to more reasonable values. Similarly we re-scale literature Salpeter based numbers4, 5, 20 to the BG03 IMF using $`M_{SP}=1.82M_{BG}`$. Our results are robust for any reasonable choice of IMF with a similar slope; to illustrate this we have also calculated masses using the popular Kennicutt IMF29.

Our approach is to find via exhaustive grid search all models consistent with the $`VIz^{}K`$ photometry in the observed frame of each galaxy. This colour set is available for all galaxies and covers rest-frame ultraviolet through near-infrared. We do not include the actual spectra (apart from the redshift information) in the fits because of variable quality and signal-to-noise, however we find the spectral classes are broadly consistent with the best photometric SED fits. The full distribution function of allowed masses were calculated by Monte-Carlo re-sampling the photometric errors. The final masses and error bars represent the mean and standard deviation of this full distribution function. Typically, we find the masses are fitted to $`\pm 0.17`$ dex in the $`K<20.6`$ sample. The stellar masses are very robust against the details of the fitting. Using the Monte-Carlo machinery to investigate the effect of different assumptions about metallicity, dust and bursts we find that the largest effect is due to bursts. If we disallowed bursts then the masses typically decrease by only 0.2 dex. Finally we note the variation in $`M/L_K`$ over the range $`1<z<2`$ is constrained by the age of the Universe (6–3 Gyr), typically the maximum range in our sample is a factor of three.

Correspondence and requests for materials should be addressed to K.G. (kgb@pha.jhu.edu)
Acknowledgments Based on observations obtained at the Gemini Observatory, which is operated by AURA under a co-operative agreement with the NSF on behalf of the Gemini partnership: NSF (U.S.), PPARC (U.K.), NRC (Canada), CONICYT (Chile), ARC (Australia), CNPq (Brazil) and CONICET (Argentina). Also based on observations made at the Las Campanas Observatory of the Carnegie Institution of Washington. Karl Glazebrook & Sandra Savaglio acknowledge generous funding from the David and Lucille Packard Foundation. Roberto Abraham gratefully acknowledges receipt of funding from NSERC and from the Government of Ontario through a Premier’s Research Excellence Award. Hsiao-Wen Chen acknowledges support by NASA through a Hubble Fellowship grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract.
Competing interests statement The authors declare that they have no competing financial interests.


text:
# Universality in two-dimensional Kardar-Parisi-Zhang growth

## I Introduction

Surface growth processes and deposition of thin films have attracted much interest from the technological point of view  and motivated the proposal of continuum and discrete models for surface and interface growth, which also play an important role in the field of non-equilibrium Statistical Mechanics. One of the most important phenomenological theories is that of Kardar, Parisi and Zhang (KPZ) , in which the time evolution of the interface described by the height function $`h`$ at position $`\stackrel{}{x}`$ and time $`t`$ is given by the KPZ equation

$$\frac{h}{t}=\nu ^2h+\frac{\lambda }{2}\left(h\right)^2+\eta (\stackrel{}{x},t).$$
(1)
Here $`\nu `$ represents a surface tension, $`\lambda `$ represents the excess velocity and $`\eta `$ is a Gaussian noise  with zero mean and variance $`\eta (\stackrel{}{x},t)\eta (\stackrel{}{x^{}},t^{})=D\delta ^d(\stackrel{}{x}\stackrel{}{x^{}})\delta \left(tt^{}\right)`$, where $`d`$ is the dimension of the substrate.

The interface width $`\xi (L,t)=\overline{h^2}\overline{h}^2^{1/2}`$ characterizes the roughness of the interface, for growth in a substrate of length $`L`$ (overbars denote spatial averages and angular brackets denote configurational averages). For short times, the interface width scales as

$$\xi t^\beta ,$$
(2)
with $`\beta `$ called the roughness exponent. For long times, a steady state is attained and the width saturates at

$$\xi _{sat}L^\alpha ,$$
(3)
with $`\alpha `$ called the roughness exponent. Eqs. (2) and (3) correspond to limits of the dynamical scaling relation of Family and Vicsek

$$\xi L^\alpha f\left(tL^z\right),$$
(4)
where the dynamical exponent $`z=\alpha /\beta `$ characterizes the crossover from the growth regime to the steady state. Galilean invariance gives $`\alpha +z=2`$ for KPZ in all dimensions . The exact scaling exponents are known in $`d=1`$, but no exact value was already obtained in two or more dimensions .

Many discrete models fall into the KPZ class, such as the restricted solid-on-solid (RSOS) model of Kim and Kosterlitz  and ballistic deposition (BD) . Numerical estimates of the scaling exponents in $`d=1`$ are consistent with the exact values  and simulations in $`d=2`$ are frequently used to estimate them. Most of the reported values of $`\alpha `$ range from $`\alpha =0.37`$ to $`\alpha =0.4`$ , which is confirmed by numerical solutions of the KPZ equation .

In 1998, assuming certain properties of the height correlation functions, Lässig obtained a quantization condition for the KPZ exponents which gave $`\alpha =2/5`$ as the only solution consistent with the range of numerical estimates . Another consequence of his work was that the moments of the steady state height distribution,

$$W_n\overline{\left(h\overline{h}\right)^n},$$
(5)
would obey power-counting, i. e. they would scale as $`W_nL^{n\alpha }`$. The second moment is the squared interface width $`\xi ^2`$. Moreover, the validity of his assumptions requires that the steady state distribution is skewed (non-zero third moment), contrary to the one-dimensional case (Gaussian distribution).

Recent numerical results of Chin and den Nijs for the RSOS and the body-centered solid-on-solid (BCSOS) models were consistent with $`\alpha =0.4`$ , but extensive simulations of Marinari et al for the RSOS model rule out that value . A recent study of the KPZ equation in the mode-coupling approximation provided an estimate $`z1.62`$ in $`d=2`$ . Despite this controversy, various works  confirm the power-counting property for models with restricted height differences (RSOS, BCSOS) and suggest that the skewness

$$S\frac{W_3}{W_{2}^{}{}_{}{}^{3/2}}$$
(6)
and the kurtosis

$$Q\frac{W_4}{W_{2}^{}{}_{}{}^{2}}3$$
(7)
are universal at the steady state. According to Chin and den Nijs, these universal values do not depend on the sign of the coefficient $`\lambda `$ of the nonlinear term of the KPZ equation (1).

Consequently, the value of KPZ exponents in $`d=2`$ is still an open question, and the universality of the values of the skewness and the kurtosis deserves to be tested in models other than those with restricted height differences. The main contribution of this paper is the analysis of simulations results of a recently proposed model for etching of crystalline solids  which belongs to the KPZ class. Additional support to some of our conclusions will be provided by simulations results of ballistic deposition and of the temperature-dependent BCSOS model . For the etching model, we will obtain $`\alpha 0.38`$ after a detailed analysis of finite-size corrections. Characteristic relaxation times will be calculated independently and will provide estimates of the dynamical exponent $`z>1.6`$, while estimates of the growth exponent give $`\beta 0.23`$. Although our error bars intercept those obtained from extensive simulations of the RSOS model , our central estimates of $`\alpha `$ are smaller and, consequently, more distant from the theoretically proposed value $`\alpha =0.4`$ . On the other hand, our estimates are near those by Colaiori and Moore  from renormalization under the mode-coupling approximation. We will also show that our data for various models confirm the universality of the steady state skewness and kurtosis. Concerning the steady state skewness, although its absolute value is universal, its sign changes with the nonlinear coefficient $`\lambda `$ of the KPZ equation.

This paper is organized as follows. In Sec. II we briefly describe the models, the simulation procedure, the methods to estimate $`W_n`$ and the method to calculate characteristic relaxation times. In Sec. III, we analyze the skewness and the kurtosis at the steady states. In Sec. IV, we analyze the finite-size estimates of the scaling exponents of the etching model, also showing some results for BD. In Sec. V we summarize our results and present our conclusions.

## II Models and simulation procedure

The model for etching of a crystalline solid of Mello et al  is illustrated in Fig. 1 in its growth version. The solids have square and simple cubic lattice structures in $`1+1`$ and $`2+1`$ dimensions, respectively. At each growth attempt, a column $`i`$, with current height $`h(i)h_0`$, is randomly chosen. Then its height is increased by one unit ($`h(i)h_0+1`$) and any neighboring column whose height is smaller than $`h_0`$ grows until its height becomes $`h_0`$. One time unit corresponds to $`L^2`$ growth attempts in $`2+1`$ dimensions. In the true etching version of this model, the columns’ heights decrease by the same quantities above. However, in this paper we will always refer to the growth version of Fig. 1 as ”the etching model”.

In the ballistic deposition (BD) model, particles are sequentially released from randomly chosen positions above the substrate, follow a trajectory perpendicular to the surface and stick upon first contact with a nearest neighbor occupied site .

In the BCSOS model defined by Chin and den Nijs , the substrate is a square lattice and the heights $`h`$ in the first (second) sublattice are restricted to assume even (odd) values. Also, the nearest neighbor columns always differ in height by $`\mathrm{\Delta }h=\pm 1`$. The energy of a given heights configuration $`\{h\}`$ is given by $`E\left(\{h\}\right)=_{i,j}\frac{1}{4}K\left(h_ih_j\right)^2`$, where $`K`$ is an inverse temperature parameter and the sum runs over all next nearest neighbor pairs. At each deposition attempt, a column $`c`$ is randomly chosen and, if the constraint of the heights difference is satisfied, then $`h(c)h(c)+2`$ with probability $`pmin(1,exp\left(\mathrm{\Delta }E\right))`$, where $`\mathrm{\Delta }E`$ is the energy change if the deposition takes place. The corresponding coefficient $`\lambda `$ (Eq. 1) changes sign at a critical point $`K_c`$ .

The etching model was simulated in lattices of lengths $`L=2^n`$ ($`n=5`$ to $`n=10`$) and lattices of lengths $`L=2^m\times 50`$ ($`m=0`$ to $`m=4`$). The maximum deposition time ranged from $`10^3`$ for the smallest lattices to $`6\times 10^4`$ for the largest ones. During half of this time or more, the systems have undoubtedly attained their steady states. For the smallest lattices, $`2\times 10^4`$ realizations were simulated, and nearly $`300`$ realizations for the largest lattices ($`L=800`$ and $`L=1024`$). BD was simulated in lattices of lengths $`L=2^n`$, from $`n=5`$ ($`L=32`$) to $`n=11`$ ($`L=2048`$). For the analysis of the skewness and of the kurtosis, it was essential to simulate this model in very large lattices, thus the number of realizations was relatively small: typically $`10^3`$ realizations for $`L512`$, $`35`$ for $`L=1024`$ and $`8`$ for $`L=2048`$. The temperature-dependent BCSOS model was simulated with $`K=0.25`$ ($`L=16`$ to $`L=128`$) and $`K=1.0`$ ($`L=16`$ to $`L=512`$). For $`L256`$, $`10^3`$ different realizations were considered, and $`40`$ realizations for $`L=512`$.

The procedure to calculate average quantities, described below, followed the same lines for all models. It was previously used in the analysis of other growth models in $`1+1`$ and $`2+1`$ dimensions .

One important point is the criterion to determine the initial time $`t_{min}`$ for estimating average quantities at the steady state. In this regime, the interface width $`\xi `$ fluctuates around an average value instead of increasing systematically, which was the case in the growth and in the crossover regions. Thus, for a fixed lattice length $`L`$, the first step was choosing a time interval $`t_{min}tt_{max}`$, with nearly constant $`\xi `$, from visual inspection of the $`\xi \times t`$ plot ($`t_{max}`$ was always the maximum simulation time). Subsequently, two tests were performed. For the first test, the interval was divided in $`5`$ subintervals and the average value of $`\xi `$ was calculated in each one, forming a sequence of estimates $`\{\xi (i)\}`$, with $`i=1,\mathrm{}5`$. If $`\xi (i)<\xi (i1)`$ at least two times along this sequence, then the average value of $`\xi `$ in the region $`t_{min}tt_{max}`$, $`\xi _{sat}^{trial}`$, was calculated. In the second test, from the plot of $`\mathrm{log}(\xi _{sat}^{trial}\xi (t))\times t`$, for $`t`$ in the crossover region, we obtained a rough estimate of the characteristic time $`\tau `$ of relaxation to the steady state, as shown in Ref. . If $`t_{min}>10\tau `$, then the interval $`t_{min}tt_{max}`$ was accepted as representative of the steady state. Otherwise, a larger value of $`t_{min}`$ was chosen and the tests were repeated (it seldom occurred).

In order to estimate the moments of the height distribution and their error bars, we used their average values within the five subintervals defined above. Final estimates of $`W_n`$ are averages of these values, and error bars were obtained from their variances.

The dynamical exponent was estimated using a recently proposed method , in which a characteristic time $`\tau _0`$, proportional to the relaxation time $`\tau `$, is calculated. For fixed $`L`$, after calculating the saturation width $`\xi _{sat}(L)`$, $`\tau _0`$ is defined through

$$\xi (L,\tau _0)=k\xi _{sat}\left(L\right),$$
(8)
with a constant $`k1`$. From relation (4), we obtain

$$\tau _0L^z.$$
(9)
Typically, the uncertainty in $`\tau _0`$ is much smaller than that of $`\tau `$, estimated from $`\mathrm{log}(\xi _{sat}\xi )\times t`$ plots - see Ref. .

We estimated $`\tau _0`$ with $`k`$ ranging from $`k=0.5`$ to $`k=0.8`$. Although these constants are not much different, the values of $`\tau _0`$ for $`k=0.5`$ and $`k=0.8`$ typically differ by a factor $`4`$, for fixed $`L`$. This method was already applied with success to calculate the dynamical exponent of other models in various universality classes .

## III Universality of skewness and kurtosis

In Fig. 2a we show the steady state skewness $`S\left(L,t\mathrm{}\right)`$ of the etching model as a function of the inverse lattice length, which provides a good linear fit of the data, with an asymptotic estimate $`S=0.26\pm 0.01`$. The absolute value of $`S`$ agrees with results of Chin and den Nijs , Shim and Landau  and Marinari et al  for the RSOS model of Kim and Kosterlitz  and for BCSOS models, but those authors obtained $`S`$ with negative sign. Here, the positive sign is related to the presence of sharp hills (see process at the left in Fig. 1) and wide valleys at the surface of the deposit, the opposite being observed in RSOS deposits.

Note that Eq. (1) is invariant under the transformations $`hh`$ and $`\lambda \lambda `$, without changing the other parameters. This transformation changes the sign of the skewness. Consequently, we expect that the sign of $`S`$ is related to that of $`\lambda `$, as previously observed in the growth regimes of $`1+1`$-dimensional KPZ systems . In fact, in the true etching version of this model, with erosion leading to decreasing heights, the sign of $`S`$ changes, corresponding to the transformation $`hh`$, $`\lambda \lambda `$.

Results for the other models contribute to this discussion. In Fig. 2b we show $`S\left(L,t\mathrm{}\right)`$ versus $`1/L`$ for the BD model. Note that $`S`$ is negative for small lattices (typically $`L<500`$), but for large lattices it becomes positive, showing that there are significant morphological differences between the steady states of small lattices and those of very large systems. This is the main reason to avoid extrapolating the data in Fig. 2b, even choosing extrapolation variables other than that $`1/L`$ (this abscissa in Fig. 2b was chosen only for illustrate the evolution of $`S`$ with $`L`$). Other consequences of this complex finite-size behavior were previously discussed in Ref. . However, it is clear from Fig. 2b that $`S`$ is asymptotically positive.

In Figs. 2c and 2d we show the steady state skewness for the temperature-dependent BCSOS models with $`K=0.25`$ and $`K=1.0`$, respectively. For $`K=0.25`$, the extrapolation of the three last data points give $`S=0.28\pm 0.015`$, reinforcing the conclusion on the universality of $`|S|`$. Due to the small number of data points, the asymptotic correction term may be other than $`1/L`$, but it would not affect that conclusion. For $`K=1.0`$, the skewness is always positive and rapidly increases with $`L`$, thus no extrapolation variable of the form $`L^\mathrm{\Delta }`$, with $`\mathrm{\Delta }>0`$, provides a reasonable linear fit. However, $`S`$ will certainly converge to a positive value as $`L\mathrm{}`$, showing that the sign of the steady state skewness changes as $`K`$ increases.

The sign of the corresponding parameter $`\lambda `$ is obtained from the size dependence of the growth velocity. The steady state growth velocity, $`v_s\left(L\right)`$, and the velocity in an infinitely large substrate at long times, $`v_{\mathrm{}}`$, obey the relation

$$v_s\left(L\right)=v_{\mathrm{}}a\lambda L^\alpha _{},$$
(10)
where $`\alpha _{}=2\left(1\alpha \right)`$ and $`a`$ is positive. Considering $`\alpha =0.38`$ (see Sec. IV), we obtain $`\alpha _{}=1.24`$. In Figs. 3a-3d, we plotted $`v_s\left(L\right)`$ versus $`x1/L^{1.24}`$ for the four models. For the etching model, the BD model and the BCSOS model with $`K=1.0`$ (Figs. 3a,b,d), $`v_s`$ decreases with $`x`$, which gives a positive $`\lambda `$, while the opposite occurs in the BCSOS model with $`K=0.25`$. It confirms that the steady state skewness of the KPZ equation has the same sign of the parameter $`\lambda `$, its absolute value being universal. The same conclusion was derived from the time dependence of the velocity in the growth regime .

At this point, we recall that Chin and den Nijs  simulated the temperature-dependent BCSOS model with $`K=0.25`$ and $`K=0.25`$, obtaining $`S0.26`$ in both cases. They concluded that $`S`$ did not depend on $`\lambda `$, but they did not estimate the values of this parameter. Since $`\lambda `$ is expected to cross zero for a positive $`K`$, we conclude that their work failed to consider the regime of positive $`\lambda `$, which is represented here by $`K=1.0`$.

Now we turn to the analysis of the kurtosis (Eq. 7). In Fig. 4a we show the steady state $`Q`$ for the etching model, as a function of $`1/L`$. The size dependence is much weaker than that of the skewness, so that extrapolation variables other than $`1/L`$ do not have a significant influence on the asymptotic estimate, $`Q=0.134\pm 0.015`$ (the large error bar is mainly a consequence of the uncertainties of the finite-size data). This estimate agrees with previous ones for models with restricted heights differences , suggesting that the steady state kurtosis is also universal. In Fig. 4b, we show the data for the BD model and for the BCSOS model with $`K=0.25`$ and $`K=1.0`$, which have significant finite-size dependence (in particular those for BD). Thus, the extrapolated values have very large error bars, but are still consistent with a universal value of $`Q`$.

## IV Roughness, dynamical and growth exponents

Our first step to estimate the roughness exponent was to calculate effective exponents $`\alpha _{(L,i)}`$ defined as

$$\alpha _{(L,i)}\frac{\mathrm{ln}\left[\xi _{sat}\left(L\right)/\xi _{sat}\left(L/i\right)\right]}{\mathrm{ln}i}.$$
(11)
It is expected that $`\alpha _{(L,i)}\alpha `$ for any value of $`i`$.

Using different values of $`i`$ in Eq. (11), we noticed that $`\alpha _{(L,i)}`$ varied with $`L`$ typically in the range $`0.33\alpha _{(L,i)}0.38`$ for $`50L1024`$, which suggests that corrections to the scaling relation (3) are relevant. Our first proposal is to assume the main scaling correction as

$$\xi _{sat}L^\alpha (a_0+a_1L^\mathrm{\Delta }),$$
(12)
where $`a_0`$ and $`a_1`$ are constants. Consequently, $`\alpha _{(L,i)}`$ is expected to vary as

$$\alpha _{(L,i)}\alpha +BL^\mathrm{\Delta },$$
(13)
where

$$B=\frac{(1i^\mathrm{\Delta })}{\mathrm{ln}(i)}\frac{a_1}{a_0}.$$
(14)
Our data with $`50L1024`$ were analyzed using four values of $`i`$ in Eq. 11: $`i=2`$, $`i=2.56`$, $`i=3.125`$ and $`i=4`$ (for noninteger $`i`$, only three or four effective exponents can be calculated). For each $`i`$, we plotted $`\alpha _{(L,i)}\times L^\mathrm{\Delta }`$ using several exponents $`\mathrm{\Delta }`$ and did least squares fits of those plots, from which the linear correlation coefficients $`r(\mathrm{\Delta },i)`$ were obtained. Since there is no argument to predict the value of $`\mathrm{\Delta }`$, we adopted the condition of maximizing the coefficient $`r`$ to extrapolate our data. In Table I, we show the exponents $`\mathrm{\Delta }`$ which gave the largest $`r`$ (best linear fits) for each $`i`$. The procedure is illustrated in Figs. 5a and 5b, in which we show $`\alpha _{(L,2)}\times L^{0.55}`$ and $`\alpha _{(L,2)}\times L^{0.65}`$, respectively, with the corresponding linear fits. The values of $`\mathrm{\Delta }`$ in Figs. 5a and 5b are those which give the best fits with $`i=2`$ and $`i=4`$, respectively (see Table I).

$`\mathrm{\Delta }`$ is expected to be independent of the particular choice of $`i`$ in Eq. (11), so the differences between the estimates in Table I are effects of the maximization of $`r`$. Moreover, other exponents $`\mathrm{\Delta }`$ near the values shown in Table I also provided reasonable linear fits of $`\alpha _{(L,i)}\times L^\mathrm{\Delta }`$ plots. In other words, large linear correlation coefficients were also obtained by considering $`0.3\mathrm{\Delta }0.8`$ for different choices of $`i`$. Consequently, it is not possible to obtain a reliable estimate of that correction exponent. On the other hand, the asymptotic $`\alpha `$ obtained from the same fits fluctuate within a narrow range. Accounting for the error bars of the data, we obtained $`\alpha =0.385\pm 0.01`$ for $`i=2`$ and $`\alpha =0.382\pm 0.01`$ for $`i=4`$.

We also checked the effect of considering a fixed correction exponent $`\mathrm{\Delta }=0.55`$ (which is near the values in Table I) for all values of $`i`$. This procedure will be called fixed $`\mathrm{\Delta }`$ method. Least squares fits of the $`\alpha _{(L,i)}\times L^{0.55}`$ plots were performed, providing the slopes $`B`$ which are shown in Table I. From Eq. (14), it is expected that $`\frac{B\mathrm{ln}i}{(1i^\mathrm{\Delta })}=\frac{a_1}{a_0}=const`$ (see also Eqs. 12 and 13), thus we also showed in Table I the corresponding estimates of $`\frac{B\mathrm{ln}i}{(1i^\mathrm{\Delta })}`$. This quantity fluctuates with $`i`$, indicating that there is no systematic trend in our results due to choice of different values of $`i`$ for calculating effective exponents. Also notice that the estimates of $`\alpha `$ from the fixed $`\mathrm{\Delta }`$ method are in same range of the those obtained from the maximization of correlation coefficients.

Since the range of lattice lengths considered here is not very large and the correction exponents $`\mathrm{\Delta }`$ estimated above are relatively small, we tried to improve our analysis with a different assumption for the scaling corrections. Contrary to the previous procedure, now we will consider a well defined form for the main scaling correction, which is an additional constant term $`\xi _I^2`$ in the dynamic scaling relation (4) for the squared interface width:

$$\xi ^2(L,t)=\xi _I^2+L^{2\alpha }g\left(tL^z\right),$$
(15)
where $`g`$ is a scaling function. $`\xi _I`$ is called intrinsic width and represents contributions of small length scales fluctuations, typical of models with large local heights differences  such as the etching model and the BD model. From Eqs. (15) and (12), the assumption of the intrinsic width as the most relevant subleading correction corresponds to a (fixed) correction exponent $`\mathrm{\Delta }=2\alpha 0.8`$. It is slightly larger than the typical values obtained in the previous analysis.

Effective exponents $`\alpha _L^{\left(I\right)}`$ which cancel the contribution of $`\xi _I^2`$ are defined as

$$\alpha _L^{\left(I\right)}\frac{1}{2}\frac{\mathrm{ln}\left[\xi _{sat}^2\left(2L\right)\xi _{sat}^2\left(L\right)\right]/\left[\xi _{sat}^2\left(L\right)\xi _{sat}^2\left(L/2\right)\right]}{\mathrm{ln}2}.$$
(16)
In Fig. 6a we show $`\alpha _L^{\left(I\right)}`$ versus $`1/L`$ for the etching model and a least squares fit of these data, which provides $`\alpha =0.383`$ asymptotically. Here, the variable $`1/L`$ in the abscissa was chosen only to illustrate the behavior of the $`\alpha _L^{\left(I\right)}`$ data. It represents a second correction term for the $`\xi `$ scaling, which is still more difficult to measure than the first correction term. Thus, we also tested other variables in the form $`L^{\mathrm{\Delta }_1}`$ to extrapolate $`\alpha _L^{\left(I\right)}`$, with $`0.5\mathrm{\Delta }_12`$ ($`\mathrm{\Delta }_1=1`$ was used in Fig. 6a). The corresponding linear fits give $`0.380<\alpha <0.387`$. The small range of the asymptotic $`\alpha `$ is a consequence of the slow variation of $`\alpha _L^{\left(I\right)}`$ with $`L`$ ($`2\%`$ from $`L=100`$ to $`L=500`$), as shown in Figs. 6a. Accounting for the error bars of the data, our final estimate is $`\alpha =0.383\pm 0.008`$.

The forms of finite-size corrections analyzed above, Eqs. (13) and (15), cannot be rigorously justified, but are based on heuristic arguments. See e. g. Ref.  and references therein. Certainly, the fact that $`\alpha _L^{\left(I\right)}`$ increases slowly with $`L`$ is a support to the assumption that the intrinsic width is the most relevant correction term to $`\xi ^2`$. Anyway, using different assumptions on the scaling corrections was essential to confirm the reliability of the above estimate of the roughness exponent.

The universality of the values of skewness and kurtosis imply that $`W_3`$ and $`W_4`$ may also be used to estimate $`\alpha `$. The effective exponents obtained from the third moment have very large error bars, but those obtained from the fourth moment behave similarly to the ones obtained from the interface width. The analogous of exponents $`\alpha _{(L,2)}`$ and $`\alpha _{(L,4)}`$ (Eq. 11) calculated with $`W_4`$ also converge to the range $`0.38\alpha 0.385`$ with strong corrections to scaling. Effective exponents which cancel the contribution of a constant additive term in the dynamic scaling relation for $`W_4`$ (analogous to $`\alpha _L^{\left(I\right)}`$) are defined as

$$\alpha _L^{(I,4)}\frac{1}{4}\frac{\mathrm{ln}\left[W_4\left(2L,t\mathrm{}\right)W_4\left(L,t\mathrm{}\right)\right]/\left[W_4\left(L,t\mathrm{}\right)W_4\left(L/2,t\mathrm{}\right)\right]}{\mathrm{ln}2}.$$
(17)
They are plotted in Fig. 6b as a function of $`1/L`$. The asymptotic estimate, obtained with the procedure described above, is $`\alpha =0.379\pm 0.012`$. It is in good agreement with the estimate from the interface width and also excludes $`\alpha =0.4`$.

The same analysis was also performed with the BD model, as shown in Fig. 7, with data for $`L512`$. Although the results are less accurate than those for the etching model, they also suggest that $`\alpha <0.4`$ asymptotically.

In order to estimate the dynamical exponent $`z`$, we calculated effective exponents $`z_L`$ defined as

$$z_L\frac{\mathrm{ln}\left[\tau _0\left(2L\right)/\tau _0\left(L/2\right)\right]}{\mathrm{ln}4},$$
(18)
using the characteristic times $`\tau _0`$ defined in Sec. II (Eqs. 8 and 9). In Figs. 8a and 8b we show $`z_L`$ versus $`1/L`$ for the etching model, obtained using $`k=0.6`$ and $`k=0.8`$ to calculate $`\tau _0`$, respectively. Here, the abscissa $`1/L`$ is also chosen to illustrate the evolution of $`z_L`$, but not to perform extrapolations of the data. Although the error bar of the data for $`L=800`$ is relatively large, those plots indicate that $`z>1.6`$. Considering the trend for large $`L`$ and different values of $`k`$, we estimate $`1.605z1.64`$. These values are consistent with the above estimates of $`\alpha `$ and the exact relation $`\alpha +z=2`$.

Finally, we estimated the growth exponent (Eq. 2) of the etching model with the same procedure previously applied with success to BD and to the Das Sarma and Tamborenea model in $`1+1`$ dimensions . The growth region for each $`L`$ begins at $`t_0=50`$ and ends at the maximum time $`\tau _{max}`$ such that the linear correlation coefficient of the data in the range $`t_0t\tau _{max}`$ exceeds a fixed value $`r_{min}`$ . Here, $`r_{min}=0.99995`$ and $`r_{min}=0.9999`$ are considered. Effective exponents $`\beta _L`$ are defined as the slopes of the linear fits of $`\mathrm{log}W\times \mathrm{log}t`$ plots using all data in the above-defined growth regions. In Fig. 9 we plot $`\beta _L`$ versus $`1/L`$ for the etching model using $`r_{min}=0.99995`$ and $`r_{min}=0.9999`$. The error bars of those effective exponents is very small. We also show in Fig. 9 the linear fits of the $`\beta _L\times 1/L`$ data for each $`r_{min}`$. Other variables in the form $`L^\mathrm{\Delta }`$ were used to extrapolate the $`\beta _L`$ data, giving asymptotic estimates $`\beta =0.229\pm 0.005`$. Within error bars, it agrees with the value $`\alpha /z=0.234\pm 0.009`$ obtained from the above estimates of $`\alpha `$ and $`z`$.

## V Conclusion

We studied four $`2+1`$ dimensional discrete growth models in the KPZ class, determining critical exponents and steady state values of the skewness $`S`$ and the kurtosis $`Q`$. Accurate estimates of the scaling exponents were obtained for the etching model proposed by Mello et al : $`\alpha =0.383\pm 0.008`$, $`1.605z1.64`$, $`\beta =0.229\pm 0.005`$. Results for the ballistic deposition model also indicate that $`\alpha <0.4`$. The presence of the intrinsic width as the main correction to the interface width scaling was considered to extrapolate the simulations data. We also obtain the estimates $`S=0.26\pm 0.01`$ and $`Q=0.134\pm 0.015`$ in the steady state regime of the etching model. Results for the BD model and of the BCSOS model, together with previous results for the RSOS model, suggest that the absolute value of $`S`$ and the value of $`Q`$ are universal, the sign of the skewness being the same of the parameter $`\lambda `$ of the corresponding KPZ equation.

The above estimate intercepts the error bar of the roughness exponent of the RSOS model given by Marinari and co-workers , $`\alpha =0.393\pm 0.003`$. However, our central estimate is significantly lower than theirs, and the result for BD confirm this trend. On the other hand, our estimate is very near that by Colaiori and Moore, $`\alpha 0.38`$, from renormalization methods . Additional support to our conclusions was given by the independent calculation of exponents $`z`$ and $`\beta `$, contrary to recent simulation works on these lines , which were limited to the calculation of exponent $`\alpha `$.

We believe that much more accurate estimates of $`\alpha `$ are difficult to be achieved with numerical simulations of this type of lattice model. However, we consider that this work provides a significant amount of numerical results indicating that the theoretically proposed value $`\alpha =0.4`$  is not valid, within the limits of the assumptions made about the form of finite-size scaling corrections. This result and the evidence of universality of the values of the skewness and the kurtosis may motivate further analytical (maybe also numerical) studies of the KPZ theory in $`2+1`$ dimensions.


text:
# The effect of morphology on the superconductor-insulator transition in 1-D nanowires

## Abstract

We study the effect of morphology on the low temperature behavior of superconducting nanowires of length $``$100 nm. A well-defined superconductor-insulator transition (SIT) is observed only in homogenous wires, in which case the transition occurs when the normal resistance is close to $`h/4e^2`$. Inhomogeneous wires, on the other hand, exhibit a mixed behavior, such that signatures of the superconducting and insulating regimes can be observed in the same sample. The resistance versus temperature curves of inhomogeneous wires show multiple steps, each corresponding to a weak link constriction (WLC) present in the wire. Similarly, each WLC generates a differential resistance peak when the bias current reaches the critical current of the WLC. Due to the presence of WLCs an inhomogeneous wire splits into a sequence of weakly interacting segments where each segment can act as a superconductor or as an insulator. Thus the entire wire then shows a mixed behavior.

Evidence for a superconductor-insulator quantum phase transition (SIT) in one-dimensional (1D) wires has been found in a number of experiments.Dynes ; Bezryadin  Yet, other studies have demonstrated a crossover, as opposed to an SIT, in thin wires where superconductivity disappears gradually, as diameter is reduced, presumably due to an increasing number of quantum phase slips (QPS).Giordano ; Lau ; Lau2  Thus the existence and possible origins of superconductor-insulator transitions in 1D remain important open problems. In particular it is not known how the SIT depends on the morphology of nanowires.

In two-dimensional system, for example, the crucial role of morphology (i.e. granularity) on the SIT is well known.Frydman ; Garno ; Finkelstein  For uniform films, as the film thickness is reduced, a reduction of the critical temperature is observed while the superconducting transition remains sharp. The SIT occurs when the square resistance of the film reaches a critical value close to the quantum resistance $`R_Q`$=$`h/4e^2`$=6.5 k$`\mathrm{\Omega }`$.Goldman  In nonhomogeneous (granular) films, on the other hand, a reduction of the film thickness results in a crossover between superconducting and insulating regimes, with very broad resistive transitions in the thinnest superconducting samples.Frydman ; Dynes2 ; Frydman2  Our goal here is to determine the morphological requirements for 1D nanowires under which an SIT can occur.

In this Communication, we present a comparative study of homogeneous and inhomogeneous nanotube-templated wires and find a qualitatively different behavior at low temperatures. Homogeneous samples (of length $``$100 nm) show an SIT which occurs when the wire’s normal resistance is close to $`R_Q`$, confirming previous results where different nanotubes were utilized as substrates.Bezryadin  Inhomogeneous wires, on the other hand, exhibit a mixed behavior displaying properties of superconductors and insulators at once. Such samples frequently show multiple steps in the resistive transitions but no resistive tails typical of QPS.Giordano ; Lau ; Lau2  We propose a model which regards the inhomogeneous wire as a sequence of weak link constrictions (WLC) connected in series. Each WLC has a certain dissipative size and corresponding normal resistance, depending on which the WLC can be either superconducting or insulating. Homogeneous (and short enough) wires do not show such mixed behavior because each wire acts as a single WLC.

The wires were fabricated using a molecular templating technique,Bezryadin  in which a single-walled carbon nanotube was suspended over a 100 nm wide trench etched into a multilayered Si/SiO2/SiN substrate. Unlike previous studies,Bezryadin ; Lau  here we use fluorinated single-wall nanotubes (FSWNT), which are known to be insulating.Margrave  The substrate with suspended nanotubes was then sputter-coatedSputtering  in one of two ways: (i) with amorphousGraybeal  Mo0.79Ge0.21 or (ii) with a slightly thinner amorphous Mo0.79Ge0.21 film followed (in the same vacuum cycle) by a 2 nm Si film. Contact pads were defined using photolithography, followed by (i) wet etching in H2O2 for MoGe or (ii) reactive ion etching followed by wet H2O2 etching for Si-coated samples. While both processes produce samples with the same geometry (Fig. 1a) and similar dimensions (Table 1), these two fabrication methods result in different wire morphologies as revealed by TEM imaging and confirmed by transport measurements (see below). The first process results in homogeneous nanowires (Fig. 1b), which we will refer to as “bare” wires. The second method yields inhomogeneous wires, which will be called “Si-coated” wires (Fig. 1c). Although these two fabrications methods give qualitatively different wires, the exact mechanism of this is not clear. The most probable explanation is that Si introduces some surface tension that leads to the creation of grains. On the other hand, if exposed to air, the outer layer of MoGe wires oxidizes unless coated by a protective film, in our case Si. To ensure that the conducting cores of bare wires have similar dimensions as the cores of the Si-coated wires, the bare wires must be made thicker to compensate for the expected surface oxidation. Since thicker wires are always considerably more homogeneous than thin ones, the bare wires exhibit a much higher degree of homogeneity.

The oxidized layer thickness (in bare wires) can be estimated as follows (Fig. 2). For each sample, the wire length, $`L`$, divided by the normal state resistance, $`R_N`$, is plotted versus SEM measured width, $`w`$ in Fig. 2a. The bare wire data (circles) can be well approximated by a linear fit (dashed line), providing evidence for their homogeneous structure. The fitting parameters give an oxidized layer thickness $``$2.3 nm (half of the fit’s x-axis intercept value) and the resistivity for MoGe $`\rho `$$``$180 $`\mu \mathrm{\Omega }\text{-cm}`$ (obtained from the fit’s slope and the wire’s average thickness, $`d_{ave}`$$``$4.7 nm, estimated under the assumption that the top 2.3 nm of the sputtered MoGe is oxidized), in agreement with published values.Bezryadin ; Graybeal  No reasonable linear fit could be obtained for the Si-coated wires (crosses), further indicating their inhomogeneous structure. In Fig. 2b we compare the $`T_c`$’s of bare MoGe films and Si-coated MoGe films (circles and crosses, respectively) plotted versus their thicknesses. By shifting the data for bare films by 2.7 nm to the left the two families of data points overlap (Fig. 2b). It is therefore concluded that the oxidized layer thickness is 2.7 nm, similar to the above estimate.

Voltage vs. current measurements, $`V(I)`$, were performed by current biasing the sample through a large ($``$1 M$`\mathrm{\Omega }`$) resistor. Zero-bias resistance, $`R(T)`$, was obtained from the slope of the linear part of the $`V(I)`$ curves as temperature was varied. Similar to Ref. Rogachev, , transport measurements were performed in 4He or 3He cryostats equipped with rf-filtered leads. The differential resistance vs. bias current, $`dV(I)/dI`$, was measured using an AC excitation on top of a DC current offset generated by a low-distortion function generator (SRS-DS360), again in series with a ($``$1 M$`\mathrm{\Omega }`$) resistor.

Zero-bias resistance vs. temperature measurements are compared for homogeneous (Figs. 3a,b) and inhomogeneous (Figs. 3c,d) nanowires. First note that the rightmost superconducting transition observed on all samples is due to the electrodes (i.e. thin MoGe films connected to the ends of the wire). Below the temperature of this transition ($`T_{c,film}`$) only the wire contributes to the measured resistance. Therefore, in what follows we will only be interested in temperatures $`T`$$`<`$$`T_{c,film}`$. We first discuss the family of $`R(T)`$ curves for bare wire samples, presented in Fig. 3a, which show a clear dichotomy. In the log-linear representation, samples A1-A8 all have negative curvature ($`d^2\text{log}(R)/dT^2`$$`<`$0) with more than half of these samples reaching immeasurably low resistances. Therefore, we refer to these samples as superconducting. Although samples A6, A7, and A8 do not reach such low resistances, their negative curvature as well as their decreasing resistance with decreasing temperature (Fig. 3b) lead us to believe they too are superconducting. Due to their smaller widths we assume that their critical temperatures are suppressedOreg  such that they do not go through their full superconducting transition within the temperature range studied. Samples A9-A12 show only the transition due to the film electrodes. Their curvature is always positive (for $`T`$$`<`$$`T_{c,film}`$) with increasing resistances as temperature is decreased. Consequently we consider these samples to be insulating. The insulating regime can be either due to a suppression of $`T_c`$ to zeroOreg  or due to a proliferation of QPSZaikin  (or both). Observation of these two regimes indicates that a superconductor-insulator transition occurs in this family of nanowires, confirming previous results.Bezryadin  As illustrated in Fig. 3b, the transition from superconducting to insulating behavior occurs when the wire’s $`R_N`$$``$$`R_Q`$. Therefore it can be suggested that the $`R_N`$ of the wire is the parameter that controls the SIT. However, since wire length was not varied over a wide range we cannot exclude the possibility that the cross-sectional area of the wire is the true control parameter.

Measurements presented in Fig. 3 demonstrate that bare and Si-coated wires are qualitatively different. In the Si-coated wires the $`R(T)`$ curves are less predictable, show multiple steps or “humps” (Figs. 3c,d), and in some cases exhibit a mixed behavior that does not always allow them to be clearly identified as either superconducting or insulating. For example, sample B4 shows three steps on the $`R(T)`$ curve and a tail at the lowest temperature, which could be the beginning of the fourth transition. The first step, as always, reflects the superconducting transition in the electrodes, but the other steps correspond to some weak links in the wire itself. Some Si-coated wires exhibit both positive and negative curvature at different temperatures. Such mixed behavior is seen in samples B5 and B6, with their resistance initially dropping as in superconducting wires (i.e. with a negative curvature) but then starting to increase as in insulating wires. The conclusion is that no clear SIT at $`R_N`$$``$$`R_Q`$ can be observed in the family of non-homogeneous wires.

A clear difference between bare and Si-coated wires is also found in the $`dV(I)/dI`$ measurements (Fig. 4). Bare MoGe wires (left column) again show one of only two distinct types of behavior. The superconducting wires exhibit a large minimum centered around zero-bias (e.g. A4 and A5) and peaks at the bias current which is equal to the critical current of the nanowire. The insulating samples show a single zero-bias resistance peak (e.g. A11 and A12). This is in agreement with previous experiments.Bezryadin

The Si-coated MoGe wires, on the other hand, frequently show various combinations of such types of behavior. For example, samples B5 has a clear resistance depression in the range $`\pm `$150 nA with the critical current peaks at the limits of this depression, and, in addition to this, a narrow zero-bias resistance peak. Such mixed behavior indicates that some parts of the wire are superconducting while some other parts act like insulating wires. Sample B6 also shows similar results. This fact again confirms that these wires are not homogeneous.

The properties of inhomogeneous wires can be understood by assuming that they contain a sequence of independent weak link constrictions, each surrounded by a dissipative region, the size of which determines the normal resistance of the WLC. Since the free energy barriers for phase slips are different in each such WLC and the constrictions are connected in series it is clear that $`R(T)`$ curves should show crossingsLau2  and multiple steps (such as for sample B4 in Fig. 3c). Moreover some links can be insulating if their normal resistance is large enough, in the same sense as homogeneous insulating wires. For example sample B5 shows a mixed behavior that is explained by assuming it contains two WLC - one superconducting and one insulating. The $`R(T)`$ curve then should show one superconducting hump and a tail with increasing resistance at lower temperatures, typical of insulating wires. Likewise, a single superconducting feature and an insulating zero-bias peak should be observed in the $`dV(I)/dI`$ curve. These characteristics are indeed observed (Fig. 3d and Fig. 4). In some wires all WLC can be superconducting, as for example in sample B1 (and possibly B2, B3 and B4). In the case of B1, two humps in the $`R(T)`$ curve (at $`T`$$``$5 K and $`T`$$``$2 K) and two superconducting peaks on $`dV(I)/dI`$, corresponding to two critical currents (at $`I`$$``$15 nA and $`I`$$``$60 nA), are observed, indicating that two superconducting WLCs are present. Finally, since each insulating WLC produces a single $`dV(I)/dI`$ peak located at zero current, the wires with a few insulating WLCs should nevertheless exhibit only one peak at zero bias, as could be the case for sample B7 in Fig. 4.

All measured homogeneous samples act as superconductors if they satisfy $`R_N`$$`<`$$`R_Q`$ and as insulators otherwise (Fig. 3b). Inhomogeneous wires frequently violate this condition, as for example sample B2 (Fig. 3d). This can be explained by the WLC model and by assuming that each WLC can be either insulating or superconducting depending on its own normal resistance. The sample B2 has shown two distinguishable peaks on its $`dV(I)/dI`$ curve (not shown), corresponding to two superconducting WLCs with slightly different critical currents. Since this sample’s normal resistance is 8.30 k$`\mathrm{\Omega }`$ and since two similar WLCs are detected in this wire, it can be estimated that the normal resistance of each WLC is $``$4.15 k$`\mathrm{\Omega }`$. Thus each WLC is superconducting and hence the entire wire is superconducting as well.

Since each WLC is manifest by a jumpwise increase of the sample resistance at the WLC’s critical current, it is reasonable to associate WLCs with phase slip centersSkocpol ; Tidecks  (which are usually positioned at weak spots along the wire). The normal resistance of a WLC is then the resistance of the dissipative region or the region which is populated by quasiparticles generated by the phase slip center. This dimension is the quasiparticle diffusion length, $`\mathrm{\Lambda }_Q`$$``$$`(D\tau _E)^{1/2}`$$``$100 nm, where $`D`$$``$$`10^4`$ m2s is the diffusion constant and $`\tau _E`$$``$$`10^{10}`$ s is the inelastic scattering time.Skocpol ; Tidecks ; Graybeal2  If constrictions are present in the wire, $`\mathrm{\Lambda }_Q`$ is further reduced and equals the distance between the constrictions. Our short homogeneous wires are not longer than $`\mathrm{\Lambda }_Q`$, so each wire acts as a single WLC. These wires therefore show a clear SIT. Inhomogeneous wires have more than one WLC (i.e. will develop more than one phase slip center under strong current bias) and consequently show a mixed behavior.

In summary, we have found that the morphology of superconducting nanowires has a strong effect on the observed superconductor-insulator transition. A clear SIT is only found in homogeneous wires whereas inhomogeneous samples show a smeared transition and mixed behavior. In the later case, the results are understood by the assumption that inhomogeneous wires are composed of weakly coupled sections connected in series. Each section exhibits either superconducting or insulating behavior while the entire wire shows a mixed behavior.

This work is supported by NSF CAREER Grant No. DMR 01-34770 and by the Alfred P. Sloan Foundation. Some fabrication was performed at CMM-UIUC supported in part by DOE grant DEFG02-96-ER45439.


text:
# Formation of giant globular cluster G1 and the origin of the M31 stellar halo

## 1 Introduction

A growing number of photometric and spectroscopic observations have suggested that G1 (=Mayall II), which is one of the brightest globular clusters belonging to M31, has very unique physical properties as a globular cluster (e.g., Meylan et al. 2001). These include a possible intrinsic metallicity dispersion among its stellar population (Meylan et al. 2001), the large central velocity dispersion of $``$ 25 km s-1 (e.g., Djorgovski et al. 1997), the very flattened shape with mean ellipticity of 0.2 and significantly high central surface brightness (Rich et al. 1996; Meylan et al. 2001), and the possible existence of a black hole with a mass of 20000 $`M_{}`$ (Gebhardt et al. 2002). In spite of its extraordinary nature, G1 is observed to be on the $`M_\mathrm{V}`$-$`\sigma _0`$ relation (where $`M_\mathrm{V}`$ and $`\sigma _0`$ are total magnitude in $`V`$band and central velocity dispersion) defined not by elliptical and dwarf galaxies but by globular clusters, which implies that G1 looks like a genuine globular cluster (Meylan et al. 2001).

One of the possible scenarios of G1 formation is that G1 is the surviving nucleus of an ancient nucleated dwarf galaxy with its outer stellar envelope almost entirely stripped by M31’s strong tidal field (Meylan et al. 1997, 2000, 2001). Such a scenario has already been suggested by Zinnecker et al. (1988) and Freeman (1993), and the viability of the scenario has been extensively discussed by many authors in terms of $`\omega `$ Cen formation (e.g., Hilker & Richtler 2000; Dinescu 2002; Gnedin et al. 2002; Zhao 2002; Bekki & Freeman 2003; Mizutani et al. 2003). However, no theoretical attempts have been made so far to investigate (1) whether M31’s tidal filed is strong enough to transform a dE,N into G1 and (2) what observable evidence of the past destruction of G1’s progenitor dE,N we can find in the M31 halo regions. The above point (2) is very important, because Ibata et al. (2001) and Ferguson et al. (2002) have recently discovered M31’s stellar halo substructures, which could have been formed from tidal destruction of M31’s satellite dwarfs.

In this paper, by using numerical simulations, we first demonstrate that G1 can be formed from a dE,N during tidal interaction between the dE,N and M31. Our fully self-consistent numerical simulations demonstrate that the stellar envelope of dE,N with $`M_\mathrm{B}`$ $``$ $`15`$ can be nearly completely stripped by the strong tidal field of M31 whereas the central nucleus can remain intact owing to its compactness. We suggest that this naked nucleus orbiting M31 is a giant globular cluster (i.e., G1). The morphological transformation from dE,Ns into very compact stellar systems was originally investigated by Bekki et al. (2001, 2003) for ultra-compact dwarfs labeled as “UCD” (Drinkwater et al. 2003) and called “galaxy threshing” (Bekki et al. 2001), though they suggested that galaxy threshing is also important for the formation of giant globular clusters such as $`\omega `$ Cen and G1. The present study confirms this early suggestion and discusses the physical relationship between UCDs, $`\omega `$ Cen, and G1.

## 2 Model

Since our numerical methods and techniques for modeling collisionless, self-gravitating systems of dE,Ns already have been described in detail by Bekki et al. (2003), we give only a brief review here. We consider that a dE,N with a mass and size similar to that observed for the dE,N types orbiting M31. The dE,N is modeled as a fully self-gravitating system and is assumed to consist of a dark matter halo, a stellar component and a nucleus: We first investigate this “three-component” dE,N model. For convenience, the stellar component (i.e., the main baryonic component) is referred to as either the “envelope” or the “stellar envelope” so that we can distinguish this component from the stellar nucleus. The density profile of the dark matter halo with the total mass of $`M_{\mathrm{dm}}`$ in the dE,N is represented by that proposed by Salucci & Burkert (2000):

$$\rho _{\mathrm{dm}}(r)=\frac{\rho _{dm,0}}{(r+a_{\mathrm{dm}})(r^2+a_{\mathrm{dm}})^2},$$
(1)
where $`\rho _{dm,0}`$ and $`a_{\mathrm{dm}}`$ are the central dark matter density and the core (scale) radius, respectively. This model profile is consistent with observations and different from the predictions of the standard CDM model (Navarro, Frenk, & White 1996, hereafter NFW). The dark matter core parameters, $`\rho _{dm,0}`$, $`a_{\mathrm{dm}}`$, and $`M_0`$ (where $`M_0`$ is the total dark matter mass within $`a_{\mathrm{dm}}`$) have a clear observed correlation, $`M_0=4.3\times 10^7(\frac{a_{\mathrm{dm}}}{\mathrm{kpc}})^{7/3}M_{}`$ (Burkert 1995).

The mass (luminosity) and the scale length of the stellar envelope of the dE,N is modeled according to the observed scaling relation of Ferguson & Binggeli (1994):

$$\mathrm{log}a_{\mathrm{dw}}[\mathrm{pc}]=0.02M_\mathrm{B}+2.6$$
(2)
for faint dwarfs ($`M_\mathrm{B}`$ $``$ $`16`$), where $`a_{\mathrm{dw}}`$ and $`M_\mathrm{B}`$ are the scale length of the exponential profile and the absolute $`B`$band magnitude, respectively. The projected density of the envelope with $`M_\mathrm{B}`$ and the total mass of $`M_{\mathrm{dw}}`$ (and $`M_{\mathrm{dw}}/L_\mathrm{B}`$ = 2) is represented by an exponential profile with a scale length $`a_{\mathrm{dw}}`$. The projected density profile of the nucleus with mass $`M_\mathrm{n}`$ is represented by a King model (King 1964) with a core radius of $`a_\mathrm{n}`$ and a central concentration parameter $`c`$ of 1.0.

The nuclei typically contribute about a few percent of the total light of dwarfs (Binggeli & Cameron 1991; Freeman 1993) and the present-day mass of G1 is estimated to be 1.5 $`\times `$ $`10^7`$ $`M_{}`$ for the King model (Meylan et al. 2001). Considering these observations, the reasonable $`M_\mathrm{B}`$ of dE,N is estimated to be $``$ $`15`$ mag. Given a value of $`M_{\mathrm{dm}}`$/$`M_{\mathrm{dw}}`$, we can determine $`a_{\mathrm{dm}}`$ from the above $`M_0`$-$`a_{\mathrm{dm}}`$ relation (Burkert 1995). For convenience, $`a_{\mathrm{dm}}`$ for $`M_\mathrm{B}`$ = $`15`$ mag and $`M_{\mathrm{dm}}`$/$`M_{\mathrm{dw}}`$ = 5 is referred to as $`a_{\mathrm{dm},0}`$ hereafter. We mainly investigate the dE,N models with $`M_\mathrm{B}`$ = $`15`$ mag, $`a_{\mathrm{dw}}`$ = 790 pc, $`M_{\mathrm{dm}}`$/$`M_{\mathrm{dw}}`$ = 5, $`a_{\mathrm{dm}}`$ = $`a_{\mathrm{dm},0}`$, $`M_\mathrm{n}`$/$`M_{\mathrm{dw}}`$ = 0.05, and $`a_\mathrm{n}`$/$`a_{\mathrm{dw}}`$ = 0.02. All of these values are reasonably consistent with observations (Binggeli & Cameron 1991; Ferguson & Binggeli 1994). We also investigate the models with $`a_{\mathrm{dm}}`$ = $`0.25a_{\mathrm{dm},0}`$ to clarify the importance of the dark matter halo structure of a dE,N in the formation processes of G1.

M31 is assumed to have the disk mass of 7.8 $`\times `$ $`10^{10}`$ $`M_{}`$, and the bulge-to-disk-ratio of 0.25, and the maximum rotational velocity of 260 km s-1, all of which are consistent with observations (e.g., van den Bergh 2000). The initial disk plane of M31 is set to be the $`x`$-$`y`$ plane of a simulation. The dE,N orbiting M31 has the initial position of ($`x`$, $`y`$, $`z`$) = ($`R_{\mathrm{apo}}\mathrm{cos}\theta `$, 0, $`R_{\mathrm{apo}}\mathrm{sin}\theta `$) and the initial velocity of ($`v_\mathrm{x}`$, $`v_\mathrm{y}`$, $`v_\mathrm{z}`$) = (0, $`\alpha V_\mathrm{c}`$, 0), where $`R_{\mathrm{apo}}`$, $`\theta `$, and $`V_\mathrm{c}`$ are the apocenter of the orbit, the inclination angle with respect to the M31’s disk (i.e., the $`x`$-$`y`$ plane), the circular velocity at the apocenter, and the parameter (0 $``$ $`\alpha `$ $``$ 1) that determines the orbital eccentricity represented by $`e_\mathrm{p}`$ (i.e, the larger $`\alpha `$ is, the more circular the orbit is). G1 is observed to have the projected distance of 40 kpc from M31 and radial velocity of $`31`$ km s-1 with respect to M31 (e.g., Meylan et al. 2001). Guided by these observations, we investigate the models with $`\theta `$ = 30 and 80, $`R_{\mathrm{apo}}`$ = 40, 80, and 160 kpc, and $`e_\mathrm{p}`$ = 0.62 ($`\alpha `$ = 0.5) and 0.18 ($`\alpha `$ = 0.9). We however describe the four representative and important models in this paper and the parameter values in these models (e.g., “fiducial model”) are given in Table 1. All the simulations have been carried out on a GRAPE board (Sugimoto et al. 1990) with the particle number of 90000.

## 3 Results

Figures 1 summarizes the morphological evolution of the envelope and the nucleus of the dE,N in the fiducial model, which shows typical behavior in morphological transformation from a dE,N into G1. As the dE,N approaches the pericenter of its orbit, the strong global tidal field of M31 stretches the envelope of the dE’N along the direction of the dwarf’s orbit and consequently tidally strips the stars of the envelope ($`T=0.6`$ Gyr). The dark matter halo, which is more widely distributed than the envelope due to its larger core radius, is also efficiently removed from the dE,N during the pericenter passage. Since the envelope (and the dark matter halo) loses a significant fraction of its mass during the passage through the pericenter, the envelope becomes more susceptible to the tidal effects of M31 after the pericenter passage. Therefore, each subsequent time the dwarf approaches the pericenter, it loses an increasingly larger fraction of its stellar envelope through tidal stripping. Finally, the envelope and the dark matter halo lose 98 % and 99 % of their initial masses, respectively, after four passages through the pericenter ($`T=4.8`$ Gyr).

The central nucleus, on the other hand, is only weakly influenced by the tidal force as a result of its compact configuration. Because of its strongly self-gravitating nature, the nucleus loses only a small amount ($`7`$%) of its mass and thus maintains its compact morphology during its tidal interaction with M31. As a result, a very compact stellar system with a negligible amount of dark matter is formed from the dE,N by $`T=4.8`$ Gyr. The total nuclear stellar mass of the remnant within $`5a_{\mathrm{dw}}`$ is $``$ 1.5 $`\times `$ $`10^7`$ $`M_{}`$, consistent with the observed mass of G1. The mass-to-light ratio, $`M/L_\mathrm{B}`$, decreases dramatically from $``$ 10 to $``$ 3 for $`r<5a_{\mathrm{dw}}`$ within 4.8 Gyr. This result clearly explains why G1 is observed to have mass-to-light ratios that are much smaller than those observed for dE,Ns ($``$ 10): galaxy threshing is most efficient in the outer regions of a dE,N where the dark matter halo dominates gravitationally. As is shown in Figure 2, both the surface density of the dark matter and that of the envelope drop by more than an order of magnitude within 4.8 Gyr. If the initial central surface brightness of the dE,N is $`\mu _B=23`$ mag arcsec-2, then the final surface brightness of the dE,N at $`r=3`$ kpc is about $`\mu _B=29.5`$ mag arcsec-2. Such a faint, low surface brightness envelope will be hard to detect, even by existing large ground-based telescopes.

The stripped stars from the envelope of the G1 progenitor dE,N form substructures and tidal tails in M31 halo region. Figure 3 shows that the tidal tail is a “rosette”, reflecting the dE,N eccentric orbit ($`e_\mathrm{p}`$ = 0.62). Most of the stellar components of the tail can be regarded as being located well outside the M31 inner halo with $`r`$ $`<`$ 40 kpc ($``$ three degrees). The $`B`$band surface brightness ($`\mu _\mathrm{B}`$) of the tail ranges from 28.7 mag arcsec-2 to 30.4 mag arcsec-2 with a mean of 30.1 mag arcsec-2. These results imply that the tidal streams formed from the G1 progenitor dE,N would be hard to detect in previous observations by Ibata et al. (2001) and Ferguson et al (2002) that mapped the M31 halo region within $`r`$ $``$ 40 kpc.

Morphological properties of tails and substructures developed during tidal destruction of a G1 progenitor dE,N depend not only on the initial orbit but also on the dark matter structure of the dE,N. Figure 4 summarizes the following three dependences. Firstly, for the model with a smaller $`e_\mathrm{p}`$ (= 0.18), in which the time scale of morphological transformation from the dE,N into G1 is much longer than that of the fiducial model, the shape of the developed tail is more like a spiral than a rosette. Secondly, the distribution of stars in the developed substructures for the model with a smaller $`R_{\mathrm{apo}}`$ (= 40 kpc) appears to have a sharp boundary within which the stars can be located in a relatively homogeneous manner. Thirdly, the surface density of stars along the tail in the model with a more compact dark matter halo (with $`a_{\mathrm{dm}}`$ = 0.25 $`a_{\mathrm{dm},0}`$) is much lower than that of the fiducial model owing to the smaller number of stripped stars. Formation of G1 from the dE,N does not occur in this model with $`a_{\mathrm{dm}}`$ = 0.25 $`a_{\mathrm{dm},0}`$ which suggests that a dE,N with the high central density of the dark matter (such as seen in NFW models) are less likely to be a G1 progenitor dE,N (See Bekki et al. 2003 for more discussions on this point).

## 4 Discussion and conclusion

If G1 originates from a dE,N with $`M_\mathrm{B}`$ = $`15`$ mag, the mean metallicity of the M31 stellar halo consisting of stars stripped from the dE,N is roughly estimated as \[Fe/H\] = $`0.96`$ (for $`BV`$ = 0.71 in the dE,N) using the observed metallicity-luminosity relation for dwarf galaxies, \[Fe/H\]=$`3.43(\pm 0.14)0.157(\pm 0.012)\times M_V`$ (e.g., Côté et al. 2000). Because of its small mass (less than $`10^9`$ $`M_{}`$) and low metallicity, the stellar halo formed from the tidal destruction of the G1 progenitor dE,N could not be the major component of the observed high density and metal-rich (\[m/H\] $``$ $`0.5`$) M31’ stellar halo (Durrell et al. 2000). Also, relatively metal-rich components (\[Fe/H\] $``$ $`0.7`$) of the stellar tail recently discovered in the M31 halo region (Ibata et al. 2001) might be less likely to be formed from the metal-poor stars stripped from G1’s progenitor dE,N.

Recently Ferguson et al. (2002) discovered a stellar halo substructure located in the proximity of G1 (referred to as “G1 clump”) and showed that total $`V`$band magnitude and $`V`$band surface brightness of the G1 clump can be estimated to be $`12.6`$ mag and 28.5 mag arcsec-2, respectively, for reasonable assumptions of dust extinction. The present study has demonstrated that metal-poor (\[Fe/H\] $``$ $`1`$) and low surface brightness ($`\mu _B30`$ mag arcsec-2) stellar halo substructures can be formed along the orbit of the G1 progenitor. We have also found that the stars in substructures close to the simulated G1 have orbital properties similar to those of the G1. Accordingly, future spectroscopic observations on radial velocities and metallicities of stars in the G1 clump will provide a new clue to the problem of whether the G1 clump can be tidal debris of the G1 progenitor dE,N.

Our simulations suggest that if a dE,N has a higher central dark matter density (i.e., more compact core), it cannot be transformed into a giant globular cluster (G1) because of the survival of its stellar envelope. Whether a dE,N can be transformed into a globular cluster by the M31 tidal field (i.e., by galaxy threshing) depends on whether the M31 tidal force is stronger than the self-gravitational force of the dark matter halo of the dE,N at the pericenter of the orbit of the dE,N. Therefore we can give the upper limit of the central dark matter density of the G1 progenitor dE,N, if we know its pericenter distance from M31. The possible central density of the dark matter halo of $`\omega `$ Cen’s progenitor dE,N is discussed in Bekki & Chiba (2003), based on the proper motion data of $`\omega `$ Cen.

Drinkwater et al. (2003) have recently discovered a new type of galactic object with effective radius of $``$ 20 pc, $`M_V`$ $``$ $`12`$ mag, and velocity dispersion of $``$ 30 km s-1 in the Fornax Cluster. These ultra-compact dwarfs (UCDs) with $`M_V`$ more than 1 mag brighter than that of G1 have been demonstrated to be formed by galaxy threshing in which dE,Ns with $`M_V`$ $`<`$ $`16`$ mag can be transformed into UCDs owing to tidal stripping of the stellar envelopes of dE,Ns (Bekki et al. 2001, 2003). We thus suggest that UCDs and giant globular clusters such as G1 and $`\omega `$ Cen can be regarded as the same class of stellar objects: The total mass or luminosity of a progenitor dE,N is the main difference between UCDs and G1. However, the difference in the location on the $`M_\mathrm{V}`$-$`\sigma _0`$ relation between UCD and G1 cannot be explained simply by the galaxy threshing scenario (Bekki et al. 2003).

The present study suggests that if G1 is not close to the pericenter of its orbit around M31, the stripped stars from the G1 progenitor dE,N can be distributed throughout the M31 outer halo region with $`R_\mathrm{p}`$ $`>`$ 40 kpc. Therefore, future observations on M31’s stellar halo need to extend to the current limit of previous surveys ($`R_\mathrm{p}`$ $``$ 40 kpc; Ferguson et al. 2002) to reveal possible tidal streams and substructures formed from dE,N. Metallicity information of the possible streams and substructures in M31’s outer halo is also useful to constrain the stellar population of the destroyed G1 progenitor dE,N. Thus, future deep, high-resolution, wide-area surveys of M31’s outer ($`R_\mathrm{p}`$ $`>`$ 40 kpc) stellar halo by wide-field cameras on large ground-based telescopes (e.g., Suprime-Cam on Subaru) will enable us to determine whether G1 originates from an ancient dE,N orbiting M31.

###### Acknowledgements.

K.B. acknowledges the Large Australian Research Council (ARC). All the simulations described here were performed with GRAPE 5 systems at the National Astronomical Observatory in Japan.