INTERNAL DLA AND THE STEFAN PROBLEM

Janko Gravner† and Jeremy Quastel‡

Abstract. The generalized internal diffusion limited aggregation (internal DLA) is a stochastic growth model on a lattice in which the dynamics has two components. Firstly, particles execute a non–degenerate zero–range dynamics (for example, independent random walks) on sites which constitute the occupied set, that is, have occupancy number exceeding α; furthermore, particles are stuck on sites with occupancy number at most α. Secondly, a finite number of sites act as sources of particles. We prove that the hydrodynamic limit of generalized internal DLA is the one–phase Stefan problem and use this to study the asymptotic behavior of the occupied set.

0. Introduction.

Internal diffusion limited aggregation is a stochastic growth model in which a set in the integer lattice grows by addition of boundary sites hit by random walks produced within the set. The model also goes by the names anti-DLA and diffusion limited erosion. It has been used in physics as a model of erosion processes, but its introduction into the mathematical literature was from the unlikely direction of an algebraic construction. Diaconis and Fulton [DF] consider a vector space whose basis elements are finite subsets of the lattice, and introduce a smash product making it a commutative algebra. The product {x} ∗ A of a singleton {x} with a set A containing it is given by ∑_{y} p(x, y)(A ∪ {y}) where p(x, y) is the probability that y is the first site in the complement of A visited by a symmetric nearest neighbor random walk starting at x. For general finite sets A and B the smash product is A ∗ B = {x_{n}} ∗ {x_{n−1}} ∗ ⋯ ∗ {x_{1}} ∗ (A ∪ B) where

1991 Mathematics Subject Classification. Primary 60K35; secondary 82A05. Key words and phrases. Interacting particle system, free boundary problem, growth model, shape

theory.

†Partially supported by U.S. National Science Foundation grant no. DMS-9703923 and Republic of Slovenia’s Ministry of Science grant J1-8542-0101-97. ‡Partially supported by U.S. National Science Foundation grant no. DMS-9504791 and Alfred P. Sloan Foundation Fellowship.

Typeset by AMS-TEX

A∩B = {x_{1}, …, x_{n}}, the key observation being that this is independent of the ordering {x_{1}, …, x_{n}}. It can now be extended to general vectors by linearity and gives a well defined product. The construction admits many generalizations, but already in this simple case there is the very interesting question of the behaviour of the random set A_{n} obtained from A_{0} = {0} by successively adding on the first point y in the complement visited by a symmetric nearest neighbor random walk starting at the origin.

In [BGL] Bramson, Griffeath and Lawler studied the asymptotic shape of A_{n} and found that in any dimension d, A_{n} is asymptotically a ball B(0, r) of radius r = a_{d}n^{1/d} around the origin where a_{d} = (dΓ(d/2)/2)^{1/d}/\sqrt{π} so that B(0, r) has volume n. More precisely, they proved that for every δ > 0, with probability one, for sufficiently large n,

B(0, (1 − δ)a_{d}n^{1/d}) ⊂ A_{n} ⊂ B(0, (1 + δ)a_{d}n^{1/d}).

In another version of the model, particles are produced at the origin as a Poisson process at rate 1 and perform continuous time random walks jumping at rate 1 to each nearest neighbour, until they hit the boundary, at which time the hitting site is added to the set. Now one has a set A_{t} and a configuration of particles on it together performing a continuous time Markov process. In dimensions d ≥ 3, [BGL] showed that the set has an asymptotic shape which is still a ball of radius a_{d}t^{1/d} around the origin. Their proof relies on the fact that t^{1/d} grows slowly enough in d ≥ 3 that particles essentially hit the boundary before the next particle is produced and therefore the model can be reduced to the discrete time case.

In dimensions 1 and 2, however, the methods of [BGL] break down and in fact computer simulations show a nontrivial density of particles in the occupied region. They conjectured that in two dimensions the set grows as K\sqrt{t} with a constant K which could, in principle, be computed with the following heuristic argument. Assume that the occupied set A_{t} is asymptotically a ball B(0, K\sqrt{t}) around the origin. The expected number of particles that hit the boundary before time T is then, after Brownian scaling,

T ∫_{0}^{1} P_{0}( sup_{0≤s≤1−t}(t + s)^{−1/2}|\sqrt{2}β_{s}| ≥ K)dt  (0.1)

where β_{s} is a standard Brownian motion starting at the origin (note that our random walks jump at rate one to each nearest neighbour and therefore the diffusion approximation is \sqrt{2}β_{s}). On the other hand the number of particles which hit by time T is clearly

T πK^{2}.  (0.2)

[BGL] conjectured that the correct growth rate is the unique K for which these expressions coincide.

In this article we prove the conjectures of [BGL] by showing that internal DLA has a hydrodynamic limit which is the one-phase Stefan problem. This is then used to study the asymptotic shape of the occupied set. For example, in the two dimensional case described above the occupied set grows as K\sqrt{t} where

e^{−K^{2}/4} = πK^{2}

which one can check is the unique K for which the (0.1) and (0.2) coincide by solving the heat equation in two dimensions with a source at the origin and Dirichlet boundary conditions on an expanding circle of size \sqrt{t}.

Our approach yields somewhat weaker shape theorems than the methods of [BGL] but has the advantage of being rather robust: It works in all dimensions and extends readily to interacting versions of the model. In particular we consider particles produced at rates of order ε^{−d} at various sites in ℝ^{d}, then performing interacting random walks on the reduced lattice εℤ^{d} at rates of order ε^{−2}. The interaction is zero-range: The jump rate of a particle depends on the number of particles at that site. We also choose some positive integer α and decree that the first α particles at a site are frozen. The set of frozen sites is the occupied set A_{ϵ}(t) and we show that it is asymptotically a ball of radius of order

