Kinetic Limits of Piecewise Deterministic Markov Processes and Grain Boundary Coarsening by Joe Klobusicky B.S., Carnegie Mellon University; Pittsburgh, PA, 2009 M.S., Carnegie Mellon University; Pittsburgh, PA, 2009 M.Sc., Brown University; Providence, RI, 2010 A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in The Division of Applied Mathematics at Brown University PROVIDENCE, RHODE ISLAND May 2014 c Copyright 2014 by Joe Klobusicky This dissertation by Joe Klobusicky is accepted in its present form by The Division of Applied Mathematics as satisfying the dissertation requirement for the degree of Doctor of Philosophy. Date Govind Menon, Ph.D., Advisor Recommended to the Graduate Council Date Robert Pego, Ph.D., Reader Date Kavita Ramanan, Ph.D., Reader Approved by the Graduate Council Date Peter Weber, Dean of the Graduate School iii Vitae Personal Information American citizen. Born January 7th , 1987 in Scranton, Pennsylvania. Education Brown University, Providence, Rhode Island USA M.S., Applied Mathematics, May 2010 Carnegie Mellon University, Pittsburgh, Pennsylvania USA M.S., Mathematics, May 2009 B.S., Mathematics, May 2009 Submitted Publications Building polyhedra by self-folding: theory and experiment, David Gracias et. al. iv Submitted, 2013. Self-assembly of mesoscale octahedral isomers, Shivendra Pandey et. al. Submitted, 2013. Papers in preparation Measure-valued solutions of a piecewise-deterministic Markov process limit, Joe Klobu- sicky, Govind Menon, and Robert Pego. In preparation, 2014. Kinetic limits of piecewise-deterministic Markov processes, Joe Klobusicky, Govind Menon, and Robert Pego. In preparation, 2014. Kinetic limits of piecewise-deterministic Markov processes applied to grain coarsen- ing, Joe Klobusicky, Govind Menon, and Robert Pego. In preparation, 2014. v Acknowledgements Our destiny is to live out what we think, because unless we live what we know, we do not even know it. -Thomas Merton Thanks to the following people who have helped bring my thoughts to action: To Govind Menon, for serving as my advisor for the past four years. Thank you for your patience as I regressed, digressed, and eventually progressed as a mathe- matician. Your commitment to excellence and clarity will serve as an example of how to approach hard problems, in mathematics and everything else. To Kavita Ramanan, for serving as a thesis committee member. Thank you for several helpful comments on queueing theory. I look forward to seeing what improvements will come from your insights. To Robert Pego: committee member, academic mentor, and overall mensch. For almost a decade, you have been a source of how to do things right. Thank you for the inordinate amount of time and concern that you gave during my years in Pitts- burgh. I’ll remember you as someone who cared a great deal about mathematics. I’ll remember you more as someone who cared enormously about passing his knowledge to others. To Johnny Guzman and Nat Trask. Thanks for letting me play numericist. I’ve vi lumped you both together because you share many things: an openness to new ideas, an academic humbleness that has substantially changed my previous rigid perceptions on mathematics, and a good sense of humor to boot. To those teachers during my upbringing who gave me encouragement. In partic- ular, I thank Mrs. Betti, the first person I can remember who admitted to enjoying mathematics. I also thank Mr. Walsh, my biology teacher and golf coach, for his timeless advice in the classroom and on the course. The term “a gentleman and a scholar” gets thrown around a lot, but you are the true embodiment of it. Most of all, to my family. Mom and Dad, your selflessness and love is so strong and unwavering that I often take it for granted. When I think about how unreason- ably blessed I am to have parents like you, all of my problems seem petty. Thanks for everything. vii Abstract of “Kinetic Limits of Piecewise Deterministic Markov Processes and Grain Boundary Coarsening” by Joe Klobusicky, Ph.D., Brown University, May 2014 The subject of this thesis is the development of a stochastic process that can be inter- preted as a model of grain coarsening, a central problem in material science. Namely, our objective is the creation of piecewise-deterministic Markov processes (PDMPs) on particles whose empirical densities converge to solutions of kinetic equations. We begin with a study of a simplified model, where particles drift toward the origin on R+ . When particles hit the origin, they are redistributed to R+ according to a fixed probability density p(x). This “dynamic shuffler” is shown to be an instance of a PDMP. From the infinitesimal generator associated with the dynamic shuffler, we then construct a martingale that we show is an approximation to a weak solution to limiting kinetic equations of the density of particles. Properties of these equations are then investigated using Laplace transforms and renewal theory. We then generalize the one-tier model by considering a k-species model, where particles travel on several tiers, and are redistributed with probability distributions that change with time. Obtaining the “correct” martingale in this “ k-species model” involves an augmented PDMP, which adds variables that keep track of certain types of jumps. The main theorem of this thesis is the existence of limiting kinetic equa- tions for densities of each species. We end by aligning the k-species model to a mean field model of grain coarsening. Several basic properties of this model which are inherent for grain systems are then shown. Finally, we run simulations of the grain coarsening PDMP and raise several conjectures on the universality of certain grain statistics. Contents Vitae iv Acknowledgments vi 1 Introduction 1 1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2 Mean-field models of grain growth . . . . . . . . . . . . . . . . . . . . 4 1.3 PDMPs and grain coarsening . . . . . . . . . . . . . . . . . . . . . . 5 1.4 Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2 An Overview of PDMPs 9 2.1 PDMPs and generators . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.1.1 A note on Skorokhod topologies . . . . . . . . . . . . . . . . . 14 3 Kinetic Limits of Piecewise-Deterministic Markov Processes: The Dynamic Shuffler 15 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 3.2 A solution via a stochastic limit . . . . . . . . . . . . . . . . . . . . . 18 3.2.1 Convergence of empirical processes . . . . . . . . . . . . . . . 21 3.3 A solution via the Laplace transform . . . . . . . . . . . . . . . . . . 28 3.4 A measure valued solution formula . . . . . . . . . . . . . . . . . . . 31 3.5 An example: a point mass under a stationary distribution . . . . . . 33 3.6 Strong solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.7 Weak solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.8 Stationarity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 4 Kinetic Limits of Piecewise-Deterministic Markov Processes: The k-Species Model 45 4.1 Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 4.2 A kinetic mixing model . . . . . . . . . . . . . . . . . . . . . . . . . . 49 4.2.1 The model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 viii 4.3 The k-species model in a PDMP framework . . . . . . . . . . . . . . 51 4.3.1 PDMP formulation of the k-species model . . . . . . . . . . . 51 4.3.2 The main theorem . . . . . . . . . . . . . . . . . . . . . . . . 53 4.3.3 Method of proof . . . . . . . . . . . . . . . . . . . . . . . . . . 57 4.4 The extended PDMP X ˜ N (t) . . . . . . . . . . . . . . . . . . . . . . . 58 4.5 Reassignment bounds and martingale equations . . . . . . . . . . . . 61 4.6 Existence of limiting measures . . . . . . . . . . . . . . . . . . . . . . 67 4.7 Convergence of boundary variables . . . . . . . . . . . . . . . . . . . 78 4.8 The existence interval revisited. . . . . . . . . . . . . . . . . . . . . . 86 4.9 Uniqueness and regularity for kinetic equations . . . . . . . . . . . . 88 5 The k-species Model Applied to Grain Coarsening 100 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 5.2 Side redistribution at singular events . . . . . . . . . . . . . . . . . . 101 5.3 PDMPs and grain coarsening . . . . . . . . . . . . . . . . . . . . . . 106 5.4 Properties of grain coarsening . . . . . . . . . . . . . . . . . . . . . . 111 5.5 Computations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 5.6 Remarks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 6 Conclusion 133 ix List of Figures 1.1 Coarsening of a grain network. Grains in material microstructure delete through annealing, causing the average grain size to increase. . 3 1.2 Grain growth via a PDMP. Particles drift according to their num- ber of sides. Each tier denotes particles of different side numbers. In this instance, a vanishing grain in tier three triggers a random reas- signment of the number of sides of grains in tiers 4,7, and 8. . . . . . 6 2.1 Topologies for cadlag functions. Left: Two functions that are “close” in the J1 topology. Note that the magnitude of their jumps are similar. Right: Two functions close in the M1 topology. Jumps are not required to be close in this case, but rather the Hausdorff distances of each function’s completed graph. . . . . . . . . . . . . . 14 3.1 A dynamic shuffler. Particles (black dots) drift left with unit speed until hitting the origin, where a particle (red dot) is redistributed according to a density p(x). . . . . . . . . . . . . . . . . . . . . . . . 17 3.2 M1 convergence of measure valued initial data. Left: Particles (small vertical segments) approximating a point mass (segment with red dot) drift toward the origin. Right: The value of F (φ)(x) before and after particles hit the origin. As the number of particles grows, their jumps form a “monotone staircase” connecting the jump of the limiting functional. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 4.1 A PDMP on four species. Particles travel on four separate copies of R+ . Velocity directions are represented by horizontal arrows, with species 3 having zero velocity. A boundary event occurs when a par- ticle (labelled by “A”) hits the origin. Here, three particles are then randomly selected (K (2) = 3), and reassigned to different species by predetermined reassignments (given by vertical arrows). In this ex- (2) (2) (2) ample, R41 = 1, R32 = 2, and R23 = 3. . . . . . . . . . . . . . . . . 51 (1) 4.2 A PDMP on two species. In this case K (1) = 1, with R21 = 1. . 57 x 5.1 Continuation through a k-degree vertex. Possible topologies of curvature flow with k-ray initial conditions, with k = 3, 4, and 5. These correspond to planar rooted trivalent tree, with the rooted vertex in each figure denoted by a dot. . . . . . . . . . . . . . . . . . 105 5.2 Reassignments of sides before and after deletions. Numbers in grains refer to side numbers. . . . . . . . . . . . . . . . . . . . . . 109 5.3 Dispersal of a fluid plug. A carrier fluid is approximated by hor- izontal tiers, and a fluid plug (the inner region of the figure) is ap- proximated by particles (denoted by dots) which travel horizontally on tiers corresponding to the vertical velocity profile V (y). . . . . . . 118 5.4 Grain number N (t). . . . . . . . . . . . . . . . . . . . . . . . . . . 119 5.5 Mean area of grains with β = .01 Mean area is plotted with line of best fit y1 = .2806 + .9214x (the lines are indistinguishable). The correlation between average area and y1 is r2 = .9998. . . . . . . . . 119 5.6 Mean area of grains with β = .1 Mean area is plotted with line of best fit y1 = .1909 + 1.0347x (the lines are almost indistinguishable). The correlation between average area and y1 is r2 = .9979. . . . . . . 119 5.7 Mean area of grains with β = 1. Mean area is plotted with line of best fit y1 = −.2796 + 1.6624x (mean area is convex). The correlation between average area and y1 is r2 = .9881. . . . . . . . . . . . . . . . 120 5.8 Comparison of side number frequency between grain coarsening PDMP (dots) at β = 1 and the Kinderlehrer model (crosses) at t = 5. . . . . 120 5.9 Frequency of side number for β = .01 for various times. . . . . . . . . 120 5.10 Frequency of side number for β = .1 for various times. . . . . . . . . . 121 5.11 Frequency of side number for β = 1 for various times. . . . . . . . . . 121 5.12 Histograms of four-sided grain densities at β = .01. . . . . . . . . . . 122 5.13 Histograms of four-sided grain densities at β = .1. . . . . . . . . . . . 123 5.14 Histograms of four-sided grain densities at β = 1. . . . . . . . . . . . 124 5.15 Histograms of six-sided grain densities at β = .01. . . . . . . . . . . . 125 5.16 Histograms of six-sided grain densities at β = .1. . . . . . . . . . . . . 126 5.17 Histograms of six-sided grain densities at β = 1. . . . . . . . . . . . . 127 5.18 Histograms of eight-sided grain densities at β = .01. . . . . . . . . . . 128 5.19 Histograms of eight-sided grain densities at β = .1. . . . . . . . . . . 129 5.20 Histograms of eight-sided grain densities at β = 1. . . . . . . . . . . . 130 5.21 Average areas of grains with sides 2-10 at β = 0. Note how for each collection of n-gons,for n = 2, . . . , M , average area increases linearly. Similar behavior holds for β = .1 and β = 1. . . . . . . . . . . . . . . 131 5.22 Side to grain deletion ratio γβ (t) at β = .01. . . . . . . . . . . . . . . 131 5.23 Side to grain deletion ratio γβ (t) at β = .1. . . . . . . . . . . . . . . . 131 5.24 Side to grain deletion ratio γβ (t) at β = 1. . . . . . . . . . . . . . . . 132 xi Chapter One Introduction 2 1.1 Introduction A well-studied phenomenon in material science concerns the coarsening of microstruc- ture found in metals, porcelains, and other common materials. Typically, these ma- terials have polycrystalline structure, where each grain has a particular orientation. Annealing induces an evolution on grain boundaries which reduces their total length. During this process, grains are deleted (and not created), causing coarsening of mi- crostructure through an increase of average grain area (see Fig. 1.1). Microstructure can determine material properties. For instance, the Hall-Petch effect gives a relation between yield stress and average grain size (see [3], section 5.2.3). More scientific and mathematical properties of grain networks can be found in [16, 36]. Grain coarsening in two dimensions can be modeled by a geometric network, or planar graph G(t) ⊂ R2 , t ∈ (0, T ), with edges that evolve in time. Faces of G are called grains, and edges are referred to as boundaries. An idealization of the coarsening process assumes that grain boundaries are isotropic, or that grain boundary energies are identical along edges of the network. Under this idealization, we define a grain network with existence interval (0, T ) as a geometric finite planar graph G(t) which satisfies the following restrictions for all t ∈ (0, T ): 1. (Herring’s Condition) G(t) is trivalent, and angles between two edges at a vertex are fixed at 120 degrees. This restriction follows from a force balance law on edges emanating from a vertex (see [20]). 2. (Mean Curvature Flow) Boundaries are smooth, and satisfy mean curvature flow. This means that the time evolution of a smooth boundary γ(x, t) ⊂ R2 3 Figure 1.1: Coarsening of a grain network. Grains in material microstructure delete through annealing, causing the average grain size to increase. is also smooth, and satisfies ∂t γ(x, t) = M σκ(x, t)~n(x, t), t ∈ (0, T ), x ∈ R2 , (1.1) with M and σ denoting grain mobility and surface energy constants, ~n(x, t) denoting the inward unit normal, and κ(x, t) denoting mean curvature. This type of motion was observed experimentally by Beck for crystallites in metals [5]. Note that grain networks retain their topological structure during their existence time. The parabolic PDE (1.1) is referred to as curve-shortening flow, and it can in fact be shown that the total edge length of a network decreases in time [4]. A theorem of curve shortening flow on grain networks, due to von Neumann and Mullins [31, 38] gives an elegant relation between a topological quantity (sides of a grain) and a geometrical quantity (area of a grain). Theorem 1. (The von Neumann-Mullins n − 6 rule). For a grain with area A and n sides, dA π = M σ (n − 6). (1.2) dt 3 One consequence of the n − 6 rule is that grains with less than six sides may shrink to a point in finite time. Such deletions, along with possible edge deletions, 4 result in networks with vertices that violate Herring’s condition. This raises the question: What types of solutions exist with initial conditions that are not grain networks? Specifically, the continuation problem for a planar network G asks the following: If we are given a planar network G with vertices that may not satisfy Herring’s conditions, is there an ε > 0 where G(t) is a grain network for t ∈ (0, ε), and G(t) → G in the Hausdorff distance as t → 0+ ? Another interesting question is whether the flow is unique, or if it can become unique under additional assumptions. Currently, there are no comprehensive answers to the above questions. Viewed as a problem in parabolic PDE theory, the current focus is local, fixing on simplified models with initial conditions of k rays meeting at a vertex. Schnurer and Shulze found a unique self-similar solution with initial conditions of three rays meeting at arbitrary angles [34]. Mazzeo and Saez generalized with the k-ray initial conditions, and discovered a bijection between self-similar solutions and k-Steiner trees of a hyperbolic metric [28]. On the computational front, techniques of flowing through a graph not satisfying the grain network conditions include level set [10] and variational [26] methods. 1.2 Mean-field models of grain growth Another approach to studying grain boundary coarsening involves using mean field assumptions to gather bulk statistics of grain behavior. A well known kinetic model is due to Fradkov [14, 17], who exploits the n − 6 rule, along with an additional deterministic rule for side redistribution, to derive kinetic equations for densities of grains with a specific number of sides. Fradkov makes the following mean field assumptions about the redistribution of sides when either a grain or edge deletes: 5 1. If a grain or side deletes, topological changes of grains will be selected in proportion to the number of sides of grains, and independent of their areas. For instance, if a singular event causes a grain to lose a side, the probability that a particular grain with j sides will drop to one with j − 1 sides is j pj = P , (1.3) k>1 kNk where Nk is the number of grains with k sides. This eliminates the notion of grains having neighbors. 2. A free parameter β is introduced which is defined as the constant ratio between side deletions and grain deletions. Such an assumption comes from the lack of a topological rule for edge growth. For a number density fn (a, t) of n sided grains with area a at time t, the Fradkov equations take the form of a transport equation with a topological source term: ∂t fn (a, t) + (n − 6)∂a fn (a, t) = Γ(f (t))(Jf )n (a, t) (a, t) ∈ (0, ∞)2 , n ≥ 2, (1.4) with a collision operator J and coupling weight Γ determined from conservation of area and polyhedral defect (networks retain an average of six sides). Well-posedness and self-similar solutions for the Fradkov model can be found in [19, 21]. 1.3 PDMPs and grain coarsening The main goal of this thesis is to understand grain statistics as a hydrodynamic limit of finitely many grains respecting deterministic drift and random topological 6 Figure 1.2: Grain growth via a PDMP. Particles drift according to their number of sides. Each tier denotes particles of different side numbers. In this instance, a vanishing grain in tier three triggers a random reassignment of the number of sides of grains in tiers 4,7, and 8. transitions. In particular, we will construct limiting kinetic equations describing the evolution of grain area densities for the class of grains with a fixed number of sides. Toward this goal, the main tool we use is the theory of piecewise deterministic Markov process (PDMP). The construction and theory of PDMPs will be discussed in Chapter 2. However, not much is lost to think of a PDMP as a generalization of a jump process for particles that drift on vector fields. For the problem of grain coarsening, we wish to build a PDMP X(t) that tracks a finite set of grains (g1 , . . . , gn ) = ((a1 , s1 ), . . . , (an , sn )) with areas ai > 0 and side numbers si ∈ N. Note that the dimension of X(t) can change with time, as grains delete from coarsening. Grain areas change at a constant rate, following the n − 6 rule until a grain shrinks to a point, or an edge deletes (we will make an assumption on side deletion rates based on total particle number). At that time, we impose mean field assumptions to change topologies of grains. 7 1.4 Overview In Chapter 2, we describe a PDMP in generality. We then state its infinitesimal generator A, and the class of functionals D(A) that give rise to martingales through the Dynkin formula. In Chapter 3, we study a simplified version of the PDMP described in Section 1.3. Here, particles drift on the positive half-line until reaching the origin, where they are reassigned according to a probability density p(x). We first prove that for empirical densities approaching nonatomic initial conditions, the densities u(x, t) of particles converge to a weak form of the limiting kinetic equation ∂t u(x, t) − ∂x u(x, t) = p(x)u(0, t) t, x ∈ R+ , (1.5) u(x, 0) = u0 (x). We then focus on constructing measure valued weak solutions of (1.5) with mild restrictions on initial data. Finally, we use a popular result of renewal theory to show the existence of an attractor under certain assumptions. Chapter 4 is a generalization of PDMPs on grain networks described in Section 1.3. Here, we show the existence of fluid limits of PDMPs that allow generalized rules and probability distributions for particle jumps between different species. As in Chapter 3, we will show the existence of a weak limit for densities uk (x, t) on Lk , 8 that satisfy the advection-reaction equations ∂t uk (x, t) + ∂x (vk (x)uk (x, t)) = (1.6)   XM− K (l) M X X (l) (l) ul (0, t)vl (0)  Wi (t)ui (x, t)1R(l) =k − K (l) Wk (t)uk (x, t) ij l=1 i=1 j=1   M K int βG(t) X X + int  wint ui (x, t)1Rijint =k − K int wkint uk (x, t) , G (t) i=1 j=1 i with initial conditions uk (x, 0) = u0k (x), k = 1, . . . , M. The solutions of (1.6) will be shown to be unique in the class of L1 ∩ L∞ (R+ ) functions. Finally, in Chapter 5, we return to the problem of grain coarsening. First, we examine the topological behavior of grain networks before and after a grain or side vanishes. After finding a finite set of rules that dictate topological changes, we then define the various free parameters in the k-species model to align with a mean field model for grain coarsening. Next, we prove some basic properties of the limiting kinetic equation, such as conservation of mass and polyhedral defect. Finally, we end the chapter with a computational investigation of grain growth, and raise several conjectures about grain statistics. Chapter Two An Overview of PDMPs 10 2.1 PDMPs and generators The main class of stochastic processes utilized in this thesis is the piecewise deter- ministic Markov Process, or PDMP, whose theory is presented in Davis [9]. For a PDMP, a Markov process is allowed deterministic drift along flow lines of a vector field before randomly jumping to a new state. In addition to covering several ex- amples from queueing theory, including the M/G/1 queue, the GI/G/1 queue, the model also encapsulates Markov chains. Our chief object of interest is the infinitesi- mal generator and martingale from Dynkin’s formula associated with the PDMP. In Chatpers 3 and 4 we plan to use these martingales to provide an approximation to limiting kinetic equations. Our goal is to describe a generalized jump process that takes place on a disjoint union of manifolds, and follows a deterministic drift between jumps. To characterize the state space, we let S be a countable set, d : S → N , where for each s ∈ S, Ms ⊂ Rd(s) is an open set. The state space is then the disjoint union a  E= Ms = (s, x) : s ∈ S, x ∈ Ms ⊂ Rd(s) . (2.1) s∈S To define the topology of E, let ιs : Ms → E be the canonical injection defined by ιs (x) = (s, x). A set A in E is open if for every s, ι−1 s (A) is open in Ms . We then define E as the Borel sets of E. This makes (E, E) a Borel space. We now define a stochastic process X(t) = (s(t), ζ(t)), with a law (X(t))t≥0 based on the following: 1. Vector fields Xs , s ∈ S. defined on Ms . 11 2. A measurable function λ : E → R+ . 3. A transition measure Q : E × (E ∪ Γ∗ ) → [0, 1]. We will formally define Γ ∗ , the exit boundary of a PDMP, shortly. To see how this fits in with a generalized jump process, points in Ms will travel according to flows defined by Xs until either a Poisson clock with intensity λ rings, or the point hits the boundary of Ms . When such an event occurs, the point jumps to a new position in E, determined by Q. The vector fields Xs are chosen so that for every z ∈ Ms , there is a unique integral curve φs (t, z) satisfying d f (φs (t, z)) = Xs f (φs (t, z)) (2.2) dt φs (0, z) = z for any smooth function f : Rd(s) → R. Here we also require that the vector fields are conservative, meaning that there exists a t > 0 where φs (r, z) is defined for r ∈ [0, t]. Now let ∂ ∗ Ms be the exit boundary of Ms , defined as  ∂ ∗ Ms = y ∈ ∂Ms : φs (t− , x) = y for some (t, x) ∈ R+ × Ms (2.3) Γ∗ is then defined as the exit boundary of our state space: a Γ∗ = ∂ ∗ Ms (2.4) s∈S 12 At a given state x ∈ E, we define the exit time as t∗ (x) = inf{t > 0 : φs (t, x) ∈ ∂ ∗ Ms }. (2.5) The stochastic process (X(t))t≥0 with initial condition X(0) = (s, z) is then defined as follows. For x = (s, z), define a survivor function F as    R   exp − t λ(s, φs (r, z))dr , t < t∗ (x), 0 Fx (t) = (2.6)   0, ∗ t ≥ t (x). The rate λ : E → R+ is a measurable function, where for every state x = (s, x) ∈ E there exists ε > 0 where the function s → λ(s, φs (s, x)) is integrable for s ∈ [0, ε). We also define Q(A; x) to be a measurable function of x for each fixed A ∈ E on x ∈ E ∪ Γ∗ , and a probability measure on (E, E) for each X ∈ E ∪ Γ∗ . Now choose a random variable T1 such that P[T1 > t] = Fx (t). Now indepen- dently choose an E-valued random variable (L, Z) with distribution Q(∙; φs (T1 , z)). The trajectory of X(t) for t ≤ T1 is then    (s, φs (t, z)), t < T1 X(t) = (2.7)   (L, Z), t = T1 . From X(T1 ), we choose the next inter-jump time T2 − T1 and X(T2 ) in a similar fashion. It can be shown that the process X(t) is Markov, and in fact, strong Markov (Section 3 of [9]). As a Markov process, the PDMP has an associated infinitesimal generator A, 13 acting on a domain D(A), defined as the set of functions f : E → R where the limit Ex (f (X(t))) − f (x) (Af )(x) = lim+ (2.8) t→0 t exists for all x ∈ E. The infinitesimal generator then takes the form Z Af (x) = X (f (x)) + λ(x) (f (y) − f (x))Q(dy; x)), f ∈ D(A). (2.9) E From here, we can use Dynkin’s formula to derive the martingale Z t Mtf := f (X(t)) − f (X(0)) − Af (X(s))ds, f ∈ D(A). (2.10) 0 The following set of sufficient conditions for membership in D(A) is given in Rolski et. al. [33] Theorem 2. A function f : E → R satisfies f ∈ D(A) if the following hold: 1. The function t → f (φ(x, t)) is absolutely continuous on [0, s∗ (x)) for every x ∈ E, where s∗ (x) = inf t {t|F (t) = 0}. 2. f (x) = limt→0 f (φ(x, t)) exists for all x ∈ Γ∗ R 3. (Boundary condition):f (x) = E f (y)Q(dy, x) for x ∈ Γ∗ 4. (Finite expectation of number of jumps): for all t ∈ R+ , x ∈ E,   m(t) X E |f (X(Ti ) − f (X(Ti−1 ))| < ∞ (2.11) i=1 14 Figure 2.1: Topologies for cadlag functions. Left: Two functions that are “close” in the J1 topology. Note that the magnitude of their jumps are similar. Right: Two functions close in the M1 topology. Jumps are not required to be close in this case, but rather the Hausdorff distances of each function’s completed graph. 2.1.1 A note on Skorokhod topologies Throughout this thesis, our focus will be on stochastic processes which are cadlag, or are right continuous with lefthand limits. Let (M, d) denote a metric space with metric d. In the following, we use the J1 Skorohod topology, denoted D([0, t], M). For our purposes, M will either be R+ or M(R+ ), the space of finite measures on R+ with the Prohorov metric [6]. The J1 topology allows for convergence of functions that “wiggle” in time as well as space, and has the following characterization (see [23]): Let R be the set of continuous reparameterizations r : [0, t] → R+ : functions that are strictly increasing, where r(0) = 0 and r(t) = t. A sequence of functions αn → α in D([0, t], M) if and only if there is a sequence rn ∈ R with 1. sups |rn (s) − s| → 0, 2. sups≤t d(αn (rn (s)), α(s)) → 0. as n → ∞. If α(t) is continuous, the functions actually converge in the local uniform topology. However, in general, the local uniform topology is strictly stronger than the Skorohod J1 topology, but this fact won’t be used, since limiting functions of our interest in this paper will have no jumps (see Remark 2). Chapter Three Kinetic Limits of Piecewise-Deterministic Markov Processes: The Dynamic Shuffler 16 3.1 Introduction This chapter examines an instance of a PDMP, called the N particle dynamic shuffler, which tracks N particles on the positive real line moving with velocity -1 until reaching the origin. Upon hitting the origin, a particle is randomly redis- tributed back to R+ , with respect to a probability density p(x) in R+ . If initial empirical densities of N particles approach a limiting density u0 (x), we will show that the random densities of the N particle dynamic shuffler at time τ > 0 approach a deterministic density u(x, τ ), which satisfies a weak form of a one-dimensional transport equation with source term, ∂τ u(x, τ ) − ∂x u(x, τ ) = p(x)u(0, τ ) τ, x ∈ R+ , (3.1) u(x, 0) = u0 (x). In this chapter, we will refer to (3.1) as the redistribution equation. Regularity conditions on u0 (x) and p(x) largely depend on the types of solutions we seek, and will be considered throughout this chapter. We will borrow methods from coagulation theory. For instance, consider the one-dimensional Allen-Cahn equation ∂t u = ∂xx u + u − u3 , which describes a coag- ulation process on the real line. Specifically, given an interval partition of R, the smallest interval merges with its two neighbors to form a single interval, and the process continues. Menon, Niethammer, and Pego used the Allen-Cahn equation as a source of motivation in [29], investigating a wide range of clustering events. The discrete nature of these clustering phenomena suggested that it was reasonable to seek measure-valued solutions. Such solutions, continuous in time in the space of probability measures, were found using an intrinsic time scale based on the number 17 Figure 3.1: A dynamic shuffler. Particles (black dots) drift left with unit speed until hitting the origin, where a particle (red dot) is redistributed according to a density p(x). of clusters in the system. The same philosophy will be applied to (1), which will use the time scale of the total number of visits to the origin by particles, Z τ t= u(0, s)ds. (3.2) 0 We should note that while the tools used in coagulation theory prove useful to this chapter, in general there is no coagulation occurring in the redistribution equation! R The total number of particles N (τ ) = R+ u(x, τ )dx is preserved, while the total R mass, or first moment, M (τ ) = R+ xu(x, τ )dx can potentially change in time. The sections of this chapter are arranged as follows. In Section 3.2, we derive (3.1) as a limit of densities of the N particle dynamic shuffler. This is done by exploit- ing the infinitesimal generator, and building a martingale approximation to (3.1). Section 3.3 derives a solution formula for (3.1) for smooth initial conditions through the Laplace transform. We then use the change of variables (3.2) in Section 3.4 to derive a measure-valued solution formula for (3.1). A specific example is presented in Section 3.5, with initial data a point mass. Sections 3.6 and 3.7 prove well-posedness of weak and strong solutions for (1), meaning that the solution formulas of Sections 3.3 and 3.4 are well-defined. In Section 3.8, we conclude the chapter by providing asymptotics for u(x, τ ) via the key renewal theorem from renewal theory. 18 3.2 A solution via a stochastic limit A discrete approximation of the behavior underlying (3.1) can be obtained through PDMPs. In the following section, we show the dynamic shuffler of N particles is a PDMP, and then take a fluid limit as N → ∞. In this limit, particle densities will converge to a weak version of (3.1). Our process is as follows. Let particles x1 , . . . , xN ∈ R+ be given. Each par- ticle moves left with unit speed until a particle hits the origin. This particle is randomly assigned a new location according to a probability distribution function p(x) ∈ C 1 (R+ ). The shuffler continues as before, with particles drifting until a new particle hits the origin. We will now show that this process can be interpreted as a PDMP. The general theory and notation of a PDMP is discussed in Section 2.1. From the viewpoint of N a dynamic shuffler with N particles, our state space is E N = (R+ ) , with a state written as x = (x1 , . . . , xN ) ∈ E N . (3.3) There is no need to consider disjoint copies of E N . This is because particles are neither created nor destroyed, and there is only one vector field, corresponding to P each particle having velocity of −1, namely X = − N ∂ i=1 ∂xi . Since all jumps occur on the boundary of E, where there is a particle at the origin, we have zero jump frequency: λ ≡ 0. For an element y ∈ ∂Γ∗ of the form y = {y1 , . . . , yi−1 , 0, yi+1 . . . , yN }, (3.4) 19 the transition kernels satisfy    p(z), x = (y1 , . . . , yi−1 , z, yi , . . . , yN ), Q(x, y) = (3.5)   0, otherwise. We thus have a well defined PMDP x(t) = (x1 (t), . . . , xN (t)) on E N . At each time t > 0, we define the random empirical measure process of particle locations: N X δxi (t) νtN (dx) = (dx). i=1 N Let φ ∈ C 1 (R+ ), we can pair an arbitrary finite measure μ ∈ P(R+ ) through inte- gration: Z hμ, φi = φ(x)μ(dx). R+ A natural choice for a functional of our PDMP is the empirical pairing F N (φ) : E → R+ , defined by N X N φ(xi ) F (φ)(x) = . (3.6) i=1 N As a Markov process, the PDMP has an associated generator and martingale determined from Dynkin’s formula. Our aim is to define a non-trivial set of test functions such that F N (φ)(x) ∈ D(A), the domain of valid functionals for the in- finitesimal generator A. We will use the sufficient conditions given in Theorem 2. The major hurdle to show F N (φ) ∈ D(A) is in boundary condition (3) of Theorem 2, which is equivalent to proving   E ΔF N (φ)(y) = 0, y ∈ Γ∗ , (3.7) 20 where, if a random variable Z(y) is distributed in E with respect to Q(∙, y), then ΔF (y) = F (Z(y)) − F (y). (3.8) In order for condition (3) to hold, we must equate the value of a particle at the origin, φ(0), with its expected value when it is redistributed back to R+ . We thus use the space of test functions  Z ∞  1 + 0 C = φ ∈ C (R ) | φ(0) = φ(x)p(x)dx, kφ k∞ < ∞ . (3.9) 0 It’s straightforward to show, using the regularity of φ ∈ C, that our empirical func- tional F N (φ), is in the domain D(A). Vector fields satisfy, by the chain rule, N X X φ0 (xi ) ∂ φ(xi ) N X (F (φ)(x)) = − =− = −hνsN , φ0 i. (3.10) i=1 ∂xi N p ∈T N i k The martingale MφN (t) associated with F N (φ) given in (2.9) and (2.10), then satisfies Z t MφN (t) = hνtN , φi − hν0N , φi − hνsN , φ0 ids. (3.11) 0 It is of note that if we allow for initial empirical data to converge to measures with jumps in its cumulative distribution function, then we should not expect con- vergence in the J1 topology. Indeed, if we initial empirical densities densities that approach a point mass, redistribution of particles concentrate during increasingly small time intervals. During this interval, the value F N (φ)(X(t)) forms an approxi- mately continuous “monotone staircase” between φ(0) and E[φ(X)] whose jumps do not approach that of the limiting function (see Fig. 3.2) The right topology to invoke in this case is the Skorokhod M1 topology, a weaker 21 Figure 3.2: M1 convergence of measure valued initial data. Left: Particles (small vertical segments) approximating a point mass (segment with red dot) drift toward the origin. Right: The value of F (φ)(x) before and after particles hit the origin. As the number of particles grows, their jumps form a “monotone staircase” connecting the jump of the limiting functional. version of the J1 topology. Convergence of functions in this case is characterized by measuring the distance between the completed graphs of cadlag functions. In particular, continuous functions can approximate a function with jumps. We refer the reader to [40] for a rigorous explanation of the M 1 topology, as well as its relation to the J1 topology. We withhold a discussion of M 1 convergence of the redistribution equation for now, and hope to address the issue in future work. Thus, in this exposition, our initial empirical measures will approach nonatomic measures. 3.2.1 Convergence of empirical processes We are now ready to state our main theorem. It may be seen as a law of large numbers limit for random empirical densities. Theorem 3. Suppose we have weak convergence of initial empirical measures ν0n → ν0 , where ν0 is nonatomic. Then the random measure process νtN converges to a deterministic limit νt under the Skorokhod topology D([0, t], P(R+ )). The limiting measure satisfies the kinetic equation Z t 0 = hνt , φi − hν0 , φi − hνs , φ0 ids (3.12) 0 for every φ ∈ C. 22 The following two lemmas will help in establishing tightness [32, 24]. Lemma 1. Let T > 0. A sequence νk (t) in D([0, T ], M(R+ )) is tight if and only if the sequence of functions hφ, νk (t)i is tight in D([0, T ], R+ ) for all φ ∈ Cb (R+ ). Theorem 4. (Aldous’ conditions) Let T > 0. A sequence of stochastic processes Xn (t) is tight in D([0, T ], R+ ) if the following two conditions hold: 1. For all rational t ∈ [0, T ) and for all ε > 0, there exists an L > 0 such that sup P(|Xn (t)| > L) ≤ ε. (3.13) n>0 2. For jump times τ1 , . . . , τm(t) , then for all ε > 0, lim lim sup sup P(|Xn ((τi + s) ∧ T ) − Xn (τi )| > ε) = 0. (3.14) r→0 n→∞ s 0: Iδ = {(a, b) ⊂ R+ , b − a < δ}. (3.15) Lemma 2. Let T > 0. Suppose ν0N → ν0 weakly, and ν0 is nonatomic. Then for all η > 0, lim sup lim sup sup P(νtN (I) > η)) = 0. (3.16) δ→0 I∈Iδ n t∈[0,T ] Proof. Let η > 0. For T = 0, we have that (3.16) follows from the fact that ν0 is nonatomic. 23 We now show (3.16) for s ∈ R+ . For S ⊂ R , we denote S+s = {x ∈ R|x−s ∈ S}. Now let T > 0. For t ∈ [0, T ], νtN (I) is equal to the sum of ν0n (I + t) and the sum of masses that are redistributed to I + t − s at time s. To estimate redistribution probabilities, we focus on the location of a single particle x(t) for t ∈ [0, T ]. Let Z x+δ Mδ = sup p(s)ds. (3.17) x∈R+ x Note that because p ∈ L1 (R+ ), Mδ → 0 as δ → 0. Now, the probability of x(t) ∈ I is bounded by the number of jumps that x(t) undergoes, multiplied by the probability that the next jump lands in the correct δ interval, or ∞ X P(x(t) ∈ I) ≤ P(x(t) has k jumps)Mδ (3.18) k=0 ∞ k ! X X ≤ Mδ P Xi < t) , k=0 i=1 where Xi , i ∈ N is a sequence of i.i.d random variables distributed with respect to p(x). A Chernoff-type bound shows that the sum over k converges. Indeed, we have k ! k ! ! X X P Xi < t = P exp − Xi ) > exp(−t) (3.19) i=1 i=1 ≤ et E[e−X1 ]k , k ∈ N, by Markov’s inequality. The series of terms (3.20) is then summable over k to a finite constant, denoted K. As particles x1 (t), . . . , xn (t) behave independently, we now have an estimate for νtN (I), given by E[νtN (I)] ≤ Ket Mδ . 24 But then we have, from Markov’s inequality, sup P(νtN (I) > η) (3.20) t∈[0,T ] E[νtN (I)] ≤ sup P(νtN (I) > η) ≤ sup t∈[0,T ] t∈[0,T ] η KeT Mδ ≤ →0 as δ → 0 η From here the proof for Aldous’ conditions is straightforward. Theorem 5. For φ ∈ C, the Aldous conditions in Theorem 4 hold for hν, φi Proof. Condition 1: This follows trivially since hφ, νtN i ≤ kφk∞ , meaning that L = kφk∞ satisfies Condition 1. Condition 2: We have Z |hφ, ντN+s∧T i − hφ, ντN i| = φ(x)(ντN+s∧T − ντN )(dx) (3.21) R+ Z ∞ ≤ (φ(x − s) − φ(x))ντN (dx) + kφk∞ ντN ([0, s]). s Since φ is continuous, the first term approaches 0 as s → 0, and by Lemma 2, we have that for any ε > 0, lim lim sup sup P(ντN ([0, s]) > ε/kφk∞ ) = 0, (3.22) s→0 N τ ∈[0,T ] which verifies Condition (2), and completes our pro of. 25 We will now take care of the martingale term, which we can show is bounded in the local uniform metric: Theorem 6. The martingale term MNφ in (3.11) is bounded by " # Bet/2 kφk∞ E sup (MNφ (s)) ≤ √ , (3.23) s∈[0,t] N where B is a constant determined by a summable series. Proof. We first note that 2kφk∞ |ΔMNφ (t)| = |MNφ (t) − MNφ (t− )| = |ΔhνtN , φi| ≤ (3.24) N It then follows, for dN (t) denoting jump times τ1 , . . . , τdN (t) , that   h i dN (t) X E MNφ (t)2 = E  ΔMNφ (τi )2  (3.25) i=1 4kφk2∞ ≤ E[dN (t)] N2 Here we again appeal to Chernoff bounds in (3.19) to help estimate the total number of jumps that occur. Let x(t) denote the stochastic process tracking a tagged particle. This gives us ∞ X E[#{jumps for p(t)}] = kP(p(t) jumps k times) (3.26) k=0 k ! ∞ X X ≤ kP Xi ≤ t ≤ et kE[e−X1 ]k ≤ Cet i=0 k=0 for some constant C. This gives us E[dN (t)] ≤ CN et , meaning 4Cet kφk2∞ E[MNφ (t)2 ] ≤ . (3.27) N 26 Now, using Doob’s inequality, we obtain  2   E sup |MNφ (s)| ≤E sup MNφ (s)2 ≤ 4E[MNφ (t)2 ] (3.28) s≤t s≤t √ Setting B = 4 C then gives our desired result. Since we have shown tightness, we can extract a subsequence νtNk → νt in the Sko- rokhod J1 topology. Furthermore, since local uniform convergence is stronger than convergence under the Skorokhod topology, we also have MNφ k (t) → 0 in D([0, T ], R). Thus our limiting measure νt satisfies equation (3.12). We now finish the proof of Theorem 11 by showing uniqueness of limits. This will be done by approximating the measure processes PDMP X N (t) by processes X N,ε (t) = (xN,ε N,ε 1 (t), . . . , xN (t)) with perturbed redistribution probabilities pε (t), given by    0 x ∈ [0, ε] pε (x) = (3.29)    R ∞p(x) x > ε. ε p(y)dy These redistributions are probabilities of a particle redistributed with respect to p, conditioned on being redistributed in [ε, ∞). Suppose now that for time t ∈ [0, ε], we have random processes νtN,ε given by N δ N,ε X x (t) νtN,ε = i . (3.30) i=1 N If for some nonatomic ν ε we have convergence of initial measures ν N,ε → ν ε in the J1 topology, then there is a unique deterministic measure valued limit νtε of νtN,ε . 27 For any open interval A, Z tZ νtε (A) = ν0ε (A + t) + pε (x + t − s)dxdFν0ε (s), (3.31) 0 A where Fν0ε (s) is the cumulative distribution function of ν0ε for time s ∈ [0, ε). That we have such a deterministic formula comes from the fact that when t ∈ [0, ε], particles only redistribute once. Thus, the total amount of redistributed mass at time t is exactly equal to Fν ε (t) . The uniqueness actually extends for all time t ∈ [0, ∞), since we can extend our unique solution to time 2ε by setting initial conditions equal to νεε , and so forth. Let LN N ε (t) and |Lε (t)| denote the support and number, respectively, of particles which at some time hit the origin and are redistributed to [0, ε] in time [0, t]. We now split the process νtN by which particles are in LN ε (t) : N − |LN X δaNi (t) X δai (t) ε (t)| νtN = + . (3.32) N N− |LN ε (t)| N aN / N i (t)∈L ε (t) aN N i (t)∈Lε (t) Now we use the martingale bounds from Theorem 6 to see that Z |LNε (t)| E[dN (t)] ε E[ ]≤ p(x) (3.33) N N 0 Z ε t ≤ Be p(x) → 0 as ε → 0 0 Therefore, the second term in (3.32) approaches 0 in L1 . Notice, however, that particles not in LN ε (t) are distributed with respect to pε . Thus, we have , for identical initial empirical measures ν0N = ν0N,ε , then for any open interval A ⊂ R+ , we have " # E sup |νsN (A) − νsN,ε (A)| → 0 as ε → 0 (3.34) s∈[0,t] 28 For subsequences n1 (N ) and n2 (N ) with ni (N ) → ∞ as N → ∞ for i = 1, 2, n (N ) n (N ) suppose that we had two converging subsequences νt 1 → νt and νt 2 → ν˜t in P(R+ ). Applying the triangle inequality then yields, " # E sup |νs (A) − ν˜s (A)| (3.35) s∈[0,t] " # " # ≤E sup |νs (A) − νsn1 (N ) (A)| + E νs (A) − νsn2 (N ) (A)| sup |˜ s∈[0,t] s∈[0,t] " # " # +E sup |νsn1 (N ) (A) − νsn1 (N ),ε (A)| + E sup |νsn2 (N ) (A) − νsn2 (N ),ε (A)| s∈[0,t] s∈[0,t] " # +E sup |νsn1 (N ),ε (A) − νsn2 (N ),ε (A)| . s∈[0,t] As n and ε were arbitrary, if we take limits of N → ∞, and then ε → 0 , we see that in fact νs = ν˜s . 3.3 A solution via the Laplace transform The goal for the rest of this chapter is to generalize (3.1) for a larger class of initial conditions, and prove a basic result about a universal attractor. Our hope is that such an analysis will highlight the use of several techniques that may be utilized for more complicated PDMPs. The existence of universal attractors is a sought after problem in coarsening systems, such as Smoluchowski’s coagulation equation [30] and grain boundary coarsening [21]. We begin our study of (3.1) by noting that the rate of redistribution in (3.1) is given by the local quantity u(0, τ ). This, if not already evident, becomes clear when we integrate along characteristics x(τ ) = x0 − τ , giving us the integral form of the 29 redistribution equation: Z τ u(x, τ ) = u0 (x + τ ) + p(x + τ − s)u(0, s)ds. (3.36) 0 Thus, well-posedness of (3.36) ultimately depends on the behavior of u(0, s). Through- out this chapter, we’ll use the Laplace transform Z v(q, τ ) := e−qx u(x, τ )dx. (3.37) R+ We use the identity for the Laplace transform of ux , Z Z −qx e ux (x, τ )dx = ∂x (e−qx u(x, τ )) − ∂x (e−qx )u(x, τ )dx (3.38) R+ Z R+ =q e−qx u(x, τ )dx − u(0, τ ) = qv − u(0, τ ) R+ Denoting L(f ) as the Laplace transform of f , (3.38) then gives us L(uτ − ux ) = v˙ − qv + u(0, τ ), and L(p(x)u(0, τ )) = Pˉ (q)u(0, τ ), where we define Pˉ (q) to be the Laplace transform of p(x). We now have (3.1) in the Laplace variable v, v˙ − qv = (Pˉ (q) − 1)u(0, τ ), (3.39) which may be easily solved for v, giving us Z τ v(q, τ ) = e v0 (q) + (Pˉ (q) − 1) qτ eq(τ −s) u(0, s)ds. (3.40) 0 This solution is similar to (3.36), but here we can use Laplace inversion to extract a formula for u(0, s). As we are taking the Laplace transform of a probability density, 30 we have that |v(q, τ )| < 1 for τ > 0, q ∈ C+ . Thus v(q, τ )e−qτ → 0 as t → ∞. We then can obtain, from the τ → ∞ limit of (3.40), Z 1 v0 (q) = e−qs u(0, s)ds. 1 − Pˉ (q) R+ Written this way, we now have a formula for u(0, τ ) based on the initial data. Notice that the left hand side is a Laplace transform in the spatial variable, whereas the right hand side is a Laplace transform in time of u(0, τ ). We now define the trace α(s) = u(0, s), and the measure describing the total number of redistribution A(ds) := α(s)ds. An application of the Fourier inversion formula now gives us Z 1 1 α(τ ) = eiξτ v0 (iξ)dξ. (3.41) 2π R 1 − Pˉ (iξ) We can apply the convolution theorem for Fourier transforms to obtain Z α(τ ) = K(τ − x)u0 (x)dx, (3.42) R+ where   −1 1 K(x) = F (x), 1 − Pˉ (q) and F −1 (f ) is the Fourier inverse of a function f . The main point here is that we now have a solution formula for α(τ ) based only on initial data, showing that well-posedness of (3.1) is equivalent to (3.41) being well-defined. 31 3.4 A measure valued solution formula We now seek measure valued solutions via a change of variables of the total amount of redistributed mass. For τ denoting the normal time coordinates, we use t to denote the change of variables given by Z τ t= α(s)ds := A(τ ). (3.43) 0 We’ll first reformulate (3.1) in terms of measures. Let Z ∞ Z ∞ Fˉτ (q) = e −qx f (x, τ )dx = e−qx Fτ (dx). (3.44) 0 0 This gives us the ODE Fˉ˙ − qF = (1 − Pˉ (q))α(τ ), (3.45) and solution formula Z τ Fˉτ (q)e−qτ − Fˉ0 (q) = (Pˉ (q) − 1) e−qs α(s)ds (3.46) 0 We have, a priori, that |Fˉτ (q)| ≤ 1 q ∈ C+ , (3.47) thus Fˉτ (q)e−qz → 0 as τ → ∞. We therefore obtain 1 α ˉ (q) = Fˉ0 (q). (3.48) 1 − Pˉ (q) Equation (3.48) is of great importance for the rest of this chapter, as it gives us a simple relation between initial data and our change of scale. 32 Again, through the convolution theorem, we can write the trace as Z ∞ α(s) = K(s − x)F0 (dx) (3.49) 0 Assuming α > 0 (we’ll give conditions for this in Corollary 1), we have that A(τ ) is strictly increasing, and that dt = α(τ ). (3.50) dτ Our ODE (3.45) is transformed as d −qτ ˉ (e F (q)) = (1 − Pˉ (q))(e−qτ α(τ )) (3.51) dτ 1 d −qτ ˉ ⇒ (e F ) = (1 − Pˉ (q))e−qτ . α(z) dτ From the chain rule, this is equivalent to d −qτ (t) ˉ (e F ) = (1 − Pˉ (q))e−qτ (t) . (3.52) dt This in turn gives us the solution formula Z t Fˉt (q)e−qτ (t) − Fˉ0 (q) = (1 − Pˉ (q)) e−qτ (s) ds. (3.53) 0 Since τ (t) → ∞ as t → ∞ , the limit of (3.53) is then Z ∞ Fˉ0 (q) = (1 − Pˉ (q)) e−qτ (s) ds (3.54) 0 and consequently Z ∞ e −qτ (t) Fˉt (q) = (1 − Pˉ (q)) e−qτ (s) ds. (3.55) t The right hand side of (3.55) is continuous in time, and as we shall see, with p(t) 33 that is positive on an interval around the origin, the left hand side is as well for any t > 0. 3.5 An example: a point mass under a stationary distribution The following example illustrates the behavior of jumps in initial data. We introduce new notation for a finite measure μ ∈ M(R+ ) by writing μ(x) := μ([0, x]) if the quantity is finite. Let the initial measure data satisfy F0 = χ[1,∞) (x) and p(x) = e−x . Using the redistribution driven change of time scale, we should expect a solution that gives extra time when a point mass hits the origin. Suppose we have u0 = p(x) = e−x . It is easy to see that the stationary solution u(x, τ ) = e−x satisfies (3.1). Then L(p(x)) = 1 q+1 , Fˉ0 (q) = e−q , and the Laplace transform for the rate function is, using (3.48) q + 1 −q −q e−q α ˉ (q) = e =e + (3.56) q q Laplace inversion then gives us α = δ0 + 1, A(s) = 1 + s (3.57) and    1 t ∈ (0, 1] τ (t) = A−1 (t) = (3.58)   t t > 1. 34 Thus, using (3.55), we then have, for t ∈ (0, 1] Z ∞ ˉ q q Ft (q) = e ( e−qτ (s) ds) (3.59) q+1 t Z ∞ q −q =e (q )((1 − t)e + e−qs ds) q+1 1 q 1 t =( )(1 − t + ) = + (1 − t). q+1 q q+1 Taking inverse Laplace transforms then gives us the solution Ft (dx) = (1−t)δ0 (dx)+ te−x dx for t ∈ (0, 1]. Similarly, we can show that Ft (dx) = e−x dx for t > 1. Thus, our change of scale immediately shifts our point mass to the origin, and then continuously redistributes a mass of one according to p(x) at rate one. After the redistribution of the point mass, the density is stationary, since u(x, τ ) = e−x is a solution of (3.1) with conditions p(x) = u0 (x) = e−x . 3.6 Strong solutions In this section, we will now work in the τ time scale, as in Section 3.3. Definition 1. A function u(x, τ ) is a strong solution to (3.1) if u(x, τ ) satisfies the integral form (3.36) and u(x, τ ) ∈ C(R+ × R+ ). Theorem 7. Let p(x) ∈ L1 (R+ ) be a probability density, and let u0 (x) be continuous. Then there exists a unique strong solution u(x, t) with u(x, 0) = u0 (x). Proof. This theorem follows immediately after showing that α(s) is integrable, since (3.42) asserts that α(s) is uniquely determined by the initial data. Since p(τ ) is a probability distribution |Pˉ (q)| < 1 for q > 0, which allows us to 35 express (3.48) as an infinite series: ∞ X α ˉ (q) = Pˉ k (q)ˉ u0 (q). (3.60) i=1 From here, we observe that from the convolution theorem, this gives us the following elegant solution for α(τ ) ∞ X α(τ ) = p∗k ∗ u0 (τ ), (3.61) k=1 where p∗k denotes k-fold self convolution. We now give a probabilistic argument to show that α(s) is well-defined and continuous. First, note that Z ∞ ∞ i ! τ X X X p∗k (s)ds = P Xj < τ , (3.62) 0 i=1 i=1 j=1 where Xj are iid random variables distributed with probability density p(x). In the sequel, we’ll define the renewal operator ∞ X Qp (τ ) = p∗k (τ ). (3.63) i=1 Using Chernoff bounds, we then have i ! X P Xj < τ ≤ eτ E[e−τ X1 ]i . (3.64) j=1 For any random variable X which is non-trivial (a point mass at 0), it’s easy to show that β(τ ) := E(e−τ X ) < 1 for τ ≥ 0. We thus have Z τ ∞ X Qp (τ )ds ≤ eτ β(τ )i < ∞. (3.65) 0 i=1 36 Since Qp (τ ) ∈ L1 ([[0, τ ]) and u0 is continuous, it follows from (3.61) that α(τ ) is also continuous. The existence and uniqueness of (3.1) then follows immediately from the solution formula (3.36). 3.7 Weak solutions That (3.1) has unique solutions for continuous initial data, and also behaves well with jumps under a suitable time scale suggests the existence of measure valued solutions. We now define weak solutions to the redistribution equation, where (3 .1) is tested against the space Cc1 (R+ ) = {b ∈ C 1 (R+ )|b(x) → 0 as x → ∞}. Definition 2. A weak solution to the redistribution equation is a map Ft : R+ → M(R+ ) where for every b ∈ Cc1 (R+ ), we have the following, R (1) The map t 7→ R+ b(x)dFt (x) is measurable. (2) Ft satisfies the following weak form of (3.1): Z ∞ Z ∞ b(x)dFt (x) − b(x)dF0 (x) = (3.66) 0 0 Z tZ ∞ Z ∞ 0 A(t)b(0) − b (x)dFs (x)ds + A(t) b(x)dP (x), 0 0 0 37 where A(t) satisfies Z tZ ∞ A(t) = K(s − x)dF0 (x)ds, (3.67) 0 0 1 K(x) = F −1 ( )(x) 1 − Pˉ (q) Our statement is then the following. Theorem 8. Suppose Fˆ ∈ M(R+ ), and P (x) is a probability measure. Assume further that Fˆ has polynomial growth, i.e. that Fˆ (x) . xc for some c < ∞. Then there exists a unique weak solution of the redistribution equation with F0 = Fˆ . If Fˆ is a probability measure, Ft is also a probability measure for all t > 0. Proof. From the Weierstrass approximation theorem, finite linear combinations of trigonometric polynomials are dense in Cc∞ (R+ ). Thus we can consider test functions b(x) = e−qx , for q > 0. This means that a weak solution exists if its Laplace transform satisfies the measure valued solution formula (3.55), which is (3.66) for b(x) = e−qx . To see this, substituting gives us Z t Fˉt (q) = Fˉ0 (q) + A(t) − q Fˉs (q)ds + A(t)Pˉ (q). (3.68) 0 As all the members of the right hand side of (3.68) are differentiable in t, so is Fˉt (q), and thus by differentiation we obtain (3.40), and therefore (3.45). We approximate Fˆ and P weakly in the space of measures by a sequence of strictly positive continuous (n) densities u0 (x) and p(n) (x), meaning that for all φ ∈ Cb (R+ ), Z Z (n) φ(x)u0 (x)dx → φ(x)dF (x), (3.69) ZR ZR + + φ(x)p(n) (x)dx → φ(x)dP (x). (3.70) R+ R+ 38 By Theorem 7, there exist strong solutions u(n) (x, τ ) of (3.1) with redistribution densities p(n) (x) that are positive as well. As shown in section 3.4, denoting Z ∞ Fˉτ(n) (q) = e−qx u(n) (x, τ )dx (3.71) 0 we then obtain Z ∞ (n) (n) (n) (s) Fˉt (q) = eqτ (t) (1 − Pˉ (n) (q)) e−qτ ds. (3.72) t (n) ˉ As n → ∞, uˉ0 (q) → Fˆ (q) and Pˉ (n) (q) → Pˉ (q) pointwise for q > 0. We also have Z t Z t (n) n A (t) = Qp(n) (s) ∗ u0 (s)ds → Q(s)dFˆ (s) = A(t). (3.73) 0 0 As shown in [29], pointwise convergence of cumulative density functions An (t) implies the pointwise convergence of CDF inverses (A(n) )−1 (t) → A−1 (t). To show converge of (3.72), we need to use the dominated convergence theorem with respect to terms involving τ (n) (t). Thus, we need to estimate τ (n) (t) = A−1(n) (t), using a more careful treatment of Chernoff’s bound. First assume p(x) has a finite moment μ. Chernoff’s bound [8] states that for any δ > 0, we can exponentially bound the tail of the sum of iid random variables by n ! X −δ 2 nμ P Xi < (1 − δ)nμ ≤e 2+δ . (3.74) i=1 Let M = b 2t μ c. Then we have for n > M , n ! n ! X X t P Xi < t =P Xi < nμ ≤ e−nμ/10 . (3.75) i=1 i=1 nμ A distribution with infinite moment can be approximated by random variables 39 (n) (n) Xi , i = 1 ∈ N, with moments μn = n. Comparing Chernoff bounds for Xi produce similar (in fact, better) bounds for tails of the sum of Xi , i ∈ N . We estimate the renewal operator by ∞ ∞ i ! X 2t X X QP (s) = P ∗k (s) ≤ + P Xj < t (3.76) i=1 μ i=M j=1 ∞ 2t X −iμ/10 2t ≤ + e ≤ + β(P ), μ i=M μ where β(P ) is some finite constant that depends only on P . This gives us a polyno- mial estimate Z t Z tZ A(t) = α(s) = Q(r)dF (s − r)ds (3.77) 0 0 R Z tZ 2r ≤ ( + β(p))dF (s − r)ds + μ Z t0 R Z 2 . c β(p)r + rdF (s − r)ds 0 μ R+ . tc+1 . It then follows A−1 (t) & t1/c+1 . (3.78) Now we can use the dominated convergence theorem along with (3.78) to give us Z ∞ lim Fˉtn (q) =e qτ (t) (1 − P (q)) e−qτ (s) ds. (3.79) n→∞ t This means that for t > 0, the measures Ftn (q) converge weakly to a measure Ft (q) that satisfies (3.55). Uniqueness follows from the uniqueness of τ (t). This follows from the uniqueness of α(t), which is evident from (3.48). We know that F0 is a probability measure, and thus satisfies Fˉ0 (0+ ) = 1. But 40 then we must also have that Fˉt (0+ ) = 1, for t > 0, since Z ∞ t Z  lim+ Fˉt = lim+ eqτ (t) (1 − Pˉ (q)) e −qτ (s) ds − e−qτ (s) ds (3.80) q→0 q→0 0 Z t 0 = Fˉ0 (0+ ) − lim+ eqτ (t) (1 − Pˉ (q)) e−qτ (s) ds = 1, q→0 0 since Pˉ (0+ ) = 1. Before we proceed with some corollaries of Theorem 8, we first show a definition and result from Tauberian theory that will be useful. First, we recall the notion of a slowly varying function, which describes those functions that are asymptotically “flat”. Definition 3. A slowly varying function at infinity L(x) : R+ → R satisfies L(xt) lim →1 for all x ∈ R+ . (3.81) t→∞ L(x) Likewise, if the above limit holds for t → 0, we say L(x) is slowly varying at 0. We now state the Hardy-Littlewood-Karamata Tauberian theorem [13] Lemma 3. Let μ ∈ M(R+ )If L is slowly varying at infinity and 0 ≤ β < ∞, then the following are equivalent: μ(x) ∼ xβ L(x) as x → ∞ and ˉ(q) ∼ q −β L(1/q)Γ(1 + β) μ as q → 0. (3.82) We proceed with some basic results: 41 Corollary 1. Let P ∈ P(R+ ) and X be a random variable distributed with respect to P . Then the following statements hold: (1) Suppose there exists ε > 0 such that P (x) > 0 for x ∈ (0, ε). Then α(t) > 0, and the measure valued solution Ft is continuous in t with respect to the weak topology of M(R+ ). (2) If E[X] = μ, then we have that A(t) ∼ t/μ as t → ∞. Proof. (1): It’s easy to show that given the assumptions for (1), we have that Q(t) > 0 for any t > 0. It then follows easily that α(t) is strictly positive, and thus A(t), and therefore A−1 is strictly increasing. From (3.55) we see that jumps in F (t) only occur when A(t) is constant, so the result follows immediately. (2) Recalling (3.54), we have Z ∞ 1 = lim Fˉ0 (q) = lim (1 − P (q)) e−qτ dA(τ ). (3.83) q→0 q→0 0 ˉ ∼ 1/μq, and so by Lemma Since (1 − P (q)) ∼ μq as q → 0, we must have that A(q) 3, A(x) ∼ q/μ. 3.8 Stationarity We now examine the convergence of initial data to a stationary solution under ap- propriate rescaling. Observe that for any probability distribution p(x), regardless of 42 moment assumptions, there is a stationary solution of (1): Z ∞ u(x, t) = u(x) = c0 ( p(y))dy), (3.84) x c0 = u0 (0). (3.85) The issue here is that for densities with no first moment, the associated stationary solution is no longer a probability density. For instance, the Cauchy density u(x) = π 2(1+x2 ) gives the stationary solution u(x) = π/2 − tan−1 (x), which is not in L1 (R+ ). The reader should also note that we shouldn’t expect solutions to converge to a stationary solution for an arbitrary redistribution measure P (x) ∈ P(R+ ). As an example, suppose we have initial data u0 (x) = 1[0,1] (x) and a redistribution measure P (x) = 1x>2 that sends all particles at the origin to x = 2. Then a solution for u(x, τ ) is then a traveling square density that is cyclic with time. Similar examples exist for any P that is an arithmetic distribution, where the support of P lies in a lattice. As we will see in the next theorem, for a large class of initial densities arithmetic distributions are the only possible instances where mixing doesn’t occur. The following result from renewal theory will help us with asymptotics (see [35]). It uses the notion of direct Riemann integrability (DRI). A function is DRI if its lower and upper Riemann sums over all of R+ converge as the mesh size of the partition approaches zero. This differs slightly from the usual Riemann integral definition, which considers Riemann sums over finite intervals [0, t], and then takes a limit as t → ∞. Theorem 9. Key Renewal Theorem. For a non-arithmetic probability measure P ∈ P(R+ ) with mean μ and u(t) that is DRI, the renewal operator QP (t) of P 43 satisfies Z 1 lim QP ∗ u(t) = u(s)ds. (3.86) t→∞ μ R+ Because of the linearity of (3.1), we can generalize slightly and prove a stability theorem for initial conditions that are Riemann integrable, or whose lower and upper Riemann sums converge as the size of the mesh approaches zero. Theorem 10. Let P ∈ P(R+ ) be nonarithmetic with finite mean μ and u0 be Rie- mann integrable probability density on R+ . Suppose Ft ∈ M(R+ ) is a weak solution of (3.1) with initial condition F ∈ M(R+ ). We then have 1 Ft → (1 − P (x)) in distribution as t → ∞. (3.87) μ Proof. Decompose the initial condition as u = u1n + u2n , where u1n = u1x≤n and u2n = u1x>n . As u1n has compact support, it is clear that it is DRI. Let Ftn,i , i = 1, 2 be the weak solution with initial data uin . We use the solution formula (3.55) and examine the limit Z ∞ lim Fˉtn,1 (q) = lim eqτ (t) e−qτ (s) ds(1 − Pˉ (q)). (3.88) t→∞ t→∞ t From the key renewal theorem and part (2) of Corollary 1, we obtain μ α(s) → R n := μn as s → ∞. (3.89) 0 u(x)dx We then change variables in (3.88) to obtain Z ∞ (1 − Pˉ (q)) lim Fˉtn,1 (q) = lim e qτ (t) e−qs α(s)ds(1 − Pˉ (q)) = . (3.90) t→∞ t→∞ A(t) μn q 44 By linearity, we can express the weak solution as Ftn = Ftn,1 (q) + Ftn,2 (q). Since (3.1) preserves total particle number, it’s evident that Ftn,2 (q) → 0 uniformly as n → ∞. However, we also have that the inverse Laplace transform of (3.90) is 1 (1 − P (x)), (3.91) μn which approaches (3.87) as n → ∞. Remark 1. A direct generalization of Theorem 10 would be to consider redistribu- tions p(x) with infinite means. Positive results would likely rely on a key renewal theorem for infinite mean random variables. Several papers have already addressed this question ([11, 1], for instance) for heavy-tailed power laws with a slowly varying factor. Chapter Four Kinetic Limits of Piecewise-Deterministic Markov Processes: The k-Species Model 46 This chapter is organized as follows. In Section 4.2, we describe dynamics involved for the k-species PDMP, including the weights and rules for determining transition probabilities. In Section 4.3, we show that the k-species PDMP is, in fact, a PDMP by defining the parameters of a general PDMP to align with the model. The main theorem involving a fluid limit for empirical particle densities is then established. Unfortunately, while the PDMP described in Section 4.3 is perfectly well-defined, the set of non-trivial functionals in the domain of its associated generator is too small to describe individual species densities. Section 4.4 addresses this issue by adding variables to the original k-species PDMP. These extra variables do not interact with the dynamics of the original stochastic process, but rather are simply meant to track certain types of jumping events that occur. With these added dimensions, we are able to produce martingales, via Dynkin’s formula, that serve as estimates of the fluid limit. The main technical details of existence of limiting densities are provided in Sec- tion 4.6. The main idea is to differentiate particles according to their paths of visiting differing species, and show empirical measures of certain paths satisfy tightness con- ditions. Section 4.7 gives a law of large numbers argument to establish convergence of terms in the martingale equation. Finally, Section 4.8 shows that existence can be established on a reasonable time interval in the sense that there solutions are well defined as each species has positive total number. 4.1 Notation Symbol Description A˜j,k N (t) Boundary running measure of subparticle process AjN (t) 47 A˜j,k int,N (t) Interior running measure of subparticle process AjN (t) β Poisson parameter for interior events φ,N Bj,k (t) Tracking dimension of particles reassignments from Lj to Lk Cb (R+ ) The space of bounded continuous functions on R+ C The space of test functions ΔX(t) The infinitesimal change of a stochastic process at time t Fl (t) total number deleted at species l gk (t) Total total number of species k limit gk0 Total initial total number of species k limit gkN (t) Total total number of species k of X N (t) G(s) Total total number of species k limit Gint (s) Weighted interior total total number of species k limit H N (t) The total deleted total number of X N (t) K (l) Number of particles reassigned from a boundary event K int Number of particles reassigned from an interior event K maxl K (l) ∨ K int κ(t) The total number of particle reassignments before time t κj,k (t) The total number of particle reassignments that send a par- ticle from Lj to Lk before time t Li Species type i |Li | Number of particles of species type i m(t) The total number of critical events before time t mint (t) The total number of interior events before time t M The number of different species (M, d) A metric space M with distance d M− The number of species with negative velocities 48 M0 The number of species with zero velocities M+ The number of species with positive velocities μN k (t) Empirical measures of species k of X N (t) μN Σ (t) Path measures of X N (t) for path Σ μk (t) Limiting empirical measures of species k μ0k Limiting initial empirical measures of species k N (t) Number of particles at time t N (l) Species selection numbers for boundary events N int Species selection numbers for interior events (l) pi Species selection probabilities for boundary events pint i Species selection probabilities for interior events ϕi (x, t) Flows corresponding to vi (x) P(x(t)) Path of a particle (l) Rij Reassignments for boundary events (l) Rij Reassignments for interior events Si (x) Support of particles in species i Σ(x(t)) Path location of a particle Te Length of existence interval τi (t) The ith jump time of a critical event. τij,k The time of the ith occurrence of a boundary event reas- signment that sends a particle from Lj to Lk . Θ(x(t)) Jump description of a particle uk (x, t) Smooth density limit of species k vi (x) The velocity at x ∈ Li (l) wi Species selection weights for boundary events wiint Species selection weights for interior events 49 (l) Wj Species selection weight fraction for boundary events of species k limit X N (t) The k-species PDMP with N initial particles ˜ N (t) X The k-species extended PDMP with N initial particles 4.2 A kinetic mixing model 4.2.1 The model Our PDMP for particle evolution tracks N (t) particles X(t) = (x1 (t), . . . xN (t)), with N (0) = N , on disjoint species Li = R+ i , i = 1, . . . , M . A species Li is equipped with a vector field vi (x) ∈ C 1 (R+ ), where vi is strictly positive, strictly negative, or identically zero. The vectors fields are also assumed to define flows ϕi (x, t) with no finite time blowup. We write M = M− + M0 + M+ , representing species with positive, stationary, and negative velocity, respectively. Specifically for i = 1, . . . , M1 we equip Li with vector fields vi < 0, while for i = M− + 1, . . . , M− + M0 , vi ≡ 0, and for i = M− + M0 + 1, . . . , M , vi > 0. Particles move according to ϕi (x, t) until one of the following two critical events occur: 1. Boundary event: A particle xj (t) ∈ R+ i , i ∈ 1, . . . M− hits the origin. 2. Interior event: A Poisson clock with a parameter rings. When a critical event occurs, particles can either be instantly deleted or transferred to another species. However, no spatial jumps occur, meaning that critical events do 50 not change the position of particles. A particle xi (t) ∈ Lj that reaches 0 is instantly deleted. We use mean field probabilities to determine where particles are reassigned at a critical event. At a boundary event triggered by a particle hitting the origin of Ll , (l) K (l) species are randomly selected according to constants wi ≥ 0, i = 1, . . . , M , (l) (l) and a probability vector p(l) (w1 , . . . , wM ), where (l) M X (l) w |Li | (l) pi = i (l) , N (l) = wj |Lj | i = 1, . . . , M. (4.1) N j=1 (l) PM (l) Note that pi ≥ 0, and i=1 pi = 1. When a species Li is chosen, one particle is then selected with equal probability 1/|Li | from particles in Li . The selected particles x1 , . . . , xK (l) will be reassigned to a new species according to a reassignment (l) Rij ∈ {1, . . . , M }, with i = 1, . . . , M, j = 1, . . . , K (l) . Explicitly, if a particle xj is in Li , then after its jump, it is assigned to LR(l) . We stress the importance of the j ij (l) parameter in Rij , since it means that the order a particle is chosen at a critical event can affect its reassignment. For interior events, we select K int particles according to a probability vector pint = (w1int , . . . , wM int ), with wiint ≥ 0 and M X (l) wiint |Li | pi = , N int = wiint |Li |. (4.2) N int i=1 int Particles are reassigned according to Rij ∈ {1, . . . , M }. After jumps occur, particles continue to drift along flow lines as before until another critical event occurs. See Fig. 4.1 for an example on four species. 51 Figure 4.1: A PDMP on four species. Particles travel on four separate copies of R+ . Velocity directions are represented by horizontal arrows, with species 3 having zero velocity. A boundary event occurs when a particle (labelled by “A”) hits the origin. Here, three particles are then randomly selected (K (2) = 3), and reassigned to different species by predetermined reassignments (2) (2) (2) (given by vertical arrows). In this example, R41 = 1, R32 = 2, and R23 = 3. 4.3 The k-species model in a PDMP framework 4.3.1 PDMP formulation of the k-species model We now show the k-species model is a specific class of PDMP, as outlined in Section 2.1. Define the possible listings of species [ S= {1, . . . , M }i , (4.3) i∈N For a species index s = (s1 , . . . , s|s| ) of length |s|, a state has the form (x1 , . . . , x|s| ) ∈ (R+ )|s| . We write the state space of particles as a  E= ((R+ )|s| )s = (s, (x1 , . . . , x|s| ) : s ∈ S, (x1 , . . . , x|s| ) ∈ (R+ )|s| . (4.4) ~s∈S A state x ∈ E can thus be described as an element of (R+ )|s| × N|s| , where we will write  x = (s1 , x1 ), . . . , (s|s| , x|s| ) . (4.5) 52 Vector fields will then be represented by N X ∂ X~s = vsi (xi ) . (4.6) i=1 ∂xi The exit boundary of Γ ∗ ⊂ E consists of particles that hit the origin. Γ∗ = {x ∈ E|there exists (si , xi ) where xi = 0, si ≤ M− } (4.7) To describe the transition kernel Q, we denote the set of states that x ∈ Γ∗ can jump to as Eb (x). This set of events is finite, and have well-defined (albeit, tedious to calculate) probabilities of a reassignment for a state x jumping to xb , denoted pb (xb , x). These probabilities are uniquely determined from the weights (l) (l) wk , numbers |Li |, and reassignments Rij , and therefore are functions of the current state x. Similarly, we can define jump probabilities for interior events as pint (xint , x). Here, xint ∈ Eint (x) denotes the set of possible states that x ∈ E can jump to. Thus, we can define      pb (y, x) x ∈ Γ∗ , y ∈ Eb (x)    Q(y, x) = pint (y, x) x ∈ E, y ∈ Eint (x) (4.8)       0 otherwise. Since each particle has a Poisson clock of parameter β, the distribution Y of the first time that a clock rings follows the distribution Y ∼ min P oisson(β) = P oisson(N β). (4.9) 1≤i≤N 53 Thus λ(x) : E → R takes the form λ(x) = β|s|. (4.10) With this, we now have a well defined PDMP which describes particle reassignments among different species. 4.3.2 The main theorem For x ∈ E, we define empirical measures μN + i (x) ∈ M(R ) by X δ xj μN i (x) = , i = 1, . . . , M. (4.11) x ∈L N j i For a process X(t), we will often write μN N i (X(t)) := μi (t). For the class of test functions C = {φ ∈ C 1 (R+ ) ∩ Cb (R+ ) : φ(0) = 0, φ0 ∈ Cb (R+ )}, (4.12) we can pair empirical measures with test functions through integration: Z X φ(xj ) hφ, μN i i = φ(x)μN i (dx) = , φ ∈ C, i = 1, . . . , M, (4.13) R+ x ∈S N j i where Si denotes the support of particles in species i, for i = 1, . . . , M . Our main result shows that empirical measures μN k (t) approximate a weak solu- tion to a transport equation, with a source term given by the flux of particles from both interior and boundary events. These measures are simplest to describe when 54 viewed as smooth limits μk (t) of particle densities μN k (t) as N → ∞ of the form uk (x, t)dx = μk (t)(dx), k = 1, . . . , M. (4.14) We also denote the total numbers Z ∞ gk (t) = uk (x, t)dx, (4.15) 0 M X M X int G(s) = gi (t), G (s) = wkint gk (t), k = 1, . . . , M. i=1 i=1 and species selection weight fractions for boundary events (l) (l) wj Wj (t) =P (l) l = 1, . . . , M− j = 1, . . . , M. (4.16) wm gm (t) The limiting number densities uk (x, t), k = 1, . . . , M , will satisfy the following system of PDE, ∂t uk (x, t) + ∂x (vk (x)uk (x, t)) = (4.17)   XM− K (l) M X X (l) (l) ul (0, t)vl (0)  Wi (t)ui (x, t)1R(l) =k − K (l) Wk (t)uk (x, t) (4.18) ij l=1 i=1 j=1   M K int βG(t) X X + int  wiint ui (x, t)1Rijint =k − K int wkint uk (x, t) (4.19) G (t) i=1 j=1 with initial conditions uk (x, 0) = u0k (x), k = 1, . . . , M. (4.20) The first line (4.17) describes advection of species densities under the velocity field vk . The next line (4.18) describes the the growth and loss of species k due to 55 boundary events. The first term of (4.18) describes rates of incoming particles to (l) species k. This is reflected in the use of indicator functions, where Rij determine species to which particles are reassigned. The next term describes loss from species k to other species. This does not depend on reassignment, so we should expect no (l) Rij terms. The final line (4.19) contain source terms due to interior events. These terms do not depend on ul (0, t), as the rate of interior events is only dependent on total total numbers gi (t). Our main theorem describes the convergence of empirical measures to a weak form of the limiting PDE (4.17)-(4.19), under the assumption of convergence of initial empirical measures to a nonatomic measure. Theorem 11. Let μ0k ∈ M(R+ ) be nonatomic. Suppose that we have weak conver- gence of initial empirical measures μ0,N k → μ0k as N → ∞ in M(R+ ). Then 1. There exists a time Te > 0 such that the μN k (t) converge to limits μk (t) along a subsequence under the Skorohod topology D([0, Te ], M(R+ )). 2. The empirical losses FlN (t) of the total number deleted in species l converge to continuous cumulative distribution functions Fl (t) along a subsequence in the mean L∞ metric:   lim E sup |FlN (s) − Fl (s)| = 0. (4.21) N →∞ s≤Te 56 3. The μk (t) satisfy the following limiting equations for all φ ∈ C : Z t hμk (t), φi − hμ0k , φi + hμk (s), φ0 vk i = (4.22) 0 M− Z t M K (l)  X XX (l) (l)  Wi (s)hμi (s), φi1R(l) =k − K (l) Wk (s)hμk (s), φi dFl (s) ij l=1 0 i=1 j=1 Z   t M K X X int βG(s)  + int −K int wkint hμk (s), φi + wiint hμi (s), φi1Rijint =k  ds, 0 G (s) i=1 j=1 for k = 1, . . . , M . Remark 2. Since individual jumps of the prelimit PDMP converge to zero as N → ∞, the limiting functions hμk (t), φi, k = 1, . . . , M are continuous in t, meaning that empirical measures also converge in the local uniform metric (see Prop. 1.17 in Chapter VI of [23]). As we will see in Section 4.7, the functions Fl (s) can be seen as an integral of the trace of measures μl (s) occurring at the boundary. In fact, for continuous solutions, we will show that dFl (s) = vl (0)u(0, s)ds. We also state a well-posedness theorem for L1 ∩ L∞ initial data, which is proved in Section 4.9. Definition 4. For k = 1, . . . , M , let uk (x, t) ∈ C([0, Te ], L1 ∩ L∞ (R+ , R+ )). We call uk (x, t) a weak solution in L1 ∩ L∞ (R+ , R+ ) of (4.17) with initial conditions u0k (x) ∈ L1 ∩ L∞ (R+ , R+ ), if for all φ ∈ C, there exist μk (t), t ∈ [0, Te ], that satisfy (4.22), with μk (t)(dx) = uk (x, t)dx, and μk (0)(dx) = u0k (x)dx. Theorem 12. For k = 1, . . . , M , there exists a unique weak solution uk (x, t) in L1 ∩ L∞ (R+ , R+ ) of (4.17) with initial conditions u0k (x) ∈ L1 ∩ L∞ (R+ , R+ ). Remark 3. Theorem 12 should not be mistaken as stating that empirical measures 57 with L1 ∩L∞ initial data converge in the uniform metric to a unique L1 ∩L∞ solution of (4.22). Rather, in the class of all possible solutions, there is one, and only one, solution that is in C([0, Te ], L1 ∩ L∞ (R+ )) Well-posedness of (4.22) for L1 data is discussed in Section 4.9. 4.3.3 Method of proof As a Markov process, the PDMP X(t)has an associated infinitesimal generator and martingale, given by 2.9 and 2.10, along with Rolski’s conditions for membership in the domain D(A)), as shown in Theorem 2. The method for proving Theorem 11 relies on constructing martingale equations that are approximations of (4.22). Unfortunately, the domain of functionals for the infinitesimal generator of X(t) is very limited. A na¨ıve approach for proving Theorem 11 might take empirical pairing functionals fkφ (x) = hφ(x), μN k i, k = 1, . . . , M. (4.23) This, however, does not satisfy the boundary condition (3) of Theorem 2 for any nontrivial class of test functions. To see this, notice that if we restrict attention to a particular species Lk , then we shouldn’t expect Δfkφ (x) to retain its expected value when a particle is deleted from or added to Lk . For an example, refer to Fig. 4.2. At the boundary event shown, for E[Δf1φ (x)] = 0 we need φ(0) = hμN N 2 , φi. But μ2 will change in time, so imposing such a restriction for any non-zero φ is unreasonable. To allow for functionals which track empirical measures on Lk , in Section 4.4 we ˜ will construct an extended PDMP X(t) from X(t) with added dimensions that track boundary events. We can then select functionals Fkφ (˜ ˜ of the x) in the domain D(A) 58 (1) Figure 4.2: A PDMP on two species. In this case K (1) = 1, with R21 = 1. ˜ These functionals satisfy condition (3) by adding a correction extended generator A. term to fkφ (x), dependent on the extra dimensions in x ˜. For each species k, we will use Fkφ (˜ x) to construct the martingales from (2.10). φ,N These martingales involve empirical pairings hμN , φi and terms Bi,j (t), which roughly describe total occurrences of boundary events. We show the limit of these equations approaches (4.22) as N → ∞. To do so, we first show the martingale term Mφk,N (t) converges to zero in the local uniform metric, a consequence of Doob’s theorem. The most technical part involves the tightness of empirical pairings . As in Chapter 2, we will use the Theorems 4 and 1 to show tightness of μk . φ,N Finally, we focus our attention on the terms Bi,j (t). Using a law of large numbers type argument, we will show these terms converge in the local uniform metric to terms involving hμk , φi and Fk (t). 4.4 ˜ N (t) The extended PDMP X To address condition (3) of Theorem 2, let us denote τij,k as the ith time that a boundary event reassigns a particle from Lj to Lk . Note that τij,k may equal τm j,k when m > i, since boundary events can reassign more than one particle. Ordering of τij,k when several particles are simultaneously reassigned from Lj to Lk is done with (l) respect to the same ordering used to determine reassignments Rjk . Suppose we have 59 κN j,k (t) particles (x1 , . . . , xκN j,k (t) ) that have been reassigned from Lj to Lk , at times τij,k . We now define boundary variables that track jumps from boundary events as κN j,k (t) φ,N X φ(xi ) Bj,k (t) = φ ∈ C. (4.24) i=1 N ˜ N (t) which tracks jumps occurring in We now define the extended PDMP X X N (t). Formally, we can write the extended state space as a n o |s|+M 2 |s|+M 2 E˜ = (R+ )s ˜ ) : s ∈ K, x = (s, x ˜∈ R+ . (4.25) ~s∈K ˜ ∈ E˜ can be written as A state x  ˜ = x1 , . . . , x|s| , b11 , . . . bM x 1 M 1 , . . . , b M− , . . . , b M− . (4.26) The quantities bji are a means of recording the number and type of boundary events ˜ which means which occur. We impose that boundary variables have no drift in E, that for vector fields in X˜ , we let X˜~s = X~s , and define the exit boundary as  ˜ ∗ (˜ x) = (y, b11 , . . . , bM ∗ Γ M− )|y ∈ Γ . (4.27) ˜ ˜ and x is seen in the extended transition kernel Q. The main difference between x For x ˜ the set of possible states E˜int (˜ ˜ ∈ E, ˜ may jump to from interior events x) that x is simply the the set Eint (x) with added unchanged boundary variables:  E˜int (˜ x) = (y, b11 , . . . , bM M− )|y ∈ Eint (x) . (4.28) ˜ ∗ (E), ˜∈Γ For x ˜ however, boundary variables will increase, corresponding to the types ˜ . Suppose a boundary event, triggered from a of boundary events that occur at x 60 particle hitting the origin of Ll , selects particles x1 , . . . , xK (l) . We then define the ˜ as new possible states of x n o x) = (y, ˜b11 , . . . , ˜bM E˜b (˜ M− )|y ∈ E b (x) , (4.29) where boundary variables change as X X M ˜bj = bj + φ(xk ) i i 1R(l) =j , i, j = 1, . . . , M. (4.30) x ∈L j=1 N ik k i Thus, boundary variables add the sum of quantities φ(xk ) to bji , where xk is reassigned from Li to Lj . We can then express the extended transition probability as      ˜∈Γ ˜) x pb (y, x ˜ ∗ , y ∈ E˜b (x)    ˜ ˜ ) = pi (y, x Q(y, x ˜ ∈ E, y ∈ E˜int (x) ˜) x (4.31)       0 otherwise. ˜ We thus have a well defined PDMP X(t) ˜ X˜ , and Q. determined from E, ˜ From here it is not hard to see that   ˜ φ,N φ,N φ,N φ,N X(t) = x1 (t), . . . , x|s| (t), B1,1 (t), . . . B1,M (t), . . . , BM− ,1 (t), . . . , BM− ,M (t) (4.32) φ,N With the parameters Bj,k (t) at our disposal, we may now define functionals X φ(x(t)) X M ˜ Fkφ,N (X(t)) = + φ,N Bj,k φ,N (t) − Bk,j (t), k = 1, . . . , M. (4.33) N j=1 x(t)∈Lk 61 4.5 Reassignment bounds and martingale equa- tions To show that Fkφ,N is a valid functional in the domain of the infinitesimal generator ˜ we give a simple bound for jumps which will be utilized throughout this paper. A, Denote the maximum number of reassigned particles at a critical event by K = K int ∨ max K (i) . (4.34) i Lemma 4. The number of particle interior events mN N int and reassignments κ (t) has bounds 1. (Expectation for fixed N ) For any N ∈ N, we have E[mN int (t)] ≤ βN t, (4.35) and E[κN (t)] ≤ KN (1 + βt). (4.36) 2. (Law of large numbers bound) mN int (t) lim sup ≤ βt a.s., (4.37) N →∞ N and κN (t) lim sup ≤ K(1 + βt) a.s.. (4.38) N →∞ N Proof. Clearly, the number of boundary events is bounded by the initial number of particles N . The expected rate of interior events mN int (t) can be bounded by the 62 Poisson parameter β(X N (t)) = β(N (t)) ≤ βN . The expected number of jumps up to time t is thus bounded by βN t, which proves (1). For part (2), we note that a Poisson process Pγ (t) of intensity γ at time t satisfies the scaling argument Pγ (t) ∼ Paγ (t/a) a > 0. (4.39) It then follows from the law of large numbers for Poisson processes, Pγ (t) →γ a.s. as t → ∞, (4.40) t that we can obtain   mN int (t) PN β (t) P(lim sup > βt) ≤ P lim sup > βt (4.41) N →∞ N N →∞ N   Pβ (N t) = P lim sup t > βt = 0, N →∞ Nt which proves part (2). ˜ corresponding to Theorem 13. For the domain of the infinitesimal generator D(A) ˜ for k = 1, . . . , M . x˜, we have Fkφ,N ∈ D(A), Proof. We need to show the four conditions of Theorem 2 hold for Fkφ,N . Conditions (1) and (2) are clear, since particles in E˜ transport along a continuous flow, and φ ∈ C 1 (R+ ). ˜ For condition (4), we need to make estimates on the frequency of jumps in E. We first note 63 h i kφk∞ φ,N E Bj,k (t) ≤ E[κN (t)] . (4.42) N Let us define the the mN (t) jump times before time t as τ1 , . . . , τmN (t) . Then from Lemma 4,   mN (t) X E ˜ i ) − F φ,N (X(τ |Fkφ,N (X(τ ˜ i−1 ))| (4.43) k i=1 N (t) h mX ≤E |hμN N k (τi ), φi − hμk (τi−1 ), φi| i=1 M X i φ,N φ,N φ,N φ,N + |Bj,k (τi ) − Bk,j (τi ) + Bj,k (τi−1 ) − Bk,j (τi−1 )| j=1 4M 4M K ≤ E[κN (t)] ≤ (βtN + N ) N N < ∞. The nontrivial part to show is the boundary condition (3). Note that for X(τ − ) ∈ Γ∗ , Z Fkφ,N (y) − F (X(t− ))Q(dy, X(t− )) = E(ΔFkφ,N (X(t))) (4.44) E   X h φ,N φ,N i = E ΔhμN k (t), φi + E ΔB j,k (t) − ΔB k,j (t) . j Assume that a boundary event at jumping time τ occurs from a particle hitting the origin in Ll , and randomly selects K (l) particles (x1 , . . . , xK (l) ). Then we have the expected change in empirical measure as a sum of the averages resulting from 64 particles leaving and entering Lk :   φ(x1 ) (l) E[ΔhμN k (τ, φ)i] = ∙E K (l) pk (τ )|x1 ∈ Lk (4.45) N XM XK (l)   (l) φ(x1 ) + pi (τ ) ∙ E |x1 ∈ Li 1R(l) =k i=1 j=1 N ij (l) M K (l) X X (l) hμN (τ ), φi K (l) wk N i =− l hμ k (τ ), φi + wi 1R(l) =k . N i=1 j=1 N (l) ij However, we see that for jumps in the boundary variables, by similar calculations, at a boundary jumping time τ X h i M X X K (l) N φ,N (l) hμi (τ ), φi E ΔBj,k (τ ) = wi 1R(l) =k (4.46) j i=1 m=1 N (l) im for expectations of particles entering species k, and X h i   φ,N (l) φ(x1 ) E ΔBk,j (τ ) = −K (l) pk (τ )E |x1 ∈ Lk (4.47) j N (l) K (l) wk = hμN k (τ ), φi Nl for particles leaving species k. In the following, define the empirical total number as M X M X N G (s) = hμN k (τ ), 1i, G int,N (s) = wkint hμN k (τ ), 1i (4.48) k=1 k=1 To each Fkφ,N , we can now construct the following martingale. 65 Theorem 14. For k ∈ {1, . . . , M }, the quantity Mkφ,N (t) = hμN N k (t), φi − hμk (0), φi+ (4.49) XM Z t φ,N φ,N 0 Bj,k (t) − Bk,j (t) − hμN k (s), φ vk i+ j=1 0   M K int N βG (s)  int int N X X K wk hμk (s), φi − wiint hμN  i (s), φi1R(l) =k ds Gint,N (s) i=1 j=1 ij is a martingale. Proof. The vector field X˜ acts only on the first N coordinates of x ˜, and satisfies N X X X ∂ φ(xi ) φ0 (xi ) X˜ (Fkφ,N (˜ x)) = vsi (xi ) = vk (xi ) (4.50) i=1 xi ∈Lk ∂xi N x ∈L N i k 0 = hμN k , φ vk i Using equation (2.9), our formula for the extended generator is then ˜ φ,N (˜ A(F x)) (4.51) k Z = hμN 0 k , φ vk i − β(˜ x) (Fkφ,N (y) − Fkφ,N (˜ x))Q(dy; x˜)) E h i 0 φ,N = hμN k , φ v k i − βN E Δ(F k (˜ x )) . The calculation of the change of the functional over interior events is similar to Theorem 13:   K int pint X X pint M K int E Δ(Fkφ,N (˜ x)) = k hμN k , φi − i hμN i , φi1R(l) =k (4.52) |Lk | i=1 j=1 |L i | ij   M K int 1 X X = K int wkint hμN k (s), φi − wiint hμN  i (s), φi1R(l) =k . Gint,N (s) i=1 j=1 ij 66 Combining terms, we see the martingale Mkφ,N (t) defined from equation (2.10) then coincides with (4.49). We now uniformly bound the martingale Mkφ,N (t), which will decay to zero in the mean L∞ metric. Theorem 15. For t < Te we have √ (4M + 2)Kkφk∞ 1 + βt E[ sup (Mkφ,N (s))] ≤ √ , (4.53) s∈[0,t] N Proof. We bound jumps in Mkφ,N (t) by noting that at most K particles are either reassigned in a critical event: M X (2M + 1)Kkφk∞ |ΔMkφ,N (t)| = |ΔhμN k (t), φi| + |Δ( φ,N Bj,k φ,N (t) + Bk,j (t))| ≤ . j=1 N (4.54) It then follows that, from Lemma 4, that   mN (t) X E[Mkφ,N (t)2 ] = E  ΔMkφ,N (τi )2  (4.55) i=1 (2M + 1)2 K 2 kφk2∞ (2M + 1)2 K 2 kφk2∞ (1 + βt) ≤ E[m N (t)] ≤ . N2 N We now use Doob’s inequality to obtain  2   E sup |Mkφ,N (s)| ≤E sup Mkφ,N (s)2 ≤ 4E[Mkφ,N (t)2 ]. (4.56) s≤t s≤t The result follows from applying (4.55) to (4.56). 67 4.6 Existence of limiting measures We base the existence interval Te > 0 off of when species will have positive total numbers almost surely. The following lemma ensures that the transition probabili- ties, and subsequently the kinetic equations (4.22), are well defined. The technical details of the the following lemma will be utilized in Lemma 6 in determining tight- ness of empirical measures. The existence time Te will eventually be extended to a reasonable time where the solution is well defined as long as each species has positive total particle number. Lemma 5. There exists a time Te where, for all k = 1, . . . , M , and 0 < t < Te , the following hold: 1. (Fraction loss of particles) 2 lim sup sup gkN (s) > gk0 a.s.. (4.57) N →∞ 0≤s≤t 3 P 2. (Total number of particles lost) The quantity H N (t) = 1 − k gkN (t) is bounded by a fraction of the total number in each species: 1 lim sup H N (t) < μ0 (R+ ) for all k ∈ {1, . . . , M } a.s.. (4.58) N →∞ 3M 2 K k 3. (Fixed time bound) The time Te is bounded by a constant determined from fixed parameters in X mink gk0 Te < . (4.59) 4βM 2 K Proof. Denote the flow determined from vector fields vi at time t on Li as ϕi (t, x), with initial condition ϕi (0, x) = x. Regardless of the choice of reassignment of par- 68 ticles, the total particle number lost from boundary deletions H N (t) can be crudely bounded by M X lim sup sup H N (s) ≤ μ0k ([0, a(t)]), (4.60) N →∞ 0≤s≤t i=1 for all N < ∞, where ϕi (s, a(t)) > 0 for all s ≤ t in Li . Since the initial measures μk are nonatomic, we can choose t1 such that M X 1 μ0k ([0, a(t1 )]) < min 2 gi0 , (4.61) k=1 1≤i≤M 3M K which shows part (2). We now show part (1). As noted in Lemma 4, the rate of interior events is bounded by βN , which means that the total particle total number lost through N interior events Hint (t) has the bound N lim sup Hint (t12 ) ≤ Kβt12 (4.62) N →∞ 1 0 < min g a.s. 1≤i≤M 6 i for a sufficiently small t12 > 0. Total loss from boundary events HbN (t) satisfies the bound M X lim sup HbN (t22 ) ≤ μ0k ([0, a(t22 )]) (4.63) N →∞ k=1 1 0 < min g , 1≤i≤M 6 i for a sufficiently small t22 > 0. Letting t2 = t12 ∧ t22 then shows part (1), and letting 1 T e = t1 ∧ t2 ∧ (4.64) 2βM K 69 gives a positive time Te > 0 that satisfies all three conditions. In the following, denote the set of intervals Iδ of length at most δ as Iδ = {(a, b) ⊂ R+ , b − a < δ}. (4.65) Our method of proof shows that the sequences μN k (t) are approximately nonatomic with probability one as N → ∞. We make this rigourous with the following defini- tion. Definition 5. A sequence of random measure processes νtn is probabilistically approximately nonatomic (p.a.n.a.) for time T if for all η > 0, lim sup lim sup sup P(νtn (I) > η) = 0. (4.66) δ→0 I∈Iδ n→∞ t∈[0,T ] We will eventually demonstrate that the p.a.n.a. property is sufficient to show the Aldous conditions of Theorem 4. We aim to express μN k (t) by a decomposition of particles though their histories as different species. With this in mind, we have the following definition, which isolates subsets of particles in species k. Definition 6. A finite measure process AkN (t) ∈ M(R+ ) × [0, Te ] is a subparticle measure process of μN k (t) if X δx AkN (t) = . (4.67) N N x∈S (t) for some support S N (t) ⊂ Lk satisfying S N (t) ⊂ Sk (X N (t)). Given a subparticle measure process AjN (t) with support S N (t) , we can then 70 define empirical measures of particles in S N (t) that jump from Lj to Lk and then flow on vk with no probability of jumping again. Definition 7. A boundary running measure A˜j,k N (t) ∈ M(R+ ) × [0, Te ] of a subparticle measure process AjN (t) with support S N (t) is a measure of the form X δϕk (ϕj (s,x),t−s) A˜j,k N (t) = (4.68) k N x∈U (s) with support U k (τ ) = {x ∈ S N (τ − )|A boundary event at τ transfers x to Lk }. (4.69) We will denote the jump times of a running measure as τ˜1 , . . . , τ˜m(t) ˜ , and the particles selected as x˜1 , . . . , x˜m(t) ˜ , where m(t) ˜ is the total number of such jumps. We can similarly define interior running measures A˜j,k int,N (t) for particles that jump from Lj to Lk due to interior events. Theorem 16. For j, k = 1, . . . , M and t ∈ [0, Te ], if the sequence of subparticle ˜ int,N (t) and measure processes RjN (t) is p.a.n.a., then so are running measures R j,k ˜ N (t). R j,k Proof. Let t ≤ Te . Suppose a particle is chosen to transfer from Lj to Lk at time s ≤ t. At a boundary event, after choosing a species Lj , particles are chosen uni- formly on Lj to reassign to another species. Thus, the (random) probability qsN (I) of reassignment of a particle in the support of RjN (s) in the correct interval at time s before traveling to I satisfies RjN (s)(ϕk (I, s − t)) qsN (I) = . (4.70) gkN (s) 71 By the p.a.n.a. property of RjN (t), for ε > 0, there exists δ > 0 where for all I ∈ Iδ lim sup sup P(qsN (I) > ε) (4.71) N →∞ s≤t 2gk0 ≤ lim sup sup P(RjN (s)(ϕk (I, s − t)) > ε) < ε N →∞ s≤t 3 The factor 2gk0 /3 arises from part 2 of Lemma 5 , which bounds the total number of a particular species. We now note that, from Markov’s inequality, and part 2 of Lemma 4, for I ∈ Iδ , ˜ N (t)(I) > η) lim sup sup P(R (4.72) k,j N →∞ t∈[0,Te ] h i E R˜ N (t)(I) j,k K(1 + βt) lim sup sup ≤ lim sup sup E(qsN (I)) N →∞ s∈[0,t] η η N →∞ s≤t K(1 + βt) h lim sup sup E(qsN (I)|qsN (I) > ε))P(qsN (I) > ε) η N →∞ s≤t i + E(qsN (I)|qsN (I) > ε))P(qsN (I) > ε) K(1 + βt) ≤ 2ε, η As ε is arbitrary, we may take δ → 0 to obtain our result. For interior particle submeasures, the proof is identical. With probability one, a particle starting in Li at time 0 visits a finite set of species before arriving at its final destination at time t. Suppose a particle x(t) at time τk jumps from Lσk to Lσk+1 , for k = 1, . . . , m. The locations of x(t) can then be symbolically described through the following path notation P(x(t)) = (Σ(x(t)), Θ(x(t))), with path location Σ(x(t)) ∈ {1, . . . , M }m , and jump description Θ(x(t)) ∈ {∂, I}m−1 . Specifically, if Σ(x(t))i = σi , and Σ(x(t))i+1 = σi+1 , then x(t) 72 jumps from Lσi to Lσi+1 at the ith jump time τi . If Θ(x(t))i = ∂, then the it h jump is a boundary event at τi , while Θ(x(t))i = I, denotes an interior event. The number of different species realized by a particle, or the length of a path, will be denoted |Σ(x(t))|. We now define path measures that keep track of particles that follow specific paths. for Σ ∈ {1, . . . , M }m , denote RΣ (t) = {x(t) ∈ Sσ|Σ| (t)|Σ(x(t)) = Σ}, (4.73) and X δx μN Σ (t) = . (4.74) N x∈RΣ (t) It’s straightforward to see that μN N Σ (t) are subparticle measures. Path measures μΣ (t) with support RΣ (t) also have the following no return property: Suppose we track a specific particle x(s) ∈ E, with initial position x0 ∈ E. If x(t) ∈ RΣ (t), and there exists t0 > t with x(s) ∈ RΣ (t0 ), then x(s) ∈ RΣ (s) for s ∈ [t, t0 ]. The no return property ensures that a particle that is a member of RΣ (t) is not permitted to exit and reenter Lk . The importance of such a property will become apparent in the following theorem. Denote the running measures of μN Σ (t) from Lσ|Σ| ˜N to Lk as μ ˜ β,N Σ,k (t) and μ Σ,k (t). We then have the following Theorem 17. μN ˜N Σ (t), μ ˜β,N Σ,k (t),and μ Σ,k (t) are p.a.n.a. for length Te . Proof. We prove the theorem by induction on the path length |Σ| . For the base case, |Σ| = 1 corresponds to particles which have not yet jumped from their initial distribution. Thus the support of μN Σ (t) is strictly contained in the support of the pullback measure (ϕσ1 )∗ μN (t), and is thus p.a.n.a., which implies that, from Theorem 73 ˜N 16, μ ˜ β,N Σ,k (t) and μ Σ,k (t) are also p.a.n.a.. For the inductive step, suppose that for |Σ| ≤ n, μN ˜N Σ (t), μ ˜ β,N Σ,k (t), and μ Σ,k (t) are ˜ = Σ∗k, where |Σ| = n and ∗ denotes concatenation. p.a.n.a.. For k = 1, . . . , M , let Σ Then μN ˜N Σ is p.a.n.a., since the support of μ N Σ,k (t) contains the support of μΣ ˜ (t), which follows from the no return property. This proves the induction step for μN Σ (t). Then ˜N we again have from Theorem 16 that for j = 1, . . . , M , μ ˜ β,N ˜ (t) and μ Σ,k ˜ (t) are Σ,k p.a.n.a., since μN ˜ (t) are subparticle measures. This completes the inductive step, Σ and our result follows. We need one more lemma to show that our empirical measures are p.a.n.a., namely that total numbers of empirical measures decrease subgeometrically with respect to the length of paths. |Σ| 3 Lemma 6. For t < Te , lim supN →∞ μN + Σ (t)(R ) ≤ c a.s., where c < 4M . ˜N Proof. If suffices to show that both lim sup N →∞ μ + N + Σ,k (t)(R ) ≤ lim sup cμΣ (t)(R ) ˜ β,N and lim sup μ + N + Σ,k (t)(R ) ≤ lim sup cμΣ (t)(R ) a.s.. We’ll first show the inequality for boundary running measures. Recall from Lemma 5, that for t < Te , we can ensure that 1. lim supN →∞ gkN (s) > 23 gk0 , and 1 2. lim supN →∞ H N (t) < g0, 3M 2 K k gk0 hold almost surely. Suppose first lim sup N →∞ μN + Σ (t)(R ) ≥ 2M . Since the maximum total number of particles that can leave a species Lk from boundary events is bounded 74 by KH(t), Condition (2) of Lemma 5 implies for k ∈ {1, . . . , M }, 1 0 ˜N lim sup μ + Σ,k (t)(R ) ≤ g . (4.75) N →∞ 3M 2 k But by the assumption on μN + Σ (t)(R ), 1 0 2 g k ≤ lim sup μN (s)(R+ ). (4.76) 3M 2 3M N →∞ Σ On the other hand, if lim sup N →∞ μN + 0 Σ (t)(R ) ≤ gk /2M , then we know from condition (1) of Lemma 5 that the probability of any particle chosen from Lk being in the support of μN 0 0 Σ (s) has probability less than (gk /2M )/(2gk /3) = 3/(4M ). Thus, 3 ˜N lim sup μ + Σ,k (t)(R ) ≤ lim sup μN + Σ (s)(R ) a.s. (4.77) N →∞ 4M N →∞ To bound interior probabilities, we again use Lemma 5, but now use condition (3), where Te satisfies KβM 2 Te < mink gk0 /4. Then, as the Poisson rate is bounded by βN , we have from Lemma 4 that the total number of species reassignments is bounded by mN int (t) mini gi0 lim sup K ≤ Kβt ≤ a.s. (4.78) N →∞ N 4M 2 Thus, we can apply an argument similar to the case of boundary events, where we again compare lim sup N →∞ μN + 0 Σ (t)(R ) with gk /2M to obtain 3 ˜β,N lim sup μ + Σ,k (t)(R ) ≤ lim sup μN (t)(R+ ) a.s. (4.79) N →∞ 4M N →∞ Σ Theorem 18. Let t < Te . The empirical measures μN k (t) are p.a.n.a. 75 Proof. Let η > 0. We can decompose μN k (t) as X ∞ X X μN k (t) = μN Σ (t) = μN Σ (t). (4.80) (Σ,Θ)∈P j=1 |Σ|=j σ|Σ| =k Let η > 0, and choose L so that μN + Σ (t)(R ) ≤ η/8 for |Σ| > L (This is possible PL−2 through Lemma 6). There are ML := M j=1 (M − 1)j−1 paths where |Σ| ≤ L and σ|Σ| = k (M initial paths and M − 1 options for each new transition). From Lemma 17, each of these path measures are p.a.n.a., meaning that lim lim lim sup sup P(μN Σ (s)(I) > η ˜) = 0 (4.81) δ→0 I∈Iδ N →∞ s∈[0,Te ] with η η˜ = , (4.82) 2ML 76 for any μN Σ (t) with |Σ| < L. But then we have, for any t ∈ [0, Te ], and I ∈ Iδ , lim sup sup P(μN k (t)(I) > η) (4.83) N →∞ t∈[0,Te ] X ∞ X  = lim sup sup P μN Σ (t)(I) > η N →∞ t∈[0,Te ] j=1 |Σ|=j σ|Σ| =k X L X ηX X N  ≤ lim sup sup P μN Σ (t)(I) + + μΣ (t)(R ) > η N →∞ t∈[0,Te ] j=1 |Σ|=j 8 j>L |Σ|=j σ|Σ| =k σ|Σ| =k X L X ηX ∞  ≤ lim sup sup P μN Σ (t)(I) + (cM )j > η N →∞ t∈[0,Te ] j=1 |Σ|=j 8 j=1 σ|Σ| =k X L X  ≤ lim sup sup P μN Σ (t)(I) > η/2 N →∞ t∈[0,Te ] j=1 |Σ|=j σ|Σ| =k ≤ 1 − lim sup sup P(μN Σ (t)(I)) ≤ η ˜) for all |Σ| < L. N →∞ t∈[0,Te ] By (4.82), this quantity approaches 0 as δ → 0. Here we show the Aldous conditions, which are a direct property of the p.a.n.a. property of μN k (t). Lemma 7. The Aldous conditions of Lemma 4 hold for μN k (s) Proof. Condition (1) follows trivially, since we have for all t ∈ [0, Te ], |hμN k , φi| < kφk∞ , (4.84) meaning that (1) holds for L = kφk∞ . For condition (2), we observe that particles 77 from time t to t + s either drift on Lk , are deleted, or are reassigned to a new species. lim sup |hφ, μN N k ((τ + s) ∧ Te )i − hφ, μk (τ )i| (4.85) N →∞ Z = lim sup φ(x)(μk ((τ + s) ∧ Te ) − μk (τ ))(dx) N N N →∞ ZR∞ + ≤ lim sup (φ(ϕ (x, s)) − φ(x))μk (τ )(dx) k N N →∞ s M− ! X N i + kφk∞ (K + 1) μi (τ )([0, ϕ (0, −s)]) + sβ a.s. i=1 The last inequality uses the fact, from Lemma 4 that almost surely at most a sβ total P − N number are reassigned from interior events, and that (K+1)( M i i=1 μi (τ )([0, ϕ (0, −s)]) particles will either be reassigned or deleted due to boundary events. We can use the continuity and boundedness of φ and ϕ to estimate (4.87). For all R > 0: Z ∞ k (φ(ϕ (x, s)) − φ(x))μN k (τ )(dx) (4.86) s ≤ kφ(ϕk (x, s)) − φ(x)kL∞ ([0,R]) + 2kφk∞ μN k (s)([R, ∞]) The first term tends to 0 as s → ∞. For the second term, we can always find R large enough where only a small total number enters [R, ∞] in any Lk . Precisely, for all ε, there is an R = R(ε) where M X M X μ0k ([R, ∞)) < ε/2, μ0k ([R/2, R]) < ε/2, (4.87) i=1 i=1 k ϕ (R/2, t) < R, t ∈ [0, Te ]. These bounds are thus independent of N and s. Thus, since R is arbitrary, (4.87) tends to 0 as s → 0. For the second term of (4.87) , notice that from the p.a.n.a. 78 property of μN k (t),we have, for η > 0,   η lim lim sup sup P μN k k (t)([0, ϕ (0, −s)]) > = 0. (4.88) s→0 N →∞ τ ∈[0,T ] e M Kkφk∞ Thus we finally have, for any η > 0, lim lim sup sup P(|hφ, μN N k (τ + s ∧ T )i − hφ, μk (τ )i| ≥ η) (4.89) s→0 N →∞ τ ∈[0,T ] e M ! X ≤ lim lim sup sup P kφk∞ K μN k k (τ )([0, ϕ (0, −s)]) ≥ η s→0 N →∞ τ ∈[0,T ] e i=1 which approaches zero from (4.89), meaning that Aldous’ conditions are satisfied. By Lemma 1, we have shown tightness of random variables hμN k (t), φi, φ ∈ Cb (R+ ). We can now use Lemma 4 to obtain existence of a limiting measure,μk ∈ M(R+ ) where for φ ∈ Cb (R+ ), there exists a subsequence Nk with hφ, μN k (s)i → hφ, μk (s)i k (4.90) in the local uniform topology for all φ ∈ Cb (R+ ). We now show the measures μk (s), k = 1, . . . , M are nonatomic almost surely. Theorem 19. Let k = 1, . . . , M and t ∈ [0, Te ]. For any point p ∈ R+ , the measures μk (t) satisfy P(μk (t)({p}) = 0) = 1. (4.91) Proof. Let η > 0. We will denote Iδp = {I ∈ Iδ |p ∈ I}. (4.92) 79 Using the p.a.n.a. property, we have, for δ > 0, P(μk (t)({p}) > η) ≤ lim sup P(μk (t)({p}) > η) (4.93) N →∞ ≤ sup lim sup P(μk (t)(I) > η) → 0 as δ → 0. I∈Iδp N →∞ As η was arbitrary, the result follows. 4.7 Convergence of boundary variables φ,N The final step in establishing Theorem (11) is to show that the quantities Bj,k (t) converge in the local uniform metric. Toward this end, we examine measures FlN , l = 1, . . . , M− which track the total number deleted from boundary events in species l. In this section, we assume limits of X N (t) are taken along subsequences where empirical measures converge local uniform metric. Throughout this section we’ll use S(a,b) the notation Y n (t) −−−→ Y (t) for convergence in the probability space S(a, b) of finite L∞ mean stochastic processes for time t ∈ [a, b], endowed with the norm   kY kS = E sup |Y (t)| . (4.94) a≤s≤b This space is complete, which is a consequence of the completeness of L∞ ([0, T ]). Also, local uniform convergence of measures μN k to μk , along with almost sure nonatomicity of μk imply that for k = 1, . . . , M , S(0,Te ) hμN −−−→ hμk , φi φ ∈ L1 (R+ ). k , φi − (4.95) 80 Let us define approximate cumulative distribution functions t bεc X  Flε (t) = μl (εk) [0, ϕl (0, −ε)) , l = 1, . . . , M− . (4.96) k=1 We can write F ε,N as the finite particle approximation of F ε : t bεc X  Flε,N (t) = μN l l (εk) [0, ϕ (0, −ε)) , l = 1, . . . , M− . (4.97) k=1 We know, through convergence of finite dimensional distributions of μk (t), and from 4.96 with φ = 1[0,ϕl (0,−ε)) , that S(0,Te ) Flε,N (t) −−−−→ F ε (t) as N → ∞ (4.98) The increasing functions Flε,N (t) assume that all particles within a ϕl (0, −ε) strip of the origin of Ll , l = 1, . . . , M , cannot be reassigned, and will thus be deleted before time ε. We then define, for t ∈ [0, Te ], Fl (t) = lim Flε (t) in S(0, Te ). (4.99) ε→0 Theorem 20. Let t ∈ [0, Te ]. The quantity Fl (t) in (4.100) is well-defined, and is equal to the total number limit in S(0, t) as N → ∞ of particles FlN (t) which hit species l. Proof. Let ε > 0. We decompose the approximated measure by the true total number lost and an error term: Flε (t) = lim FlN (t) + El,N ε (t) in S(0, t). (4.100) N →∞ 81 Here, FlN (t) denotes the empirical CDF of the total number that actually hit the ε origin of Ll . The error El,N denotes the sum over k = 1, . . . , bt/εc of the total number not in [0, ϕl (0, −ε)) ∩ Ll at time εk that hit the origin of Ll during time [ε, ε(k + 1)], minus the sum of particles in [0, ϕl (0, −ε)) ∩ Ll at time εk that do not hit the origin of Ll during time [ε, ε(k + 1)]. ε A simple bound for El,N (t) is simply the total number of particles that are re- assigned in [0, ε) before time t. The following argument for bounding redistribution probabilities in [0, ε] will mirror the bound in Theorem 16. Specifically, denote qsε,N as the (random) probability that a critical event selects a particle in [− M  [0, ϕl (0, −ε)) ∩ Ll . (4.101) l=1 Probabilities can be expressed as M− X μN l j (s)([0, ϕ (0, −ε)) ∩ Ll ) qsε,N = . (4.102) i=1 gjN (s) Then, fixing any η > 0, by the p.a.n.a. property, we can find an ε0 > 0 where any ε < ε0 satisfies lim sup sup P(qsε,N > η) (4.103) N →∞ s≤t M− X 2gk0 ≤ lim sup sup P([0, ϕ (0, −ε)) ∩ Ll ) > l η) < η. i=1 N →∞ s≤t 3 From the p.a.n.a. property, this gives us the bound lim sup E(sup |El,N ε (s)|) ≤ (K + βt)(2η). (4.104) N →∞ s≤t 82 As η, this means ε S lim sup El,N (t) − →0 as ε → 0. (4.105) N →∞ We thus have that the ε limit of Flε (t), is Cauchy, meaning that there exists a limit Fl (t) in the S(0, t) norm. The limit satisfies lim Flε (t) = lim FlN (t) S(0, t). (4.106) ε→0 N →∞ In particular, for a continuous deterministic solution dμl (s) = ul (x, s)dx, we have that Fl (t) has the form bεc Z t X ϕ(0,−ε) Fl (t) = lim ul (x, εk)dx (4.107) ε→0 0 k=1 Z t = ul (0, s)vl (0)ds, l = 1, . . . , M− . 0 We aim to show convergence of the tracking parameters to quantities involving limiting particle densities. Namely, we wish to show, for k = 1, . . . , M , that M X Z tX M− K (l) M XX φ,N φ,N (l) Bj,k (t) − Bk,j (t) − Wj (s)hμj (s), φi1R(l) =k (4.108) ji j=1 0 l=1 i=1 j=1 M− X (l) S(0,t) + K (l) Wk (t)hμk (s), φidFl (s) −−−→ 0. l=1 By linearity of expectation, this is equivalent to focusing on one tracking parameter, and showing that as N → ∞, Z tX M− K (l) X φ,N S(0,Te ) (l) Bj,k (t) −−−−→ Wj (s)hμj (s), φidFl (s)1R(l) =k . (4.109) ji 0 l=1 i=1 83 S(0,Te ) This is done to showing two facts. The first is that FlN (t) −−−−→ Fl (t), which we have shown in Lemma (20). We then must show that, as N → ∞, Z tX M− K (l) X φ,N S(0,Te ) (l) Bj,k (t) −−−−→ Wj (s)hμj (s), φidFl (s)1R(l) =k . (4.110) ji 0 l=1 i=1 Proving (4.111) may seem tedious, but it follows quickly from a bound that is similar to a proof of the law of large total numbers for random variables with moment generating functions (see [37]). Lemma 8. Suppose we are given independent (thought not necessarily identical) ran- N (i) dom variables Dj , for j = 1, . . . , N (i), and N (i) → ∞ as i → ∞. Suppose further that the variables all have well-defined moment generating functions Mj,N (i) (t), and that for each i ∈ N, j = 1, . . . , N (i), we have −∞ < μ := lim inf E(Dij ), lim sup E(Dij ) := μ < ∞. (4.111) i,j i,j We then have N (i) N (i) N (i) N (i) X Dj X Dj lim inf , lim sup ∈ [μ, μ] a.s. (4.112) i j=1 N (i) i j=1 N (i) Proof. We utilize the most basic version of the Chernoff bound, which states that for a random variable X with moment generating function M (t), P(X > a) < min e−at M (t), P(X < a) < min e−at M (t) (4.113) t>0 t<0 N (i) This would imply that, for fixed i, and Dj with moment generating function 84 Mj,N (i) , N (i) X N (i)   Dj −(μ+η)tN (i) P( > μ + η) ≤ min e max(Mj,N (i) (t)) N (i) (4.114) j=1 N (i) t>0 j N (i) ≤ min max e−(μ+η)t Mj,N (i) (t) . t>0 j Note that the function fj,N (i) (t) = e−(μ+η)t Mj,N (i) (t) satisfies fj,N (i) (0) = 1. Also, we have that for a sufficiently large n1 ∈ N, fj,N 0 (i) (0) < −η/2 for i > n1 . Thus there exists λ ∈ [0, 1) such that   N (i) N (i) X Dj P > μ + η  < λN (i) . (4.115) j=1 N (i) Similarly, we can find n2 ∈ N and γ ∈ [0, 1), where for sufficiently large i > n2 ,   N (i) N (i) X Dj P < μ − η  < γ N (i) . (4.116) j=1 N (i) Now, define the event   XN (i) D N (i)  j Dηi = / [μ − η, μ + η] ∈ , (4.117)  N (i)  j=1 P∞ It follows from (4.116) and (4.117) that i=1 P(Dηi ) < ∞, which from the Borel- Cantelli lemma implies that P(Dηi occurs infinitely often in i) = 0. The sets Dηi are monotonically decreasing in probability as η → 0. Thus   N (i) N (i) X Dj 0 = lim P(Dni ) = P(lim Dni ) = P  / [μ, μ] , ∈ (4.118) η→0 η→0 j=1 N (i) which is what we sought to prove. 85 To apply Lemma 8, we approximate the integral in (4.111) by a Riemann sum so that if we define the following continuous process K (l) X (l) Ql (t) = Wj (t)hμj (t), φi1R(l) =k , s ∈ [0, Te ] (4.119) ij i=1 S(0,Te ) then we have, for some δ > 0, since FlN (t) −−−−→ Fl (t), there exists a partition P = Pδη with mesh a1 < ∙ ∙ ∙ < a|P| that satisfies Z t X Ql (s)dFl (s) − Ql (ai )ΔFlN (ai ) < δ. (4.120) 0 ai ∈P S(0,Te )) where we impose that the mesh of the partition P is defined so that for an arbitrary η > 0, E[|Ql (s) − Ql (t)|] < η for s, t ∈ [ai , ai+1 ] i = 1, . . . , |P |. (4.121) The reason for defining the partition based on a modulus of continuity for Ql is to approximate expectations of reassignments in a short period of time. More specifi- cally, we can define boundary variables in terms of their behavior in each segment of P , as φ,N X X ZlN (tk ) Bi,j (t) = . (4.122) P N tk ∈[ak ,ak+1 ] Here, assuming a boundary event triggered at Ll , the random variables ZlN (tk ) at P boundary jump time tk assign 0 if Li is unaffected, and m φ(xm ) for particles xm ∈ Li that are reassigned to Lj . Using the mean-field species probabilities, we can 86 see that for t ∈ [0, Te ], E[ZlN (t)] = E[QN l (t)], (4.123) where K(l) X (l),N QN l (t) = Wj (t)hμN j (t), φi1R(l) =k , s ∈ [0, Te ]. (4.124) ij i=1 However, since QN l (t) → Ql (t) in the local uniform metric, it follows that that for some boundary jump time t ∈ [ai , ai+1 ], lim sup E[ZlN (t)] ≤ Ql (ai ) + η, lim inf E[ZlN (t)] ≥ Q(ai ) − η, (4.125) N →∞ N But then from Lemma 8, as the total number of jump times t1 , . . . , tml (i) in [ai , ai+1 ] tends to infinity, we have with probability 1 that N (i) X N (i) N (i) Zl (tj ) Zl (tj ) lim inf , ∈ [Ql (ai ) − η, Ql (ai ) + η] , a.s. (4.126) i,j j=1 m l (i) m l (i) from which we can infer that over the partition interval [ai , ai+1 ],   M− N (i) X ml (i) X N (i) Zl (tj ) lim sup  − Ql (pi )ΔFl (pi )  (4.127) N →∞ N j=1 ml (i) l=1 S(ai ,ai+1 ) < ηM1 ΔFl (pi ). which summing over the partition gives us Z tX M− φ,N lim sup Bj,k (t) − Ql (s)dFl (s) < ηM1 Fl (t). (4.128) N →∞ 0 l=1 S(0,Te ) As η is arbitrarily small and Fl (t) < 1, (4.111) follows. 87 4.8 The existence interval revisited. Our interval for existence [0, Te ] is solely based on three conditions of Lemma 5. Running the equation on its maximal interval, it then follows that there exists a time t1 and constants 0 < c1 , c2 < 1, c3 < ∞ where one of the following three events occurs: 1. There exists k ∈ {1, . . . , M } where gk (t1 ) = (1 − c1 )gk0 . 2. There exists k ∈ {1, . . . , M } where H(t1 ) = c2 gk0 . 3. t1 = c3 minM 0 i=1 gi . At time t1 we can run the equation afresh by setting the initial conditions as μnk,0 = μnk (t1 ). The limiting distribution μk,0 is nonatomic, so we can find a time t2 > 0 as before. Continuing in this fashion, we can obtain a sequence of extension times ti , i ∈ N where the weak form of the kinetic equation (4.22) is well defined for each [0, ti ]. This turns out to be the right method for determining a reasonable maximum interval for (4.22). Theorem 21. Given a sequence of extension times tj % Tˉ < ∞, there exists k ∈ {1, . . . , M } where lim gk (tj ) = 0. (4.129) tj →Tˉ− (l) This is a reasonable time parameter, since species selection probabilities Wi (s) (l) may be undefined based on the selection of weights wi . Of course, specific weights may allow for a species to be empty with well defined selection probabilities, and in this case nothing would prevent a solution for existing until enough species exhaust 88 their particles to make the transition probabilities ill-posed. Determining exact con- ditions for which species can delete would be tedious and not very informative, so we’ll only consider extending solutions to with strictly positive species total numbers. Proof. Assume, to the contrary, that there exist α, ε > 0 where for gk (t) > α > 0 for some t ∈ [Tˉ − ε, Tˉ ). As the extension times are converging to some finite Tˉ < ∞, condition (3) can only be satisfied finitely many times. Similarly, condition (2) cannot occur infinitely often, as that would imply the deletion of infinite particle total number. Thus, condition (1) occurs infinitely often. This means that there exists a species Lk and arbitrarily small intervals tj+1 − tj where gk (tj+1 ) − gk (tj ) ≥ c1 α. This, however, implies an infinite total number of particle density transported in a finite time, which cannot occur by Lemma 4. 4.9 Uniqueness and regularity for kinetic equa- tions Uniqueness of the weak solutions satisfying (4.22) with continuous bounded initial data will result from the following steps. 1. Prove uniqueness for a small time for mild solutions of kinetic densities, ob- tained through the variation of parameters. This draws heavily from the tech- niques used in [18]. 2. Show solutions change continuously with respect to initial data, and that dif- ferentiable initial data implies a differentiable bounded solution. 89 3. Take a limiting sequence for differentiable initial data that converges to arbi- trary continuous bounded initial data. 4. Extend the existence interval to a reasonable time, as described in Theorem 21. Step 1 Our goal is to set up a contraction for limiting equations for the k-species system. The following exposition is very technical, but the main idea is simple: severely restricting time, we can find an interval of uniqueness for positive continuous initial data by bounding L1 and L∞ norms. We begin with distributional equations of tier densities, ∂t uk (x, t) + ∂x (vk (x)uk (x, t)) = (4.130)   XM− K (l) M X X (l) (l) ul (0, t)vl (0)  Wi (t)ui (x, t)1R(l) − K (l) Wk (t)uk (x, t) ij l=1 i=1 j=1   M K int βG(t)  X X + wiint ui (x, t)1Rijint =k − K int wkint uk (x, t) Gint (t) i=1 j=1 with initial conditions uk (x, 0) = u0k (x), k = 1, . . . , M. (4.131) Integral equations may be written, from Duhamel’s formula (or variation of pa- 90 rameters), as uk (x, t) = u0k (ϕk (x, −t))+ (4.132) Z tX M1 M X X K (l) (l) ul (0, s)vl (0) Wi (s)ui (ϕk (x, s − t), s)1R(l) ij 0 l=1 i=1 j=1 M1 X (l) − K (l) ul (0, s)vl (0)Wk (s)uk (ϕk (x, s − t), s) l=1 βG(s) + − K int wkint uk (ϕk (x, s − t), s) Go (s) M K int X X  + wiint ui (ϕk (x, s − t), s) 1Rijint =k ds. i=1 j=1 We will be working in the Banach space M n=1 ∈ C(R , R ) : kf k1 + kf k∞ < ∞, X M = {f = (fn )M + M min kfn k1 > 0}, (4.133) n=1 with norms defined as the maxima of the traditional norms of functions fi : R → R: M M kf k∞ = max kfi k∞ , kf k1 = max kfi k1 (4.134) i=1 i=1 The relevant Banach space for our equations in time and space is M Y M = {h ∈ C([0, T ], X M ) : sup khk1 + sup khk∞ < ∞, sup min khn k1 > 0}. t t t n=1 (4.135) We also desire for our initial conditions to be positive, or in the set M G = {f ∈ X M : min ≥ 0 ∀k ∈ {1, . . . , M }, x ∈ R+ }. (4.136) i=1 91 Written with respect to flux terms (l) (l) wj ul (0, s)vl (0) Vj = P (l) , (4.137) m w m g m we now construct the iteration map I, where a fixed point of I gives a mild solution: I(f ) = (In (f ))M i=1 (4.138) Ik (f )(x, t) = u0k (ϕk (x, −t)) Z tX M1 X M X K (l) (l) + Vi (s)fi (ϕk (x, s − t), s)1R(l) =k ij 0 l=1 i=1 j=1 M1 X (l) − K (l) Vk (s)fk (ϕk (x, s − t), s) l=1 βG(s) + − K int wkint fk (ϕk (x, s − t), s) Go (s) M K int X X  + wiint fi (ϕk (x, s − t), s) 1Rijint =k ds i=1 j=1 Define the weighted total number for the j th species as M X Z (j) (j) g (u(t)) = wm um (x, t)dx, j = 1, . . . , M− . (4.139) m=1 R+ Our closed subspace, where we will show a contraction, is then B = {h ∈ Y M , kh(t)k1 ≤ A1 , kh(t)k∞ ≤ A∞ , (4.140) hj (u(t)) ≥ A ∀t ∈ [0, T ]}, Where we have constants gm (u0 ) A1 > ku0 k1 , A∞ > ku0 k∞ , A < min . (4.141) m 2 92 Throughout the next , we will also use the notation w˜ = max wij v˜ = max vi (0). (4.142) i,j j We may now state our main uniqueness theorem: Theorem 22. For initial data u0 ∈ G, there exists a unique continuous mild solution u(x, t) to the integral equation (4.167) for a small time t ∈ [0, T ] where T is given by constants T < min{T1 , T , Tv , Tf }, (4.143) where A(A1 − ku0 k1 ) A(A∞ − ku0 k∞ ) T1 = min{ , }, (4.144) KM− wA ˜ 1 (˜ ˜ − wA v A∞ + 2β) KM ˜ ∞ (˜ v A∞ + 2β) ku0 k1 /2 − A T = , and (4.145) KM− w˜ ˜ v A1 (A∞ + 2β)/A  A∞ w˜ ˜ v + q˜β Tf = 2M− K (4.146) A  2M w(˜ ˜ v + 1)  −1 2 ∙ 1+ (A 1 + A ∞ (1 + β w)) ˜ . A2 Proof. We must show the map I is iterative and contractive. To see it retains L1 and L∞ bounds, note KM− wA ˜ 1 (˜ v A∞ + 2β) sup k(Iu)(t)k1 ≤ ku0 k1 + T < A1 , (4.147) t A KM− wA ˜ ∞ (˜ v A∞ + 2β) sup k(Iu)(t)k∞ ≤ ku0 k∞ + T < A∞ . t A 93 For weighted total numbers, we can show gj (u0 ) ˜ v K l M − A1 w˜ gj (I(u(t)) > −T (A∞ + 2β) > A. (4.148) 2 A We now bound the difference of flux terms, using the fact that the denominator of each of these terms satisfies D(f ) ≥ A, and |Wlj (u) − Wlj (˜ u)| (4.149) (l) (l) qj ul (0, s)vl (0) qj u˜l (0, s)vl (0) ≤ P (j) − P (j) | wm gm (u) wm gm (˜ u) ˜ v M− (A1 ku − u˜k∞ + A∞ ku − u˜k1 ) w˜ ≤ A2 Similarly, bounding interior event terms gives us G(u) G(˜u) 2wA ˜ 1 k˜ u − uk∞ − int ≤ . (4.150) int G (u) G (˜ u) A2 94 Our main estimates on the difference of the iterative map I are then sup kI(u)(t) − I(˜ u)(t)k1 (4.151) t Z M− M K (l) ∞X XX (l) (l) ≤ T sup Vj (u(s))uj (x, s) − Vj (˜ u(s))˜ uj (x, s) 1R(l) =k t ij 0 l=1 i=1 j=1 M− X (l) (l) (l) + K Vk (u(s))uk (x, s) − K (l) Vk (˜ u(s))˜ uk (x, s) l=1 βG(u(s)) βG(˜ u(s)) int int + int K int wkint uk (x, s) − int K wk u˜k (x, s) G (u(s)) G (˜ u(s)) X βG(u(s)) M K int X βG(˜u(s)) int int + int wj uk (x, s)ds − int wk u˜k (x, s) 1R(l) =k dx i=1 j=1 G (u(s)) G (˜ u(s)) ij h i (l) (l) (l) ≤ 2T M− K sup(sup |Vj (u(t))|ku(t) − u˜(t)k1 + A1 sup |Vj (u(t)) − Vj (˜ u(t)))| t l,j l,j h β w˜ G(u(t)) G(˜u(t)) i +2T M− K sup(ku(t) − u˜(t)k1 + A1 β w˜ int − int | A t G (u(t)) G (˜ u(t)) h A w˜ ∞ ˜v + β w ˜ (l) (l) ≤ 2T M− K ku(t) − u˜(t)k1 + A1 (sup |Vj (u(t)) − Vj (˜ u(t)))| A l,j G(u(t)) G(˜ u(t)) i +A1 β w˜ int − int G (u(t)) G (˜ u(t)) Similarly, sup kI(u)(t) − I(˜ u)(t)k∞ (4.152) t h A w˜ ∞ ˜v + β w ˜ (l) (l) ≤ 2T M− K ku(t) − u˜(t)k∞ + A∞ (sup |Vj (u(t)) − Vj (˜ u(t)))| A l,j G(u(t)) G(˜u(t)) i + A∞ β w˜ int − int . G (u(t)) G (˜ u(t)) We now finally show our contraction. Using bounds for differences of flux terms (l) Vj (u(t)) and G(u(t))/Gint (u(t)), we then have 95 sup kI(u)(t) − I(˜ u)(t)k∞ + sup kI(u)(t) − I(˜ u)(t)k1 (4.153) t t ≤ T Tf−1 (sup ku(t) − u˜(t)k1 + sup ku(t) − u˜(t)k∞ ) t t < (sup ku(t) − u˜(t)k1 + sup ku(t) − u˜(t)k∞ ) t t by our choice of T < Tf−1 . This shows that we have a short time contraction, and thus uniqueness of the mild solution on [0, T ]. Step 2 We will now show continuous dependence upon initial conditions. Theorem 23. Mild solutions depend continuously upon initial data. Proof. Suppose u1 (x, t) and u2 (x, t) are solutions with respective initial data u1,0 (x) and u2,0 (x). Define L(t) = ku1 (x, t) − u2 (x, t)k∞ + ku1 (x, t) − u2 (x, t)k1 (4.154) We can use the estimates from (4.152) and (4.153) to obtain, for some positive constant C, Z t 1,0 2,0 1,0 2,0 L(t) ≤ ku (x) − u (x)k∞ + ku (x) − u (x)k1 + C L(s)ds. (4.155) 0 Taking Gronwall’s inequality from 4.156 , we obtain L(t) ≤ (ku1,0 (x) − u2,0 (x)k∞ + ku1,0 (x) − u2,0 (x)k1 ) exp(Ct), (4.156) 96 which gives us our required continuous dependence on initial conditions. Our next task is to show that in the case of continuous functions, mild solutions are equivalent to weak solutions. This is actually done by first analyzing solutions with differentiable initial data, whose equivalence to weak and mild data can be easily shown. Theorem 24. The mild solution uk (x, t) of (4.167) with differentiable initial data u0k (x) is differentiable in R+ × [0, T ], where T is the interval of existence defined in Theorem 22. Assuming differentiable initial data u0 (x) = ∂x U0 (x), we know from Theorem 22 the existence of a unique continuous solution U (x, t). Define the difference quotient dh (U (x, t)) as U (x + h, t) − U (x, t) dh (U (x, t)) = (4.157) h From the definition of weak solution, this gives us the bounds kdh (U (t))k∞ ≤ kdh (U0 (x))k∞ + (4.158) Z t j N (s) K sup(sup(Vi (s)uj (0, s)) + ) kdh (U (s))k∞ ds (4.159) s≤t i,j Ns (0) 0 where K < ∞ is a fixed constant based on the parameters of the system of equations. Applying Gronwall’s inequality then gives us kdh (U (t, x))k∞ ≤ kdh (U0 (x))k∞ exp(C(g(t))) (4.160) with C(g(t)) < ∞ is based off of fixed parameters and total numbers g(t). This means that dh (U (x, t)) has a derivative u(x, t) := ∂x U (x, t). Showing that this is 97 actually the solution to 4.167 with initial conditions u0 (x), we perform the same estimates to obtain kdh (U (t, x))−u(x, t)k∞ ≤ kdh (U0 (x))−u0 (x)k∞ exp(C(g(t)) → 0 as h → 0. (4.161) Step 3 Our next step is to show the equivalence of mild and weak solutions for continuous initial conditions. This is done through approximating continuous initial data with differentiable initial data. Theorem 25. For differentiable initial conditions, there exists a unique differentiable solution u(x, t) that satisfies both weak and mild equations. Proof. Suppose u(x, t) ∈ C 1 (R+ × [0, T ]) is a weak solution to (4.22). Integration by parts and differentiation of (4.22) then yields 0 = h−∂t uk (x, t) + (uk (x, s)vk )0 , φi (4.162) DX Kl M1 X X M1 X j + Vl (s)uj (x, s)vl (0)ul (0, s) − Kl Vlk (s)uk (x, s) l=1 i=1 Ri (j)=k l=1 l βN (s) h E Ks X X k + − Ks qs uk (x, s) + qsj Cj (s), φ ds. No (s) i=1 i Rs (j)=k from which we deduce that u(x, t) is also a strong solution, and thus a mild solution. Similarly, assuming u(x, t) is a mild solution with differentiable initial data, we know that it is unique and solves the strong equation, as it is differentiable, and therefore solves the weak equation as well. 98 We now show equivalence of weak and mild continuous solutions. Theorem 26. For continuous initial conditions, there exists a unique continuous solution u(x, t) that satisfies both weak and mild equations. Proof. We have already shown the uniqueness of mild solutions. We are left to show that weak solutions are unique, and that the solutions also satisfies the mild formulation. Let u(x, t) ∈ C(R+ , [0, T ]) be a weak solution of (4.167) with initial conditions u0 (x), and let u ˜(x, t) ∈ C(R+ , [0, T ]) be a mild solution of (4.22) also with initial conditions u0 (x). Now let an approximating sequence un,0 (x) ∈ C 1 (R+ , [0, T ]) of u0 be given, with respect to the Z = L∞ ∩ L1 norm given by ku − vkZ = ku − vk∞ + ku − vk1 . (4.163) By Theorem 24, there exist unique differentiable weak solutions un (x, t) for initial data un,0 . Also, un (x, t) approximate the solution u ˜ (x, t) ∈ C(R+ , [0, T ]) in the Z norm. As they are differentiable, they solve the strong, and thus mild equations. But by Theorem 23, the un (x, t) also approximate u(x, t), the unique mild solution in the Z norm. This establishes that u ˜ (x, t) = u(x, t). Step 4 For the final step, we extend the uniqueness parameter T . This is, however, straight- forward, as we have already devised existence on an interval Tˉ as described in The- orem 21. Theorem 27. The uniqueness interval for a continuous weak solution u(x, t) may be extended to t ∈ [0, Tˉ ]. 99 Proof. Let S = sup {t|u(x, t) has a unique solution on [0, t)} (4.164) t∈[0,Tˉ ) From step (1), we know that S > 0. Now assume that S < Tˉ . From Theorem 21, solutions u(x, t) are strictly positive and bounded in [0, Tˉ ). As solutions are continuous, we may define uˆ0 (x) = lim u(x, t). (4.165) t→S But then we have uniqueness of initial conditions for some small time T > 0 with initial conditions u ˆ 0 (x). But patching together solutions establishes uniqueness on [0, S + T ), in contradiction of the definition of S. Proof of Theorem 12 Proof. We may write mild equations for continuous u(x, t) in the form of Stieltjes integrals as uk (x, t) = u0k (ϕk (x, −t)) + (4.166) Z tX M1 X M X K (l) (l) Wi (s)ui (ϕk (x, s − t), s)1R(l) ij 0 l=1 i=1 j=1 M1 X (l) − K (l) ul (0, s)vl (0)Wk (s)uk (ϕk (x, s − t), s)dFl (s) l=1 βG(s) + − K int wkint uk (ϕk (x, s − t), s) Go (s) M K int X X  + wiint ui (ϕk (x, s − t), s) 1Rijint =k ds, i=1 j=1 where dFl (s) = ul (0, s)vl (0)ds. Now, let u0k ∈ Z, and take a sequence un,0 → u0 in Z, where un,0 1 + k ∈ C (R ). By continuous dependence on initial conditions, we know that 100 there exists a limiting function u(x, t) where un (x, t) → u(x, t) in Z. As functions in L1 it makes no sense to write dFl (s) = vl (0)ul (0, s)ds since we may not define ul pointwise. However, from section 4.7 we have the non-pointwise characterization bc Z t X ϕl (0,−ε) Fl (t) = lim ul (x, εk)ds, l = 1, . . . , M− . (4.167) ε→0 0 k=1 As unl (x, t) → ul (x, t) in Z, for dFln (t) = unl (0, t)vl (0)dt, this gives kFl (t) − Fln (t)k∞ (4.168) b t c Z ϕl (0,−ε) X ≤ lim kul (x, εk) − unl (x, εk)k∞ ds ε→0 0 k=1 ≤ sup kul (x, t) − unl (x, t)kZ ∙ tvl (0) → 0 as n → ∞ t∈[0,Te ] Thus, it is clear to see that the limit in Z of un (x, t) satisfies (4.167), as well as the weak form (4.22). Uniqueness of the solution is manifest from the continuous dependence of parameters of continuous solutions. Chapter Five The k-species Model Applied to Grain Coarsening 102 5.1 Introduction The purpose of this chapter is to align the k-species PDMP given in Chapter 3 with a mean-field model of grain coarsening. We will assume many of the mean field assumptions from the Fradkov model, given in Section 1.2. In fact, the only significant point of departure is in the selection of the flipping parameter β. After defining the appropriate parameters, we will proceed to show some basic properties of the fluid limit, such as conservation of polyhedral defect and area, that are present in instances of the finite particle model. Finally, we will examine simulations of grain PDMP. 5.2 Side redistribution at singular events We now take a closer look at particular singular events of grain and edge deletion. Flow through a singular point in a grain network is highly nonunique (see [12] for some illuminating examples). For the grain coarsening problem, we are interested in what happens when a face or side in a grain network shrinks to a point. In this case, we look for grain networks that solve the continuation problem outlined in Section 1.1 with initial conditions of G that satisfy the k-ray condition. This is defined by a finite planar graph G, with smooth edges αi (x), that is trivalent at all vertices except for one vertex, which is of degree k. To rule out possible pathologies, we’ll make the following assumptions of evolving networks G(t) ⊂ R2 , for t ∈ (0, ε) under curvature flow which enforce topological changes to take place at vertices: Continuation assumptions: Suppose G consists of n smooth boundaries αi (x) : (0, 1) → R2 , i = 1, . . . , n. Then there exist evolving boundaries αi (x, s) : (0, 1) × 103 (0, ε) → R2 and neighborhoods Uv (s), s ∈ (0, ) for each vertex v ∈ V where 1. Uv (s) → {v} in the Hausdorff distance. S Sn 2. G(s) \ v∈V Uv (s) = i=1 αi ((c1i (s), c2i (s)), s), where c1i (s) → 0 and c2i (s) → 1 as s → 0. 3. limt→0+ α(x, t) = α(x), x ∈ (0, 1). 4. Uv (s) ∩ G(s) is connected for s ∈ (0, ε), v ∈ V . Assumptions 1-3 forbid topological changes at points in G that are not vertices. The connectedness assumption is a reasonable one when we consider the underlying physical motivation of the problem, where each grain is understood as having a different material property. Solutions that disconnect at a degree k vertex without the connectedness assumption have have several different grains of positive area immediately merge regions. It may seem permissible, however, for new edges to arise from vertices in G. The following theorem shows with the continuation assumptions imposed on G, face creation is impossible, and that there are a limited number of topologies for G(t). Theorem 28. Suppose G has k-ray conditions. For a grain network G(s) with initial conditions G satisfying the continuation assumptions, there are at most 5 possible topologies of G(s) when k = 5, 2 topologies when k = 4, and 1 topology when k = 3 or 2. All solutions in a sufficiently small neighborhood of the initial degree k vertex are trees. Proof. Let G be a network with k-ray conditions, with a degree k vertex at v˜, and boundaries αi (x), x ∈ (0, 1) emanating from v, for i = 1, . . . , l. For some ε > 0, 104 suppose G(s), s ∈ (0, ε) is a grain network satisfying the continuation assumptions, with corresponding evolving boundaries αi (x, s) and neighborhoods Uv for v ∈ V . We now fix a time s ∈ (0, ε) and focus on U := Uv˜ . For a network with V vertices, E edges, and F faces, the Euler characteristic χ(GU (s)) = V − E + F is 1 − k. This characteristic includes the edges αi that cross ∂U , but neither the vertices of these edges that lie outside of U , nor faces that intersect the compliment of U . This can be seen by explicitly calculating the Euler characteristic for simple examples of subnetworks (k edges meeting at a point in U , for instance), and recalling from Euler’s theorem that χ is invariant for any planar graph contained in U with k edges that intersect ∂U . We claim that no grains in GU (s) can have less than 7 sides. First note that no topological changes occur in GU (t) during time t ∈ (0, ε). However, by the n − 6 rule, this means that any grain with n < 7 sides and area a > 0 at time s ∈ (t, t + ] would have area a + (s − t)(n − 6) > 0 at time t, which is a contradiction to the initial configuration of edges. The main claim, then, is that no network GU (s) can exist which contains only faces of more than 6 sides. The first step is an identity relating the number of edges and vertices. Since GU (s) is trivalent, we can count edges by double counting. For every v ∈ GU (s), there are three edges. We sum over the vertices to obtain the edge count, noting that every edge, except those which intersect ∂U , is counted twice. This gives us 3V + k 2E − k E= ⇒V = . (5.1) 2 3 We can also estimate E with relation to F . Supposing each face of F has at least 7 sides, the number of edges is at least 7F/2, since in the worst case every edge is 105 shared between two faces. We also have an extra k edges (those intersecting ∂U ) which are not part of any face in GU (s). This gives us 7F E≥ + k. (5.2) 2 This estimate is sufficient for the 2 and 3 ray cases. To see this, we use (5.2), (5.1), and the Euler formula to give us V −E+F =1−k (5.3) 2E − k ⇒ −E+F =1−k 3 E 2k ⇒F =1+ − 3 3 7F k 2k ≥1+ + − 6 3 3 ⇒ F ≤ 2(k − 3) For the case of 4 and 5 rays, however, we need a slightly stronger estimate. Toward this end, we can assume, without loss of generality, that U is a circle, and G(s) ∩ ∂U consists of k points p1 , . . . , pk that are ordered clockwise. If any edges corresponding to pi meet at a vertex, then we can draw a new neighborhood U 0 with k − 1 points where GU 0 (s) is a subnetwork. This process may possibly be continued until we are left with three or less points in ∂U ∩ G(s), in which case we are done. Otherwise, for k > 3, between any two points pk and pk+1 with vertices vk , vk+1 ∈ GU (s), mod k, there exists edges E1k , . . . , Elkk that connect pk and pk+1 , and are also edges of a face F k in the complement of GU (s). The sequence of edges E11 , . . . , El11 , . . . , E1k , Elkk forms a circuit whose edges do not cross that connects the vertices vi , i = 1, . . . k. Thus each edge Eji , i = 1, . . . , k, j = 1, . . . , li occurs in only one face of GU (s). This implies that we have undercounted (5.2) by at least k/2 106 Figure 5.1: Continuation through a k-degree vertex. Possible topologies of curvature flow with k-ray initial conditions, with k = 3, 4, and 5. These correspond to planar rooted trivalent tree, with the rooted vertex in each figure denoted by a dot. edges, which gives us the new estimate 7F + 3k E≥ . (5.4) 2 We now combine (5.1), (5.4), and the Euler characteristic formula to obtain V −E+F =1−k (5.5) E 2k ⇒F =1+ − 3 3 7F + 3k 2k ≥1+ − 6 3 ⇒ F ≤ k − 6. 107 We therefore search for solutions containing no faces. Such networks are precisely trees. Specifically, we seek trivalent trees with k leaves corresponding to edges con- taining pi , i = 1, . . . k. Since the endpoints of the leaves are fixed, the number Mn of such trees is equal to the planar trivalent rooted trees of n leaves, which is equal to the n − 1st Catalan number Cn−1 [22]. Specifically, M2 = 1, M3 = 1, M4 = 2, and M5 = 5. The explicit graph types are shown in Fig. 5.1. However, regardless of topological type, flowing through a four and five sided grain deletions gives the same topological transitions, as shown in Fig. 5.2. 5.3 PDMPs and grain coarsening We now demonstrate that a mean-field model of grain coarsening is a specific example of the k-species PDMP explored in [27] . The one parameter constants are as follows, where we index our species i = 2, . . . , M + 1 according to the n − gons (grains with n-sides) they describe: General Variable Value for grain coarsening M M > 6 (free parameter) M− 4 β β > 0 (free parameter) vi i−6 K (l) K 2 = 2, K (3) = 3, K (4) = 2, K (5) = 3, K int = 4 Describing our parameters in more detail: 108 • A particle of size a of species Ln = R+ , n = 2, . . . , M , corresponds to an n-gon with area a. • By the n − 6 rule, the velocity vn (a) of an n-gon is the constant vn (a) ≡ n − 6. Thus, M− = 4. • Constants K (i) , K int correspond to the number of side redistributions for each type of grain and side deletion, e.g. five sided grain deletions affect three grains, K (5) = 3. For n-gon probability weights, we impose the mean-field assumption that n-gon selection is proportional to n. To obtain a closed system, we forbid two-sided grains to lose sides, three-sided grains to lose two sides, and M −sided grains to gain a side. We will not consider one-sided grains, as they are relatively rare in actual metal networks (about .1%, see [15]). We now give the explicit weights for n-gon selection in grain coarsening. These are based on the Fradkov assumption that grains are selected based only on their side number. The weights are then    j, j ∈ {4, . . . , M }, (2) wj = (5.6)   0, j ∈ {2, 3},    j, j ∈ {3, . . . , M }, (3) wj = (5.7)   0, j = 2,    j, j ∈ {3, . . . , M }, (4) wj = (5.8)   0, j = 2, 109    j, j ∈ {3, . . . , M − 1}, (5) wj = (5.9)   0, j ∈ {2, M },    j, j ∈ {3, . . . , M − 1}, wjint = (5.10)   0, j ∈ {2, M }. Reassignments are in accordance with side redistribution under grain and side (l) deletion. Recall the upper index l = 1, . . . M− in Rkj refers to the side number of the deleted grain, k = 2, . . . , M refers to the side number of a neighboring grain before undergoing side reassignment, and j = 1 . . . , K (l) refers to the specific reassignment of the j t h grain that is reassigned sides. Refer to Fig. 5.2 for a pictorial description distribution of sides before and after grain and side deletions. Note that we disallow reassignments that create 1 and M + 1 sided grains. The explicit reassignments are then    k − 2, k ∈ {4, . . . , M }, j ∈ {1, 2}, (2) Rkj = (5.11)   0, k ∈ {2, 3},    k − 1, k ∈ {3, . . . , M }, j ∈ {1, 2, 3}, (3) Rkj = (5.12)   0, k = 2,    k − 1, k ∈ {3, . . . , M }, j ∈ {1, 2}, (4) Rkj = (5.13)   0, k = 2,      k − 1, k ∈ {3, . . . , M − 1}, j ∈ {1, 2},    (5) Rkj = k + 1, k ∈ {3, . . . , M − 1}, j = 3, (5.14)       0, k ∈ {2, M }, 110 Figure 5.2: Reassignments of sides before and after deletions. Numbers in grains refer to side numbers.      k − 1, k ∈ {3, . . . , M − 1}, j ∈ {1, 2},    int Rkj = k + 1, k ∈ {3, . . . , M − 1}, j ∈ {3, 4}, (5.15)       0, k ∈ {2, M }. We now give an explicit form to the limiting kinetic equations of [27] for grain coarsening. Assuming a continuous density of grain areas uk (x, t), we can write a system of kinetic equations for grains with k sides: ∂t uk (x, t)+(k−6)∂x uk (x, t) = hk+ k− k+ k− grain (u, t)−hgrain (u, t)+hside (u, t)−hside (u, t). (5.16) The four source terms of the right hand side of (5.16) describe, in order, addition and deletion of grains from Lk due to grain deletions, and addition and deletion from Lk due to side deletions. Explicitly, they are, for k = 2, . . . , M , (2) (3) hk,+ grain (u, t) = 8u2 (0, t)Wk+2 (t)uk+2 (x, t) + 9u3 (0, t)Wk+1 uk+1 (x, t) (5.17) (4) (5) +4u4 (0, t)Wk+1 (t)uk+1 (x, t) + 2u5 (0, t)Wk+1 (t)uk+1 (x, t) (5) +u5 (0, t)Wk−1 (t)uk−1 (x, t), 111 (2) (3) hk,− grain (u, t) = uk (x, t)[8u2 (0, t)Wk (t) + 9u3 (0, t)Wk (t) (5.18) (4) (5) +4u4 (0, t)Wk (t) + 3u5 (0, t)Wk (t)], 2βG(t) int hk,− side (u, t) = int [w uk−1 (x, t) + wk+1 uk+1 (x, t)], (5.19) Gint (t) k−1 4βG(t) int hk,+ side (u, t) = w uk (x, t), (5.20) Gint (t) k (l) for tier weights wi , wiint , total numbers G(t), Gint (t), and species selection weight (l) fraction Wi , Wiint defined in (4.22). Our main result of the previous chapter is the convergence of empirical PDMPs to the weak form of (5.16) for absolutely continuous measures μk (t) ∈ M(R+ ), obtained by pairing solutions μk with a test function φ ∈ C and integrating along R+ . The weak form of (5.16) is then, for k = 2, . . . , M , Z t hCk (t), φi − hCk (0), φi + (k − 6) hCk (s), φ0 i (5.21) 0 k+ k− k+ k− = Hgrain (u, t) − Hgrain (u, t) + Hside (u, t) − Hside (u, t), with Z t k+ (2) Hgrain (u, t) =2 Wk+2 (s)hCk+2 (s), φidF2 (s) + (5.22) 0 Z t Z t (3) (4) 3 Wk+1 (s)hCk+1 (s), φidF3 (s) +2 Wk+1 (s)hCk+1 (s), φidF4 (s) + (5.23) 0 0 Z ∞ (5) (5) 2 Wk+1 (s)hCk+1 (s), φi + Wk−1 (s)hCk−1 (s), φidF5 (s), 0 Z t k,− (2) Hgrain (u, t) =2 hCk (s), φiWk (s)dF2 (s) (5.24) 0 Z t Z t (3) (4) +3 hCk (s), φiWk (s)dF3 (s) + 2 hCk (s), φiWk (s)dF4 (s) 0 Z 0t (5) +3 hCk (s), φiWk (s)dF5 (s), 0 112 Z t k,− 2βN (s) int int Hside (u, t) = [wk−1 hCk−1 (s), φi + wk+1 hCk+1 (s), φi]ds, (5.25) 0 No (s) Z t k,+ 4βN (s) int Hside (u, t) = w hCk (s), φids. (5.26) 0 No (s) k A point of departure from the Fradkov model is the choice of side deletion parameter β. For the PDMP model, each grain is equipped with a Poisson clock that activates with average time β. When such a clock is activated, four grains are selected to change side number in accordance with rules for side flipping. This differs from Fradkov’s model, which assumes side deletion over grain deletion occurs at a constant ratio β. Graphs depicting the ratio of side to grain deletions are depicted in Fig. 5.22-5.24. . From the previous chapter, we have seen that there exists limiting subsequence that converges to a solution of (5.21). We now show that the solution retains several properties from the pre-limit PDMP. 5.4 Properties of grain coarsening One advantage to the stochastic method of proving existence of solutions to (5.21) is that we can use properties from the finite particle PDMP model to prove, with little difficulty, that the same properties hold for the hydrodynamic limit. Our main l.u. tool is the observation that since hμN −→ hμk (t), φi for all φ ∈ Cb (R+ ), we can k (t), φi − l.u. choose φ ≡ 1 to obtain the convergence of total numbers gkN (t) −−→ gk (t). Theorem 29. The following properties hold for the limiting distributions μk (t) 1. (Conservation of Polyhedral Defect). If the initial average side number 113 for grains is 6, it remains so: M X M X (k − 6)gk0 = 0 ⇒ P (t) := (k − 6)gk (t) = 0. (5.27) k=2 k=2 2. (No runoff at infinity). All loss of total number occurs from boundary events, or M X M X 5 X N (t) := gk (t) = gk0 − Fi (t). (5.28) k=2 k=2 i=2 3. (Decreasing densities). The total number is decreasing in time: N (t) ≤ N (s) for t ≥ s. (5.29) Proof. Suppose we have a limiting initial measure μ0k with zero polyhedral defect. For i ∈ N, let μkNi ,0 ∈ M(R+ ) be a sequence of initial empirical distributions of particles with zero polyhedral defect, satisfying μkNi ,0 → μ0k in M(R+ ). Constructing such a sequence is always possible. To see this, choose integers Ni → ∞ ,where PM P j Ni = Ni2 + ∙ ∙ ∙ + NiM , k k=2 (k − 6)Ni = 0, and Ni / k j Ni → gi,0 . For each n-gon in Ln , randomly select Nik particles with respect to the distribution μ0k . By the Glivenko-Cantelli theorem, with probability one the empirical distributions μkNi ,0 will converge to μ0k in M(R+ ). For the empirical measures μNi ,0 (t), it is straightforward to see conservation of polyhedral defect. This is because transport of grains does not change a grain’s side number, and one can check that any critical event has a net zero change in polyhedral defect. Thus M X P Ni (t) := (k − 6)gkNi (t) = 0. (5.30) k=2 Since gk is the local uniform limit of gkNi (t), we obtain property (1). 114 N0 For property (2), we again choose a sequence of empirical initial measures μk , → μ0k in M(R+ ), this time with no further qualifications. For the empirical densities, we have M X M X 5 X gkN (t) = N gk,0 − FiN (t). (5.31) k=2 k=2 i=2 l.u. In the previous chapter we’ve shown that FiN (t) −−→ Fi (t), so taking the local uni- form limit proves (2). Statement (3) then follows immediately from (5.31), as the cumulative distributions Fi are increasing. Our next task is to show that our system conserves area, which is a consequence of conservation of polyhedral defect. In this case, we’ll use the weak form of the limiting equations. Here, we denote the identity function id : x 7→ x. Theorem 30. (Conservation of area). For initial distributions which satisfy M X hμ0k , idi < ∞, P (0) = 0, (5.32) k=2 we have M X M X hμk (t), idi = hμ0k , idi (5.33) k=2 k=2 Proof. For r, j > 0, we use functions φr,j : R+ → R+ defined by      x x r + 1j , where ψ j : [0, 1j ] → R+ is a C 1 (R+ ) function with ψ 0 (0) = 1 and ψ 0 ( 1j ) = 0, so that φr,j ∈ C. If we sum (5.16) over grain classes, a simple (though lengthy) calculation 115 shows that M X k+ k− k+ k− Hgrain (u, t) − Hgrain (u, t) + Hside (u, t) − Hside (u, t) = 0. (5.35) k=2 Intuitively, this identity simply states that critical events change side numbers of grains, but not their areas. We are left with the relation M X M X M X Z t r,j hμk (t), φ i = r,j hμk (0), φ i − (k − 6) hμk (s), (φr,j )0 i. (5.36) k=2 k=2 k=2 0 Notice, however, that Z 1 j r,j 0 hμk (s), (φ ) i = μk (s)([0, r)) + ψ j dμk (s) → gk (s) as r, j → ∞. (5.37) 0 This implies that taking the limit as j, r → ∞ in (5.36) gives M X M X Z t hμk (t), idi = hμk (0), idi − P (s)ds. (5.38) k=2 k=2 0 Our theorem then follows from conservation of polyhedral defect. 5.5 Computations Simulations for PDMP grain coarsening were performed in Matlab. Grains over several different initial configurations and β parameters were tested. Illustrations, unless otherwise noted, are with N = 9 × 104 initial grains with a maximum of 20 sides. Uniform initial distributions have 104 n-gons, for n = 2, . . . , 10, dis- tributed uniformly in [0, 1]. At time t = 5 the number of grains drops an order of magnitude (see figure 5.4). We observe the following: 116 Mean area growth Simulations of grain networks which directly approximate mean curvature flow sug- gest that the average grain area grows linearly. (see, [10, 2, 39], for instance). In our simulations, the plot the average area hAβ (t)i at time t for the grain coarsening PDMP with varying Poisson parameters β. For small values of β, hAβ (t)i is approx- imately linear (see Fig. 5.5,5.6). The same behavior holds for average area of grains with k sides, k = 1, . . . , M . For larger values of β, coarsening becomes convex (Fig. 5.7). Distribution of n-gons For low values of β, the distribution of n-gons is approximately stationary (see Fig. 5.9) . However, as β increases, distributions of n-gons n = 3, . . . , 9 become more uniform as time increases (Fig. 5.11). This is similar to the numerical experiments of Fradkov [17]. Plots of n-gon densities for β = 0 and the Kinderlehrer model described in [26] for stationary densities (time t = 5 for both models) shows a similarity of statistics (see Fig. 5.8). Grain densities and diffusion In Fig. 5.12-5.20, we plot histograms for grain densities of 4, 6, and 8 sides under β values of .01, .1, and 1. Each bar represent a bin of 40 grains. From grain deletions, the number of bars decreases, corresponding to a decrease in number for increasing time. The densities for all n-gons drift to the right with time (see Fig. 5.21). The parameter β appears to play a diffusing role, as larger β also correspond to increased 117 variance. Side to grain deletion ratio As observed in Section 5.3, the PDMP model differs from the Fradkov model in its choice of parameter β. Thus, we will denote γβ (t) as the (time dependent) ratio of side deletions to grain deletions. For a uniform initial distribution, under differing β, γβ (t) exhibits similar concave profiles (see Fig. 5.22-5.24). 5.6 Remarks A major advantage to viewing grain growth statistics as approximations of a hydro- dynamic limit is the relative ease of implementing simulations. As we have seen, we can immediately raise several conjectures about grain behavior. For example, for β = 0, does there exist a stable attractor for the distribution of n-gon densities? Several details of the limiting densities make this question difficult. First, individual n-gon densities do not have stable attractors. Also, we still have yet to establish ex- istence for an infinite time, which would be central in an asymptotic analysis. Such long time existence was established for the Fradkov system in [19], but the existence of universal attractors in this model remains unanswered. The role of the flipping parameter β plays an interesting role in altering grain statistics. The general trend noted in simulations is that β tends to act as a diffusion agent between n-gons. The rule for side flipping is similar to the particle transfer rule in the following model, as described in [7]. Here, we are approximating the 118 dispersal of a fluid plug with initial density u0 (x, y) in a carrier fluid. If the fluid has a vertical velocity profile of V (y), the governing equations for the evolving density u(x, y, t) have the form α ut − V (y)ux = uyy (5.39) 2 u(x, y, 0) = u0 (x, y). The factor α is the diffusivity factor for the plug. Taken as a k-species PDMP, we can approximate the domain carrier fluid by a large number of horizontal tiers, and then approximately the fluid plug as particles traveling on these tiers with velocity V (y ∗ ) for a tier located at y = y ∗ . This model has no particle deletions (the fluid is conserved), but each particle is equipped with a Poisson clock of parameter α, that when set, transfers a particle to the next highest species with probability 1/2, and to the next lowest species with probability 1/2. This, rule, however, is roughly the same as side flipping. In both models, for n particles that are selected in an intermediate species, roughly n/2 will transfer upwards by one species, and n/2 will lower by a species. Thus, it is reasonable to expect the β parameter in grain coarsening to act in a similar manner as the diffusion parameter α for the fluid plug. This is, of course, not precise, as we have not actually shown the right scaling of vertical species interspacing for the fluid plug, nor have we considered how grain deletions may interact with diffusion due to side deletion. 119 Figure 5.3: Dispersal of a fluid plug. A carrier fluid is approximated by horizontal tiers, and a fluid plug (the inner region of the figure) is approximated by particles (denoted by dots) which travel horizontally on tiers corresponding to the vertical velocity profile V (y). 120 Figure 5.4: Grain number N (t). Figure 5.5: Mean area of grains with β = .01 Mean area is plotted with line of best fit y1 = .2806 + .9214x (the lines are indistinguishable). The correlation between average area and y1 is r2 = .9998. Figure 5.6: Mean area of grains with β = .1 Mean area is plotted with line of best fit y1 = .1909 + 1.0347x (the lines are almost indistinguishable). The correlation between average area and y1 is r2 = .9979. 121 Figure 5.7: Mean area of grains with β = 1. Mean area is plotted with line of best fit y1 = −.2796 + 1.6624x (mean area is convex). The correlation between average area and y1 is r2 = .9881. Figure 5.8: Comparison of side number frequency between grain coarsening PDMP (dots) at β = 1 and the Kinderlehrer model (crosses) at t = 5. Figure 5.9: Frequency of side number for β = .01 for various times. 122 Figure 5.10: Frequency of side number for β = .1 for various times. Figure 5.11: Frequency of side number for β = 1 for various times. 123 300 t=1 t=2 250 t=3 t=4 t=5 200 Frequency 150 100 50 0 0 5 10 15 Area Figure 5.12: Histograms of four-sided grain densities at β = .01. 124 300 t=1 t=2 t=3 250 t=4 t=5 200 Frequency 150 100 50 0 0 5 10 15 Area Figure 5.13: Histograms of four-sided grain densities at β = .1. 125 500 t=1 450 t=2 t=3 400 t=4 t=5 350 300 Frequency 250 200 150 100 50 0 0 5 10 15 Area Figure 5.14: Histograms of four-sided grain densities at β = 1. 126 500 t=1 450 t=2 t=3 400 t=4 t=5 350 300 Frequency 250 200 150 100 50 0 0 5 10 15 Area Figure 5.15: Histograms of six-sided grain densities at β = .01. 127 500 t=1 450 t=2 t=3 400 t=4 t=5 350 300 Frequency 250 200 150 100 50 0 0 5 10 15 Area Figure 5.16: Histograms of six-sided grain densities at β = .1. 128 350 t=1 t=2 300 t=3 t=4 t=5 250 Frequency 200 150 100 50 0 0 5 10 15 Area Figure 5.17: Histograms of six-sided grain densities at β = 1. 129 t=1 t=2 200 t=3 t=4 t=5 150 Frequency 100 50 0 0 2 4 6 8 10 12 14 16 18 20 Area Figure 5.18: Histograms of eight-sided grain densities at β = .01. 130 t=1 t=2 200 t=3 t=4 t=5 150 Frequency 100 50 0 0 2 4 6 8 10 12 14 16 18 20 Area Figure 5.19: Histograms of eight-sided grain densities at β = .1. 131 180 t=1 160 t=2 t=3 140 t=4 t=5 120 Frequency 100 80 60 40 20 0 0 2 4 6 8 10 12 14 16 18 20 Area Figure 5.20: Histograms of eight-sided grain densities at β = 1. 132 Figure 5.21: Average areas of grains with sides 2-10 at β = 0. Note how for each collection of n-gons,for n = 2, . . . , M , average area increases linearly. Similar behavior holds for β = .1 and β = 1. Figure 5.22: Side to grain deletion ratio γβ (t) at β = .01. Figure 5.23: Side to grain deletion ratio γβ (t) at β = .1. 133 Figure 5.24: Side to grain deletion ratio γβ (t) at β = 1. Chapter Six Conclusion 135 In the previous chapters, we constructed two examples of PDMPs, the dynamic shuffler and the k-species model, and proved hydrodynamic limit theorems for their particle densities. While the dynamic shuffler can be seen as an illuminating example for illustrating PDMPs, it would be inaccurate to call it a simplification of the k- species model. Unlike the k-species model, the dynamic shuffler conserves total number, and critical events change particle sizes. This mixing behavior allows us to analyze attractors. The fluid limit of the k-species model can be interpreted as a result in queueing theory for aggressive crowds. In our model, particles can be seen as customers who are served upon reaching the origin. A served customer triggers other customers to “cut” into other lines. This carries similarities to a queueing model of call centers from Kang and Ramanan [25] that analyzes fluid limits of queues with “impatient” customers, who after waiting a certain time without service, remove themselves from the model. Customers are served on a First-Come-First-Serve basis, and service times are not immediate. The weak formulations of the limiting kinetic equations (3.5) in [25] and (4.22) include Stieltjes integrals with cumulative distributions of served customers. In the k-species model, if the cumulative densities Fl (t), l = 1, . . . , M− are sufficiently smooth, then their derivatives dFl (t) = vl (0)ul (0, t)dt are described by behavior at a single point of limiting densities ul (x, t). This is why we require both L∞ and L1 estimates in the well-posedness analysis of Section 4.9. The underlying motivation for analyzing such queues, however, is an explanation of grain boundary coarsening. Direct simulation of coarsening through calculat- ing discrete curvatures involves a great deal of computation. Also, such methods inevitably face the fundamental question of flowing through a singularity, which re- mains an unresolved issue. On the other hand, mean field models assume, a priori, the existence of continuous fluid limits, and therefore close the question of measure 136 valued behavior. Interpreting the k-species model as a mean field description of grain growth carries several advantages. For instance, prelimiting empirical distributions are, by definition, measure-valued. This allows for an investigation of coarsening in a much larger class of initial data. This interpretation also simplifies proofs of certain grain statistics, such as the conservation of area and polyhedral defect, by allow us to pass limits in the uniform metric. Finally, it should be mentioned that simulation of the grain coarsening PDMP is essentially a sorting problem, and thus is not difficult to implement. The generality of the k-species model allows for a wide variety of interpretation. However, there are certainly more places to generalize. For instance, we may want to consider countable tiers, or particles that advect in domains of Rn . Also, a better understanding of well-posedness for measure-valued initial data would exploit the advantages for working with measure-valued empirical distributions in the k-species PDMP. Regularity results, along improved results on existence times, can illuminate universal grain statistics. There are several open questions that deserve attention. For instance, grain densities do not converge to stationary solutions, but certain statistics, such as the distribution of n-gons, seem to reach a steady state. For the dynamic shuffler, we have addressed universal attractors for non-arithmetic redis- tributions, but understanding limiting cyclic behavior of point masses on a lattice would give a more complete cataloguing of long time behavior for measure-valued so- lutions. These questions, among others, all share a common goal: describe complex phenomena as a limit of finite systems that obey simple rules. Bibliography [1] K. K. Anderson and K. B. Athreya, A strong renewal theorem for gen- eralized renewal functions in the infinite mean case, Probability Theory and Related Fields, 77 (1988), pp. 471–479. [2] M. Anderson, D. Srolovitz, G. Grest, and P. Sahni, Computer simu- lation of grain growth I. kinetics, Acta Metallurgica, 32 (1984), pp. 783–791. [3] T. L. Anderson, Fracture mechanics: fundamentals and applications, CRC press, 2005. [4] K. Barmak, E. Eggeling, M. Emelianenko, Y. Epshteyn, D. Kinder- lehrer, R. Sharp, and S. Taasan, An entropy based theory of the grain boundary character distribution, Discrete Contin. Dyn. Syst, 30 (2011), pp. 427– 454. [5] P. Beck, Metal Interfaces, (1952), p. 208. [6] P. Billingsley, Convergence of probability measures, vol. 493, John Wiley & Sons, 2009. [7] H. Bruus, Theoretical microfluidics, vol. 18, Oxford University Press, 2008. [8] H. Chernoff, A measure of asymptotic efficiency for tests of a hypothesis based on the sum of observations, The Annals of Mathematical Statistics, 23 (1952), pp. 493–507. [9] M. H. Davis, Piecewise-deterministic Markov processes: A general class of non-diffusion stochastic models, Journal of the Royal Statistical Society. Series B (Methodological), (1984), pp. 353–388. [10] M. Elsey et al., Diffusion generated motion for grain growth in two and three dimensions, Journal of Computational Physics, 228 (2009), pp. 8015–8033. [11] K. B. Erickson, Strong renewal theorems with infinite mean, Transactions of the American Mathematical Society, 151 (1970), pp. 263–291. [12] L. C. Evans and J. Spruck, Motion of level sets by mean curvature I, J. Diff. Geom, 33 (1991), pp. 635–681. 137 138 [13] W. Feller, Introduction to probability theory and its applications, Vol. II POD, (1974). [14] V. Fradkov, A theoretical investigation of two-dimensional grain growth in the gas approximation, Philosophical Magazine Letters, 58 (1988), pp. 271–275. [15] V. Fradkov, A. Kravchenko, and L. Shvindlerman, Experimental in- vestigation of normal grain growth in terms of area and topological class, Scripta metallurgica, 19 (1985), pp. 1291–1296. [16] V. Fradkov and D. Udler, Two-dimensional normal grain growth: topolog- ical aspects, Advances in Physics, 43 (1994), pp. 739–789. [17] V. Fradkov, D. Udler, and R. Kris, Computer simulation of two- dimensional normal grain growth (the gas approximation), Philosophical Mag- azine Letters, 58 (1988), pp. 277–283. [18] R. Henseler, A kinetic model for grain growth, doctoral thesis. [19] R. Henseler, M. Herrmann, B. Niethammer, and J. J. Vel a ´ zquez, A kinetic model for grain growth, arXiv preprint arXiv:0807.3529, (2008). [20] C. Herring, Surface tension as a motivation for sintering, in Fundamental Contributions to the Continuum Theory of Evolving Phase Interfaces in Solids, Springer, 1999, pp. 33–69. [21] M. Herrmann, P. Laurenc ¸ ot, and B. Niethammer, Self-similar solutions to a kinetic model for grain growth, Journal of Nonlinear Science, 22 (2012), pp. 399–427. [22] N. Hungerbu ¨ hler, The isomorphism problem for Catalan families, J. Com- bin. Inform. System Sci, 20 (1995), pp. 129–139. [23] J. Jacod and A. N. Shiryaev, Limit theorems for stochastic processes, vol. 288, Springer-Verlag Berlin, 1987. [24] A. Joffe and M. Me ´tivier, Weak convergence of sequences of semimartin- gales with applications to multitype branching processes, Advances in Applied Probability, (1986), pp. 20–65. [25] W. Kang et al., Fluid limits of many-server queues with reneging, The Annals of Applied Probability, 20 (2010), pp. 2204–2260. [26] D. Kinderlehrer, I. Livshits, and S. Ta’asan, A variational approach to modeling and simulation of grain growth, SIAM Journal on Scientific Comput- ing, 28 (2006), pp. 1694–1715. [27] J. Klobusicky et al., Kinetic limits of piecewise deterministic Markov pro- cesses (preprint), (2013). [28] R. Mazzeo and M. Saez, Self similar expanding solutions of the planar net- work flow, arXiv preprint arXiv:0704.3113, (2007). 139 [29] G. Menon, B. Niethammer, and R. Pego, Dynamics and self-similarity in min-driven clustering, Transactions of the American Mathematical Society, 362 (2010), pp. 6591–6618. [30] G. Menon and R. L. Pego, Approach to self-similarity in Smoluchowski’s coagulation equations, Communications on Pure and Applied mathematics, 57 (2004), pp. 1197–1232. [31] W. W. Mullins, Two-dimensional motion of idealized grain boundaries, Jour- nal of Applied Physics, 27 (1956), pp. 900–904. [32] S. Roelly-Coppoletta, A criterion of convergence of measure-valued pro- cesses: application to measure branching processes, Stochastics: An Interna- tional Journal of Probability and Stochastic Processes, 17 (1986), pp. 43–65. [33] T. Rolski, H. Schmidli, V. Schmidt, and J. Teugels, Stochastic pro- cesses for insurance and finance, vol. 505, Wiley, 2009. [34] O. C. Schnu ¨ rer and F. Schulze, Self-similarly expanding networks to curve shortening flow, arXiv preprint math/0702698, (2007). [35] R. Serfozo, Renewal and regenerative processes, in Basics of Applied Stochas- tic Processes, Springer, 2009, pp. 99–167. [36] C. V. Thompson, Grain growth and evolution of other cellular structures, Solid State Physics, 55 (2001), pp. 269–314. [37] H. Tijms, Understanding probability, Cambridge University Press, 2012. [38] J. Von Neumann, Discussion: grain shapes and other metallurgical applica- tions of topology, Metal Interfaces, (1952). [39] F. Wakai, N. Enomoto, and H. Ogawa, Three-dimensional microstructural evolution in ideal grain growthgeneral statistics, Acta Materialia, 48 (2000), pp. 1297–1311. [40] W. Whitt, Stochastic-process limits: an introduction to stochastic-process lim- its and their application to queues, Springer, 2002.