\{ \begin{matrix}t^{1/d} & d ≥ 2; \\ (t log t)^{1/2} & d = 1. \end{matrix}

We comment briefly on the relation between our scaling and that in [BGL]. If we consider time t = 1 in our model and some small ε > 0, then we will have created O(ε^{−d}) particles which have each taken O(ε^{−2}) steps of size ε producing an occupied set which is a ball of order 1. This corresponds to taking times of order ε^{−2} in [BGL], at which time order ε^{−2} particles will have been created and the occupied set will be a ball of radius O(1) only if we rescale the lattice width to be O(ε^{2/d}). In d ≥ 3 therefore, the model in [BGL] is, in some sense, a renormalized version of our model, with each of their particles corresponding to ε^{2−d} of our particles. The surprise is that the size of the occupied set is the same in both models, which means that the dynamics is essentially unaffected by this renormalization.

The convergence to the Stefan problem is proved using a fairly standard method in hydrodynamics, the H_{−1} method, with technical modifications to deal with the lousy ergodicity properties of the model, and the singular behaviour of the Stefan problem at creation points. The H_{−1} method was developed in [Y] and [CY] to study the metastability and nonequilibrium properties of conservative Ginzburg-Landau models. It is probably the easiest method in hydrodynamic scaling limits; however, it is restricted to gradient systems. In particular, it completely avoids the use of entropy arguments and the two block estimate, which in the internal DLA models seem to be hopeless. In fact the only type of ergodicity result needed as input is the characterization of the translation invariant, invariant measures for the process in infinite volume which is proved using coupling methods of [A]. This unfortunately forces us to assume the attractiveness (monontonicity) of the dynamics.

1.Notatation and results.

We now describe the generalized internal DLA model. Particles are moving on the lattice εℤ^{d}. We let η_{x} denote the number of particles at x ∈ εℤ^{d}. The rate of jumping of particles from x to nearest neighbor site y is ε^{−2}a(η_{x}) where a(n) = g(n − α) for n > α and a(n) = 0 for n ≤ α corresponding to the freezing of the first α particles. g(n) are the rates of an attractive zero range process and are assumed to satisfy for some C < ∞

C^{−1} ≤ g(n + 1) − g(n) ≤ C

for all n. We will also assume that g(n) is approximately linear, ie. |g(n) − Kn| ≤ C for all n = 0, 1, … for some finite C and K. Of course by rescaling time we can and will assume K = 1 so the assumption is that for all n = 0, 1, …, for some finite C,

|g(n) − n| ≤ C.

Particles are created at a finite number of sites x_{i} = ⌊x_{i}⌋ at rates c_{i}. i = 1, …, n. We use x to denote locations on ℝ^{d}, x to denote locations on εℤ^{d} and ⌊x⌋ to denote the lattice part of x. The generator of the process is

L = ε^{−2}L_{0} + ε^{−d}L_{c}

where

L_{0}f(η) = ∑_{x,e} a(η_{x}) (f(η^{x,x+e}) − f(η))  (1.1)

where e denote the basis vectors of length ε in our lattice and η_{z}^{x,y} = η_{y} + 1 if z = y, η_{x} − 1 if z = x and η_{z} otherwise, and

L_{c}f(η) = ∑_{i=1}^{n} c_{i} (f(η^{x_{i}}) − f(η))

where η_{z}^{x} = η_{z} + 1 if z = x and η_{z}^{x} = η_{z} otherwise. The corresponding zero range process, with generator

L_{1}f(η) = ∑_{x,e} g(η_{x}) (f(η^{x,x+e}) − f(η))

has as invariant measures the product measures with marginals

Z^{−1}(γ)γ^{n}/[g(n) ⋯ g(1)].

The factor Z(γ) is the normalization to make it into a probability measure. The parameter γ is related to the density ρ by the formula

ρ = γZ′(γ)/Z(γ).

Define

λ(u) = \{ \begin{matrix}γ(u − α) & if u ≥ α; \\ 0 & if 0 ≤ u < α. \end{matrix}

Note that from the assumptions on g(n) we have for some finite C,

C^{−1}|u_{1} − u_{2}| ≤ |λ(u_{1}) − λ(u_{2})| ≤ C|u_{1} − u_{2}|,  u_{1}, u_{2} ∈ [α, ∞).  (1.2)

The upper bound actually holds trivially for all u_{1}, u_{2}. Note also that the assumptions imply that |γ(ρ) − ρ| ≤ C for all ρ and therefore

|λ(u) − u| ≤ C  (1.3)

for all u.

The Stefan problem corresponding to our microscopic system is to find functions ρ(x, t) and s(x) on ℝ^{d} satisfying

\{ \begin{matrix}\frac{∂ρ}{∂t} = Δγ(ρ) + ∑_{i=1}^{n} c_{i}δ_{x_{i}} & in the region \{ρ > 0\}; \\ ρ = 0 & on s(x) ≥ t; \\ ∇_{0}ρ · ∇s = −α/γ′(0) & on s(x) = t, \end{matrix}  (1.4)

where ∇_{0}ρ denotes the gradient taken from the inside of the region s(x) < t. Note the last equation can be rewritten ∇_{0}γ(ρ) · ∇s = −α. The set s(x) ≤ t is the occupied set at time t and the density ρ(x, t) evolves according to the nonlinear heat equation with Dirichlet boundary conditions on the occupied set, whose boundary in addition moves in the outward normal direction at a rate proportional to ∇γ(ρ).

We make the following entropy transformation to obtain a weak formulation of the problem. Let ρ be a classical solution. Let φ be a smooth test function with compact support in [0, T) × ℝ^{d}. Multiplying by φ and using Green’s identity,

∫_{0}^{T} ∫ [ρ\frac{∂φ}{∂t} + γ(ρ)Δφ] dxdt+∫_{0}^{T} ∫_{s(x)=t} φ∇γ(ρ)· \frac{∇s}{|∇s|} dSdt+∫_{0}^{T} ∑_{i=1}^{n} c_{i}φ(t, x_{i})dt = 0

Using ∇γ(ρ) · ∇s = −α on s(x) = t the middle term gives

−α ∫_{0}^{T} ∫_{s(x)=t} φ\frac{dS}{|∇s|}dt = −α ∫ φ(s(x), x)dx = ∫_{0}^{T} ∫ α1_{ρ>0} \frac{∂φ}{∂t}dxdt.

Let

σ(ρ) = \{ \begin{matrix}ρ + α & if ρ > 0; \\ 0 & if ρ = 0. \end{matrix}

We conclude that

∫_{0}^{T} ∫ [γ(ρ)Δφ + σ(ρ)\frac{∂φ}{∂t}] dxdt + ∫_{0}^{T} ∑_{i=1}^{n} c_{i}φ(t, x_{i})dt = 0,  (1.5)

for all smooth test functions φ with compact support in [0, T) × ℝ^{d}. Now let ρ satisfy (1.5) and define

u = σ(ρ).

Note that γ(ρ) = λ(u). Thus u is a weak solution of

\frac{∂u}{∂t} = Δλ(u) + ∑_{i=1}^{n} c_{i}δ_{x_{i}}  (1.6)

in the sense that for each smooth test function φ with compact support in [0, T) × ℝ^{d},

∫_{0}^{T} ∫ [λ(u)Δφ + u\frac{∂φ}{∂t}] dxdt + ∫_{0}^{T} ∑_{i=1}^{n} c_{i}φ(t, x_{i})dt = 0.  (1.7)

We will call such a ρ a weak solution of the Stefan problem. We will prove below an existence and uniqueness theorem for such weak solutions. Note that even in the simple case of two creation sites, at some time the clusters will meet and the weak solution will fail to satisfy (1.4) in a classical sense. We will not pursue the issue here of in what sense weak solutions satisfy (1.4) but refer the reader to the extensive literature (see, for example, [M]).

The main results of this paper are:

Theorem 1.1. For each t > 0, as ε → 0, the configuration η^{ε}(t) of the generalized internal DLA process converges weakly in probability to the unique weak solution u(t) of the Stefan problem (1.4).

To state the shape theorems, let us introduce some notation. Let A_{t}^{ε} be the occupied set in the generalized internal DLA process with creation at rate ε^{−d} at the origin. Note that A_{t}^{1} is in the occupied set when creation and diffusion occur at the same rate, hence it corresponds to the occupied set for the model in [BGL]. In all dimensions we denote by B(x, r) the Euclidean ball of radius r centered at x.

Definition. The total variation distance dist(A, A) between a set A ⊂ εℤ^{d} and a Borel set A ⊂ ℝ^{d} is defined as follows: Let μ_{A} be the the measure on ℝ^{d}, given by

μ_{A}(·) = ε^{d}[# of points in A ∩ ·],

where A ⊂ ℝ^{d} by the natural embedding with every point having mass ε^{d}. Also, let μ_{A} be the measure on ℝ^{d} with density 1_{A} dx. Then

dist(A, A) = ∥μ_{A} − μ_{A}∥_{T V} ,

where ∥ · ∥_{T V} is the total variation norm.

Theorem 1.2. 1. In any dimension d and for each t > 0,

dist(A_{t}^{ε}, A_{t}) → 0

as ε → 0 in probability where A_{t} = {x : s(x) ≤ t} is the occupied set in the Stefan problem (1.4).

2. Assume d = 1. There exist constants K_{1} and K_{2} such that for sufficiently large t,

B(0, K_{1}\sqrt{2t log t}) ⊂ A_{t} ⊂ B(0, K_{2}\sqrt{2t log t}).  (1.8)

If g(n) = n (the case of independent random walks), then, for arbitrary δ > 0 and K_{1} = 1 − δ and K_{2} = 1 + δ, (1.8) holds for sufficiently large t and it continues to hold with probability one if A_{t} is replaced with A_{t}^{1}.

3. Assume d = 2. There exists a constant K such that

A_{t} = B(0, K\sqrt{t})

and

dist(A_{t}^{1}/\sqrt{t}, B(0, K)) → 0

in probability, as t → ∞. In the case g(n) = n, K is the solution of

e^{−K^{2}/4} = παK^{2}.

4.Assume d ≥ 3. There exist constants K_{1} and K_{2} such that for sufficiently large t,

B(0, K_{1}t^{1/d}) ⊂ A_{t} ⊂ B(0, K_{2}t^{1/d}).  (1.9)

If g(n) = n, then, for arbitrary δ > 0 and K_{1} = (1 − δ)a_{d} and K_{2} = (1 + δ)a_{d}, (1.9) holds for sufficiently large t and it continues to hold with probability one if A_{t} is replaced with A_{t}^{1}. Here, a_{d} = (dΓ(d/2)/2)^{1/d}/\sqrt{π}.

When g(n) = n, (1.9) for A_{t}^{1} was proved already in [BGL] (the random walks in this paper have different rates, but the argument applies unchanged). In the case of zero–range dynamics, strengthening (1.8) or (1.9) to make K_{1} and K_{2} arbitrarily close remains an open problem, and so does proving (1.8) or (1.9) for A_{t}^{1} instead of A_{t}.

2. Stefan problem.

In this section we prove some preliminary results about the Stefan problem which we will need. Although these results are fairly standard, we were not able to find any references in the vast literature on the Stefan problem (see references in [M]) dealing with seemingly natural situation of delta function sources.

The idea behind the uniqueness result we will prove is that if u and v are two solutions of (1.6) the H_{−1} norm contracts at rate ∫ (u − v)(λ(u) − λ(v)) which is positive because λ is nondecreasing. Note however, that in dimensions d ≥ 4, the solutions to our Stefan problem are not even L^{2} functions due to the Green’s function like singularities at the creation points. So we have to identify and subtract the singularities.

In the proceeding it may help to keep in mind the asymptotic behaviour of the Stefan problem with delta function sources: As will be descibed in detail in Section 4, in dimensions 1 and 2 the solution is always in L^{2} and grows, at rate t^{1/2} in d = 1 and logarithmically in d = 2. In dimensions 3 and higher one can check easily in the linear case γ(u) = u that the solution grows to the solution of Δλ(u_{∞}) = − ∑_{i=1}^{n} c_{i}δ_{x_{i}} and in the nonlinear case one expects the same behaviour.

Let P_{x}^{ε}(t) be the number of particles created at x up to time t in the internal DLA process. Of course P_{x}^{ε}(t) is just a Poisson process with rate ε^{−d}c_{i} if x = x_{i} and vanishes otherwise.

Let Δ_{ε} be the lattice Laplacian on εℤ^{d},

[Δ_{ε}ϕ]_{x} = ε^{−2} ∑_{|e|=ε}[ϕ_{x+e} − ϕ_{x}]

and consider the lattice approximation to the Stefan problem,

\frac{∂u^{ε}}{∂t} = Δ_{ε}λ(u^{ε}) + dP^{ε}  (2.1)

on εℤ^{d} with u^{ε}(0) = 0. We compare it with the ‘free’ problem

\frac{∂ρ^{ε}}{∂t} = Δ_{ε}ρ^{ε} + dP^{ε}  (2.2)

also on εℤ^{d} with ρ^{ε}(0) = 0. Note that ρ^{ε}(x, t) is given explicitly by ∑_{i=1}^{n} c_{i} ∫_{0}^{t} p_{x−x_{i}}^{ε}(t− s)dP_{x_{i}}^{ε}(s) where p_{x}^{ε}(t) is the heat kernel on εℤ^{d}. An elementary computation shows that for any 1 ≤ p < d/(d − 2) if d ≥ 3, any 1 ≤ p < ∞ if d = 2 and any 1 ≤ p ≤ ∞ if d = 1, and any T ≥ 0 there exists a finite constant C_{p}(T) so that for all ε > 0 and 0 ≤ t ≤ T,

E[ε^{d} ∑_{x}[ρ_{x}^{ε}(t)]^{p}] ≤ C_{p}(T).  (2.3)

It is also clear that for any δ > 0 and any T > 0 there exists a finite constant C(δ, T) such that for all ε > 0 and 0 ≤ t ≤ T,

E[ε^{d} ∑_{|x−x_{i}|≥δε^{−1}}[ρ_{x}^{ε}(t)]^{2}] ≤ C(δ, T).  (2.4)

Note that lim_{δ→0} C(δ, T) < ∞ only in d ≤ 3.

Lemma 2.1. Let u^{ε} be the solution of (2.1) and ρ^{ε} be the solution of (2.2). Then for each T > 0 there exists a finite constant C independent of ε such that,

∫_{0}^{T} ∥u^{ε}(t) − ρ^{ε}(t)∥_{ℓ^{2}(εℤ^{d})}dt ≤ C.  (2.5)

Proof. Let g_{ε} be the kernel of −Δ_{ε}^{−1}, ie.

[−Δ_{ε}^{−1}f]_{x} = ε^{d} ∑_{y} g_{ε}(x − y)f_{y}  (2.6)

and let

∥f∥_{−1}^{2} = ε^{2d} ∑_{x,y} g_{ε}(y − x)f_{y}f_{x}.  (2.7)

We have

∥u^{ε}(t) − ρ^{ε}(t))∥_{−1}^{2} = −2 ∫_{0}^{t} ε^{d} ∑_{x}(u_{x}^{ε} − ρ_{x}^{ε})(λ(u_{x}^{ε}) − ρ_{x}^{ε})ds.

We use the bounds |λ(u) − u| ≤ C and ε^{d} ∑_{x} u_{x}^{ε}, ε^{d} ∑_{x} ρ_{x}^{ε} ≤ C to bound this by

−2 ∫_{0}^{t} ε^{d} ∑_{x}(u_{x}^{ε} − ρ_{x}^{ε})^{2}ds + C.

We conclude that the absolute value of the first term is bounded by a constant independent of ε. □

Remark. In dimensions d ≤ 3, both terms in (2.5) are in ℓ^{2} uniformly in ε so the lemma is easier. This is proved as follows. If g_{x} = ∫_{0}^{t} p_{x}^{ε}(t − s)ds then by Schwarz’s inequality g_{x}^{2} ≤ ∫_{0}^{t}[p_{x}^{ε}(s)]^{2}s^{−c}ds ∫_{0}^{t} s^{c}ds. Since ∥p^{ε}(s)∥_{ℓ^{2}} ≤ Cs^{−d/2} we have that ∥g∥_{ℓ^{2}} is finite if we can choose −1 < c < (2 − d)/2, ie. if d ≤ 3. For equation (2.2) one has the integral against a Poisson process instead of Lebesgue measure, but the same result holds with probability one.

For equation (2.1) one has the added complication of the nonlinearity. However if ρ^{ε} is a function on the lattice we can write Δ_{ε}γ(ρ^{ε}) = ∑_{e} a_{x,x+e}[u_{x+e}^{ε} −u_{x}^{ε}] where a_{x,x+e} = \frac{γ(u_{x+e}^{ε})−γ(u_{x}^{ε})}{u_{x+e}^{ε}−u_{x}^{ε}} are bounded above and below. So \frac{∂ρ^{ε}}{∂t} = Δ_{ε}γ(ρ^{ε}) is the forward equation for a reversible random walk and therefore one knows that the transition probabilities are bounded above and below by Gaussians, by standard parabolic estimates for divergence form diffusion equations [D]. Therefore the above argument shows that the solution of the inhomogeneous equation \frac{∂ρ^{ε}}{∂t} = Δ_{ε}γ(ρ^{ε}) + dP^{ε} is square integrable in dimensions 3 and lower. Finally the same estimate holds for solutions of (2.1) by comparison with this case, because u^{ε} ≤ ρ^{ε} + α.

Define Ψ on ℝ_{+} by

Ψ(u) = \{ \begin{matrix}u − α & if u ≥ α; \\ 0 & otherwise. \end{matrix}

Lemma 2.2. Let u^{ε} be the solution of (2.1). Then for each δ > 0 there exists a constant C(δ) < ∞ depending only on δ and u_{0} such that

E[ ∫_{0}^{T} ε^{d} ∑_{|x−x_{i}|≥δε^{−2}} ∑_{|e|=ε} ε^{−2}{Ψ(u_{x+e}^{ε}) − Ψ(u_{x}^{ε})}^{2}dt] ≤ C(δ).  (2.8)

Remark. The idea of the lemma is that if one has a solution of \frac{∂u}{∂t} = Δλ(u) then the rate of contraction of the L^{2} norm is 2 ∫ λ′(u)|∇u|^{2}dx so the time integral of this term can be controlled. Note that λ′(u) = 0 for 0 ≤ u ≤ α which is the origin of the cutoff Ψ in the lemma.

Proof. Let ϕ(x) be a smooth function vanishing in a neighbourhood of the x_{i}. We have

ε^{d} ∑_{x}[u_{x}^{ε}(t)]^{2}ϕ(x) = ∫_{0}^{t} ε^{d} ∑_{x} u^{ε}Δ_{ε}λ(u^{ε})ϕds.

After a summation by parts the last term becomes

\begin{matrix}− ∫_{0}^{t} ε^{d−2} ∑_{x,e}(u_{x+e}^{ε} − u_{x}^{ε})(λ(u_{x+e}^{ε}) − λ(u_{x}^{ε}))ϕ(x)ds \\ − ∫_{0}^{t} ε^{d−2} ∑_{x,e} u_{x+e}^{ε}(λ(u_{x+e}^{ε}) − λ(u_{x}^{ε}))(ϕ(x + e) − ϕ(x))ds. \end{matrix}

Applying Schwarz’s inequality and using the assumptions on λ we see that for some 0 ≤ c and C < ∞ this is bounded above by

\begin{matrix}− c ∫_{0}^{t} ε^{d−2} ∑_{x,e}(Ψ(u_{x+e}^{ε}) − Ψ(u_{x}^{ε}))^{2}ϕ(x)ds \\ + C ∫_{0}^{t} ε^{d−2} ∑_{x,e}[u_{x+e}^{ε}]^{2}(ϕ(x + e) − ϕ(x))^{2}ϕ^{−1}(x)ds. \end{matrix}

Now we choose a smooth ϕ which is unity on |x − x_{i}| ≥ δ, vanishes on a neighborhood of the x_{i}, and has (ϕ(x + e) − ϕ(x))^{2}ϕ^{−1}(x) uniformly bounded in ε. (Note this can be done since f(r) = r^{2} satisfies |f′|^{2} ≤ Cf.) By the previous lemma u^{ε} have a uniform bound in L^{2} in any region not containing the creation points x_{i}. Hence the last term is bounded independent of ε. We conclude that the first term is also bounded independent of ε which is what was to be proved. □

Now we are able to state the main existence and uniqueness theorem for the Stefan problem (1.6). Let ρ be the solution of

\frac{∂ρ}{∂t} = Δρ + ∑_{i=1}^{n} c_{i}δ_{x_{i}}.  (2.9)

with ρ(0) = 0. Note that (2.9) can be solved explicitly and the solution is ρ(x, t) = ∑_{i=1}^{n} c_{i} ∫_{0}^{t} p(x − x_{i}, t − s)ds where p(x, t) is the heat kernel on ℝ^{d}.

Theorem 2.3. Let u^{ε} be solutions of the lattice version of the Stefan problem (2.1). There exist a subsequence converging weakly to the unique solution u of the Stefan problem (1.6) with u(0) = 0 satisfying

∫_{0}^{T} ∥u(t) − ρ(t)∥_{L^{2}(ℝ^{d})}^{2}dt < ∞  (2.10)

Furthermore the corresponding subsequence of Ψ(u^{ε}) converge strongly to Ψ(u) away from the creation points x_{1}, …, x_{n}.

Proof. From (2.5) we obtain

∫_{0}^{T} ∥u^{ε}(t) − ρ^{ε}(t)∥_{ℓ^{2}(εℤ^{d})}^{2}dt ≤ C.

Since ρ^{ε} converge to ρ weakly (and strongly away from creation points) the bound suffices to let us take a weakly convergent subsequence converging to u which by lower semicontinuity satisfies (2.10). The strong convergence of Ψ(u^{ε}) to Ψ(u) can be deduced from (2.8). This allows us to pass to the limit in the weak formulation of the discrete equation (2.1) to see that u is a weak solution of the Stefan problem (1.6). Now suppose u and v are two such solutions. From (2.10) we obtain

∫_{0}^{T} ∥u(t) − v(t)∥_{L^{2}(ℝ^{d})}^{2}dt < ∞.

By a well known theorem [LM], if V ⊂ H ⊂ V′ are separable Hilbert spaces with V and V′ being in duality relative to the inner product of H and we have a curve u(t) : 0 ≤ t ≤ T satisfying u ∈ L^{2}(0, T; V) and u′ ∈ L^{2}(0, T; V′), then u is almost everywhere equal to a function continuous from [0, T] into H and for 0 ≤ t_{1} ≤ t_{2} ≤ T, ∥u(t_{2})∥_{H}^{2} − ∥u(t_{1})∥_{H}^{2} = 2 ∫_{t_{1}}^{t_{2}}⟨u(t), u′(t)⟩dt. Therefore

∥u(T) − v(T)∥_{H_{−1}(ℝ^{d})}^{2} ≤ −2 ∫_{0}^{T} ∫ (u − v)(λ(u) − λ(v))dxdt.

Since λ is nondecreasing we conclude that ∥u(T) − v(T)∥_{H_{−1}(ℝ^{d})} = 0. This is true for every T > 0 and therefore u ≡ v. □

One of the nice consequences of having the lattice approximation to the Stefan problem is that it allows us to prove very easily comparison principles for solutions even though the creation rates are singular. For example, passing to the limit after an easy computation on the lattice level using the maximum principle we obtain the following result which will be useful in Section 4.

Lemma 2.4. Let u be the weak solution of (1.6)with bounded initial condition u_{0} and v be a subsolution (respectively supersolution) satisfying

\frac{∂v}{∂t} − Δλ(v) ≤ ∑_{i=1}^{n} c_{i}δ_{x_{i}} v(0, x) = v_{0}(x)  (resp. ≥).  (2.11)

If v_{0} ≤ u_{0} then v ≤ u (resp. ≥).

3. Hydrodynamic limit.

Recall that the H_{−1} norm ∥ · ∥_{−1} on εℤ^{d} is given by ∥f∥_{−1}^{2} = ε^{2d} ∑_{x,y} G_{y−x}^{ε}f_{x}f_{y} where G^{ε} is the kernel of −Δ_{ε}^{−1} on εℤ^{d} (Green function). The main result of this section is

Theorem 3.1. Let η^{ε} be the configuration of the generalized internal DLA process and u^{ε} the solution of the discretized Stefan problem (2.1) with the same initial condition. Then

lim_{ε→0} E [∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2}] = 0.

Here and in the rest of the section, E denotes expectation with respect to the generalized internal DLA process η^{ε}(·). The theorem will be proved in Lemmas 3.2-3.9.

Lemma 3.2. Let η^{ε} be the configuration of the generalized internal DLA process and u^{ε} the solution of the discretized Stefan problem (2.1) with the same initial condition. Then

∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2} = 2 ∫_{0}^{T} ε^{d} ∑_{x} V (η_{x}^{ε}, u_{x}^{ε})dt + M^{ε}(T)  (3.1)

where

V (n, u) = −(n − u)(a(n) − λ(u)) + a(n)  (3.2)

and M^{ε}(t) is a martingale.

Proof. This is an explicit computation. We compute

d [(η_{x} − u_{x})(η_{y} − u_{y})] ={ [Δ_{ε}(a(η_{x}) − λ(u_{x}))] (η_{y} − u_{y})

+ (η_{x} − u_{x}) [Δ_{ε}(a(η_{y}) − λ(u_{y}))] − ε^{−2}1(|x − y| = 1)(a(η_{x}) + a(η_{y}))

+ ε^{−2}1(x = y) ∑_{e}(a(η_{x}) + a(η_{x+e}))}dt + dM_{xy}

where M_{xy} are martingales and then use the fact that [−Δ_{ε}G^{ε}]_{0} = ε^{−d}. □

Now consider the generalized IDLA process running in a finite box Λ with reflecting boundary conditions and no creation. The set of invariant measures will be denoted μ_{Λ}^{β} running over the parameter set β ∈ B. They are described as follows: For each N ≥ α|Λ| there is exactly one invariant measure with N particles. Under that measure each η_{x} ≥ α and the configuration η_{x}−α, x ∈ Λ is distributed according to the canonical invariant measure for our reference zero range process with N − α|Λ| particles. On the other hand, if N < α|Λ| then for each configuration of N particles in Λ with no more than α particles at each site, the Dirac mass at that configuration is invariant. So the parameter set B consists really of each fixed configuration with no more than α particles per site, together with each N ≥ α|Λ|.

Lemma 3.3. The marginal on Λ of any translation invariant, invariant measure for L_{0} on ℤ^{d} is a mixture of the μ_{Λ}^{β}.

The proof of this key lemma is postponed until Section 5. We now continue on with the proof of Theorem 3.1.

Lemma 3.4. Let V be as in (3.2) and μ be a translation invariant, invariant measure for the process. Then for any u,

E_{μ}[V (η_{0}, u)] ≤ 0.

Proof. Let Λ be a box around the origin and μ_{Λ} be the marginal of μ on Λ. By the previous lemma μ_{Λ} is a mixture of the extremal invariant measures μ_{Λ}^{β} on Λ. Let f(n) be any function of whole numbers which vanishes if n ≤ α and μ_{Λ}^{β} be an extremal invariant measure on Λ with density v_{β}. Let ν_{γ} be an invariant measure for the zero-range dynamics, ie. a product measure with marginals ν_{γ}(η_{x} = n) = Z^{−1}(γ)γ^{n}[g(n) ⋯ g(1)]^{−1}. We have

E_{μ_{Λ}^{β}}[f(η_{0})] = \{ \begin{matrix}0 & if v_{β} ≤ α; \\ E_{ν_{γ(v_{β} −α)}}[f(η_{0} − α)] + O(|Λ|^{−1}) & if v_{β} > α. \end{matrix}

In the second case the error O(|Λ|^{−1}) comes from changing the measure μ_{Λ}^{β} with fixed density to the product measure. Note that both a(n) and na(n) vanish if n ≤ α. Furthermore if v > α,

\begin{matrix}E_{ν_{γ(v−α)}}[a(η_{0})] = λ(v), \\ E_{ν_{γ(v−α)}}[η_{0}a(η_{0})] = (v + 1)λ(v). \end{matrix}

Letting |Λ| → ∞ we have

E_{μ}[V (η_{0}, u)] = − ∫ (v_{β} − u)(λ(v_{β}) − λ(u))dΥ_{β}

where Υ_{β} is a probability measure on the parameter set for the extremal invariant measures. This is nonpositive since λ in nondecreasing. Note the apparent lack of uniformity in the argument is easily resolved: Except for small densities, the extremal invariant measures on Λ are parametrized by the density itself, and because of the stationarity of μ one has an easy cutoff of large densities uniform in Λ. □

Consider the right hand side of (3.1). If we take the expectation, the martingale disappears and because of the time averaging we will be evaluating V (η_{x}, u) under an invariant measure. By Lemma 3.4 the result is nonpositive, giving the hydrodynamic limit. This argument is made precise in Lemma 3.9. However, there are several cutoffs needed first. Lemma 3.5 provides a basic estimate that will be used repeatedly. Lemma 3.6 cuts off the bad points near the creation sites where the density is unbounded. Lemma 3.7 cuts off large values of V, and Lemma 3.8 provides for some averaging in space and time in order to obtain a translation invariant, invariant measure when we look at the the limit of the distribution around a macroscopic site as we let ε → 0.

Lemma 3.5. 1. Let p < d/(d − 2). For each T > 0 there is a constant C_{p}(T) < ∞ such that for all ε > 0,

E [∫_{0}^{T} ε^{d} ∑_{x}[η_{x}^{ε}(t)]^{p}dt] ≤ C_{p}(T).

2. There exists a function C(δ) < ∞ for each δ > 0 such that for all ε > 0,

E [∫_{0}^{T} ε^{d} ∑_{|x−x_{i}|≥δ}[η_{x}^{ε}(t)]^{2}dt] ≤ C(δ).

Proof. Follows from the similar bound on u and Lemma 3.2 which implies that the expectation of the L^{2} norm of the difference of u and η can be controlled. □

Lemma 3.6. Let V be as in (3.2). Then

E [∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2}] ≤ 2E[ ∫_{0}^{T} ε^{d} ∑_{|x−x_{i}|≥δε^{−1}} V (η_{x}^{ε}, u_{x}^{ε})dt] + Ω_{1}(δ, ε, T) (3.3)

where for each T > 0,

lim sup_{δ→0} lim sup_{ε→0} Ω_{1}(δ, ε, T) = 0.

Proof. Let Ω_{1}(δ, ε, T) = C ∫_{0}^{T} ε^{d} ∑_{{x:|x−x_{i}|<δ}} η_{x}dt where C is large enough that a(n) ≤ Cn. Since we chucked a negative term from (3.1) the inequality holds. We can estimate Ω_{1} by part 1 of Lemma 3.5 and Čebyšev’s inequality. □

Lemma 3.7. Let ϕ_{ℓ}(x) = x if x ≤ ℓ and ϕ_{ℓ}(x) = ℓ if x ≥ ℓ be the cutoff at level ℓ. Let

V_{ℓ}(n, x) = −ϕ_{ℓ}((n − x)(a(n) − λ(x))) + ϕ_{ℓ}(a(n)).

Then

\begin{matrix}E [∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2}] \\ ≤ 2ε^{d} ∑_{|y−x_{i}|≥δ} ∫_{0}^{T} Av_{(x,s)∈B_{a}(y,t)}E [V_{ℓ}(η_{x}^{ε}(s), u_{y}^{ε}(t))] dt + Ω_{2}(ℓ, ε, a, δ, T) \end{matrix}  (3.4)

where

lim sup_{δ→0} lim sup_{a→0} lim sup_{ℓ→∞} lim sup_{ε→0} Ω_{2}(ℓ, ε, a, δ, T) = 0.

Proof. Since one of the error terms from (3.3) is negative all we have to show is that

lim sup_{ℓ→0} lim sup_{ε→0} E [∫_{0}^{T} ε^{d} ∑_{x} |ϕ_{ℓ}(a(η_{x})) − a(η_{x})|dt] = 0

which follows from Lemma 3.5 and Čebyšev’s inequality. □

Lemma 3.8. Let V be as in (3.2). For each σ > 0, let B_{σ}(y, t) = {(x, s) : |x − y| ≤ σ, |s − t| ≤ σ} ⊂ εℤ^{d} × [0, T]. Then

\begin{matrix}E [∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2}] \\ ≤ 2ε^{d} ∑_{|y−x_{i}|≥δ} ∫_{0}^{T} Av_{(x,s)∈B_{σ}(y,t)}E [V_{ℓ}(η_{x}^{ε}(s), u_{y}^{ε}(t))] dt + Ω_{3}(σ, δ, ε, T) \end{matrix}  (3.5)

where for each T > 0,

lim sup_{δ→0} lim sup_{σ→0} lim sup_{ε→0} Ω_{3}(σ, δ, ε, T) = 0.

Here Av denotes the average.

Proof. Note that from the special form of V, and using the bounds a(n) ≤ Cn and λ′ ≤ C, and Schwarz’s inequality, for each ζ > 0,

[V_{ℓ}(η_{x}, u_{y}) − V_{ℓ}(η_{x}, u_{x})] ≤ Cζ^{−1}η_{x}^{2} + ζ(Ψ(u_{y}) − Ψ(u_{x}))^{2}.

By Lemma 3.5 we have

lim sup_{ε→0} ε^{d} ∑_{|y−x_{i}|≥δ} ∫_{0}^{T} Av_{(x,s)∈B_{σ}(y,t)}E [V (η_{x}(s), u_{y}(t)) − V (η_{x}(s), u_{x}(s))] dt

≤ Cζ^{−1} + ζ ∫_{0}^{T} ∫ Av_{(x,s)∈B_{σ}(y,t)}(Ψ(u(x, s)) − Ψ(u(y, t)))^{2}dydt.

By the strong convergence of Ψ(u^{ε}) to Ψ(u) this last term vanishes in the limit of small σ for each fixed δ. Sending ζ → ∞ gives the lemma. □

In the final lemma it will be convenient to think of our configuration η on ℤ^{d} instead of εℤ^{d} because we want to study the limiting measure. Note however that time is still sped up by a factor of ε^{−2}. We also need some notation. Let μ^{ε}(t, dη) be the distribution of the η at time t and let

μ̄_{ε,σ}(y, t, dη) = Av_{(x,s)∈B_{σ}(y,t)}τ_{x}μ^{ε}(s, dη)

so that Av_{(x,s)∈B_{σ}(y,t)}E [V_{ℓ}(η_{x}(s), u_{y}(t))] can be rewritten ∫ V_{ℓ}(η_{0}, u_{y}(t))μ̄_{ε,σ}(y, t, dη). If we further define a measure ν_{ε} on ℝ^{d} × [0, T] by ν_{ε}(dy, dt) = ε^{d} ∑_{|y−x_{i}|≥δ} δ_{y}(dy)dt and

β_{ε,σ}(dy, dt, dη) = μ̄_{ε,s}(y, t, dη)ν_{ε}(dy, dt)

then what we have on the right hand side of (3.5) is ∫ V_{ℓ}(η_{0}, u_{y}(t))β_{ε,σ}(dy, dt, dη).

Lemma 3.9. For each fixed σ > 0 the family {β_{ε,σ}(dy, dt, dη)}_{ε>0} is weakly compact and any weak limit β̄(dy, dt, dη) is of the form μ(y, t, dη)ν(dy, dt) where ν is Lebesgue measure (on {|y−x_{i}| ≥ δ}×[0, T]) and for each fixed y and t, μ(y, t, dη) is a translation invariant, invariant measure of the internal DLA process (without creation).

Proof. Let c_{x} ∈ ℓ^{2}(ℤ^{d}) with c_{x} > 0. Note that β_{ε,σ} has no mass on {|y| ≥ CT + 1} × [0, T] for some constant C since CT + 1 is an easy bound for the occupied set of the process up to time T. Let K_{M} = {|y| ≤ CT + 1, 0 ≤ t ≤ T, η_{z} ≤ Mc_{z}^{−1}, z ∈ ℤ^{d}} ⊂ ℝ^{d} × [0, T] × ℕ^{ℤ^{d}} be weakly compact sets increasing to the support of β_{ε,a} as M → ∞. By Čebyšev’s inequality it is straighforward to see that we have the following bound on the measure of the complement

β_{ε,σ}(K_{M}^{C}) ≤ M^{−2}∥c∥_{ℓ^{2}}^{2}E [∫_{0}^{T} Av_{|y−x_{i}|≥(δ−σ)ε^{−1}}η_{y}^{2}dt]

which vanishes in the limit M → ∞ by Lemma 3.5. This establishes the weak compactness. Since it is constructed as the limit of large averages of bounded measures the translation invariance of μ is immediate.

We now prove μ is invariant. Let f be any bounded local function and L_{0} be the generator of the process without creation (1.1).

Eβ̄[L_{0}f(η)] = lim_{ε→0} ε^{d} ∑_{|y−x_{i}|≥δ} ∫_{0}^{T} Av_{(x,s)∈B_{σ}(y,t)}E [L_{0}f(τ_{x}η)] dt.

By the translation invariance of the dynamics,

Av_{|x−y|≤σ}E [L_{0}τ_{x}f] = E [L_{0} {Av_{|x−y|≤σ}τ_{x}f}] .

By the definition of the generator,

\begin{matrix}∫_{t−σ}^{t+σ} E [L_{0} \{Av_{|x−y|≤σ}τ_{x}f\} (s)] ds \\ = ε^{2}E [Av_{|x−y|≤σ}τ_{x}f] |_{t−σ}^{t+σ} − ε^{2−d} ∫_{t−σ}^{t+σ} E [L_{c} \{Av_{|x−y|≤σ}τ_{x}f\}] ds. \end{matrix}

The right hand side is bounded by Cε^{2} where C is independent of y and t since the creation operator only hits finitely many shifts of the local function f. Therefore Eβ̄[L_{0}f] = 0 and we conclude that μ is invariant. □

Proof of Theorem 3.1. From the previous results we have, since V_{ℓ} is continuous with respect to weak convergence in η,

lim sup_{ε→0} E [∥η^{ε}(T) − u^{ε}(T)∥_{−1}^{2}]

≤ 2 ∫_{0}^{T} ∫_{|y−x_{i}|≥δ} E_{μ(y,t)} [V_{ℓ}(η_{0}, u(y, t))] dydt + Ω_{4}(ℓ, ε, σ, δ, T)

where for each T > 0,

lim sup_{δ→0} lim sup_{σ→0} lim sup_{ℓ→∞} lim sup_{ε→0} Ω_{4}(ℓ, ε, σ, δ, T) = 0.

Here for each y and t, μ(y, t) is a translation invariant, invariant for L. By Lemma 3.4 and the monotone convergence theorem we have

lim_{ℓ→∞} E_{μ}[V_{ℓ}(η_{0}, x)] ≤ 0.

Theorem 3.1 follows. □

4. Shape theory.

This section is devoted to convergence results for the occupied sets in the internal DLA and Stefan problems. Throughout, we will assume that the creation is at the origin only at constant rate c_{1} = 1. We will denote by A_{t}^{ε} ⊂ εℤ^{d} the occupied set {x : η_{x} > 0} and by A_{t} = {x : ρ(x, t) > 0} the occupied set of the Stefan problem (1.4).

Theorem 4.1. For every t > 0,

lim_{ε→0} dist(A_{t}^{ε}, A_{t}) = 0

in probability.

Remark. What follows immediately from the hydrodynamic limit is that macroscopic densities on the occupied set A_{t}^{ε} and the occupied set in the Stefan problem A_{t} agree. On the other hand, convergence in total variation means (1) that on the complement of a slightly larger set than the occupied set in the Stefan problem the density of particles is going to zero and, (2) that on a slightly smaller set than the occupied set in the Stefan problem the density of “holes” is going to zero. Note that the first statement follows immediately from the weak convergence, but the second does not because the density could be macroscopically correct while on the microscopic scale there could be clusters of free particles surrounded by seas of holes. The proof of Theorem 4.1 given below utilizes a bootstrapping argument to demonstrate that this does not happen. Actually, computer simulations leave no doubt that A_{t}^{1}/\sqrt{t} → A_{1} a.s. in the even stronger Hausdorff metric.

A free particle here means one at a site where there are more than α particles and the number of holes at a site x is α − η_{x} if η_{x} ≤ α.

Proof. Let t be a fixed time and γ > 0. Let 1/2 > δ > 0 and τ > 0. Let O_{δ,τ} ⊂ ℝ^{d} be the set where the solution of the Stefan problem (1.6) is larger than α + 2^{d+1}δ at time t − τ. Note that O_{δ,τ} → A_{t} at τ, δ → 0. Divide the set O_{δ,τ} into boxes of side length a = \sqrt{τ}. Let n be an integer large enough that α(1 − c_{0}δ/2)^{n} < γ where 0 < c_{0} < (2α)^{−1} is a small physical constant which will be described at the end of the proof. From the hydrodynamic limit we know that there are at least δ(aε^{−1})^{d} free particles in a smaller box of side length a/2 with the same center in each little box at each time t_{0}, ..., t_{n−1} where t_{i} = t − \frac{n−i}{n}τ, with probability which goes to one as ε → 0. We will now show that given this many free particles, the density of holes in each box is smaller than γ at time t with probability approaching one as ε → 0. Note that we can, without loss of generality, assume that any free particle which leaves the box in which it started is annihilated immediately, for the true number of holes which is destroyed is only decreased in this way.

Now we fix attention on one of our little lattice boxes on a time interval [t_{i}, t_{i+1}] and study how many holes are removed. Our box has side length L = aϵ^{−1} sites and there are N = δL^{d} free particles in a cocentric box of side length L/2. Let M denote the number of holes in our box (counted according to multiplicity). The free particles are symmetrically labelled 1 through N (ie. all permutations are equally likely) and they move by random walks until either 1. they hit a hole, in which case the particle dies and the hole multiplicity at that site is decreased by one, or 2. they hit the boundary where they are killed. Let X_{i} be the indicator function that the ith particle hits a hole before time t_{i+1} and let

β = E[X_{i}].

Note that β ≤ M/N. We are now going to show that Nβ > c_{0}δM also. We have

\begin{matrix}E[X_{i}] = & E[X_{i} | \frac{1}{N − 1} ∑_{j≠i} X_{i} > 2β]P (\frac{1}{N − 1} ∑_{j≠i} X_{i} > 2β) \\ + E[X_{i} | \frac{1}{N − 1} ∑_{j≠i} X_{i} ≤ 2β]P (\frac{1}{N − 1} ∑_{j≠i} X_{i} ≤ 2β). \end{matrix}

Note that E[(N − 1)^{−1} ∑_{j≠i} X_{i}] = β and that the X_{i} are negatively correlated so V ar[(N − 1)^{−1} ∑_{j≠i} X_{i}] ≤ β(1 − β)(N − 1)^{−1}. Now either β > 1/2 in which case Nβ > N/2 > c_{0}δM or β ≤ \sqrt{β(1 − β)} in which case by Čebyšev’s inequality,

P (\frac{1}{N − 1} ∑_{j≠i} X_{i} > 2β) ≤ \frac{1}{N − 1}.

Hence we have for large enough N,

E[X_{i}] ≥ \frac{1}{2}E[X_{i} | \frac{1}{N − 1} ∑_{j≠i} X_{j} ≤ 2β].  (4.1)

If \frac{1}{N−1} ∑_{j≠i} X_{j} ≤ 2β we have at least R = M − 2β(N − 1) holes remaining at time t. There are two possibilities: If R ≤ 0 then (as long as N ≥ 2),

Nβ ≥ \frac{1}{2}M ≥ c_{0}δM.

If R > 0 then, the probability for a particle starting inside the smaller box to hit a hole, of which there are at least R throughout, before time τ without hitting the boundary is at least the probability that at time τ it is within the set of holes, and has not hit the boundary, which by an elementary Brownian scaling argument is seen to be at least a fixed small constant c times the relative volume of the set of holes remaining vs. the size of the box,

E[X_{i} | \frac{1}{N − 1} ∑_{j≠i} X_{j} ≤ 2β] ≥ cR/L^{d}

which, together with (4.1) gives

Nβ ≥ c_{0}δM

for a fixed small constant c_{0}. This is the promised definition of c_{0}. Thus we have shown that always Nβ > c_{0}δM.

Now suppose that β > N^{−1/2}. By Čebyšev’s inequality again,

P (∑_{i=1}^{N} X_{i} ≤ c_{0}δM/2) ≤ Cε^{d/2}

where C is a finite constant (which does depend on a, δ and c_{0}). On the other hand if β ≤ N^{−1/2} we have that M = O(N^{1/2}), ie the density of holes has become much less than γ so we can stop. We repeat the argument in each time interval until we have either stopped, or used all n steps. In the worst case the probability that at least c_{0}δM/2 are killed in every single interval is at least (1 − Cε^{d/2})^{n} which goes to zero as ε → 0 so even in this case we are left with less than (1 − c_{0}δ/2)^{n} of the original holes, so the density of holes is less than γ. Since γ > 0 was an arbitrary positive constant, we have shown that the density of holes at time t is arbitrarily small. □

When d = 2, the creation and diffusion have the same scaling factor, yielding the following result.

Corollary 4.2. Assume that d = 2. As t → ∞, dist(A_{t}^{1}/\sqrt{t}, A_{1}) → 0 in probability as t → ∞.

We now show that in d = 2, A_{1} is a circular disc with the conjectured radius by constructing an explicit radially symmetric solution of (1.4).

Putting the equation is polar coordinates, the problem is to show that there exists a constant K ∈ (0, ∞), and a function ρ(r, t) which is smooth in the interior of {(r, t) : r ≤ K\sqrt{t}\} and satisfies

\frac{∂ρ}{∂t} = ( \frac{∂^{2}}{∂r^{2}} + \frac{1}{r} \frac{∂}{∂r} ) γ(ρ) + \frac{1}{2πr} δ_{0}

with boundary conditions

ρ(K\sqrt{t}, t) = 0,

\frac{∂ρ}{∂r}(K\sqrt{t}, t) = −\frac{Kα}{2\sqrt{t}}.

We look for self–similar solutions of the form ρ(r, t) = v(rt^{−1/2}). Writing s = rt^{−1/2}, F = γ^{−1}, and u = γ(v), one gets

\begin{matrix}u′′(s) + ( \frac{1}{s} + \frac{1}{2}sF′(u)(s) ) u′(s) + \frac{1}{2πs}δ_{0} = 0 \\ u(K) = 0, u′(K) = −\frac{1}{2}αK. \end{matrix}  (4.2)

It follows that, off s = 0, (su′)′ = −\frac{1}{2}s^{2}F′(u)u′ and therefore

su′(s) = −\frac{1}{2}αK^{2} + ∫_{s}^{K} \frac{1}{2}σ^{2}F′(u)(σ)u′(σ) dσ.  (4.3)

In addition, u(K) = 0 and su′(s) → −1/(2π) as s → 0. Recall first that F′ is bounded above and below by positive constants and note that (4.2) implies that u′ < 0. Moreover, (4.3) and Gronwall’s inequality imply that, for some constant C (which depends only on γ), −\frac{1}{2}αK^{2}e^{CK^{2}} ≤ su′(s) ≤ −\frac{1}{2}αK^{2}. As su′(s) is also increasing on (0, K), it has a limit as s → 0, which varies continuously between 0 and −∞ as K varies between 0 and ∞. Therefore, there exists a K such that lim_{s→0} su′(s) = −1/(2π), and therefore K can be chosen so that a solution of (4.2) exists (such a solution must be unique, since it produces a solution of (1.4)).

When γ(u) = u, and therefore F′(u) = 1, (4.2) can be explicitly solved, giving

u(s) = ∫_{s}^{K} \frac{1}{2πs}e^{−s^{2}/4} ds

which forces K = K(α) to be the solution of

e^{−K^{2}/4} = παK^{2}.  (4.4)

For example, to 5 significant digits, K = 0.54372 if α = 1, K = 0.39138 if α = 2, and K ∼ 1/\sqrt{πα} as α → ∞.

In particular, Corollary (4.3) now implies the following shape theorem.

Theorem 4.2. In dimension d = 2, there exists a constant K, such that

dist(A_{1}(t)/\sqrt{t}, B(0, K)) → 0

as t → ∞. In the linear case γ(u) = u, K is given by (4.4).

Now we address the long term behaviour of both the particle system and the Stefan problem in dimensions other that 2. For constant creations at a single site there is now no self-similar solutions. One way to circumvent this problem is to make the creation rate variable with time. Assume that the creation rate at the origin c_{1} is a function of t and satisfies c_{1}(t)t^{−(d−2)/2} → 1 as t → ∞. Then, after speeding up time by factor ε^{−2}, and sending ε → 0, the limiting PDE becomes

\frac{∂ρ}{∂t} = Δγ(ρ) + t^{(d−2)/d}δ_{0},

or in polar coordinates,

\frac{∂ρ}{∂t} = ( \frac{∂^{2}}{∂r^{2}} + \frac{1}{r} \frac{∂}{∂r} ) γ(ρ) + \frac{Γ(d/2)}{2π^{d/2}}r^{1−d}t^{(d−2)/2}δ_{0}

where Γ is the gamma function. The convergence is proved in the same way as for constant creation rate. A self–similar solution ρ = v(rt^{−1/2}) exists for any γ (with the same proof as before). In the case γ(u) = u,

v(s) = \frac{Γ(d/2)}{2π^{d/2}} ∫_{s}^{K} s^{1−d}e^{−s^{2}/4} ds , where \frac{Γ(d/2)}{π^{d/2}}e^{−K^{2}/4} = αK^{d}.

Now we study the problem of long time behaviour in dimensions different from 2 in the case of constant creation rate at the origin. It is proved in [BGL] that when d ≥ 3, γ(u) = u and α = 1, the asymptotic behavior of A_{t}^{1} is the same as when the source is discrete; that is, with probability 1, for every δ > 0,

((1 − δ)t^{1/d} · B(0, K)) ∩ ℤ^{d} ⊂ A_{t}^{1} ⊂ (1 + δ)t^{1/d} · B(0, K),  (4.5)

eventually in t. Here, K^{d} = \frac{dΓ(d/2)}{2π^{d/2}}, so that B(0, K) has unit volume. Our first theorem shows that the occupied set in the Stefan problem expands at the same rate.

Proposition 4.3. Assume that d ≥ 3, γ(ρ) = ρ, and α = 1 and δ > 0 is fixed. For sufficiently large t,

(1 − δ)t^{1/d} · B(0, K) ⊂ A_{t} ⊂ (1 + δ)t^{1/d} · B(0, K),

with the same K as in 4.3. In the case of general γ, there exist two constants 0 < K_{1} < K_{2} < ∞, so that

t^{1/d} · B(0, K_{1}) ⊂ A_{t} ⊂ t^{1/d} · B(0, K_{2})

.

Proof. The idea is to construct sub and supersolutions of the Stefan problem (1.4) out of the explicit solution

v(r, t) = ∫_{0}^{t}(4πs)^{−d/2}e^{−r^{2}/4s} ds = \frac{1}{4π^{d/2}} r^{2−d} ∫ ∞_{\frac{r^{2}}{4t}} s^{\frac{d}{2} −2}e^{−s}ds

of the heat equation ∂v/∂t = Δv + δ_{0}. Consider first the case γ(ρ) = ρ. We have

\frac{∂v}{∂r} (r, t) = −\frac{r^{1−d}}{2π^{d/2}} ∫ ∞_{\frac{r^{2}}{4t}} s^{(d−2)/2}e^{−s}ds = −\frac{r^{1−d}}{2π^{d/2}} (Γ(d/2) + O(r^{d}t^{−d/2})).  (4.6)

Let us take a ball of radius Ct^{1/d}. The rate of expansion of the boundary is Ct^{\frac{1}{d} −1}/d. Note that this coincides with −\frac{∂v}{∂r} (Ct^{1/d}, t) up to order t^{1/d−d/2} exactly when C = K where K^{d} = \frac{dΓ(d/2)}{2π^{d/2}} as in (4.5), ie. exactly when the ball has volume t. Define

ρ(r, t) = v(r, t) − v(Ct^{1/d}, t)

for r ≤ Ct^{1/d} and v = 0 outside the ball. If C > K then, after a certain time, ρ is a supersolution of the Stefan problem, and if C < K then it is eventually a subsolution. The proposition then follows from Lemma 2.4.

For the case of general γ, note that the spatial coupling term Δγ(ρ) can be written as ∇ · D∇ where from the assumptions γ′ = D is uniformly bounded above and below away from zero. Aronson’s estimates [D] tell us that there are positive finite constants c_{1}, c_{2}, b_{1}, b_{2} such that

c_{1}(4πs)^{−d/2}e^{−b_{1}|x−y|^{2}/4(t−s)} ≤ p(s, x, t, y) ≤ c_{2}(4πs)^{−d/2}e^{−b_{2}|x−y|^{2}/4(t−s)}

where p(s, x, t, y) is the solution at (y, t) of the corresponding heat equation, starting from a Dirac mass at x at time s. Therefore the flux across the boundary of the ball of radius Ct^{1/d} is bounded above and below in terms of the corresponding flux in the linear case (4.6). But now we can repeat the same argument as in the linear case to obtain the upper and lower bounds. □

In one dimension the proportion of created particles which are killed converges to 0 as time progresses. Therefore, as pointed out in [BGL], the rate of boundary expansion is of much slower order than t. Again, [BGL] have a conjecture in the linear case which we now state and prove as a theorem.

Theorem 4.4. Assume that d = 1, γ(u) = u and α = 1. Then for every ε > 0, with probability 1, for sufficiently large t,

(1 − ε)\sqrt{2t log t} · [−1, 1] ∩ ℤ ⊂ A_{t}^{1} ⊂ (1 + ε)\sqrt{2t log t} · [−1, 1].  (4.7)

Proof. Choose 1 < c < 1 + ε. Let N_{t} be the number of particles created during [0, t] which do not exit [−c\sqrt{2t log t}, c\sqrt{2t log t}] before time t. Let p_{s,t} the probability of the event that a random walk started at the origin at time s ∈ [0, t] visits a site outside [−c\sqrt{2t log t}, c\sqrt{2t log t}] during the time interval [s, t]. Then, by Skorohod embedding,

p_{s,t} ≤ p_{0,t} ≤ C · P (W_{t} ≥ c\sqrt{2t log t}) ≤ C · t^{−c^{2}/2}.

Then P (N_{t} ≥ \sqrt{t}) ≤ E(exp(−\sqrt{t} + N_{t})), where E(exp(N_{t})) can be estimated using the upper bound on p_{s,t}, giving

P (N_{t} ≥ \sqrt{t}) ≤ exp ( −\sqrt{t} + Ct^{1−c^{2}/2} ) .

Hence, P (N_{t} ≤ \sqrt{t} eventually) = 1. The upper bound in (4.7) now follows since (1 + ε − c)\sqrt{t log t} > \sqrt{t} for large t.

To prove the lower bound, let A_{t}^{1} = [−L_{t}, R_{t}]. We will first show that, with large probability, at least one of L_{t} and R_{t} must be at least c\sqrt{t log t}, for some small c. Then we will show that both L_{t} and R_{t} must satisfy this lover bound; in the final step we will bootstrap c close to \sqrt{2}.

To this end, fix c < 1 − 2ε and let c_{1} = c/(1 − 2ε). Let us momentarily consider the free particle system with creation only, without killing. In this case, let N_{t} (resp. N′_{t}) be the number of particles created in [0, εt], whose position is at least c\sqrt{2t log t} (resp. at most −c\sqrt{2t log t}) at time t. Again, P (N_{t} ≤ c\sqrt{2t log t}) ≤ E(exp(c\sqrt{2t log t} − N_{t}), and normal approximation can be used to show that

E(exp(−N_{t})) ≤ exp(−Cϵt^{1−c_{1}^{2}/2}).  (4.8)

so that P (N_{t} ≥ c\sqrt{t log t} eventually) = 1 and by symmetry the same is true for N_{t} replaced by N′_{t}. The natural coupling between the internal DLA and the free particle system is the one in which the particles in the internal DLA are tagged but allowed to continue to diffuse. Under this coupling, if N_{t} and N′_{t} are both larger than c\sqrt{t log t}, then R_{t} ∨ L_{t} ≥ c\sqrt{t log t}.

By the previous paragraph, it follows that for large enough t, R_{εt} ∨ L_{εt} ≥ cε\sqrt{t log t}. Assume, for example, that L_{εt} ≥ cε\sqrt{t log t}. Now make particles created after time εt execute free random walks, and define N_{t}^{1} to be the number of particles created in [εt, 2εt], whose position is at least c\sqrt{2t log t} at time t and which do not go below −L_{εt} during the time interval [εt, t]. A single particle will satisfy these two conditions with probability at least εt^{−c_{1}^{2}/2}, and so (4.8) holds with N_{t} replaced by N_{t}^{1} and ε replaced by ε^{2}. Thus it follows that, with probability 1, for large t, R_{t} ∧ L_{t} ≥ cε\sqrt{t log t}.

For the final bootstrapping step, notice that what we proved so far implies that R_{εt} ∧ L_{εt} ≥ cε^{2}\sqrt{t log t} for t large enough. Then repeat the argument in the previous paragraph with ε replaced by ε^{2} to show that R_{t} ∧ L_{t} ≥ c\sqrt{t log t} with probability at least exp(−Cε^{3}t^{1−c_{1}^{2}/2}). □

Again, it turns out that A_{t} has exactly the same rate of expansion.

Proposition 4.5. Assume that d = 1, γ(v) = v and α = 1, and fix an ε > 0. Then, for large t,

(1 − ε)\sqrt{2t log t} · [−1, 1] ⊂ A_{t} ⊂ (1 + ε)\sqrt{2t log t} · [−1, 1].

In the case of general γ, there exist two constants 0 < K_{1} < K_{2} < ∞, so that K_{1}\sqrt{t log t}· [−1, 1] ⊂ A_{t} ⊂ K_{2}\sqrt{t log t} · [−1, 1].

Proof. With notation as in the proof of Proposition 4.3, we need to prove that for large t,

−\frac{∂v}{∂r} (C\sqrt{2t log t}, t) ≥ C\sqrt{\frac{log t}{2t}},  (4.9)

when C < 1, while the inequality is reversed when C > 1. This follows from the fact that the integral on the left of (4.9) is asymptotic to a constant multiple of t^{−C/2}(log t)^{−1/2}. The case of general γ again follows as in the proof of Proposition 4.3 by using Aronson’s estimates. □

Remark. The upper and lower bounds obtained in Propositions 4.3 and 4.5 can easily be extended to the case of multiple creation sites or α > 1, thereby proving the shape theorems in this cases as well.

5. Invariant measures.

In this section we prove Lemma 3.3 which identifies the infinite volume invariant measures of our degenerate zero–range processes. All the key ideas are contained in [A], but our process is sufficiently different (and [A] contains enough misprints) that we are compelled to provide a fair amount of details.

The family of measures μ^{β} on configurations η ∈ {0, 1, … }^{ℤ^{d}} comprises the following two sets of measures:

(1) for every u ≥ α, the translation invariant product measures in which, for every x, η_{x} ≥ α and η_{x} − α is distributed as the marginal for the canonical invariant measure for the zero-range dynamics L_{1} with density ρ = u − α and

(2) all Dirac masses on configurations in {0, 1, …, α}^{ℤ^{d}}.

It is clear that any such μ^{β} is invariant. What we need is the converse.

Theorem 5.1. Any invariant, translation invariant measure for L_{0} defined in (1.1) is a mixture of μ_{β}. More precisely, any extremal invariant, translation invariant measure is one of μ_{β}.

The only methods available at the present time to prove such a result rely on either attractiveness or duality. Since there is no useful duality available in our context we have to assume attractiveness and use a method based on coupling. The coupled process makes the usual attempt to “move the two processes in unison as much as possible.” Formally, the generator of the coupled process is given by

L̄f(η, ξ) = ∑_{x∼y}1_{η_{x}>ξ_{x}}(a(η_{x}) − a(ξ_{x}))(f(η^{xy}, ξ) − f(η, ξ))
+1_{η_{x}<ξ_{x}}(a(ξ_{x}) − a(η_{x}))(f(η, ξ^{xy}) − f(η, ξ))
+a(η_{x} ∧ ξ_{x})(f(η^{xy}, ξ^{xy}) − f(η, ξ)).

In the statements of Lemmas 5.2 -5.5, we assume that ν̄ is an invariant, translation invariant measure for the coupled process.

Lemma 5.2. For any two neighboring sites x_{0} and y_{0}:

ν̄{(η, ξ) : η_{x_{0}} < ξ_{x_{0}}, ξ_{y_{0}} < η_{y_{0}}, α < η_{y_{0}}} = 0.  (5.1)

Proof. Let f(η, ξ) = (ξ_{x_{0}} − η_{x_{0}})_{+}. Then

L̄f = ∑_{x∼x_{0}}1_{η_{x_{0}} >ξ_{x_{0}}} (a(η_{x_{0}}) − a(ξ_{x_{0}}))(f(η^{x_{0}x}, ξ) − f(η, ξ))
+1_{η_{x_{0}} <ξ_{x_{0}}} (a(ξ_{x_{0}}) − a(η_{x_{0}}))(f(η, ξ^{x_{0}x}) − f(η, ξ))
+a(η_{x_{0}} ∧ ξ_{x_{0}})(f(η^{x_{0}x}, ξ^{x_{0}x}) − f(η, ξ))
+1_{η_{x}>ξ_{x}}(a(η_{x}) − a(ξ_{x}))(f(η^{xx_{0}}, ξ) − f(η, ξ))
+1_{η_{x}<ξ_{x}}(a(ξ_{x}) − a(η_{x}))(f(η, ξ^{xx_{0}}) − f(η, ξ))
+a(η_{x} ∧ ξ_{x})(f(η^{xx_{0}}, ξ^{xx_{0}}) − f(η, ξ))
= ∑_{x∼x_{0}}0
(I)  +1_{ξ_{x_{0}} >η_{x_{0}}} (a(ξ_{x_{0}}) − a(η_{x_{0}}))(−1)
+0
(II)  +1_{η_{x}>ξ_{x},ξ_{x_{0}} >η_{x_{0}}} (a(η_{x}) − a(ξ_{x}))(−1)
(III)  +1_{η_{x}<ξ_{x},ξ_{x_{0}} ≥η_{x_{0}}} (a(ξ_{x}) − a(η_{x}))(+1)
+0.

Then

Eν̄( ∑_{x∼x_{0}} expression (I)) = − ∑_{x∼x_{0}} ∑_{k>ℓ}(a(k) − a(ℓ))ν̄(ξ_{x_{0}} = k, η_{x_{0}} = ℓ)

while

Eν̄( ∑_{x∼x_{0}} expression (II)) = ∑_{x∼x_{0}} ∑_{k>ℓ}(a(k) − a(ℓ))ν̄(ξ_{x} = k, η_{x} = ℓ, ξ_{x_{0}} ≥ η_{x_{0}}).

Therefore, by translation invariance,

Eν̄(L̄f) ≤ Eν̄( ∑_{x∼x_{0}} expression (III))

≤ − ∑_{k<ℓ}(a(ℓ) − a(k))ν̄(ξ_{y_{0}} = k, η_{y_{0}} = ℓ, ξ_{x_{0}} > η_{x_{0}}).

However, as ν̄ is invariant, Eν̄(L̄f) = 0. Since the rates a(ℓ) are strictly incresing for ℓ > α, ν(ξ_{y_{0}} = k, η_{y_{0}} = ℓ, ξ_{x_{0}} > η_{x_{0}}) = 0 whenever k < ℓ and α < ℓ, which is equivalent to the statement of the lemma. □

Lemma 5.3. For any two neighboring sites x_{0} and y_{0}:

ν̄{(η, ξ) : η_{y_{0}} > α, η_{x_{0}} < α} = 0.  (5.2)

Proof. This time, let f = (η_{x_{0}} − α)_{+}. Then

L̄f = ∑_{x∼x_{0}} ( a(η_{x_{0}})1_{η_{x_{0}} >α}(−1) + a(η_{x})1_{η_{x_{0}} ≥α}(+1) ) .

Therefore

Eν̄(L̄f) = ∑_{x∼x_{0}} ∑_{k>α} a(k) ( −ν̄(η_{x_{0}} = k) + ν̄(η_{x} = k, η_{x_{0}} ≥ α) )

= − ∑_{x∼x_{0}} ∑_{k>α} a(k)ν̄(η_{x} = k, η_{x_{0}} < α).

This implies that, for each k > α, ν̄(η_{y_{0}} = k, η_{x_{0}} < α) = 0. □

We will shorten max η = max{η_{x} : x ∈ ℤ^{d}} and min η = min{η_{x} : x ∈ ℤ^{d}}.

Lemma 5.4. ν̄(max η > α, min η < α) = 0.

Proof. It suffices to show that (5.2) holds for arbitrary sites x_{0} and y_{0}. To this end we will prove, by induction on k ≥ 1, that (5.2) holds for all x_{0}, y_{0} for which ||x_{0}−y_{0}||_{1} = k. To prove the k − 1 → k step, take a y_{1} with ||x_{0} − y_{1}||_{1} = k − 1 and ||y_{1} − y_{0}|| = 1. Choose the function f = 1_{η_{y_{1}} >α,η_{x_{0}} <α}. Then, by the induction hypothesis, f = 0 ν̄-a.s., and therefore L̄f ≥ 0 ν̄-a.s. However, Eν̄(L̄f) = 0, thus so is Eν̄ of the part of L̄f which corresponds to a particle jumping from y_{0} to y_{1} in η. That is,

0 = Eν̄ ( a(η_{y_{0}})1_{η_{y_{1}} ≥α,η_{x_{0}} <α} )

= ∑_{k>α} a(k)ν̄(η_{y_{0}} = k, η_{y_{1}} ≥ α, η_{x_{0}} < α)

= ∑_{k>α} a(k)ν̄(η_{y_{0}} = k, η_{x_{0}} < α),

the last inequality being the consequence of Lemma 5.2 applied to y_{0} and y_{1}. Thus ν̄(η_{y_{0}} = k, η_{x_{0}} < α) = 0 for any k > α. □

Lemma 5.5. The two marginals of ν̄ are monotonely coupled, that is

ν̄ ( ( ∩_{x∈ℤ^{d}}{(η_{x} − α)_{+} ≥ (ξ_{x} − α)_{+}} ) ∪ ( ∩_{x∈ℤ^{d}}{(η_{x} − α)_{+} ≤ (ξ_{x} − α)_{+}} ) ) = 1.

Proof. This proof is very similar to the previous one, so we just point out the main steps. It suffices to prove that equation (5.1) holds for arbitrary sites x_{0} and y_{0}, that is, for ||x_{0}−y_{0}||_{1} = k for k ≥ 1. Again, take a y_{1} with ||x_{0}−y_{1}||_{1} = k−1 and ||y_{1}−y_{0}|| = 1, and choose the function f = 1_{η_{y_{1}} >ξ_{y_{1}} ,ξ_{x_{0}} >η_{x_{0}} ,α<η_{y_{1}}} . As before, the induction hypothesis implies that Eν̄ of the part of L̄f which corresponds to a particle jumping from y_{0} to y_{1} in η, but not in ξ, is 0:

0 = Eν̄ ( (a(η_{y_{0}}) − a(ξ_{y_{0}}))1_{η_{y_{0}} >ξ_{y_{0}}} f(η^{y_{0}y_{1}}, ξ) )

= ∑_{k>ℓ}(a(k) − a(ℓ))ν̄(η_{y_{0}} = k, ξ_{y_{0}} = ℓ, η_{y_{1}} ≥ ξ_{y_{1}}, ξ_{x_{0}} > η_{x_{0}}, α ≤ η_{y_{1}})

= ∑_{k>ℓ}(a(k) − a(ℓ))ν̄(η_{y_{0}} = k, ξ_{y_{0}} = ℓ, ξ_{x_{0}} > η_{x_{0}}, α ≤ η_{y_{1}}),

by Lemma 5.2. It follows that

0 = ν̄(η_{y_{0}} > ξ_{y_{0}}, η_{y_{0}} > α, ξ_{x_{0}} > η_{x_{0}}, α ≤ η_{y_{1}})
= ν̄(η_{y_{0}} > ξ_{y_{0}}, η_{y_{0}} > α, ξ_{x_{0}} > η_{x_{0}}),

by Lemma 5.3. This completes the proof. □

Proof of Theorem 5.1. Assume that μ_{1} is an extremal invariant translation invariant measure for the process. Pick another extremal invariant measure μ_{2} (to be specified later). Then there exists an extremal invariant measure ν̄ for the coupled process, with marginals μ_{1} and μ_{2} (Lemma 4.3 in [A]).

By Lemma 5.4, one can decompose the measure μ_{1} as follows:

μ_{1}(η ∈ ·) =μ_{1}(η ∈ · | min η < α)μ_{1}(min η < α)
+μ_{1}(η ∈ · | max η > α)μ_{1}(max η > α)
+μ_{1}(η ∈ · |η ≡ α)μ_{1}(η ≡ α).

The three conditional measures above are invariant and translation invariant. Since μ_{1} is extremal, one of the three non–conditional probabilities above must be 1. If either μ_{1}(min η < α) = 1 or μ_{1}(η ≡ α) = 1, then μ_{1}(max η ≤ α) = 1 so there is nothing to prove. If μ_{1}(max η > α) = 1, then η_{x} ≥ α for every x μ_{1}-a.s. Assume this for the rest of the proof.

Now let μ_{1} be an extremal invariant measure with u = E_{μ_{1}}(η_{0}) ≥ α. Moreover, let μ_{2} to be the extremal measure μ_{u} from subset (1) of the family μ^{β}; that is, μ_{u} is the product measure in which, for every x, η_{x} ≥ α and η_{x} − α has the marginal distribution of the canonical invariant measure for the non-degenerate zero-range process L_{1}. Then, by Lemma 5.5 and the fact that ν̄ is extremal (see the proof of Lemma 4.5 in [A]),

ν̄ ( ∩_{x∈ℤ^{d}}{η_{x} ≥ ξ_{x}} ) = 1 or ν̄ ( ∩_{x∈ℤ^{d}}{η_{x} ≤ ξ_{x}} ) = 1.  (5.3)

As μ_{1} and μ_{2} have the same density, Eν̄(η_{x} − ξ_{x}) = 0 for every x. This forces η_{x} = ξ_{x} ν̄-a.s. in either case of (5.3), that is, ν̄({(η, ξ) : η ≡ ξ}) = 1. In particular, μ_{1} = μ_{u}. □

References.

[A] Andjel E., Invariant measures for the zero-range processes. Ann. Prob. 10 525–547 (1982).

[BGL] Bramson, M., Griffeath, D., Lawler, G., Internal diffusion limited aggregation, Ann. Prob.20 (1992), no. 4, 2117–2140.

[CY] Chang, C.C., Yau, H.-T.,Fluctuations of one-dimensional Ginzburg-Landau models in nonequilibrium. Comm. Math.. Phys. 145 (1992), no. 2, 209–234.

[D] Davies, E.B., Heat Kernels and Spectral Theory, Cambridge University Press, 1989.

[DF] Diaconis, P., Fulton, W., A growth model, game, an algebra, Lagrange inversion, and characteristic classes, Ren. Sem. Mat. Univers. Politecn. Torino, 49 95–119 (1991).

[F] Friedman, A., Variational principles and free boundary problems, Wiley-Interscience, New York, 1982.

[L] Lawler, G.F., Subdiffusive fluctuations for internal diffusion limited aggregation, Ann. Prob. 23 (1995), 71–86.

[LM] Lions, J.-L., and Magenes, E., Non-homogeneous Boundary Value Problems and Applications, Grundlehren der Mathematischen Wissenschaften Vols. 181-182, Springer, New York-Heidelberg, 1972.

[M] Meirmanov, A.M., The Stefan problem, Walter de Gruyter, Berlin, 1992.

[S] Spohn, H., Interface motion in models with stochastic dynamics, J. Stat. Phys. 71 (1993).

[Y] Yau, H.-T., Metastability of Ginzburg-Landau model with a conservation law. J. Statist. Phys. 74 (1994), no. 3-4, 705–742.

Janko Gravner, Department of Mathematics, University of California, Davis, CA 95616, USA

INTERNAL DLA 35

E-mail address: gravner@math.ucdavis.edu

Jeremy Quastel, Departments of Mathematics and Statistics, University of Toronto, 100 St. George Street, Toronto, Ontario M5S 3G3, Canada

E-mail address: quastel@math.toronto.edu