# Asymptotic Analysis of a bi-monomeric nonlinear
Becker-Döring system

Marie Doumic, Klemens Fellner, Mathieu Mezache, Juan J L Velázquez

<figure><img src="image_1.png" /></figure>

## To cite this version:

Marie Doumic, Klemens Fellner, Mathieu Mezache, Juan J L Velázquez. Asymptotic Analysis of a bi-monomeric nonlinear Becker-Döring system. Nonlinearity, 2025, 10.1088/1361-6544/adc3e5. hal-04635659v2

HAL Id: hal-04635659

https://hal.science/hal-04635659v2

Submitted on 1 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>
PAPER • OPEN ACCESS

# Asymptotic analysis of a bi-monomeric nonlinear Becker-Döring system

To cite this article: Marie Doumic et al 2025 *Nonlinearity* 38 065010

View the [article online](#) for updates and enhancements.

## You may also like

- Enhanced dissipation for the Poiseuille-Couette flow in a finite channel
Qionglei Chen and Zhen Li

- Global dynamics of an SIS epidemic model with cross-diffusion: applications to quarantine measures
Jiawei Chu and Zhi-An Wang

- Almost global existence for some Hamiltonian PDEs on manifolds with globally integrable geodesic flow
D Bambusi, R Feola, B Langella et al.
# Asymptotic analysis of a bi-monomeric nonlinear Becker–Döring system

Marie Doumic<fnref n="1" />, Klemens Fellner<fnref n="2" />,*, Mathieu Mezache<fnref n="3" /> and Juan J L Velázquez<fnref n="4" />

E-mail: klemens.fellner@uni-graz.at, marie.doumic@inria.fr, mathieu.mezache@inrae.fr and velazquez@iam.uni-bonn.de

Received 29 October 2024; revised 3 March 2025
Accepted for publication 21 March 2025
Published 21 May 2025

<figure><img src="image_3.png" /></figure>

CrossMark

Recommended by Professor Arik Yochelis

## Abstract

To provide a mechanistic explanation of sustained then damped oscillations observed in a depolymerisation experiment, a bi-monomeric variant of the seminal Becker–Döring system has been proposed in Doumic et al (2019 J. Theor. Biol. 480 241–61). When all reaction rates are constant, the equations are the following:

$$ \frac{dv}{dt} = -vw + v \sum_{j=2}^{\infty} c_j, \quad \frac{dw}{dt} = vw - w \sum_{j=1}^{\infty} c_j, $$

$$ \frac{dc_j}{dt} = J_{j-1} - J_j, \quad j \ge 1, \quad J_j = wc_j - vc_{j+1}, \quad j \ge 1, \quad J_0 = 0, $$

where v and w are two distinct unit species, and $c_i$ represents the concentration of clusters containing i units.

<figure><img src="image_4.png" /></figure>

Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

<footnote n="1">¹ Ecole Polytechnique, Inria, CNRS, Institut Polytechnique de Paris, route de Saclay, 91128 Palaiseau Cedex, France</footnote>

<footnote n="2">² Department of Mathematics and Scientific Computing, University of Graz, Heinrichstrasse 36, Graz 8010, Austria</footnote>

<footnote n="3">³ Université Paris-Saclay, INRAE, MaIAGE (UR 1404), 78350 Jouy-en-Josas, France</footnote>

<footnote n="4">⁴ University of Bonn, Institute for Applied Mathematics, Endenicher Allee 60, D-53115 Bonn, Germany</footnote>

<footnote>* Author to whom any correspondence should be addressed.</footnote>
We study in detail the mechanisms leading to such oscillations and characterise the different phases of the dynamics, from the initial high-amplitude oscillations to the progressive damping leading to the convergence towards the unique positive stationary solution. We give quantitative approximations for the main quantities of interest: period of the oscillations, size of the damping (corresponding to a loss of energy), number of oscillations characterising each phase. We illustrate these results by numerical simulation, in line with the theoretical results, and provide numerical methods to solve the system.

Keywords: polymerisation-depolymerisation reactions, oscillations, asymptotic analysis, Lotka–Volterra system, drift-diffusion

Mathematics Subject Classification numbers: 34E10, 34C41, 92C42, 92D25, 92E20

Contents

<table><tr><td>1. Introduction</td><td>3</td></tr><tr><td>2. Coupling LV and Becker–Döring systems</td><td>4</td></tr><tr><td>2.1. A first overview of the dynamics</td><td>4</td></tr><tr><td>2.2. Study of the unperturbed LV system</td><td>7</td></tr><tr><td>2.3. Cluster size distribution dynamics along one LV cycle, initial time</td><td>9</td></tr><tr><td>2.4. Main results and sketch of the dynamics phases</td><td>11</td></tr><tr><td>3. Phases I and II: from one LV cycle to the next</td><td>12</td></tr><tr><td>3.1. An approximate transport and diffusion equation</td><td>17</td></tr><tr><td>3.2. From the initial cycle to the next</td><td>18</td></tr><tr><td>3.3. Energy decay and iterative mapping</td><td>22</td></tr><tr><td>4. Phases I and II: overall dynamics</td><td>25</td></tr><tr><td>4.1. Early Phase I: initial increase of the characteristic length</td><td>25</td></tr><tr><td>4.2. Continuation of Phase I and Phase II: size distribution evolution</td><td>25</td></tr><tr><td>4.3. Overall dynamics during phases I and II</td><td>29</td></tr><tr><td>5. Phase III: energy reduction from order ε to order ε<sup>2</sup></td><td>31</td></tr><tr><td>5.1. Early Phase III: matching condition</td><td>31</td></tr><tr><td>5.2. Dynamics of Phase III: a new scaling</td><td>32</td></tr><tr><td>5.3. End of Phase III: energy decay</td><td>36</td></tr><tr><td>6. Phase IV: oscillations decay and stabilization</td><td>39</td></tr><tr><td>6.1. End of Phase III and early Phase IV: decay of the oscillations</td><td>39</td></tr><tr><td>6.2. Phase IV: stabilization by means of a parabolic equation</td><td>43</td></tr><tr><td>7. Discussion</td><td>45</td></tr><tr><td>7.1. Other choices of initial concentrations</td><td>45</td></tr><tr><td>7.2. Concluding remarks</td><td>46</td></tr><tr><td>Data availability statement</td><td>48</td></tr><tr><td>Appendix A. Extra computations for the LV system</td><td>48</td></tr><tr><td>Appendix B. Linearised stability around the positive steady state</td><td>52</td></tr><tr><td>Appendix C. Numerical computations</td><td>55</td></tr><tr><td>C.1. The continuous model</td><td>55</td></tr><tr><td>C.2. The discrete model</td><td>58</td></tr><tr><td>References</td><td>58</td></tr></table>
# 1. Introduction

We are interested in describing the damped oscillations of the following chemical reaction system, first introduced in [10]:

$$ 
\begin{aligned}
\mathcal{V} + \mathcal{W} & \rightarrow 2\mathcal{W}, \\
\mathcal{W} + C_j & \rightarrow C_{j+1}, \quad i \ge 1, \\
C_j + \mathcal{V} & \rightarrow C_{j-1} + 2\mathcal{V}, \quad i \ge 2.
\end{aligned} 
\quad (1) $$

In this chemical reaction network, clusters containing $i$ monomers, represented by the species $C_i$, with $i \ge 1$, grow by polymerisation events adding a $\mathcal{W}$-monomer and shrink by catalytic depolymerisation induced by $\mathcal{V}$-monomers.

The goal of introducing this model is to explain sustained, though damped, oscillations, experimentally observed during the time-course of protein fibrils depolymerisation experiments [10, 12], which are also displayed when simulating (1), see [10]. These experiments raised significant interest since the classical Becker-Döring model for the polymerisation/de-polymerisation of polymers features a well-known Ljapunov functional, which conflicts with sustained oscillations. The need to propose a new reaction system like (1) was therefore pointed out by the experimentalists in the hope to gather new insights into the mechanics of prion dynamics. The question of how variants of standard aggregation-fragmentation systems can give rise to persistent or persistent then damped oscillations has attracted increasing interest in recent years. For example, another modification of the Becker-Döring system has been proposed in [19, 20]. This system leads to sustained oscillations due to a very different mechanism, namely the atomisation of large polymers into monomers, see also [7, 13]. Such rapid shortening events also appear to be responsible for damped oscillations observed numerically in a model of microtubule dynamics [16].

Notice that there are two remarkable differences between the model (1) and the classical Becker-Döring model [4], as it can be found for instance in [2, 17, 23, 24]. On the one hand, there are two different types of monomers $\mathcal{V}$ and $\mathcal{W}$, whereas $C_1$ represents the smallest cluster, which is not necessarily a monomeric species. This leads to the conservation of the total number of clusters, a quantity which is not preserved in the Becker-Döring system.

As a second difference, the reaction yielding the depolymerisation of a cluster of size $j$ into a $\mathcal{V}$-monomer and a cluster of size $(j-1)$ is not modelled as a linear process, but as a nonlinear reaction catalyzed by a $\mathcal{V}$-monomer itself. We emphasise that the presence of several types of oligomers [1] and monomers as well as the possibility of the uncommon catalytic depolymerisation reaction have been then confirmed by other observations [22].

When assuming constant reaction rates, we obtain the following model:

$$ \frac{dv}{dt} = -vw + v \sum_{j=2}^{\infty} c_j, \quad v(0) = v^0, \quad (2) $$

$$ \frac{dw}{dt} = vw - w \sum_{j=1}^{\infty} c_j, \quad w(0) = w^0, \quad (3) $$

$$ \frac{dc_j}{dt} = J_{j-1} - J_j, \quad J_j = wc_j - vc_{j+1}, \quad c_j(0) = c_j^0, \quad j \ge 1, \quad J_0 = 0. \quad (4) $$

In this system, $v(t)$ and $w(t)$ denote the concentrations of the two monomeric species $\mathcal{V}$ and $\mathcal{W}$ and $c_j(t)$ represents the concentrations of $C_j$, the polymers, clusters or aggregates containing $j$ units/monomers at time $t$. System (2)-(4) is the particular case of the system (1) obtained when assuming that the reaction rates are size-independent.
A first numerical and theoretical study has been carried out in [10]. It gave necessary and sufficient conditions for the existence of a positive steady state and studied the stability of boundary steady states (proposition 4 in [10]). Numerical simulations and comparison with well-known oscillatory systems (subsection 5.3. in [10]) provided evidence for sustained damped oscillations converging towards the positive steady state, when it exists.

The aim of the present study is to characterise in full detail the mechanisms driving the system from pronounced oscillations to equilibrium. The results will be justified using formal asymptotic expansions (see for instance [14, 18] or [5]), complemented with numerical simulations yielding results analogous to those obtained in the asymptotic formulae. Some rigorous mathematical results that support the scenario for the solutions described in this paper are proved in the companion paper [11].

The plan of this paper is the following. In section 2, we describe in detail the model studied and provide some important notations and lemmas which are used throughout the paper. In section 2.2, we summarize the form of the solutions of a rescaled version of the Lotka–Volterra (LV) model, which we will identify as approximation to the behaviour of the monomeric concentrations $(v, w)$, and which will be shown to drive the oscillatory behaviour of the solutions of (2)–(4) (section 2.3). Finally, in section 2.4, we sum-up the main results and describe the processes that evolve the initial distribution of clusters to the final equilibrium distribution. We have found it convenient to decompose this evolution into four different Phases I, II, III, and IV, which are characterised by the range of values of the LV energy associated to the monomer concentrations $(v, w)$ as well as the range of a typical cluster size $k$ containing most of the mass of the cluster distribution $\{c_k\}_{k \in \mathbb{N}}$.

Section 3 describes in detail the evolution of the cluster concentrations $\{c_k\}$ during one LV cycle of the concentrations $(v, w)$. This analysis is relevant throughout the phases that we denote as Phases I and II. The evolution of the concentrations $\{c_k\}$ along several LV cycles, and the slow changes occurring from one cycle to the next during Phases I and II, is then described in section 4. Phases III and IV are described in sections 5 and 6 respectively. Section 7 summarises and discusses the main conclusions of this article as well as possible generalisations. Finally, three appendices at the end of the paper some detailed calculations about the form of the solutions of the LV problem, on analysis of the linearised problem around the steady states of (2)–(4) and a description of the numerical methods we used.

**Notation 1.1.** We will use extensively the classical asymptotic notation $f \sim g$ as $x \to x_0$ to indicate that $\lim_{x \to x_0} \frac{f(x)}{g(x)} = 1$. When $x_0$ is not specified it means that $f \sim g$ when $\varepsilon \to 0$. We will use also the notation $f \ll g$ as $x \to x_0$ to indicate that $\lim_{x \to x_0} \frac{f(x)}{g(x)} = 0$. We will use the notation $f \simeq g$ to indicate informally that the two functions $f$ and $g$ can be expected to be approximately equal for a suitable range of values. We will use the notation $f \approx g$ to indicate that the functions $f$ and $g$ are of the same order of magnitude for some range of values of their argument. More precisely, there exists a constant $C > 1$ such that $\frac{1}{C}f \le g \le Cf$. We will use the notation $f \lesssim g$ to indicate that there is a constant $C > 0$ such that $f \le Cg$.

## 2. Coupling LV and Becker–Döring systems

### 2.1. A first overview of the dynamics

When looking at (2)–(4), we first observe that there are two conserved quantities, namely

$$M = v + w + \sum_{j=1}^{\infty} j c_j, \quad \varepsilon = \sum_{j=1}^{\infty} c_j. \qquad (5)$$
The first conserved quantity represents the total mass $M > 0$ whereas the second is the number of clusters. This second conservation law is an important difference between our system and the classical Becker-Döring system, in which the smallest clusters and monomers are identical. Our study will show its importance to the dynamical behaviour: Since our system preserves the number of clusters, the question is how its size distribution evolves with time.

We remark that we can assume $M = 1$ without loss of generality because the equations (2)–(4) are invariant under the change of variables $v \to Mv$, $w \to Mw$, $c_j \to Mc_j$, and $\varepsilon \to M\varepsilon$, which transforms the equations (5) in

$$1 = v + w + \sum_{j=1}^{\infty} j c_j, \quad \varepsilon = \sum_{j=1}^{\infty} c_j. \qquad (6)$$

Notice that since $v \ge 0$ and $w \ge 0$ these equations imply that $\varepsilon \le 1$ and the identity $\varepsilon = 1$ takes place only if $v = w = 0$ and $c_j = \delta_{j,1}$.

The aim of our study is to quantitatively analyse the asymptotic behaviour of the solutions of the problem (2)–(4), (6) for $\varepsilon \to 0$. The smallness of $\varepsilon$ may have two biological interpretations: either only few polymers exist initially and the initial mass lies mainly in the two monomeric species $v$ and $w$, or the mass is distributed among monomers and clusters with a very large average cluster size $i_M$ so that $i_M(0) := \sum_{j \ge 1} i_{c_j}(0) = O(1/\varepsilon)$. This last case is typical of amyloid fibrils depolymerisation or fragmentation experiments [3, 26].

Using the conservation of the number of clusters, we can rewrite (2)–(4) in a simpler form:

$$\frac{dv}{dt} = -vw + v(\varepsilon - c_1), \qquad (7)$$

$$\frac{dw}{dt} = vw - \varepsilon w, \qquad (8)$$

$$\frac{dc_j}{dt} = J_{j-1} - J_j, \quad \forall j \ge 1, \quad J_0 = 0, \quad J_j = wc_j - vc_{j+1}, \quad j \ge 1. \qquad (9)$$

The global-in-time existence and uniqueness of a nonnegative solution to this system has been proved in [10] (theorem 2), under the assumption that the initial nonnegative state $(v^0, w^0, c_j^0)$ satisfies

$$v^0 + w^0 + \sum_{j=1}^{\infty} j^2 c_j^0 < \infty.$$

There is a family of boundary steady states for which $w=0$. We can have either $v=0$ and an arbitrary distribution of concentrations $c_j$ satisfying the constraint $\sum_{j \ge 1} j c_j = 1$, or alternatively, $v = 1 - \varepsilon$, $c_1 = \varepsilon$, and $c_j = 0$ for $j \ge 2$. All the boundary steady states happen to be unstable as soon as $\varepsilon < 1/2$, what is assumed from now on (since we even assume $\varepsilon \ll 1$), see proposition 4 in [10]. Existence and uniqueness of a unique positive steady state (corollary 4 in [10]) for $\varepsilon < 1/2$ is given by

$$\theta := 1 - \frac{1}{2\varepsilon} \left( 1 - \sqrt{(1 - 2\varepsilon)^2 + 4\varepsilon^2} \right) \sim \varepsilon, \qquad (10)$$
$$\bar{c}_1 := \varepsilon\theta \sim \varepsilon^2, \quad \bar{c}_i := (1 - \theta)^{i-1}\bar{c}_1, \quad \bar{v} := \varepsilon, \quad \bar{w} := \varepsilon - \bar{c}_1.$$

In appendix B, we prove that this positive steady state is locally linearly stable.
The system (7)–(8) for $(v, w)$ may be viewed as a $c_1$-perturbation of a LV system, which would permanently oscillate in case $c_1 \equiv 0$. The unperturbed LV system, i.e. (7)–(8) with $c_1 \equiv 0$, is governed by the Hamiltonian

$$E = v + w - 2\varepsilon - \varepsilon \log \left( \frac{vw}{\varepsilon^2} \right). \qquad (11)$$

The hamiltonian structure of the (unperturbed) LV system satisfied by $(v, w)$ may be seen by the transformation $(p, q) = (\log(\frac{v}{\varepsilon}), \log(\frac{w}{\varepsilon}))$, which satisfies the system $(p', q') = (-\partial_q E, +\partial_p E)$. It features the steady state $v = w = \varepsilon$, which corresponds to the minimum of the Hamiltonian $E=0$. All other trajectories of the unperturbed LV system (7)–(8) with $c_1 \equiv 0$ are closed periodic orbits for positive values of $E$ and we define $T(E)$ as the associated time period.

For the full system (7)–(9), in which $v, w$ interact with the clusters $c_j$ and vice versa, the function $E(t)$ is neither a conserved quantity nor a Ljapunov functional. The absence of a suitable entropy is one difficulty in the study of (7)–(9). More precisely, the following formula for the change of the energy holds

$$\frac{dE}{dt} = (\varepsilon - v)c_1. \qquad (12)$$

As we shall demonstrate, solutions to the full system (7)–(9) feature approximate LV cycles, where $v$ and $w$ oscillate around the steady state values of the full system $\bar{v} = \varepsilon$ and $\bar{w} = \varepsilon - \bar{c}_1$; see lemmas 2.2 to 2.6 for a full description of the dynamics of one LV cycle. Hence we see from (12) that the sign of $\frac{dE}{dt}$ changes over every time period. The energy level $\bar{E}$ of the steady state of the full system is computed from (10):

$$\bar{E} = \bar{v} + \bar{w} - 2\varepsilon - \varepsilon \log \left( \frac{\bar{v}\bar{w}}{\varepsilon^2} \right) = -\bar{c}_1 - \varepsilon \log \left( 1 - \frac{\bar{c}_1}{\varepsilon} \right) \sim \varepsilon \frac{\theta^2}{2} \sim \frac{\varepsilon^3}{2}. \qquad (13)$$

We expect that solutions to our model system (7)–(9), which depart from $E \approx 1$, will slowly decrease the mean value of the energy $E$ along successive approximate LV cycle, each of these cycles being characterised by a given level of energy $E$. Quantifying precisely this trend is the aim of this article. We will show that the decay-on-average of the values of $E$ takes place by means of an intricate mechanism, in which the concentrations $c_j$ oscillate in cluster sizes $j$ with a transport speed given by $w - v$.

For small perturbations of the steady state (10), we expect $c_1$ to remain in the order of $\bar{c}_1 \sim \varepsilon^2$. Else, we can expect $c_1$ to be very small, except for the times in which $\langle j \rangle$, the mean value of the concentrations $c_j$, reaches its minimum values. This happens when depolymerisation dominates, namely for $v \gg \varepsilon$. During these times, even if $c_1$ still remains small, it modifies $E$ the most. This explains why we expect that in average $\Delta E = \int_t^{t+T(t)} \frac{d}{ds}E(s)ds$, where $T(t)$ is the time (pseudo) period, to be negative.

In fact, the values of $c_1$ remain small during most times, except for a possible initial transient state (for instance if the system departs from $c_j = \varepsilon\delta_1(j)$). Therefore, the change of $E$ will be small in each cycle of the variables $v, w$, but eventually its value will decrease and approach the minimum value $E \sim \frac{\varepsilon^3}{2}$ after sufficiently long times.

A full description of this slow approach is the aim of this article. To do so, we first focus on the dynamics during one (approximately unperturbed) LV cycle. In section 2.2 we begin by analysing the different stages of the dynamics of $(v, w)$ over one LV cycle (lemmas 2.2 to 2.6). This done, in section 2.3 we can focus on the size distribution dynamics, which appears close to a drift-diffusion equation with coefficients given in terms of the previously studied $(v, w)$.
This strategy is the basis of our asymptotic long-time study and explains how the perturbation slowly modifies the LV cycles; we sketch it in section 2.4.

## 2.2. Study of the unperturbed LV system

According to our above argumentation, $c_1$ is expected to be very small, i.e. $c_1 \ll 1$, during most of the evolution; it will be shown a posteriori that this ansatz is self-consistent. In this section, we thus study in detail the dynamics over one cycle of $(v, w)$ solutions to the unperturbed LV system, i.e. (7)–(8) with $c_1 \equiv 0$, namely

$$ \frac{dv}{dt} = -vw + v\epsilon, \quad \frac{dw}{dt} = vw - \epsilon w, \quad w(0) = v(0) \ge \epsilon, \qquad (14) $$

where the last condition defines a convenient starting- and end-point of a LV cycle. Notice that the steady state of (14) is at

$$ v = w = \epsilon, \quad E = 0. \qquad (15) $$

We recall that the energy $E$ defined by (11) is constant along the solutions of (14): the trajectories associated to $(v, w)$ solution of (14) are the level lines of the function $E$. The solutions $(v, w)$ are periodic, oscillating around the equilibrium value (15).

In order to describe the behaviour of the solutions of (14) we first use some rescaling properties of the solutions, summed-up in the following lemma.

**Lemma 2.1 (Scaling properties).** Denote by $(v, w)(t; E, \epsilon)$ the unique solution to (14) with $v(0) = w(0) > \epsilon$ and $E$ defined by (11), that is

$$ E = 2(v(0) - \epsilon) - 2\epsilon \log \left( \frac{v(0)}{\epsilon} \right). \qquad (16) $$

Note that (16) is strictly monotone increasing in $v(0) > \epsilon$ and that trajectories of (14) are uniquely defined in terms of values $E > 0$. Let $T(E, \epsilon)$ denote the period of the corresponding LV cycle. By a scaling argument we have the following identities:

$$ (w(t; E, \epsilon), v(t; E, \epsilon)) = \left( Ew \left( Et; 1, \frac{\epsilon}{E} \right), Ev \left( Et; 1, \frac{\epsilon}{E} \right) \right) \qquad (17) $$

and

$$ T(E, \epsilon) = \frac{T\left(1, \frac{\epsilon}{E}\right)}{E}. \qquad (18) $$

A large part of this paper will be concerned with the dynamics of the solutions during the range of times in which $\frac{\epsilon}{E} \ll 1$. Due to (17) and (18), in order to compute the asymptotic behaviour of $w$, $v$, $T$ for this range of values of $\epsilon$ and $E$, it is enough to compute their asymptotic behaviour for $E=1$ and $\epsilon \to 0$. Given $E=1$, we first note that (16) implies

$$ v(0) = w(0) \sim \frac{1}{2} + \epsilon \log \left( \frac{1}{\epsilon} \right) + \epsilon (1 - \log(2)) + O(\epsilon^2 \log(\epsilon)) \quad \text{as } \epsilon \to 0. \qquad (19) $$

Figure 1 depicts the $(v, w)$ phase portrait of a LV cycle and the associated time evolution of the concentrations. We divide a LV cycle into five stage subject to the following five lemmas, each describing the corresponding contribution to the period $T$. The proofs of these lemmas by standard asymptotic expansions are carried out in appendix A.
<figure><img src="image_5.png" /><figcaption>Figure 1. Numerical solution of the Lotka–Volterra system (14) for $\varepsilon = 0.02$ and $v(0) = w(0) = \frac{1}{2}$. To avoid numerical dissipation we rather solved the system for $(\log(v), \log(w))$, see appendix C.1, System (144). Left: phase portrait for $(v, w)$ illustrating the different stages of the trajectory over one period $T$. The dashed purple line is the subset where $v = w$, the dotted blue lines are $w = \varepsilon$ or $v = \varepsilon$. Right: time-evolution of the solution $(v, w)$. The interval $t_2 - t_3$ is the slowest timescale, all the other time intervals are faster, the fastest being $t_0 - t_1$ and $t_4 - t_5$.</figcaption></figure>

**Lemma 2.2 (First stage: fast decay of $v$ from $v(0)$ to $\varepsilon$, increase of $w$).** Let $E$ be given by (11) and $\varepsilon \ll 1$. Let $(v, w)$ be the unique periodic solution to (14) corresponding to $E = I$, i.e. subject to initial data at time $t_0 = 0$ asymptotically characterised in (19). We define $t_1$ by

$$t_1 := \min_{t \ge 0} \{t, v(t) = \varepsilon\}.$$

On $(0, t_1)$, $w$ increases, $v$ decreases, and as $\varepsilon \to 0$ the solution is approximated by

$$\begin{align}
v(t) &= 1 - \frac{1}{1 + e^{-t}} + O(\varepsilon e^{-t} \log(\varepsilon) t), \nonumber \\
w(t) &= \frac{1}{1 + e^{-t}} (1 + O(\varepsilon \log(\varepsilon) t)), \quad t_1 \sim -\log(\varepsilon), \nonumber \\
w(t_1) &= 1 - \varepsilon \log(\varepsilon) + \varepsilon + o(\varepsilon), \quad v(t_1) = \varepsilon. \tag{20}
\end{align}$$

**Lemma 2.3 (Second stage: decay of $v$ and $w$).** Under the assumptions of lemma 2.2, we define

$$t_2 := \min_{t \ge t_0} \{t, w(t) = \varepsilon\}.$$

On $(t_1, t_2)$, $w$ and $v$ decrease, and as $\varepsilon \to 0$ we have

$$t_2 - t_1 \sim \frac{1}{\varepsilon} \log \left( \frac{1}{\varepsilon} \right). \tag{21}$$

**Lemma 2.4 (Third stage: $v$ increases, $w$ decreases, minimum of $v + w$).** Under the assumptions of lemma 2.2, we define

$$t_3 := \min_{t \ge t_2} \{t, v(t) = \varepsilon\}.$$
On $(t_2, t_3)$, $w$ decreases, $v$ increases to $\varepsilon$, and as $\varepsilon \to 0$ we have

$$t_3 - t_2 \sim \frac{1}{\varepsilon^2}, \quad w \sim \varepsilon \exp(-\varepsilon(t-t_2)), \quad v \sim v(t_2) e^{-1} \exp(\varepsilon(t-t_2)). \quad (22)$$

**Lemma 2.5 (Fourth stage: v and w increase, symmetric to Stage 2).** Under the assumptions of lemma 2.2, we define

$$t_4 := \min_{t \ge t_3} \{t, w(t) = \varepsilon\}.$$

On $(t_3, t_4)$, $w$ and $v$ increase, and as $\varepsilon \to 0$ the time interval is approximated by

$$v \sim \varepsilon e^{\varepsilon(t-t_3)}, \quad t_4 - t_3 \sim \frac{1}{\varepsilon} \log \left( \frac{1}{\varepsilon} \right). \quad (23)$$

**Lemma 2.6 (Fifth stage: v decreases, w increases, end of the cycle).** Under the assumptions of lemma 2.2, we define

$$t_5 := \min_{t \ge t_4} \{t, v(t) = w(t) > \varepsilon\}.$$

On $(t_4, t_5)$, $w$ increases, $v$ decreases, and as $\varepsilon \to 0$ the time interval is approximated by

$$t_5 - t_4 \sim \log \left( \frac{1}{\varepsilon} \right). \quad (24)$$

We remark that the period $T$ is given as $T = t_5 - t_0 = t_5$. We also notice that the third stage is the longest, giving the total cycle its asymptotic length. During this stage almost nothing visible happens since both $v$ and $w$ are exponentially small. By contrast, the first and last stages, where $v$ and $w$ are largest, are extremely fast. These considerations drive all the following analysis, where we distinguish between a transport part (first and last stages of the cycle) and a transport-diffusion part (second to fourth stages). This is made clearer in the following section where we analyse the dynamics of the cluster size distribution along one LV cycle.

## 2.3. Cluster size distribution dynamics along one LV cycle, initial time

In section 2.2, we have described the dynamics of $(v, w)$ along one LV cycle, neglecting the perturbation given by $c_1$. Let us now reason similarly for the cluster size distribution: assuming that $(v, w)$ is given by the unperturbed LV system, how does the cluster distribution evolve?

At this point, we remark that the equations (9) for the evolution of the clusters can be rewritten as

$$\frac{dc_j}{dt} = \left(\frac{w-v}{2}\right)(c_{j-1} - c_{j+1}) + \left(\frac{w+v}{2}\right)(c_{j-1} + c_{j+1} - 2c_j), \quad j \ge 2, \quad (25)$$

$$\frac{dc_1}{dt} = vc_2 - wc_1. \quad (26)$$

For given $v$ and $w$, the equations (25) for $j \ge 2$ constitute discrete-in-size convection-diffusion equations while equation (26) for the smallest cluster concentration $c_1$ forms a kind of dynamical boundary condition. The approximation of discrete-in-size polymerisation/depolymerisation by a continuous-in-size convection-diffusion PDE is well-known in the literature, see e.g.
[24]. If we assume that the concentrations $c_j$ change slowly in the variable $j$ (in a sense that will be precised later), we can read (25)–(26) as a discretized version of the PDE

$$ \frac{\partial c}{\partial t}(j,t) + \tilde{V}(t) \frac{\partial c}{\partial j}(j,t) = \frac{d(t)}{2} \frac{\partial^2 c}{\partial j^2}(j,t). \quad (27) $$

where we introduce

$$ \tilde{V}(t) = (w - v)(t), \quad d(t) = (w + v)(t). \quad (28) $$

In particular, we approximate the discrete-in-size cluster concentrations $c_j(t)$ by a smooth continuous-in-size distribution function $c(j,t)$ where $j \in (0, \infty)$ denotes now—with a slight abuse of notation—a continuum size variable. Hence, we assume

$$ c_j(t) \approx c(j,t), \quad j \in \mathbb{N}, \quad |c_{j+1} - c_j| \ll 1, \quad \left| \frac{\partial c}{\partial j} \right| \ll 1. \quad (29) $$

The PDE (27) constitutes a valid approximation for sufficiently large $j$, [24]. For cluster sizes $j$ of order one the approximation (27) breaks down and in order to calculate the concentrations $c_j$ we need to use the system of equations (25)–(26). Hence, we denote by the outer layer the range in size where the approximation (27) holds and by the boundary layer the range where we have to consider the system (25)–(26).

In the following, we assume that there exists initially a characteristic cluster size, denoted as $L_0 := \sum_{j=1}^{\infty} j c_j = \frac{1-v-w}{\varepsilon} \gg 1$ (and $L_n$ for the nth following LV cycle), around which the concentrations $c_j$ give a relevant contribution to the cluster distribution $\{c_j\}_{j \in \mathbb{N}}$. More precisely, we assume at $t=0$ when $v=w > \varepsilon$ that the initial concentrations are given approximately by means of

$$ c_j(0) \approx c(j,0) \approx \frac{\varepsilon}{L_0} \varphi_0 \left( \frac{j}{L_0} \right), \quad j \in (0, \infty). \quad (30) $$

We assume that $\varphi_0$ is a measure defined in $[0, \infty)$ which can contain a Dirac delta in $x=0$, modelling the fact that the mass contained in $c_1$ may be of the same order of magnitude as the mass contained in large clusters: we thus define

$$ \varphi_0(x) = \psi_0(x) + m_0 \delta(x), \quad (31) $$

where $\psi_0$ is a smooth function and supported in $(0, \infty)$ and $m_0 \ge 0$ the mass contained in the Dirac delta. Thanks to (6) and w.l.o.g. we consider $M=1$, which implies $\sum_{j=1}^{\infty} j c_j = O(1)$ as well as $\sum_{j=1}^{\infty} c_j = \varepsilon$. For concentrations of the form (30) it follows that we have $L_0 = O(\frac{1}{\varepsilon})$. Assuming $\varepsilon \ll 1$ in (6), we typically have $L_0, L_n \gg 1$. We also impose

$$ \int_0^\infty \varphi_0(dx) = \int_0^\infty x \varphi_0(dx) = 1, \quad (32) $$

which defines uniquely $L_0$ by $L_0 = \frac{\sum j c_j}{\varepsilon} \in [1, \frac{1}{\varepsilon})$. Note that we understand $\int_0^\infty = \int_{[0, \infty)}$, i.e. we include the Dirac mass at 0.

The heuristic considerations in this section provides us with the qualitative behaviour of the clusters along one LV cycle: their size distribution is transported at the speed $\tilde{V}(t)$, and is dif-fused at the diffusion rate $\frac{d(t)}{2}$. We remark for instance that if $w > v$, we have that $\tilde{V}(t) > 0$ (see figures 4(c) and 3(a) below as part of a more detailed discussion). In particular the convective
term in (27) shifts the distribution of clusters towards larger values of j (see figures 4(d) and 3(b)). Similar to lemma 2.1 in section 2.2, we obtain the following lemma 2.7, which collects scaling properties concerning the transport and the total diffusion along one LV cycle.

**Lemma 2.7.** With the notations and assumptions of lemma 2.1, let us define the characteristic curve $Y(t; E, \varepsilon)$ and the total diffusion $D(E, \varepsilon)$ by

$$Y(t; E, \varepsilon) := \int_{0}^{t} (w - v)(s; E, \varepsilon) \, ds, \quad D(E, \varepsilon) := \frac{1}{2} \int_{0}^{T(E, \varepsilon)} d(s) \, ds. \quad (33)$$

We have the relations

$$Y(t; E, \varepsilon) = Y(Et; 1, \frac{\varepsilon}{E}), \quad D(E, \varepsilon) = D(1, \frac{\varepsilon}{E}). \quad (34)$$

The next lemma gathers results concerning $D$, $Y$, $T$ and the average values of $v$ and $w$.

**Lemma 2.8.** Under the assumptions of lemma 2.2, we have

$$T(1, \varepsilon) \sim \frac{1}{\varepsilon^2}, \quad \int_0^T v(s) \, ds = \int_0^T w(s) \, ds = \varepsilon T \sim \frac{1}{\varepsilon}, \quad D(1, \varepsilon) \sim \frac{1}{\varepsilon}. \quad (35)$$

Let us denote $Y(t) = Y(t; 1, \varepsilon)$ and $D = D(1, \varepsilon)$. We have the following behaviour of $Y$ and $D$ along one LV cycle.

* On $(0, t_2)$, $Y$ is increasing.

* The main contribution to $D$ is due to $w$ in the interval $(t_1, t_2)$ (see figure 3(a)) and $v$ in $(t_3, t_4)$, (see figure 3(e)). Elsewhere, the contribution is negligible.

* On $(t_2, t_3)$, $Y$ is almost constant, its largest value $Y_{\max} \sim \frac{1}{\varepsilon}$ being reached at $t_2 + \frac{t_3 - t_2}{2}$,

* On $(t_3, t_5)$, $Y$ is decreasing.

## 2.4. Main results and sketch of the dynamics phases

A main goal of this paper is to determine the evolution of the characteristic cluster length $L$ and its approximate cluster size distribution $\varphi$ over time; In particular for the nth LV cycle, we will denote by $(L_n, \varphi_n)$ the corresponding typical length and size distribution. It will turn out that after the variables $v$, $w$ have approximately followed several LV cycles (14), the energy $E$ changes and the concentrations $\{c_j\}_{j \in \mathbb{N}}$ can be approximated by means of a formula of the form (30), for new choices of the length $L_n$ and the function $\varphi_n$.

We now assume that the initial length $L_0$ is much smaller than $\frac{1}{\varepsilon}$. Otherwise, if $L_0$ is of order $\frac{1}{\varepsilon}$ for the initial distribution of clusters, then the first of the four phases described below is skipped by the dynamics. We also remark that $L_0 \ll \frac{1}{\varepsilon}$ in (6) implies that the sum $\sum_{j=1}^{\infty} j c_j$ is much smaller than 1 for small $\varepsilon$. Thus, $(v+w) \simeq 1$ for small $\varepsilon$ and we can assume $E=1$ and use the analysis carried out in section 2.2, lemmas 2.2 to 2.5. Actually, this approximation holds, not only for the initial time, but whenever we have $L \ll \frac{1}{\varepsilon}$.

In our case that initially $L_0 \ll \frac{1}{\varepsilon}$, the evolution of the concentrations towards the equilibrium distribution undergoes four different phases. In each phase a different mechanisms governs the evolution of the concentrations $\{c_j\}_{j \in \mathbb{N}}$, $v$, $w$ and the energy $E$ and the characteristic length $L$. We describe the leading magnitudes associated to these phases as well as the mechanism yielding the dissipation of energy $E$. At the end of these four phases the concentrations $\{c_j\}_{j \in \mathbb{N}}$, $v$, $w$ approach to an equilibrium. Notice that the approach to equilibrium takes place by means of
an involved procedure and we do not see (at least currently) that it can be captured by means
of a Ljapunov function.

*   **Phase I: The energy *E* remains near a constant of order 1.** The period $T(E, \varepsilon)$ of the LV oscillations is of order $\frac{1}{\varepsilon^2}$ (see section 2.2). The evolution of the cluster distributions takes place by means of a combination of an oscillatory motion of the cluster distribution towards larger values of $j$ and backwards to clusters with $j$ of order one, combined with a spreading of the cluster concentration distribution in the space of cluster sizes (see section 2.3). The distribution spreads an amount of order $\frac{1}{\sqrt{\varepsilon}}$ in the space of cluster sizes in each LV cycle, and if we assume that initially $L = L_0 \approx 1$ we require a number of cycles $n$ of order $\frac{1}{\varepsilon}$ to arrive to $L = L_{\text{fin}} \approx \frac{1}{\varepsilon}$ (see section 4.3, proposition 4.3). From this point on, the order of magnitude of $L$ remains in the order of $\frac{1}{\varepsilon}$.

*   **Phase II: decay of the energy *E* from order 1 to order *ε*.** The function $\varphi$ approximating the concentrations $c_j$ has a nontrivial evolution (corollary 3.3 and proposition 4.1). The energy *E* decreases from a value close to 1 to $E_{\text{fin}} \approx \varepsilon$. The distribution of clusters evolves by means of a combination of an oscillatory displacement towards large cluster sizes combined with their spread in the space of cluster size as in Phase I (lemma. 3.4). The period *T* of the LV oscillations is of order $\frac{E}{\varepsilon^2}$ (lemmas 2.1 and 2.8). The number of LV oscillations taking place during this phase is of order $n \approx \frac{-\log(\varepsilon)}{\varepsilon}$ (see section 4.3, proposition 4.3).

*   **Phase III: Decay of the energy *E* from order *ε* to order *ε*³.** The concentrations of clusters $c_j$ remain nearly constant for large sizes $j$, but the concentrations $c_j$ with $j$ of order one oscillate during the LV oscillations of the concentrations $v, w$ (see section 5, proposition 5.2). The period of the oscillations $T$ remains in the order of $\frac{1}{\varepsilon}$. The number of LV oscillations $n$ taking place during this phase is of order $n \approx \frac{-\log(\varepsilon^2)}{\varepsilon}$ (see section 5.3). The oscillations for the cluster concentrations $j$ of order one are progressively damped and finally become negligible (see section 5.3, proposition 5.4).

*   **Phase IV: Final trend to the equilibrium.** The energy *E* remains of order *ε*³. The concentrations $c_j$ evolve towards their equilibrium value and their values can be approximated using a nonlinear second order parabolic equation (see section 6, proposition 6.2). During this phase, all the concentrations $c_j$ with $j$ of order one can be approximated by a single function $C(0^+, t)$ which changes in the same time scale in which the whole distribution of concentrations $\{c_j\}_{j \in \mathbb{N}}$ approach to equilibrium. The oscillations being negligible, the total duration for this last phase is in the order of $\frac{1}{\varepsilon^3}$.

An sketch overviewing these four phases is provided by figure 2.

**3. Phases I and II: from one LV cycle to the next**

As pointed out in section 2.3, the dynamics of the {$c_j$}$_{j \in \mathbb{N}}$ is a discrete-in-size convection-diffusion process. We now describe in detail the evolution of the cluster concentrations {$c_j$}$_{j \in \mathbb{N}}$ during one LV cycle of the monomer concentrations (v, w). We will illustrate how this mechanism leads to an approximation of the cluster distribution in terms of the ansatz (31) consisting of an aggregated mass for the smallest cluster and a smooth profile function for large cluster sizes during the Phases I and II. For this reason, we assume in all of this section that the energy satisfies $\varepsilon \ll E \lesssim 1$, hence all the lemmas of sections 2.2 and 2.3 are valid. We use them in order to find out what is the decay of energy along one LV cycle, and how to go from a cycle $n$ to a cycle $n+1$.
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Overview sketch of the four asymptotic phases. Phase I: nearly constant energy $E \sim 1$, nearly unperturbed LV-cycles with period $T \sim 1/\epsilon^2$ causing oscillatory advection of the cluster distribution, diffusive spreading of the cluster distribution by about $1/\sqrt{\epsilon}$ per period; Phase II: decay of $E \searrow \epsilon$, damped LV-cycles with $T \sim E/\epsilon^2$, oscillatory advection and diffusive spreading reaching a balance for large clusters; Phase III: different decay of $E \searrow \epsilon^3$, damped LV-cycles with $T \sim 1/\epsilon$, only small clusters are left to feature oscillatory behaviour; Phase IV: nearly constant energy $E \sim \epsilon^3$, no more oscillations, final exponential equilibration of monomer- and cluster concentrations.</figcaption></figure>

Figure 3 plots phase portraits of the monomeric concentrations ($v$, $w$) and snapshots of the approximate cluster size distribution $c(j, t)$ on time intervals $[t_1, t_2]$, $[t_2, t_3]$ and $[t_3, t_4]$ according to lemmas 2.3–2.5. It is for didactic reasons that figure 3 starts with lemma 2.3 describing stage two: The corresponding evolution of the size distribution in figure 3(b) shows very nicely the convective transport to larger cluster sizes ($\tilde{V} \simeq w$ in (27)–(28) due to $v$ being small on $[t_1, t_2]$) combined with the diffusive spreading occurring in this stage according to $d \simeq w$. Convection and diffusion slow down as $w$ decays in the phase portrait figure 3(a).

Next in figure 3(c), the monomer concentrations ($v$, $w$) are small compared to $\epsilon$ on the interval $[t_2, t_3]$ of lemma 2.4. Hence, the cluster size distribution undergoes only negligible convection and diffusion in figure 3(d), making $[t_2, t_3]$ the longest time interval within one LV cycle.

The dynamics picks up again in stage four described in lemma 2.5 and plotted in figure 3(e), which sees the growth of $v$ on the time interval $[t_3, t_4]$ while $w$ remains small. Accordingly, figure 3(f) shows convection of the size distribution towards smaller clusters combined with diffusive spreading. We highlight that the resizing of cluster towards smaller sizes is naturally limited by the smallest cluster size $c_1$ where we observe the formation of an aggregate.

Figure 4 shows in subfigure 4(b) that the aggregation at the smallest cluster is further enhanced on the time interval $[t_4, t_5]$ of stage five, lemma 2.6, during which the phase portrait of ($v$, $w$) reaches the end time $T$ of a LV cycle in figure 4(a). Finally, figures 4(d) and (d) show the first stage of a LV cycle as in lemma 2.2. We remark again that we decided to plots
<figure><img src="image_7.png" /><figcaption>(a) Phase portrait of $(v, w)$ for
$t \in [t_1, t_2] = [3, 191]$, cf. Lemma 2.3</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(b) Evolution of cluster distribution $c(j, t)$ for $t \in [t_1, t_2] = [3, 191]$</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>(c) Phase portrait of $(v, w)$ (log scale) for $t \in [t_2, t_3] = [191, 2541]$, cf. Lemma 2.4</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>(d) Evolution of cluster distribution $c(j, t)$ for $t \in [t_2, t_3] = [191, 2541]$</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>(e) Phase portrait of $(v, w)$ for $t \in [t_3, t_4] = [2541, 2729]$, cf. Lemma 2.5</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>(f) Evolution of cluster distribution $c(j, t)$ for $t \in [t_3, t_4] = [2541, 2729]$</figcaption></figure>

**Figure 3.** Phase portraits and snapshots of the cluster size distribution for one cycle of the LV system during Phase I in the case $\epsilon = 0.0212$. See appendix C for the details about the numerical simulations.
<figure><img src="image_13.png" /><figcaption>(a) Phase portrait of $(v, w)$ for $t \in [t_4, t_5 = [2729, 2732]$, cf. Lemma 2.6</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>(b) Evolution of cluster distribution $c(j, t)$ for $t \in [t_4, t_5 = T] = [2729, 2732]$</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(c) Phase portrait of $(v, w)$ for $t \in [T, T + t_1] = [2733, 2736]$, cf. Lemma 2.2</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>(d) Evolution of cluster distribution $c(j, t)$ for $t \in [T, T + t_1] = [2733, 2736]$</figcaption></figure>

<caption>**Figure 4.** Phase portraits and snapshots of the cluster size distribution during the fastest region of the LV phase portrait in Phase I. See appendix C for the details about the numerical simulations.</caption>

this first time interval last, since the evolution of the cluster distribution during the first stage when the aggregate of the smallest cluster begins to convect towards larger clusters is difficult to read without knowing the previous discussion.

Combining the time scales characterised in lemma 2.8 with the above discussion of figures 3 and 4 suggests to separate the dynamics of each LV cycle into two principal steps: 1) a 'pure transport' step, during which the diffusion is negligible—this is the interval $[0, t_1]$ and symmetrically $[t_4, t_5]$; we could even extend these intervals to the intermediate time points $[0, t_{1,2}]$ and $[t_{3,4}, t_5]$ (see the proof of lemma 2.3 in the appendix where these times are defined); 2) a 'transport and diffusion' step, where both transport and diffusion are important, which is given by the symmetric intervals $[t_{1,2}, t_2]$ and $[t_3, t_{3,4}]$. The remaining time interval $[t_2, t_3]$ is a kind of lag time where almost nothing happens, neither transport nor diffusion, due to the exponential smallness of both $v$ and $w$.
The following proposition characterises one iteration step of the size distribution of clusters as plotted in figures 4(a) and (b) (by taking the final depicted time points) from a first appearance (which we shall call $t=0$ rather than $t=T$ for the sake of simplicity) to its second appearance after one period $T$.

**Proposition 3.1 (Initial cycle and iterative formula).** Let $M=1$, $\varepsilon \ll E_0 \ll 1$, and $v(0) = w(0) > \varepsilon$ such that (11) is satisfied. Let $1 \ll L_0 \lesssim \frac{1}{\varepsilon}$, and $c_j(0)$ defined by (30) for some measure $\varphi_0$ satisfying (31)–(32). Let $T_0 := T(E_0, \varepsilon)$ be the time period of the unperturbed LV system (14), $D_0 := D(E_0, \varepsilon)$ the diffusion defined by (33), and $\sigma_0 := \frac{D_0}{(L_0)^2}$.

Defining $(\varphi_1, L_1)$ such that $c_j(T_0) = \frac{\varepsilon}{L_1}\varphi_1(\frac{j}{L_1})$ and $\varphi_1$ satisfies (31)–(32), we have

$$ \varphi_1(x) \approx \frac{L_1}{L_0} \chi_{(0,\infty)}(x) \int_0^\infty G\left(\frac{L_1}{L_0}x - \eta; \sigma_0^2\right) \varphi_0(\eta) \, d\eta \\ + \delta(x) \int_{-\infty}^0 d\zeta \int_0^\infty G\left(\zeta - \eta; \sigma_0^2\right) \varphi_0(\eta) \, d\eta, \quad (36) $$

where $G$ is the fundamental solution for the heat equation:

$$ G(\xi; s) = \frac{1}{\sqrt{4\pi s}} \exp\left(-\frac{\xi^2}{4s}\right) \quad (37) $$

and $L_1$ is given by

$$ \frac{L_1}{L_0} = \int_0^\infty \varphi_0(\eta) \, d\eta \int_0^\infty xG(x - \eta; \sigma_0^2) \, dx \\ = 1 + \int_0^\infty xdx \int_0^\infty G(x + \eta; \sigma_0^2) \varphi_0(\eta) \, d\eta. \quad (38) $$

**Remark 3.2.** If we had that $L_0$ is of order one (i.e. the initial size concentration is concentrated in cluster sizes $j$ of order one), then the solution during a LV cycle would spread the cluster concentrations to regions of size $\sqrt{\frac{E}{\varepsilon}} \gg 1$. Therefore, in a single cycle the width of the region in which the clusters are concentrated would become much larger than one; this explains why we have assumed in proposition 3.1 that $L_0 \gg 1$. Notice that with this assumption, if there is a fraction of clusters in the concentrations $c_j$ with $j$ of order one, we should include a Dirac mass at in $\varphi_0$ at $x=0$ containing this fraction of clusters.

**Heuristic proof of proposition 3.1.** We provide a heuristic proof of proposition 3.1 divided into the following steps: In section 3.1, we analyse the ODE system as being close to a continuous convection–diffusion equation: its dynamics along one LV cycle is thus described in lemma 3.4 below. This result is used in section 3.2 to obtain (36)–(38) by solving the diffusion equation during one cycle with the help of the fundamental solution to the heat equation. $\square$

It is clear that as long as the assumptions of the proposition 3.1 remain valid, we can apply it to successive time periods $T_n$ themselves defined by successive energies $E_n$. Hence, the following corollary holds as long as $E_n \gg \varepsilon$.

**Corollary 3.3.** Under the assumptions and notations of proposition 3.1, the same result applies between two successive cycles $n$ and $n+1$ as long as $E_n \gg \varepsilon$. By denoting accordingly $\varphi_n$, $L_n$, $m_n$, $\psi_n$ and $\sigma_n^2$, we have the following iterative formulae:
$$
\begin{align}
\psi_{n+1}(x) &= \frac{L_{n+1}}{L_n} \int_0^\infty G\left(\frac{L_{n+1}}{L_n}x - \eta, \sigma_n^2\right) \psi_n(\eta) \,d\eta \nonumber \\
&\quad + m_n \frac{L_{n+1}}{L_n} G\left(\frac{L_{n+1}}{L_n}x, \sigma_n^2\right), \\
m_{n+1} &= \int_{-\infty}^0 d\zeta \int_0^\infty G(\zeta - \eta, \sigma_n^2) \psi_n(\eta) \,d\eta + \frac{m_n}{2}, \\
\frac{L_{n+1}}{L_n} &= 1 + \int_0^\infty xdx \int_0^\infty G(x + \eta, \sigma_n^2) \psi_n(\eta) \,d\eta \nonumber \\
&\quad + m_n \int_0^\infty xG(x, \sigma_n^2) \,dx. \tag{39}
\end{align}
$$

3.1. An approximate transport and diffusion equation

As a first step to establish proposition 3.1, let us study the dynamics of the cluster size dis-
tribution during one cycle. We use equations (25)–(26) as written in section 2.3 and interpret
them as a discretised drift-diffusion equation (27). Moreover, it will be useful to rescale time
in (27) and consider the following convection–diffusion equation

$$
\frac{\partial c}{\partial \tau}(j, \tau) + \frac{2(w-v)}{(w+v)} \frac{\partial c}{\partial j}(j, \tau) = \frac{\partial^2 c}{\partial j^2}(j, \tau), \quad \tau(t) = \int_0^t \frac{d(s)}{2} ds, \quad (40)
$$

where we recall that $d = w + v$ and lemmas 2.7, 2.8, i.e.

$$
D(E, \varepsilon) := \frac{1}{2} \int_{0}^{T(E, \varepsilon)} d(s) ds \sim \frac{E}{\varepsilon},
$$

which is a consequence of $d = O(\varepsilon)$ for those part of the LV cycle which contribute to the
integral the most. Hence, one can image the time-rescaling as roughly $\tau = O(\varepsilon t)$. The time-
rescaled equation (40) is a key step in understanding these dynamics.

**Lemma 3.4.** Under the assumptions and notations of proposition 3.1, we define $\Delta t := \sqrt{\frac{1}{E_0\varepsilon}}$.
During the time interval $[0, T_0]$, the size distribution $c_j(t)$ approximately evolves according to
the following dynamics:

1. Over the time interval $[0, \Delta t]$, transport dominates diffusion: i.e. $c_j(t) \approx c_{j-Y(t)}(0)$ is transported along $Y(t) > 0$ defined by (33), and moves towards larger sizes.

2. Over the time interval $t \in [\Delta t, T_0 - \Delta t]$, $c_j(t)$ is close to the solution to (40), hence $c_j(t) \approx u(x, \tau(t))$ with $j = x + Y(t)$, $\tau$ and $d$ defined by (28), and $u$ solution to the pure diffusion equation

$$
\frac{\partial u}{\partial \tau}(x, \tau) = \frac{\partial^2 u}{\partial x^2}(x, \tau), \quad E_0 \Delta t \le \tau \le D_0 - E_0 \Delta t. \tag{41}
$$

During this time interval, the total diffusion is in the order of $D_0 \approx \frac{E_0}{\epsilon}$: the size distribution
diffuses in a range $\Delta x \approx \sqrt{\frac{E_0}{\epsilon}}$.

3. In the time interval $t \in [T_0 - \Delta t, T_0]$, transport dominates diffusion again. We have $c_j(t) \approx c_{j-Y(t)+Y(T_0-\Delta t)}(T_0 - \Delta t)$ and $Y(t) < 0$. The mass moving again towards smaller sizes, for $j < \Delta x$, the mass accumulates in $c_1$.
**Proof.** W.l.o.g. and for the sake of simplicity, we set $E_0 = 1$ as the general case follows by the change of variables detailed in lemmas 2.1 and 2.7. Let us begin by proving the central behaviour of Step 2 given by (41). To do so, we use (40) and identify—with some abuse of notations—its solution $c(j, \tau)$ with $c(x, \tau)$. We have $\tau(T(E_0, \varepsilon)) = D_0$ defined in (33), and it follows that (40) becomes (41).

The approximation (41) is valid on the whole size space $(0, +\infty)$ (and even $(-\infty, +\infty)$ in the sense that we do not need any boundary condition to solve the heat equation) if the mass is contained in large enough clusters, i.e. $Y(t)$ large enough, which is true, thanks to lemma 2.8, on an interval $[\Delta t, T_0 - \Delta t]$ for some $\Delta t \in (t_1, t_2)$ to be specified. In this case, we expect to have $c_1$ very small, so that the equations for $v$ and $w$ (7)–(8) may be approximated by the unperturbed LV system (14), and the approximation of $Y$ done in lemma 2.8 is a posteriori validated.

During this time interval, the range of values of cluster sizes in which the concentrations $c$ are relevant spreads according to the equation (41). Notice that this spreading process can be thought as diffusion in the space of cluster sizes, although this diffusion is not related to any physical diffusion in the physical space, but it is just a mathematical description that allows to compute the evolution of the concentrations in the space of cluster sizes. Let us denote the heat equation semigroup $\exp(\tau \partial_x^2)$ as $S(\tau)$. Assuming $\Delta t \ll D_0$ (this will be confirmed a posteriori), the total time for diffusion during $[\Delta t, T_0 - \Delta t]$ is thus in the order of $\tau(T_0) = D_0 \sim \frac{E_0}{\varepsilon}$, see lemmas 2.7 and 2.8. Due to the properties of the heat semigroup, the characteristic length for diffusion is thus of order $\sqrt{\frac{E_0}{\varepsilon}} \ll D_0$.

This provides us with the order of magnitude for $\Delta t$: at the final stage $[T_0 - \Delta t, T_0]$ of the LV cycle, when the concentrations front returns to values of $j$ of order one, the concentrations $c_j$ for $j = O(1)$ got diffused over a region of order $\sqrt{\frac{E_0}{\varepsilon}}$ cluster sizes (see figures 3(f) and 4(b)). The transport speed $\frac{dY}{dt} = w(t) - v(t)$ in this region is of order $-E_0$ by lemmas 2.1 and 2.5, thus the time interval during which the size region $[1, \sqrt{\frac{E_0}{\varepsilon}}]$ hits $c_1$ is of order

$$ \Delta t = \frac{1}{E_0} \sqrt{\frac{E_0}{\varepsilon}} = \sqrt{\frac{1}{E_0 \varepsilon}} $$

Concerning the dynamics of the clusters in the boundary layer during the arrival of the concentrations front, the mass diffused accumulates in $j = 1$ (see figures 3(f) and 5(a)). For clusters of order one, we can approximate the concentrations $c_j$ by means of the solutions of

$$ \frac{dc_j}{d\tau} = (c_{j+1} - c_{j-1}) + (c_{j-1} + c_{j+1} - 2c_j), \quad j \ge 2 \qquad (42) $$

$$ \frac{dc_1}{d\tau} = 2 \left( c_2 - \frac{\varepsilon}{E_0} c_1 \right). \qquad (43) $$

To conclude, it remains to prove that the diffusion may be neglected in the intervals $[0, \Delta t]$ and $[T_0 - \Delta t, T_0]$; this is due to the fact that the time interval is much smaller than the total diffusion $D_0$. The approximations in Steps 1 and 2 by a pure transport equation are thus valid. $\square$

### 3.2. From the initial cycle to the next

In this section, we (heuristically) obtain (36) and (38) of proposition 3.1. These two formulae may be viewed as linking a prototypical size distribution at a first time point $T_0$ (or in general
<figure><img src="image_17.png" /><figcaption>**Figure 5.** Scheme of the evolution of the cluster sizes of order one.</figcaption></figure>

$T_n$, for the $n$th cycle) to the size distribution after one more cycle at time $T_1$ (respectively $T_{n+1}$) provided that the initial energy $E_0$ (or $E_n$) of this period satisfies $\varepsilon \ll E_0, E_n \ll 1$.

We remark that as a first prototypical size distribution at $T_0$, we consider the size distribution obtained after finishing an initial LV cycle, at the time point when the values of $w$ and $v$ cross the half-line $\{w = v > \varepsilon\}$ and we already have aggregated a positive amount of mass in clusters with size $j$ of order one, and more precisely in the smallest cluster size $j=1$ as shown in figure 4(b) (see also figure 5(b)) below. This mass comes from two sources, first, a possible contribution of the original size-distribution $\varphi_0(1/L_0)$ which has been transported back and forth—and diffused, hence lowered—and secondly, to the diffusion of small clusters $\varphi_0(j/L_0)$
with $j$ within a width $\sqrt{\frac{E}{\epsilon}} \gg 1$, which will accumulate to some fraction on $j=1$ during the transport back to the origin.

We remark that the conditions in (6) imply that

$$1 = v + w + \varepsilon L_0 \int_0^\infty x \varphi_0(x) \, dx, \quad \int_0^\infty \varphi_0(x) \, dx = 1$$

where we approximate the sums in (6) by means of integrals. In order to define $L_0$ in a unique manner, we assume (32) and obtain

$$1 = v + w + \varepsilon L_0, \quad \int_0^\infty \varphi_0(x) \, dx = 1 \qquad (44)$$

First, we assume that we can neglect the effect of the boundary conditions arising for clusters with size $j=1$ - we verify that this assumption is valid a posteriori. Then, lemma 3.4 allows us to compute the concentrations after one LV cycle from the concentrations at the beginning of the LV cycle (with energy $E$) by evaluating the heat semigroup $S(\tau)$ with time step $\tau = D(E, \varepsilon)$. This would give the approximation

$$\frac{\varepsilon}{L_1} \varphi_1 \left( \frac{j}{L_1} \right) = \frac{\varepsilon}{L_0} S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)} (\cdot) \right] (j) \quad \text{for } j \gg \sqrt{\frac{E_0}{\varepsilon}}, \qquad (45)$$

where $\chi_{(0,\infty)}(\cdot)$ is the characteristic function in the half-line $\{j > 0\}$. This approximation can be expected to be valid for cluster sizes away of the boundary layer, i.e. for $j \gg \sqrt{\frac{E_0}{\varepsilon}}$. We notice that if we were to take negative values of $j$, the right-hand side would still be positive, which is impossible with the discrete system. As already stated in lemma 3.4, this mass, instead of being diffused in the region of negative clusters, accumulates in clusters with $j=1$ (see the last frames in figures 3(f), 4(b) and 5(b)).

Then, to define $\varphi_1$ in terms of $\varphi_0$ we must include all this mass in a Dirac mass at the origin. This gives

$$\frac{\varepsilon}{L_1} \varphi_1 \left( \frac{j}{L_1} \right) = \frac{\varepsilon}{L_0} \chi_{(0,\infty)}(j) S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)}(\cdot) \right] (j) \\ + \frac{\varepsilon}{L_0} \delta(j) \int_{-\infty}^0 \left( S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)}(\cdot) \right] \right) (\xi) \, d\xi$$

or after some simplifications

$$\varphi_1 \left( \frac{j}{L_1} \right) = \frac{L_1}{L_0} \chi_{(0,\infty)}(j) S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)}(\cdot) \right] (j) \\ + \frac{L_1}{L_0} \delta(j) \int_{-\infty}^0 \left( S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)}(\cdot) \right] \right) (\xi) \, d\xi.$$
We rewrite this equation using the fundamental solution for the heat Equation $G(\xi; s)$ defined by (37), and we obtain

$$ \varphi_1 \left( \frac{j}{L_1} \right) = \frac{L_1}{L_0} \chi_{(0,\infty)}(j) \int_0^\infty G(j-\xi; D_0) \varphi_0 \left( \frac{\xi}{L_0} \right) d\xi \\ + \frac{L_1}{L_0} \delta(j) \int_{-\infty}^0 dk \int_0^\infty G(k-\xi; D_0) \varphi_0 \left( \frac{\xi}{L_0} \right) d\xi. $$

We write $x = \frac{j}{L_1}$ and use the change of variables $\xi = L_0\eta$, and $k = L_0\zeta$:

$$ \varphi_1(x) = L_1 \chi_{(0,\infty)}(x) \int_0^\infty G(L_1 x - L_0 \eta; D_0) \varphi_0(\eta) d\eta \\ + \delta(x) L_0 \int_{-\infty}^0 d\zeta \int_0^\infty G(L_0 \zeta - L_0 \eta; D_0) \varphi_0(\eta) d\eta, $$

where we use that $L_1\delta(L_1x) = \delta(x)$.

We now use that $G(\xi, s) = \frac{1}{\lambda}G(\frac{\xi}{\lambda}, \frac{s}{\lambda^2})$ for any $\lambda > 0$, and taking $\lambda = L_0$ we obtain

$$ \varphi_1(x) = \frac{L_1}{L_0} \chi_{(0,\infty)}(x) \int_0^\infty G\left(\frac{L_1}{L_0}x - \eta; \frac{D_0}{(L_0)^2}\right) \varphi_0(\eta) d\eta \\ + \delta(x) \int_{-\infty}^0 d\zeta \int_0^\infty G\left(\zeta - \eta; \frac{D_0}{(L_0)^2}\right) \varphi_0(\eta) d\eta, $$

which is (36).

In order to conclude the computation of $\varphi_1$ we need to determine $L_1$. To this end we impose the condition (32), namely $\int_0^\infty x\varphi_1(x)dx = 1$. Notice that since we assume that $\int_0^\infty \varphi_0(x)dx = 1$ we automatically have $\int_0^\infty \varphi_1(x)dx = 1$. We thus have

$$ 1 = \int_0^\infty x \varphi_1(x) dx = \frac{L_1}{L_0} \int_0^\infty x dx \int_0^\infty G\left(\frac{L_1}{L_0}x - \eta; \sigma_0^2\right) \varphi_0(\eta) d\eta. $$

Then

$$ \frac{L_1}{L_0} = \int_0^\infty \int_0^\infty G(x - \eta; \sigma_0^2) x dx \varphi_0(\eta) d\eta, $$

which is the first equality of (38). The second one has the advantage of highlighting that $L_1 > L_0$. To obtain it, we use the symmetry of $G$ to write

$$ \frac{L_1}{L_0} = \int_0^\infty \left( \int_{-\infty}^\infty G(x-\eta; \sigma_0^2) x dx - \int_{-\infty}^0 G(x-\eta; \sigma_0^2) x dx \right) \varphi_0(\eta) d\eta \\ = \int_0^\infty \eta \varphi_0(\eta) d\eta - \int_0^\infty \int_{-\infty}^0 G(x-\eta; \sigma_0^2) x dx \varphi_0(\eta) d\eta \\ = 1 + \int_0^\infty \varphi_0(\eta) \int_0^\infty x G(x+\eta; \sigma_0^2) dx d\eta $$
and this is the second equality of (38). It ends the proofs of proposition 3.1. This also provides us with the the formulae in (39). Taking as the initial step $\psi_n$ and $m_n$ with $\text{Supp}(\psi_n) \subset (0, \infty)$ as in (31) and denoting accordingly $D_n$ and $\sigma_n^2 = \frac{D_n}{(L_n)^2}$ using $E = E_n$ in (33), we obtain the formulae in (39) and finish the proof of corollary 3.3 by using the formulae in proposition 3.1, replacing the corresponding terms and iterating over the LV cycles.

### 3.3. Energy decay and iterative mapping

As a first result in proposition 3.5, we compute the energy decay over a first LV cycle and we obtain a semi-explicit representation formula and provide an estimation in terms of $\sqrt{\varepsilon}$. As a second step, proposition 3.6 quantifies the change of energy $E$ and the change of typical cluster length scale $L$ over following LV cycle iterations as long as the energy remains relatively big in comparison with $\varepsilon$. Please note that some cases in the proof of proposition 3.6 require a more detailed analysis compared to the present section and are therefore postponed to section 4.2. While this might be sub-optimal in terms of presentation, we felt that as a result, proposition 3.6 belongs to this section and should not be stated any later.

**Proposition 3.5.** Under the assumptions and notations of proposition 3.1, the main change of the energy with a LV cycle occurs in the interval $[T_0 - \Delta t, t_4]$ and may be approximated as follows:

$$E_1 - E_0 \sim - \int_{-\infty}^{0} |s| W_0(s) ds, \qquad (46)$$

with $W_0$ defined by

$$W_0(\cdot) = \frac{\varepsilon}{L_0} S(D_0) \left[ \varphi_0 \left( \frac{\cdot}{L_0} \right) \chi_{(0,\infty)}(\cdot) \right]. \qquad (47)$$

Moreover, since $\sigma_0 := \frac{D_0}{(L_0)^2} \approx 1$, the following estimation holds

$$|\Delta_{E_0}| := \int_{-\infty}^{0} |s| W_0(s) ds \approx \varepsilon L_0 \lesssim \sqrt{\varepsilon}. \qquad (48)$$

With propositions 3.1 and 3.5, we can now obtain formulae quantifying the change of energy and the change of the cluster length scale over instances of following LV cycles.

**Proposition 3.6.** Under the assumptions and notations of propositions 3.1 and 3.5 and corollary 3.3, as long as $E_n \gg \varepsilon$, the following estimation holds

$$|\Delta_{E_n}| := \int_{-\infty}^{0} |s| W_n(s) ds \approx m_n \sqrt{E_n \varepsilon} + \frac{E_n}{L_n} \psi_n(0) \lesssim \sqrt{\varepsilon} \qquad (49)$$

if $\sigma_n^2 \ll 1$. Moreover, iteration formulae for the energy $E_{n+1}$ and the cluster length scale $L_{n+1}$ are given by

$$|E_{n+1} - E_n| \approx \frac{E_n}{L_n}, \qquad \frac{L_{n+1}}{L_n} = 1 + C_n, \qquad C_n = \frac{|\Delta_{E_n}|}{\varepsilon L_n} \approx \frac{E_n}{\varepsilon L_n^2}. \qquad (50)$$
**Proof of proposition 3.5.** As a first step, we prove the expression (46). We recall that the energy dynamics is given by (12): $\frac{dE}{dt} = (\varepsilon - v)c_1$, so that we need to compute where the contribution to $\int (\varepsilon - v)c_1 dt$ is dominant. We also recall the partition of the LV phase space and the corresponding partition of the time period introduced in section 2.2 and in lemma 3.4.

* During $[0, t_1]$, by lemmas 2.1 and 2.2 we have $c_1 \le \varepsilon$ and $v \in [\varepsilon, \frac{E_0}{2}]$ hence $|\int_0^{t_1} (\varepsilon - v)c_1 dt| \lesssim \varepsilon \frac{E_0}{2} t_1 \sim -\frac{\varepsilon}{2} \log(\varepsilon)$.

* During $[t_1, \Delta t]$, we have $v \le \varepsilon$ hence $|\int_{t_1}^{\Delta t} (\varepsilon - v)c_1 dt| \lesssim \varepsilon^2 \Delta t \approx \frac{\varepsilon^{\frac{3}{2}}}{\sqrt{E_0}}$.

* During $[\Delta t, T_0 - \Delta t]$, $c_1 \lesssim e^{-\sqrt{\frac{E_0}{\varepsilon}}} \ll \varepsilon^4$ hence $\int_{\Delta t}^{T_0 - \Delta t} (\varepsilon - v)c_1 dt \ll \varepsilon^2$ is negligible.

* During $[T_0 - \Delta t, t_4]$, we compute below the contribution $\int_{T_0 - \Delta t}^{t_4} (\varepsilon - v)c_1 dt$.

* During $[t_4, T_0]$, analogously to the time $[0, t_1]$, we have $|\int_{t_4}^{T_0} (\varepsilon - v)c_1 dt| \lesssim -\frac{\varepsilon}{2} \log(\varepsilon)$.

It remains to compute the contribution to the change of energy during the time $[T_0 - \Delta t, t_4]$. We have seen in lemma 3.4 that we can approximate the size distribution dynamics by the pure transport equation along $Y(t)$, and moreover in $[T_0 - \Delta t, t_4]$ we have $v \sim E_0$ and $w \lesssim \varepsilon$, hence in the variable $\tau$ we have $\frac{dY}{d\tau} \sim -2$. The discontinuity of the characteristic function $\chi_{(0,\infty)}(\cdot)$ results in the formation of a concentration front, that can be described when it approaches to clusters of order one by means of the function

$$c(\cdot, \tau) = \frac{\varepsilon}{L_0} S(D_0) \left[ \varphi_0 \left( \frac{-2(\tau_* - \tau)}{L_0} \right) \chi_{[0, \infty)}(\cdot) \right] \quad (51)$$

where $\tau_* = D_0$ is the period of the LV cycle. It is convenient to rewrite this formula defining the function $W_0$ by (47), then

$$c(j, \tau) = W_0(j - 2(\tau_* - \tau)).$$

Here, we use $j$ as continuum variable. The flux of monomers towards clusters of order one can thus be approximated as $c(0^+, t) = W_0$ leading (43) to be written

$$\frac{dc_1}{d\tau} = 2W_0(-2(\tau_* - \tau)).$$

Then, using that the function describing the concentrations front $W_0$ decays fast for large negative values of its argument so that integrating between $\tau(T_0 - \Delta t)$ and $\tau$ as the same order of magnitude as integrating between 0 or even $-\infty$ and $\tau$, we obtain

$$c_1(\tau) \sim 2 \int_0^\tau W_0(-2(\tau_* - s)) ds \sim \int_{-\infty}^{2(\tau - \tau_*)} W_0(s) ds \quad \text{for } \tau < \tau_* \quad (52)$$

In order to estimate the order of magnitude of the change of the energy due to the interaction of the concentrations waves with the regions with cluster sizes $j$ of order one, we use the fact that $W_0$ is of order $\frac{\varepsilon}{L_0}$ and it has a width of order $\sqrt{\frac{E_0}{\varepsilon}}$ (cf (47)). We remark that this scaling properties for $W_0$ are based in the assumption that the contribution of the Dirac in $\varphi_0$ gives a contribution to $W_0$ of the same order of magnitude that the part of $\varphi_0$ in $\{x > 0\}$.

We can then compute the change of the energy using (12). We recall that $\frac{dE}{dt} = (\varepsilon - v)c_1$, hence $\frac{dE}{d\tau} = \frac{2(\varepsilon - v)}{(w+v)}c_1 \approx -2c_1$ in $[T_0 - \Delta t, t_4]$ since $v \simeq E_0$ and $w \lesssim \varepsilon$. Then, as long as $E \gg$
$\epsilon$, we have the following approximation during the range of times in which $c_1$ contributes significantly to the change of $E$

$$ \frac{dE}{d\tau} \approx -\frac{2E}{E}c_1 = -2c_1 \qquad (53) $$

Combining (52) and (53) and claiming again that the integrals are negligible away from $[T_0 - \Delta t, t_4]$, we obtain

$$ \begin{aligned} E_1 - E_0 &\approx -4 \int_0^{\tau_*} \int_0^{\tau} W_0(-2(\tau_* - s)) \, ds \, d\tau = \int_{-\tau_*}^0 W_0(s) \, ds \\ &\approx - \int_{-\infty}^0 W_0(s) |s| \, ds. \end{aligned} $$

In order to prove the estimation (48), we compute

$$ \begin{aligned} \Delta_{E_0} &= \frac{\epsilon}{L_0} \int_{-\infty}^{0} |s| \int_{0}^{\infty} \varphi_0 \left( \frac{y}{L_0} \right) G(s-y; D_0) \, dy \, ds \\ &= \frac{\epsilon}{L_0^2} \int_{0}^{\infty} s \int_{0}^{\infty} \varphi_0 \left( \frac{y}{L_0} \right) G \left( \frac{s+y}{L_0}; \sigma_0^2 \right) \, dy \, ds \\ &= \epsilon L_0 \int_{0}^{\infty} s \int_{0}^{\infty} \varphi_0(y) G(s+y; \sigma_0^2) \, dy \, ds \\ &= \epsilon L_0 \int_{0}^{\infty} s \left( m_0 G(s; \sigma_0^2) + \int_{0}^{\infty} \psi_0(y) G(s+y; \sigma_0^2) \, dy \right) \, ds. \end{aligned} $$

Using $\int_0^\infty sG(s; \sigma^2)ds = \frac{\sigma}{\sqrt{\pi}}$, we obtain

$$ \Delta_{E_0} = \epsilon L_0 \left( \frac{m_0 \sigma_0}{\sqrt{\pi}} + \int_0^\infty s \int_0^\infty \psi_0(y) G(s+y; \sigma_0^2) \, dy \, ds \right). \qquad (54) $$

Since $\sigma_0 \approx 1$, the sum of the first and second terms is in the order of $\epsilon L_0$ since $m_0 + \int_0^\infty \psi_0(x)dx = 1$, which shows (48). $\square$

**Proof of proposition 3.6.** The formula of the energy decay (46) as well as its estimate (54) holds over the iterations of the LV cycles as long as $E_n \gg \epsilon$. The only change in comparison to the first cycle is that if $\sigma_n \ll 1$, we have that

$$ \Delta_{E_n} = \epsilon L_n \left( \frac{m_n \sigma_n}{\sqrt{\pi}} + \int_0^\infty s \int_0^\infty \psi_n(y) G(s+y; \sigma_n^2) \, dy \, ds \right), $$

where the first term is in the order of $\epsilon L_n m_n \sigma = \epsilon m_n \sqrt{D_n} = m_n \sqrt{E_n \epsilon}$ (cf lemma 2.8), whereas the second term may be approximated by

$$ \epsilon L_n \psi_n(0) \sigma_n^2 \int_0^\infty \int_y^\infty (z-y) G(z; 1) \, dz \, dy \approx \epsilon \psi_n(0) \frac{D_n}{L_n} = \frac{E_n}{L_n} \psi_n(0). $$
Concerning the changes in the cluster distribution over the LV cycles, for $\frac{L_{n+1}}{L_n}$, we have

$$ \frac{L_{n+1}}{L_n} = 1 + C_n, $$

with

$$ C_n = \int_0^\infty \int_0^\infty xG(x+\eta, \sigma_n^2) \varphi_n(\eta) d\eta dx = \frac{1}{\varepsilon L_n} \int_{-\infty}^0 |s| W_n(s) ds = \frac{|E_{n+1} - E_n|}{\varepsilon L_n}, $$

so that it only remains to prove either the estimation of $E_{n+1} - E_n$ or the one for $C_n$ to get the other.

In the cases that $L_n \approx \frac{1}{\sqrt{\varepsilon}}$, we easily obtain (50) from (48) since then $E_n \approx 1$, so $\varepsilon L_n \approx \sqrt{\varepsilon} \approx \frac{E_n}{L_n}$, which finishes the proof of proposition 3.6 in this case. However, the proof is much more involved in case $L_n \gg \frac{1}{\sqrt{\varepsilon}}$, or equivalently, for $\sigma_n \ll 1$. It then requires the analysis carried out in section 4.2 below. □

# 4. Phases I and II: overall dynamics

## 4.1. Early Phase I: initial increase of the characteristic length

With an initial energy $E_0 \approx 1$, we have noticed that even if $L_0 \ll \frac{1}{\sqrt{\varepsilon}}$, the diffusive effects within one LV cycle yield $L_1 \approx \frac{1}{\sqrt{\varepsilon}}$, whereas the energy remains $E_1 \approx 1$ since $\Delta E_0 \approx \varepsilon L_0 \ll 1$ by proposition 3.5. Over following LV cycles, by (50) we have that $L_n$ increases much faster than $E_n$ decreases. Departing from $L_1 \approx \frac{1}{\sqrt{\varepsilon}}$, we have $C_1 \approx 1$. Hence, after a few cycles (in the order of $-\log(\varepsilon)$) we have $L_n \gg \frac{1}{\sqrt{\varepsilon}}$, so that $\sigma_n^2 = \frac{E_n}{\varepsilon L_n^2} \ll 1$ (see proposition 3.6). The change in the size distribution for the early Phase I is illustrated in the figure 6. Figure 7 plots the change in the energy E for one iteration of the LV cycle due to the evolution of $c_1$.

## 4.2. Continuation of Phase I and Phase II: size distribution evolution

The change of energy in one LV cycle being of order at most $\sqrt{\varepsilon}$ (see proposition 3.5 and 3.6), the energy $E_n$ first remains of order one while $\sigma_n \ll 1$ (Phase I) then decreases to order $\varepsilon$ (Phase II).

To describe the dynamics after the early Phase I, we are in the case of $\sigma_n \ll 1$ (see proposition 3.6) so that we need to evaluate both $m_n$ and $\psi_n(0)$ to prove (50). It is thus convenient to separate the size domain in three parts: the Dirac delta for $c_1$ (providing $m_n$), a boundary layer for sizes $1 \ll i \lesssim \frac{1}{\sqrt{\varepsilon}}$, and large sizes $i \gg \frac{1}{\sqrt{\varepsilon}}$, so that $\psi_n(0)$ may be viewed as the limit for $x \to \infty$ of the boundary layer as well as the limit $x \to 0$ of the large-sizes part of the domain. This is expressed in the following proposition.

**Proposition 4.1.** Let $\varepsilon \ll E_n$ and $\frac{1}{\sqrt{\varepsilon}} \ll L_n$, so that $\sigma_n^2 = \frac{E_n}{\varepsilon L_n^2} \ll 1$. Set $x = \frac{j}{L_n}$. Following the notations of corollary 3.3, we have the following approximations for the size distribution $\psi_n$.

* For $x > 0$ not small, more precisely $x = \frac{j}{L_n}$ with $j \gg \sqrt{\frac{E_n}{\varepsilon}}$, we have

$$ \psi_n(x) \approx \psi(x) := \frac{2}{\pi} e^{-\frac{x^2}{\pi}}, \qquad (55) $$

from which we also achieve the proof of (50).
<figure><img src="image_18.png" /><figcaption>(a) Size distribution $c(j, t)$ when $v = w > \epsilon$, i.e. for $t = 0, T_0, T_1, T_2$.</figcaption></figure>

<caption>(b) Size distribution $c(j, t)$ when $v = w < \epsilon$, i.e. for $t \approx \Delta t, T_0+\Delta t, T_1+\Delta t, T_2+\Delta t$.</caption>

<figure><img src="image_19.png" /><figcaption>(a) Phase portrait of the monomers</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>(c) Evolution of the energy $E$</figcaption></figure>

<caption>**Figure 7.** Change of the energy $E$ (11) for one iteration of the LV cycle when $c_1 \neq 0$ in equations (7) and (8).</caption>
* For $x = \frac{j}{L_n}$ with $j \lesssim \sqrt{\frac{E_n}{\epsilon}}$, defining $x = \sigma_{n+1}\xi$, $m_n = \sigma_n M_n$ and $\psi_n(\sigma_n\xi) = U_n(\xi)$, we have that $(U_n, M_n) \approx (U, M)$ solutions of the following system:

$$ U(\xi) = \int_{0}^{\infty} G(\xi - \zeta; 1) U(\zeta) d\zeta + MG(\xi; 1) , \quad \xi > 0 \quad (56) $$

$$ \frac{M}{2} = \int_{-\infty}^{0} d\zeta \int_{0}^{\infty} G(\zeta - \xi; 1) U(\xi) d\xi \quad (57) $$

$$ U(\infty) = \frac{2}{\pi} \quad (58) $$

**Proof.** In all the following dynamics, we have that $|\frac{\Delta E_n}{\epsilon L_n}| \ll 1$ so that (50) implies $\frac{L_{n+1}}{L_n} \sim 1$ in the limit $\sigma_n \to 0^+$. For convenience, we recall the first two iteration formulas (39) previously stated in corollary 3.3:

$$ \begin{aligned} \psi_{n+1}(x) &= \frac{L_{n+1}}{L_n} \int_0^\infty G\left(\frac{L_{n+1}}{L_n}x - \eta, \sigma_n^2\right) \psi_n(\eta) d\eta + m_n \frac{L_{n+1}}{L_n} G\left(\frac{L_{n+1}}{L_n}x, \sigma_n^2\right), \\ m_{n+1} &= \int_{-\infty}^0 d\zeta \int_0^\infty G(\zeta - \eta, \sigma_n^2) \psi_n(\eta) d\eta + \frac{m_n}{2}. \end{aligned} \quad (59) $$

In the following we approximate (39) resp. (59) using the expansion formula

$$ \int_{0}^{\infty} G(x-y; \sigma^2) \psi(y) dy = \psi(x) + \sigma^2 \partial_{xx} \psi(x) + o(\sigma^2) \quad \text{as} \quad \sigma \to 0^+ \quad (60) $$

which holds not uniformly but provided that $x > 0$ is fixed while and $\sigma \to 0^+$. On the other hand, the term $m_n \frac{L_{n+1}}{L_n} G\left(\frac{L_{n+1}}{L_n}x; \sigma_n^2\right)$ in the first equation of (39) resp. (59) is exponentially small in $\sigma_n^2$ if $x$ is of order one. We then obtain using (60) the following approximation

$$ \psi_{n+1}(x) = \frac{L_{n+1}}{L_n} \left[ \psi_n + \sigma_n^2 \partial_{xx} \psi_n \right] \left( \frac{L_{n+1}}{L_n} x \right), \quad x > 0, \quad (61) $$

if $\sigma_n^2$ tends to zero. Notice that the iterative formula (61) does not depend on $m_n$. This suggests in particular that the ratio $\frac{L_{n+1}}{L_n}$ can be determined independently on $m_n$ if $\sigma_n^2 \ll 1$. We will assume that $\psi_n$ stabilizes to a steady state. In addition, we use the fact that $\sigma_n^2$ tends to zero. We will assume that $\frac{L_{n+1}}{L_n} \sim 1 + a\sigma_n^2$ as $\sigma_n^2 \to 0$, where $a > 0$ must be determined. Then, looking for steady states $\psi$ of (61) and neglecting lower order terms, as well as using the approximation $\psi((1+a\sigma_n^2)x) \approx \psi(x) + a\sigma_n^2 x \partial_x \psi(x)$, we obtain that $\psi$ solves

$$ \sigma^2 (a\psi(x) + ax\partial_x\psi(x) + \partial_{xx}\psi(x)) = 0 $$

where we keep only the lower order terms in $\sigma^2$. Then

$$ a\psi(x) + ax\partial_x\psi(x) + \partial_{xx}\psi(x) = 0, \quad x > 0. \quad (62) $$

We can expect to have the approximation $\psi_n \approx \psi$ as $n \to \infty$. Moreover, we remark that we can expect to have $m_n \to 0$ as $n$ becomes larger enough. Indeed, this is a consequence of the fact
that as $\sigma_n^2 \to 0$, also the amount of mass transferred to the region $x < 0$ by the heat semigroup tends to zero. Therefore, the normalisation condition (32) implies that approximately

$$ \int_0^\infty \psi(x) dx = \int_0^\infty x\psi(x) dx = 1. \qquad (63) $$

The integrable solutions of (62) have the form $\psi(x) = Ce^{-\frac{ax^2}{2}}$ where $C$ is an arbitrary real constant. The normalisation conditions (63) imply $C = a = \frac{2}{\pi}$. Then we obtain (55), and also

$$ \frac{L_{n+1}}{L_n} = 1 + \frac{2}{\pi}\sigma_n^2, \quad \sigma_n^2 = \frac{D(E_n, \varepsilon)}{(L_n)^2}. $$

as $n$ becomes large enough. This also achieves the heuristic proof of (50).

We now consider the description of the functions $\psi_n(x)$ for small $x$, i.e. for $x$ in the order of $\sigma_n$. In the first equation of (39) resp. (59), we use $\frac{L_{n+1}}{L_n} \sim 1$ and introduce the new variables $x = \sigma_{n+1}\xi$, $\psi_n(\sigma_n\xi) = U_n(\xi)$. Then

$$ U_{n+1}(\xi) \sim \int_0^\infty G(\sigma_{n+1}\xi - \eta; \sigma_n^2) \psi_n(\eta) d\eta + m_n G(\sigma_{n+1}\xi; \sigma_n^2). $$

Defining $\eta = \sigma_n\zeta$ and using $G(y; \sigma^2) = \frac{1}{\sigma}G(\frac{y}{\sigma}; 1)$, we obtain

$$ U_{n+1}(\xi) \sim \int_0^\infty G\left(\frac{\sigma_{n+1}}{\sigma_n}\xi - \zeta; 1\right) U_n(\zeta) d\zeta + \frac{m_n}{\sigma_n} G\left(\frac{\sigma_{n+1}}{\sigma_n}\xi; 1\right). $$

We now remark that $\frac{\sigma_{n+1}}{\sigma_n} \sim 1$ as $\sigma_n \to 0$. Moreover, using the rescaling $m_n = \sigma_n M_n$, we introduce the following equation as the definition for the rescaled profile $U_n$ of boundary layer describing the concentrations near $x = 0$:

$$ U_{n+1}(\xi) = \int_{0}^{\infty} G(\xi - \zeta; 1) U_{n}(\zeta) d\zeta + M_{n}G(\xi; 1), \quad \xi > 0. \qquad (64) $$

In order to obtain a closed system for both $U_n$ and $M_n$, we use the second equation in (39) resp. (59). With the definitions of $U_n$ and $M_n$, we obtain the equation

$$ M_{n+1} = \frac{L_{n+1}}{L_n} \int_{-\infty}^{0} d\zeta \int_{0}^{\infty} G(\zeta - \xi; 1) U_n(\xi) d\xi + \frac{L_{n+1}}{L_n} \frac{M_n}{2}. $$

Approximating then $\frac{L_{n+1}}{L_n}$ by 1, we obtain

$$ M_{n+1} = \int_{-\infty}^{0} d\zeta \int_{0}^{\infty} G(\zeta - \xi; 1) U_{n}(\xi) d\xi + \frac{M_{n}}{2}. \qquad (65) $$

The equations (64)–(65) describe the iterative dynamics of the boundary layer near $x = 0$ yielding the cluster concentrations as well as the mass of the peaks appearing there. These equations must be solved combined with the matching condition that results from the fact that $\psi_n(\sigma_n\xi) = U_n(\xi)$. Then, since $\sigma_n \ll 1$, it is natural to impose the matching condition
$U_n(\infty) = \psi_n(0^+)$. Taking into account (55), we must then impose the following matching condition for the system (64)–(65)

$$U_n(\infty) = \frac{2}{\pi}. \qquad (66)$$

It is natural to assume that the solutions of (64)–(66) approach a stationary solution for large $n$. Therefore, they become close to the solutions of the problem (56)–(58), and it ends the (heuristic) proof. □

In the companion paper [11], it will be proved that there exists a unique solution of (56)–(57) satisfying the boundary condition (58). The resulting function $U(\xi)$ and the rescaled mass $M$ describe the concentration of clusters $c_j$ for $j$ in the order of $\frac{1}{\sigma_n}$, with $1 \ll \frac{1}{\sigma_n} \ll L_n$.

**Remark 4.2.** We can use the third equation in (39) to estimate $\frac{L_{n+1}}{L_n}$ to check the correctness of the approximation (50) that has been obtained with a different approach, approximating the evolution of the cluster sizes with $x$ of order one by means of a differential operator. We have

$$\frac{L_{n+1}}{L_n} = 1 + \sigma_n^2 \int_0^\infty ydy \int_0^\infty G(y+\xi; 1) U_n(\xi) d\xi + M_n \sigma_n^2 \int_0^\infty yG(y; 1) dy.$$

Using then the approximations $U_n \approx U$ and $M_n \approx M$, we obtain the following approximation

$$\frac{L_{n+1}}{L_n} = 1 + \sigma_n^2 \left( \int_0^\infty ydy \int_0^\infty G(y+\xi; 1) U(\xi) d\xi + M \int_0^\infty yG(y; 1) dy \right).$$

Therefore, both terms in the formula of $\frac{L_{n+1}}{L_n}$ give a comparable contribution. Another way to see this result is to go back to (49) in proposition 3.6: replacing $\psi(0)$ by $\frac{2}{\pi}$ and $m$ by $\sigma_n M$ we have

$$|\Delta_{E_n}| \approx M\sigma_n\sqrt{E_n\varepsilon} + \frac{E_n}{L_n}\frac{2}{\pi} \approx \frac{E_n}{L_n}\left(M + \frac{2}{\pi}\right),$$

and here again both terms are of similar order.

### 4.3. Overall dynamics during phases I and II

From (50) we can now detail the overall dynamics during Phases I and II, that we summarize in proposition 4.3. Let us first define the constant

$$A := \int_0^\infty ydy \int_0^\infty G(y+\xi; 1) U(\xi) d\xi + M \int_0^\infty yG(y; 1) dy. \qquad (67)$$

**Proposition 4.3.** Departing from $E_0 \approx 1$ and $L_0 \approx \frac{1}{\sqrt{\varepsilon}}$, with $A$ defined by (67), we have the following dynamics.

* Phase I: as long as $n \ll \frac{1}{\varepsilon}$, we have

$$L_n \sim \sqrt{\frac{2An}{\varepsilon}}, \quad E_n \sim 1. \qquad (68)$$
* Phase II: $L_n$ increases to $\frac{1}{\epsilon}$ and $E_n$ decreases, and at the end of Phase II, i.e. for $\frac{1}{\epsilon} \ll n \lesssim -\frac{\log(\epsilon)}{\epsilon}$, we have

$$L_n \sim \frac{1}{\epsilon} (1 - e^{-1}e^{-A\epsilon n}), \quad E_n \sim e^{-1}e^{-A\epsilon n}. \qquad (69)$$

**Proof.** The whole proof is based on studying the sequence $(L_n, E_n)$ defined by (50). We have already described the early Phase I in section 4.1: after a few cycles we have $L_n \gg \frac{1}{\sqrt{\epsilon}}$, and then by proposition 4.1 $\psi_n$ approaches $\frac{2}{\pi}e^{-\frac{x^2}{\pi}}$ and we can approximate $M_n$ and $U_n$ by the solutions of (56)–(58). We can then approximate the evolution of $L_n$ using that $C_n \sim A\sigma_n^2$ so that (50) may be written in the following more precise way:

$$\frac{L_{n+1}}{L_n} - 1 \sim \frac{AE_n}{\epsilon(L_n)^2}, \quad E_{n+1} - E_n \sim -A\frac{E_n}{L_n}. \qquad (70)$$

### Phase I asymptotics

As long as $E_n \sim 1$, we can simplify (70) by writing

$$L_{n+1} \sim L_n + \frac{A}{\epsilon L_n}.$$

Hence

$$\frac{(L_n)^2}{2} \sim \frac{(L_0)^2}{2} + \frac{An}{\epsilon}$$

for $1 \ll n \ll \frac{1}{\epsilon}$, and we obtain the asymptotics (68).

### Phase II asymptotics

We define a new set of variables, namely

$$L_n = \frac{1}{\epsilon} \ell(s), \quad E_n = e(s), \quad s = \epsilon n \qquad (71)$$

Then (70) becomes

$$\frac{\ell(s + \epsilon) - \ell(s)}{\epsilon} \sim \frac{Ae(s)}{\ell(s)}, \quad \frac{e(s + \epsilon) - e(s)}{\epsilon} \sim -\frac{Ae(s)}{\ell(s)}.$$

It is natural to approach the left-hand side of these equations by derivatives. We then obtain the following system of ODEs

$$\frac{d}{ds}\ell(s) = \frac{Ae(s)}{\ell(s)}, \quad \frac{d}{ds}e(s) = -\frac{Ae(s)}{\ell(s)}, \qquad (72)$$

and since at the beginning of Phase II we still have $L_n \ll \frac{1}{\epsilon}$ and $E_n \sim 1$, the matching conditions for $L_n$, $E_n$ yield the initial conditions

$$\ell(0) = 0, \quad e(0) = 1. \qquad (73)$$
The initial conditions stated in (73) makes the system (72) singular, however the initial conditions are compatible since the solution $\ell$ is bounded away from 0 exponentially fast. Equations (72)–(73) yield the evolution of $L_n$, $E_n$ during Phase II. These equations can be solved in the following implicit form

$$ -\frac{1}{A}\ell - \frac{1}{A}\log(1-\ell) = s, \quad e = 1-\ell. \qquad (74) $$

Notice that this equation (or directly (72)) implies the asymptotic behaviour

$$ \ell(s) \sim \sqrt{2As}, \quad e(s) \sim 1 \quad \text{as } s \to 0. $$

We obtain that $\ell$ increases from $\ell=0$ for $s=0$ to $\ell=1$ as $s$ tends to infinity. On the other hand $e$ decreases from $e=1$ for $s=0$ to $e=0$ as $s$ tends to infinity. We have the following asymptotic behaviour for $\ell$ and $s$

$$ \ell(s) \sim (1 - e^{-1}e^{-As}), \quad e(s) \sim e^{-1}e^{-As} \quad \text{as } s \to \infty. \qquad (75) $$

Combining (71) and (75) yields the evolution of $L_n$, $E_n$ during Phase II. The formula for $\ell(s)$ in (75) implies that $L_n$ remains of order $\frac{1}{\epsilon}$ during the Phase II. $\square$

**Remark 4.4.** The approximations (70) have been computed under the assumption that this displacement $\frac{E_n}{\epsilon}$ is very large. However, the formula for $e(s)$ in (75) as well as the fact that $E_n = e(s)$ imply that for $n$ of order $\frac{1}{A}\frac{1}{\epsilon}\log(\frac{1}{\epsilon})$ the energy $E_n$ becomes of order $\epsilon$ and then the approximations (70) are not any longer valid. This marks the beginning of Phase III. We also notice that during this phase the concentrations $c_k$ experience a displacement of order $\frac{E_n}{\epsilon}$ in the space of cluster sizes $k$ (see lemmas (2.7) and (2.8)).

# 5. Phase III: energy reduction from order $\epsilon$ to order $\epsilon^2$

## 5.1. Early Phase III: matching condition

Phase III begins when $E_n$ becomes of order $\epsilon$, so that all the analysis made for Phases I and II, which rely on the assumption $\frac{E_n}{\epsilon} \gg 1$, is no more valid. However, we can infer the initial state at the beginning of Phase III as the limit value obtained at the end of Phase II, which is illustrated in figure 8: this is expressed in the following lemma.

**Lemma 5.1.** With the notations of the previous results, at the end of Phase II and beginning of Phase III, we have

$$ E_n \approx \epsilon, \quad L_n \approx \frac{1}{\epsilon}, \quad T_n \approx \frac{1}{\epsilon}, \quad Y_{\max}^n \approx 1, \quad c_j(T_n) \sim \epsilon^2 \psi(\epsilon j), \quad j \gg 1, \qquad (76) $$

with $\psi$ defined by (55). Moreover, we have $v(T_n) = w(T_n) = O(\epsilon)$ and $c_j(T_n) \approx U(j)$ for $j = O(1)$.

**Proof.** The energy $E_n \approx \epsilon$ is the definition of the beginning of Phase III, so that we deduce from proposition 4.3 that $n \approx -\frac{\log(\epsilon)}{\epsilon}$. Hence, by (69) we have $L_n \approx \frac{1}{\epsilon}$. Thanks to lemmas 2.7 and 2.8, the maximal displacement for clusters $Y_{\max}^n \approx \frac{E_n}{\epsilon}$ tends to become of order 1, and the time period for one cycle $T_n \approx \frac{E_n}{\epsilon}$ tends to be of order $\frac{1}{\epsilon}$. Since $E_n$ is of order $\epsilon$ we obtain that $w$ and $v$ are of order $\epsilon$ too (cf (11)).
<figure><img src="image_21.png" /><figcaption>**Figure 8.** Numerical simulation of the evolution of the monomers and the energy in Phase II (cf section 4.3). The colour scheme indicates the time evolution in figures (a) and (b). (Left) The phase portrait of $(v, w)$ illustrates the trajectories of the monomers' concentrations when the Energy $E_n$ starts to have an exponential decay at the end of Phase II (see proposition 4.3). (Right) The y-axis vertical axis scale is logarithmic and the horizontal line is $E_n \approx \varepsilon$ where $\varepsilon = 0.02$ in the numerical simulation. The decreasing affine line in red is an illustration of the exponential decay of the energy that takes place in Phase II (see proposition 4.3).</figcaption></figure>

For the size distribution, we recall its description done in proposition 4.1: a Dirac mass of weight $m_n \approx \sigma_n M$, a boundary layer of width $\sigma_n$ in the variable $x = \frac{j}{L_n}$ and the function $\frac{\varepsilon}{L_n} \psi(\frac{j}{L_n})$ for $\frac{j}{L_n} \gg \sigma_n$. We now have $\sigma_n = \frac{\sqrt{D_n}}{L_n} = \frac{\sqrt{E_n}}{\sqrt{\varepsilon L_n}} \approx \varepsilon$ so that the boundary layer width vanishes: the approximation by $\frac{\varepsilon}{L_n} \psi(\frac{j}{L_n}) \sim \varepsilon^2 \psi(\varepsilon j)$ is valid for any $j \gg 1$, and for the boundary layer we have $c_j \approx U(j)$. $\square$

## 5.2. Dynamics of Phase III: a new scaling

From this initial state, it is natural to introduce the following change of variables in order to describe Phase III:

$$w = \varepsilon W, \quad v = \varepsilon V, \quad c_k = \varepsilon^2 C_k, \quad \tau = \varepsilon t. \qquad (77)$$

Notice that the variable $\tau$ here is different from the variable $\tau$ defined for Phases I and II. Then the equations (7)–(9) become

$$\frac{dV}{d\tau} = -VW + V(1 - \varepsilon C_1), \qquad (78)$$

$$\frac{dW}{d\tau} = VW - W, \qquad (79)$$

$$\frac{dC_j}{d\tau} = \bar{J}_{j-1} - \bar{J}_j, \quad \forall j \ge 1, \quad \bar{J}_0 = 0, \quad \bar{J}_j = WC_j - VC_{j+1}, \quad j \ge 1. \qquad (80)$$

It is also natural to define a rescaled LV energy associated to the problem (78)–(80). Notice that if we set $\varepsilon = 0$ in (78)–(79), we obtain the conservation of the rescaled LV energy $\tilde{E}$, i.e.

$$E = \varepsilon \tilde{E}, \quad \tilde{E} = V + W - 2 - \log(VW). \qquad (81)$$
During the Phase III we expect the rescaled energy $\tilde{E}$ to be reduced from values of order one to values of order $\varepsilon^2$, which is the order of magnitude of its equilibrium value. Notice that the solutions of (78)–(80) are small perturbations of the LV equation, though the perturbation is different from the one in Phases I and II: here, for $\tilde{E} \approx 1$, we have that $V$ and $W$ always remain in the order of one, and each phase of the LV cycle lasts for a time of order one in the variable $\tau$. The perturbative term $-\varepsilon VC_1$ yields changes in the rescaled energy $\tilde{E}$ that is illustrated numerically in figures 9 and 10.

The main difference between Phase III, described by means of (78)–(80), and Phases I and II is that the change of energy of the LV oscillations does not take place in the specific interval $[T - \Delta t, t_4]$ of the LV cycle, but during the whole cycle. Moreover, in the analysis of Phases I and II, the main change of the energy during a LV cycle takes place when the concentration wave is closer to cluster sizes $j$ of order one and during those times $c_1$ is much larger than the concentrations $c_j$ with $j \neq 1$. On the contrary, during Phase III, the rescaled concentrations $C_j$ with $j$ of order one have the same order of magnitude as $C_1$. Due to this, we cannot use a perturbative argument to approximate the values of $C_1$ as it was made in the analysis of Phases I and II, but we need to study a problem which involves the whole sequence $\{C_j\}_{j \in \mathbb{N}}$. Let us gather the results for Phase III in the following proposition.

**Proposition 5.2.** Let $T_3^{\text{in}}$ denote the beginning of Phase III. Let us depart from $E_n, L_n, T_n$ and $c_j(T_n)$ described by (76) and define $V, W, \tau$ and $C_j$ by (77), $\tilde{E}_n := \frac{1}{\varepsilon}E_n$ and $\tau_n^{\text{in}} = \varepsilon T_3^{\text{in}}$. On each cycle $\tau_n > \tau_n^{\text{in}}$, we can approximate the dynamics with the following system taken on $[\tau_n, \tau_{n+1}]$:

$$ \frac{dV}{d\tau} = -VW + V, \qquad (82) $$

$$ \frac{dW}{d\tau} = VW - W, \qquad (83) $$

$$ \frac{dC_j}{d\tau} = \frac{W-V}{2}(C_{j-1}-C_{j+1}) + \frac{V+W}{2}(C_{j-1}-2C_j+C_{j+1}), \quad j \ge 2, \qquad (84) $$

$$ \frac{dC_1}{d\tau} = VC_2 - WC_1, \quad \tilde{E}_k = V + W - 2 - \log(VW), \quad V(0) = W(0) > 1, \qquad (85) $$

and with the additional matching condition

$$ \lim_{j \to \infty} C_j(\tau) = \frac{2}{\pi} \quad \forall \tau \in [\tau_n, \tau_{n+1}]. \qquad (86) $$

The change of energy between two successive cycles is then given by

$$ \tilde{E}_{n+1} - \tilde{E}_n \approx \varepsilon \int_{\tau_n}^{\tau_{n+1}} (1 - V(\tau)) C_1(\tau) d\tau, \qquad (87) $$

with $(C_1, V)$ solutions to (82)–(86).

**Proof.** Let us first prove that $C_j \approx \psi(\varepsilon j)$ for $j \approx \frac{1}{\varepsilon}$: this implies the matching condition (86), the system (82)–(85) being written for $j$ large but $\varepsilon j$ small. This is obtained through a transport-diffusion equation as for Phases I and II. We depart from the initial condition given by lemma 5.1
<figure><img src="image_22.png" /><figcaption>(a) Monomers</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>(c) Size distribution</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>(e) Energy</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>(b) Monomers</figcaption></figure>

<figure><img src="image_26.png" /><figcaption>(d) Size distribution</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>(f) Energy</figcaption></figure>

**Figure 9.** Numerical illustration of the evolution of the size distribution of the clusters, the monomers and the energy (81) in Phase III over one cycle (c.f the rescaled system described in the equations (78)–(80)). The time evolution indicated by the colours is the same for the three figures (a), (c) and (e) (resp. (b), (d) and (f)). The figures (c) and (d) are truncated and only show the evolution of the size distribution for sizes smaller than 50. The time $t=0$ corresponds to the end of Phase II and the beginning of Phase III in this simulation.
<figure><img src="image_28.png" /><figcaption>Figure 10. Numerical simulation of the time evolution of the Energy $\tilde{E}$ defined in (81) during Phase III and the beginning of Phase IV. The time $t=0$ corresponds to the end of Phase II and the beginning of Phase III in this simulation. The y-axis vertical axis is logarithmic and $\varepsilon = 0.02$. The decay of the energy is exponential-like until it reaches an order of magnitude of $\varepsilon^2 \approx 10^{-4}$; the order of magnitude of the energy oscillations then become of the same order of magnitude as the energy itself, and then decay (end of Phase III).</figcaption></figure>

$$C_j(\tau_n^{\text{in}}) = \psi(\varepsilon j), \quad \tau_n^{\text{in}} = \varepsilon T_3^{\text{in}} \quad \text{for } j \text{ large,} \qquad (88)$$

We rewrite (80) as a discrete transport-diffusion equation, as we did in (25)

$$\frac{dC_j}{d\tau} = \left(\frac{W-V}{2}\right)(C_{j-1} - C_{j+1}) + \frac{V+W}{2}(C_{j-1} - 2C_j + C_{j+1}). \qquad (89)$$

Denoting as $x$ the variable $\varepsilon j$ and approximating the discrete derivatives in the previous formula by means of continuous derivatives as we did in (27) (writing $C_j(\tau)$ as $C(x, \tau)$ with $x = \varepsilon j$), we obtain

$$\frac{\partial C(x, \tau)}{\partial \tau} = \varepsilon (V - W) \frac{\partial C(x, \tau)}{\partial x} + \frac{V+W}{2} \varepsilon^2 \frac{\partial^2 C(x, \tau)}{\partial x^2}. \qquad (90)$$

Equation (90) allows to estimate the variation of the concentrations $C_j$ with large $j$:

* The periods $\tau_n > \tau_n^{\text{in}}$ of each LV cycle is of order one during Phase III.

* $(V-W)$ is periodic in each LV cycle, hence the transport rate $\varepsilon (V-W)$ cause an oscillation of $C(x, \tau)$ in the variable $x$ of order $\varepsilon$ or equivalently, in the variable $j$ the amplitude of the oscillations is of order one.

* Since the second order term $\frac{V+W}{2} \varepsilon^2 \frac{\partial^2 C(x, \tau)}{\partial x^2}$ is $O(\varepsilon^2)$, significant changes in the shape of the concentrations $C(x, \tau)$ only happen over times $\tau$ of order $\frac{1}{\varepsilon^2}$.

To sum-up, over times $\tau$ much smaller than $\frac{1}{\varepsilon^2}$, we can neglect the modifications of the concentrations $C(x, \tau_n)$, and only small oscillations during each period. Notice that the equations (78)–(80) suggest that there is a change of the rescaled LV energy $\tilde{E}$ of order $\varepsilon$ in each LV cycle as long as $\tilde{E}$ is of order one, so that the number of LV cycles needed for $\tilde{E}$ to become small is of order $\frac{1}{\varepsilon}$. The end of Phase III, characterised by $\tilde{E}$ small, is studied below in proposition 5.4, and concerns a number of cycles in the order $-\frac{\log(\varepsilon^2)}{\varepsilon}$. Since one LV period is of order one
in the variable $\tau$ during the whole Phase III, we deduce that the total duration of Phase III is of order $-\frac{\log(\varepsilon^2)}{\varepsilon} \ll \frac{1}{\varepsilon^2}$ in the variable $\tau$, hence the external concentrations $C(x, \tau) \approx \psi(x)$ remain frozen during Phase III.

This allows to obtain the matching condition (86) (cf prop 4.1). We then obtain (82)–(85) to approximate $C_1$ to the leading order. We expect to have a stable periodic solution of (82)–(86) for each value of the rescaled energy $\tilde{E}$. Notice that the solution of (82)–(83) can be computed for each value of $\tilde{E} \ge 0$ independently on the values of $C_j$. Thus, (84)–(86) become a system of infinitely many ODEs with prescribed functions $W, V$ for each value of $\tilde{E}$.

We end the proof of proposition 5.2 by computing the change of rescaled energy $\tilde{E}$ during each cycle: we use that

$$ \frac{d\tilde{E}}{d\tau} = \varepsilon(1 - V)C_1, $$

which yields (87).

<figure><img src="image_29.png" /></figure>

**Remark 5.3.** Contrarily to Phase II, where (70) provided an explicit formula to estimate the decay of energy, we have in Phase III no way to approximate (82)–(86) in order to estimate the right-hand side of (87): in general this has to be done numerically, except at the end of Phase III, for which we have the following result.

### 5.3. End of Phase III: energy decay

As already said, we expect to have a number of cycles in the order of $\frac{1}{\varepsilon}$ before reaching $\tilde{E}_n \ll 1$. We now compute the integral on the right-hand side of (87) at the end of Phase III, i.e. for $\tilde{E}_n \to 0$. We obtain the following result.

**Proposition 5.4.** Under the assumptions and notations of proposition 5.2, for the range $\varepsilon^2 \ll \tilde{E}_n \ll 1$, we can approximate the decay of energy by

$$ \tilde{E}_{n+1} - \tilde{E}_n \approx -\varepsilon a \tilde{E}_n, \qquad (91) $$

for a given constant $a > 0$ of order one. Hence the energy decays exponentially fast at the end of Phase III. The assumptions of proposition 5.2 are valid until $\tilde{E}_n \approx \varepsilon^2$, which happens after a number of cycles in the order of $\frac{1}{a\varepsilon}\log(\frac{1}{\varepsilon^2})$. At the end of this period, we thus have $\tilde{E}_n = O(\varepsilon^2)$, $V-1 = O(\varepsilon)$, $W-1 = O(\varepsilon)$, $C_j \approx \psi(\varepsilon j)$ for $j \approx \frac{1}{\varepsilon}$.

**Remark 5.5.** We notice that, as at the end of Phase II and Phase III, and despite the fact that the scalings differ, the energy decays exponentially, with a rate in the order of $\varepsilon$ in both cases (cf figures 8(b) and 10). However the constants ($A$ and $a$) are different.

**Remark 5.6.** The energy decay at the end of Phase III is also characterizing the damping of the oscillations of the trajectories of the monomers’ concentrations. proposition 5.4 proves that the trajectories $(V, W)$ are enclosed in a ball of centre $(1, 1)$ and radius $\varepsilon$. However, one can note that the steady-state (10) for the monomers is shifted from the point $(1, 1)$, hence more precise estimate can be found for $V-1$ after $\frac{1}{a\varepsilon}\log(\frac{1}{\varepsilon^2})$ cycles. The damping of the oscillations are studied more precisely in section 6 and the results will be refined in lemma 6.1 using a more precise asymptotic expansion when linearizing around the steady-state.
**Proof.** Let us linearise the set of equations (82)–(86), writing

$$V = 1 + \alpha, \quad W = 1 + \beta, \quad C_j = \frac{2}{\pi}(1 + \eta_j)$$

where $|\alpha|, |\beta|$ and $|\eta_j|$ are small. Then, keeping only the linear terms in $\alpha, \beta$ and $\eta_j$, we obtain the following equations

$$\frac{d\alpha}{d\tau} = -\beta, \qquad (92)$$

$$\frac{d\beta}{d\tau} = \alpha, \qquad (93)$$

$$\frac{d\eta_j}{d\tau} = (\eta_{j-1} - 2\eta_j + \eta_{j+1}), \quad j \ge 2, \qquad (94)$$

$$\frac{d\eta_1}{d\tau} = (\alpha - \beta) + \eta_2 - \eta_1, \quad \lim_{j \to \infty} \eta_j(\tau) = 0 \text{ for all } \tau. \qquad (95)$$

The rescaled energy $\tilde{E}$ can be approximated for small values as

$$\tilde{E} = \frac{1}{2}(\alpha^2 + \beta^2).$$

Then, the solution of (92)–(93) can be written, up to a translation of the origin of time $\tau$ as

$$\alpha = \sqrt{2\tilde{E}\cos(\tau)}, \quad \beta = \sqrt{2\tilde{E}\sin(\tau)}.$$

We notice that we now have periods of the LV cycle of $\tau_n \sim 2\pi$. In fact, it is enough to solve the reference problem

$$\frac{d\varphi_j}{d\tau} = (\varphi_{j-1} - 2\varphi_j + \varphi_{j+1}), \quad j \ge 2, \qquad (96)$$

$$\frac{d\varphi_1}{d\tau} = e^{i\tau} + \varphi_2 - \varphi_1, \qquad (97)$$

$$\lim_{j \to \infty} \varphi_j(\tau) = 0 \text{ for all } \tau. \qquad (98)$$

Since solutions $\{\eta_j\}_{j \in \mathbb{N}}$ to (94)–(95) are recovered from solutions $\{\varphi_j\}_{j \in \mathbb{N}}$ to (96)–(98) by means of

$$\eta_j = \sqrt{2\tilde{E}[\text{Re}(\varphi_j) - \text{Im}(\varphi_j)]}, \quad j \ge 1.$$

We look for solutions of the problem (96)–(98) in the form

$$\varphi_j = A_j e^{i\tau}$$

where the coefficients $A_j$ solve

$$iA_j = A_{j-1} - 2A_j + A_{j+1}, \quad j \ge 2 \qquad (99)$$

$$iA_1 = A_2 - A_1 + 1. \qquad (100)$$

We can look for particular solutions of (99) in the form

$$A_j = (r)^{j-1}, \quad j \in \mathbb{N},$$
which yields the two roots

$$r_{\pm} = \frac{1}{2} \left[ (2+i) \pm \sqrt{4i-1} \right].$$

We compute numerically $|r_{-}| \approx 0.48053 < 1$ while $|r_{+}| > 1$. Therefore, in order to obtain a solution $\varphi_{j}$ satisfying (98), we must have

$$A_j = K_0 (r_-)^{j-1}, \quad j \in \mathbb{N}$$

for some $K_{0} \in \mathbb{C}$. In order to determine $K_{0}$ we use (100). Then

$$K_0 [(1 + i) - r_{-}] = 1.$$

Therefore,

$$\varphi_j = \frac{(r_-)^{j-1}}{[(1+i)-r_-]} e^{i\tau} = \frac{(r_-)^{j-1}}{\left[\frac{i}{2} + \sqrt{i-\frac{1}{4}}\right]} e^{i\tau}, \quad j \ge 1.$$

Thus,

$$\eta_j = \sqrt{2\tilde{E}} \left[ \text{Re} \left( \frac{(r_-)^{j-1}}{\left[ \frac{i}{2} + \sqrt{i-\frac{1}{4}} \right]} e^{i\tau} \right) - \text{Im} \left( \frac{(r_-)^{j-1}}{\left[ \frac{i}{2} + \sqrt{i-\frac{1}{4}} \right]} e^{i\tau} \right) \right], \quad j \ge 1.$$

or, equivalently

$$\eta_j = \sqrt{2\tilde{E}\text{Re}} \left( \frac{(r_-)^{j-1}(1+i)}{\left[ \frac{i}{2} + \sqrt{i-\frac{1}{4}} \right]} e^{i\tau} \right) = \sqrt{2\tilde{E}\text{Re}} ((1+i)\varphi_j), \quad j \ge 1.$$

We can now approximate the integral on the right-hand side of (87) as follows:

$$\begin{aligned}
\tilde{E}_{n+1} - \tilde{E}_n &= \varepsilon \int_0^{\tau_n} (1 - V(\tau; \tilde{E}_n)) C_1(\tau; \tilde{E}_n) d\tau \\
&= -\frac{2\varepsilon}{\pi} \int_0^{\tau_n} \alpha(\tau; \tilde{E}_n) (1 + \eta_1(\tau; \tilde{E}_n)) d\tau \\
&= -\frac{2\sqrt{2\tilde{E}\varepsilon}}{\pi} \int_0^{\tau_n} \cos(\tau) \eta_1(\tau; \tilde{E}_n) d\tau \\
&= -\frac{4\varepsilon\tilde{E}}{\pi} \text{Re} \left( (1+i) K_0 \int_0^{2\pi} \cos(\tau) e^{i\tau} d\tau \right) \\
&= -\frac{4\varepsilon\tilde{E}}{\pi} \text{Re} \left( (1+i) K_0 \int_0^{2\pi} \cos^2(\tau) d\tau \right) = -4\varepsilon\tilde{E}\text{Re}((1+i)K_0)
\end{aligned}$$

with

$$\text{Re}((1+i)K_0) = \text{Re} \left( \frac{(1+i)}{\frac{i}{2} + \sqrt{i-\frac{1}{4}}} \right) = 0.92505 > 0.$$
Then

$$ \tilde{E}_{n+1} - \tilde{E}_n = -4\varepsilon\tilde{E}_n \operatorname{Re}((1+i)K_0) = -a\varepsilon\tilde{E}_n \quad \text{where } a \approx 3.6922 > 0. \quad (101) $$

Phase III continues until the time in which the contribution of the term $-\varepsilon VC_1$ in (78) becomes comparable to the term $-VW + V = -V\beta$, so that the approximation of (78) by (82) becomes invalid. Therefore, we need to compare the terms $\varepsilon C_1$ and $\beta$. Since $C_1$ remains of order one during the whole phase, Phase III ends when $\beta$ becomes of order $\varepsilon$. The order of magnitude of $\beta$ is of order $\sqrt{\tilde{E}_n}$. Thus, Phase III is valid as long as $\tilde{E}_n \gtrsim \varepsilon^2$. Notice that (101) yields an exponential decay for $\tilde{E}_n$ with the form $\exp(-an\varepsilon)$. It then follows that $\tilde{E}_n$ becomes of order $\varepsilon^2$ after $\frac{1}{a\varepsilon}\log(\frac{1}{\varepsilon^2})$ LV cycles. $\square$

# 6. Phase IV: oscillations decay and stabilization

As outlined in the previous section, the final Phase IV starts once the energy level is in the order of $\varepsilon^3$ (or $\varepsilon^2$ for the rescaled variables of Phase III). We recall that this is also the order of magnitude of our system's equilibrium energy value $\tilde{E}$, see (13). We decompose Phase IV into two successive stages: first, a phase similar to the end of Phase III, where now the oscillations (yet no longer the energy) decay and finally become negligible. Secondly, the trend to equilibrium occurs through an approximate parabolic equation, with non-oscillatory coefficients.

## 6.1. End of Phase III and early Phase IV: decay of the oscillations

In the previous section, proposition 5.4 gives the decay of the energy until the end of Phase III as well as the characterisation of the monomers and the cluster distribution until the end of Phase III, i.e. until $\tilde{E}_n \approx \varepsilon^2$, resp. $E_n \approx \varepsilon^3$. The following lemma also characterizes the end of Phase III, but allows us to describe more precisely the damping of the oscillations occurring at the very end of Phase III and the beginning of Phase IV.

**Lemma 6.1.** Under the assumptions of proposition 5.4, the oscillations become negligible after a number of cycles in the order of $\frac{\log(\varepsilon^{-1})}{\varepsilon}$ and we have

$$ V = 1 + O(\varepsilon^{3/2}), \quad W = 1 - \varepsilon C_{\infty} + O(\varepsilon^{3/2}), \quad C_1 = C_{\infty} + O(\varepsilon). $$

During these cycles, the changes in the size distribution $C_j$ have been negligible, so that $C_{\infty} = \frac{2}{\pi}$ and $C_j \approx \psi(\varepsilon j)$.

**Proof.** In (78), when linearizing and taking $W = 1 + \beta$ (cf Proof of proposition 5.4), the term $-\varepsilon VC_1$ becomes comparable to $-V\beta$ and it cannot be any longer ignored to the leading order. In order to understand the evolution of $W$, $V$, $C_j$ during this phase, it is natural to introduce new variables $\tilde{\alpha}$, $\tilde{\beta}$ and $\tilde{\eta}_j$ defined by

$$ V = 1 + \varepsilon \tilde{\alpha}, \quad W = 1 - \varepsilon C_{\infty} + \varepsilon \tilde{\beta}, \quad C_j = C_{\infty} (1 + \varepsilon \tilde{\eta}_j), $$

where $C_{\infty} = \frac{2}{\pi}$ at the 'beginning of the very end' of Phase III but is expected to slowly change during Phase IV. Then, keeping the leading terms in (78)–(80) we obtain the following approximate model

$$ \frac{d\tilde{\alpha}}{d\tau} = -\tilde{\beta} - \varepsilon\tilde{\alpha}\tilde{\beta} - \varepsilon C_{\infty}\tilde{\eta}_1 \quad (102) $$
$$ \frac{d\tilde{\beta}}{d\tau} = \tilde{\alpha}(1 - \varepsilon C_{\infty}) + \varepsilon\tilde{\alpha}\tilde{\beta} \qquad (103) $$

$$ \frac{d\tilde{\eta}_j}{d\tau} = \varepsilon \left( \frac{\tilde{\beta} - \tilde{\alpha}}{2} \right) (\tilde{\eta}_{j-1} - \tilde{\eta}_{j+1}) + (\tilde{\eta}_{j-1} - 2\tilde{\eta}_j + \tilde{\eta}_{j+1}) \qquad (104) $$

$$ -\varepsilon C_{\infty} (\tilde{\eta}_{j-1} - \tilde{\eta}_{j}) + \varepsilon \frac{\tilde{\beta} + \tilde{\alpha}}{2} (\tilde{\eta}_{j-1} - 2\tilde{\eta}_{j} + \tilde{\eta}_{j+1}), \quad j \ge 2 $$

$$ \frac{d\tilde{\eta}_1}{d\tau} = \tilde{\eta}_2 - \tilde{\eta}_1 + \tilde{\alpha} - \tilde{\beta} + C_{\infty} + \varepsilon (\tilde{\alpha}\tilde{\eta}_2 - \tilde{\beta}\tilde{\eta}_1 + C_{\infty}\tilde{\eta}_1). \qquad (105) $$

We keep only the leading order terms to get

$$ \frac{d\tilde{\alpha}}{d\tau} = -\tilde{\beta} \qquad (106) $$

$$ \frac{d\tilde{\beta}}{d\tau} = \tilde{\alpha} \qquad (107) $$

$$ \frac{d\tilde{\eta}_j}{d\tau} = \tilde{\eta}_{j-1} - 2\tilde{\eta}_j + \tilde{\eta}_{j+1}, \quad j \ge 2 \qquad (108) $$

$$ \frac{d\tilde{\eta}_1}{d\tau} = \tilde{\eta}_2 - \tilde{\eta}_1 + \tilde{\alpha} - \tilde{\beta} + C_{\infty}. \qquad (109) $$

We notice that the system obtained only differs from (92)–(95) by the source term $+C_{\infty}$ in the equation for $\tilde{\eta}_1$, which expresses the fact that there is a constant non-negligible difference between the influx of monomers (due to the depolymerisation rate $V$) and the outflux due to the polymerisation $W \approx V - C_{\infty}$.

We first solve (106) and (107) as we did for (92) and (93): we notice that solutions to (106) and (107) are $2\pi$ periodic (as for (92) and (93)), which implies that periods $\tau_n$ of (102) and (103) are $2\pi + O(\varepsilon)$. With respect to the energy of (102)–(103), we have

$$ \tilde{E} \approx \varepsilon^2 \frac{\tilde{\alpha}^2 + (\tilde{\beta} - C_{\infty})^2}{2} =: \varepsilon^2 \hat{E} + \varepsilon^2 \frac{C_{\infty}^2 - 2\tilde{\beta}C_{\infty}}{2}, $$

where we have defined $\hat{E} := \frac{\tilde{\alpha}^2 + \tilde{\beta}^2}{2}$, so that up to choosing a proper initial time we have $\tilde{\alpha} = \sqrt{2\hat{E}}\cos(\tau)$ and $\tilde{\beta} = \sqrt{2\hat{E}}\sin(\tau)$. The oscillations of the energy $\tilde{E}$ are now of the same order of magnitude as its value (see figure 10, values around $E = 10^{-4}$, for an illustration), as they are expressed by the term $-\varepsilon^2\tilde{\beta}C_{\infty}$, so that we also define its average over a time period:

$$ \langle \tilde{E} \rangle = \varepsilon^2 \Delta_{\hat{E}} + \varepsilon^2 \frac{C_{\infty}^2}{2}, $$

so that it tends to its equilibrium $\varepsilon^2 \frac{C_{\infty}^2}{2}$ when $\Delta_{\hat{E}}$ vanishes. We also compute

$$ \frac{d\hat{E}}{d\tau} = \varepsilon^{-2} \frac{d\tilde{E}}{d\tau} + C_{\infty} \frac{d\tilde{\beta}}{d\tau} = -\varepsilon\tilde{\alpha}C_{\infty}\tilde{\eta}_1 - \varepsilon\tilde{\alpha}\tilde{\beta}(C_{\infty} + \tilde{\alpha} - \tilde{\beta}), \qquad (110) $$

so that over a period $\tau_n = 2\pi + O(\varepsilon)$, we have

$$ \Delta_{\hat{E}} \sim -\varepsilon C_{\infty} \frac{1}{2\pi} \int_{\tau_n}^{\tau_n + 2\pi} \tilde{\alpha}(s) \tilde{\eta}_1(s) ds, \qquad (111) $$

and we are left to compute $\tilde{\eta}_1$ in order to estimate the decay of energy.
**Solution to (108) and (109).** Recognizing in (108) and (109) a linear system with two sources $\tilde{\alpha} - \tilde{\beta}$ and $C_\infty$, we superimpose the solution $\tilde{\eta}_j^{(1)}$ to (106)–(109) with $C_\infty = 0$—this is exactly the solution to (94)–(95) computed above for $\eta_j$, where we simply replace $\alpha, \beta, \tilde{E}$ by $\tilde{\alpha}, \tilde{\beta}$ and $\tilde{E}$ respectively—with the solution $\tilde{\eta}_j^{(2)}$ of the system with a constant source $C_\infty$, namely

$$ \frac{d\tilde{\eta}_j^{(2)}}{d\tau} = \tilde{\eta}_{j-1}^{(2)} - 2\tilde{\eta}_j^{(2)} + \tilde{\eta}_{j+1}^{(2)}, \quad j \ge 2 \qquad (112) $$

$$ \frac{d\tilde{\eta}_1^{(2)}}{d\tau} = \tilde{\eta}_2^{(2)} - \tilde{\eta}_1^{(2)} + C_{\infty}. \qquad (113) $$

To solve (112)–(113), let us compute the system satisfied by its Laplace transform $L_j(z) := \int_0^\infty e^{-z\tau} \tilde{\eta}_j^{(2)}(\tau)d\tau$:

$$ \begin{aligned} zL_1(z) &= L_2(z) - L_1(z) + \frac{C_\infty}{z}, \\ zL_j(z) &= L_{j-1}(z) - 2L_j(z) + L_{j+1}(z). \end{aligned} $$

We look for solutions $L_j(z) = B(z)\theta(z)^{j-1}$, to find

$$ \theta^2 - \theta(2+z) + 1 = 0, $$

so that

$$ \theta_{\pm}(z) = 1 \pm \frac{z}{2} \pm \sqrt{z + \frac{z^2}{4}}. $$

Since we want $|\theta| < 1$ to ensure that the solutions vanish when $j \to \infty$, we have to choose $\theta_-(z)$ which satisfies $|\theta_-(z)| < 1$ for $z \notin \{0, 4\}$. We then compute $B(z)$ with the equation for $L_1$, and find

$$ zB(z) = B(z)(\theta_{-}(z) - 1) + \frac{C_{\infty}}{z}, $$

so that finally

$$ L_j(z) = \frac{C_\infty \theta_-(z)^{j-1}}{z(z - \theta_-(z) + 1)} = \frac{C_\infty \left(1 + \frac{z}{2} - \sqrt{z + \frac{z^2}{4}}\right)^{j-1}}{z\left(\frac{z}{2} + \sqrt{z + \frac{z^2}{4}}\right)}. $$

We can now compute the inverse Laplace transform, and obtain, taking $\delta > 0$,

$$ \tilde{\eta}_j^{(2)}(z)(\tau) = \frac{C_\infty}{2i\pi} \int_{\delta-i\infty}^{\delta+i\infty} e^{z\tau} \frac{\left(1 + \frac{z}{2} - \sqrt{z + \frac{z^2}{4}}\right)^{j-1}}{z\left(\frac{z}{2} + \sqrt{z + \frac{z^2}{4}}\right)} dz. $$

We can deform the contour of integration $(\delta - i\infty, \delta + i\infty)$ by making $\delta \to 0$ and taking the contour $(-\infty - i\delta, -i\delta) \cup \{|z| = \delta, \text{Re}(z) > 0\} \cup (i\delta, -\infty + i\delta)$ so that the main contribution
for the inverse Laplace transform comes from the half circle $C := \{|z| = \delta, \text{Re}(z) > 0\}$, on which we compute

$$ \tilde{\eta}_1^{(2)} \approx \frac{C_\infty}{2i\pi} \int_C e^{z\tau} z^{-\frac{3}{2}} dz \approx \frac{C_\infty \sqrt{\tau}}{2i\pi} \int_C \frac{e^\eta d\eta}{\eta^{3/2}} \approx \frac{C_\infty \sqrt{\tau}}{i\pi} \int_C \frac{e^\eta d\eta}{\eta^{1/2}} = C \frac{C_\infty \sqrt{\tau}}{i\pi}, $$

for a given constant $C$. This allows us to compute the change of energy due to the contribution of $\tilde{\eta}_1^{(2)}$, see below.

Let us also verify that for $j \gg 1$ the influence of $\tilde{\eta}_j^{(2)}$ remains negligible on the timecourse of this phase, which is in the order of $1/\epsilon$, so that the large-sizes distribution remains approximately unchanged. For $z \ll 1$, we compute

$$ L_j(z) \sim C_\infty \frac{(1 - \sqrt{z})^{j-1}}{z^{\frac{3}{2}}} \sim C_\infty \frac{e^{-j\sqrt{z}}}{z^{\frac{3}{2}}}. $$

Thus, computing the inverse Laplace transform for $j \gg 1$—more precisely $j = \sqrt{\tau}x$, since the typical timescale for the discrete diffusion equation is $O(\sqrt{\tau})$—and differentiating it in time, we obtain

$$ \frac{\partial}{\partial \tau} \tilde{\eta}_j^{(2)} \approx \frac{C_\infty}{2i\pi} \int_C e^{z\tau} \frac{e^{-j\sqrt{z}}}{\sqrt{z}} dz = \frac{C_\infty}{2i\pi\sqrt{\tau}} \int_C e^{\eta} \frac{e^{-x\sqrt{\eta}}}{\sqrt{\eta}} d\eta, $$

where we have used the change of variables $\eta = z\tau$. We can see by integrating by parts that this expression decreases faster than any polynomial when $x \to \infty$, which ensures that the contribution of $\tilde{\eta}_j^{(2)}$ vanishes for $j \to \infty$. Moreover, since the time derivative of $\tilde{\eta}_j^{(2)}$ is in the order of $1/\sqrt{\tau}$, the time needed for the clusters distribution departing from $j \approx 1$ to come to $j \approx \frac{1}{\epsilon}$ is of order $\frac{1}{\epsilon^2}$, which appears much larger than the total time of this period in the $\tau$ variable, which is in the order of $\frac{\log(\epsilon^{-1})}{\epsilon}$, see below.

**Computation of the change of energy.** Writing $\tilde{\eta}_1 = \tilde{\eta}_1^{(1)} + \tilde{\eta}_1^{(2)}$, we can now gather the two contributions for the change of energy (cf (101) and (111)), writing

$$ \Delta_{\hat{E}} = \Delta_{\hat{E}}^{(1)} + \Delta_{\hat{E}}^{(2)}, \\ \approx -\epsilon a \hat{E} - \epsilon \sqrt{2\hat{E}} \frac{C_\infty^2}{2\pi} \frac{C}{i\pi} \int_{\tau}^{\tau+2\pi} \sqrt{\tau} \cos(\tau) d\tau. $$

As already seen above, we have $a > 0$ so that the first term ensures an exponential decay for $\hat{E}$. However, the sign of the second term changes according to the exact time-point $\tau \in [\tau_n, \tau_n + 2\pi]$, and since it is weighted by $\sqrt{\hat{E}}$ it may after a while reveal dominant. However, at the beginning of this stage, we have $\hat{E} \gg 1$ so that the second term is negligible compared to the first one. We can thus assume $\tau > \tau_n \gg 1$ without loss of generality. We now estimate the integral in $\Delta_{\hat{E}}^{(2)}$ for large $\tau$,
**Remark 6.3.** The system (114)–(115) has been proposed in [15] in the context of a chemo-taxis model, when the chemical is produced at the boundary, and studied in [8, 9]. Our case corresponds to the critical case $M = 1$, where the long-time convergence towards $e^{-x}$ has been proved in [8, Theorem 1.1.], under the assumption that the second moment is initially finite.

**Proof.** At the end of Phase III, we have seen that $\tilde{\alpha}$ and $\tilde{\beta}$ became negligible compared to $C_\infty$ and to $\tilde{\eta}_j$. We can approximate system (108)–(109) by

$$ \frac{d\tilde{\eta}_j}{d\tau} = \tilde{\eta}_{j-1} - 2\tilde{\eta}_j + \tilde{\eta}_{j+1}, \quad j \ge 2 \qquad (117) $$

$$ \frac{d\tilde{\eta}_1}{d\tau} = \tilde{\eta}_2 - \tilde{\eta}_1 + C_\infty, \qquad (118) $$

or equivalently, by recalling $C_j = C_\infty(1 + \varepsilon\tilde{\eta}_j)$ and taking the equations for $C_j$ for $j$ of order one, we can approximate (108)–(109) as

$$ \frac{dC_j}{d\tau} = (C_{j-1} - 2C_j + C_{j+1}), \quad j \ge 2 \qquad (119) $$

$$ \frac{dC_1}{d\tau} = C_2 - C_1, \qquad (120) $$

where, at the beginning of Phase IV, $C_j$ is approximately constant equal to $C_\infty = \frac{\pi}{2}$ for $j$ large. We notice that equations (119)–(120) yield an evolution for $C_j$ independent of $\tilde{\alpha}, \tilde{\beta}$. equation (119) is a discrete diffusion equation, so that solutions of (119)–(120) converge, in times $\tau$ of order one, to $C_j = C_\infty$.

To understand what happens for larger times, we now examine the evolution of the concentrations $C_j$ with $j$ of order $\frac{1}{\varepsilon}$. Let us recall (89), where we keep only the leading order term for the diffusion and the $\varepsilon$ order for the transport:

$$ \frac{dC_j}{d\tau} = -\frac{\varepsilon C_\infty}{2} (C_{j-1} - C_{j+1}) + (C_{j-1} - 2C_j + C_{j+1}). $$

As done for the previous phases (we recall that since the end of Phase II we have $L_n \sim \frac{1}{\varepsilon}$), we use the approximation of $C_j(\tau)$ as $C(x, \tau)$ with $x = \varepsilon j$. Then

$$ \frac{\partial C(x, \tau)}{\partial \tau} = -\varepsilon^2 C_\infty \frac{\partial C(x, \tau)}{\partial x} + \varepsilon^2 \frac{\partial^2 C(x, \tau)}{\partial x^2}. \qquad (121) $$

This equation must be solved with the initial condition (116) obtained in lemma (6.1) at the beginning of Phase IV. The value of $C_\infty$ can be approximated identifying it with the value obtained for the outer concentrations $C(x, \tau)$ as $x \to 0^+$. Introducing also the time scale $\bar{\tau} = \varepsilon^2(\tau - \tau_{in})$ where $\tau_{in} = \varepsilon T_4^{in}$ is the time when we assume the Phase IV to begin, we obtain (114).

In order to determine the boundary condition to be imposed at $x=0$ for the solutions to (114), we use the second condition in (6), that is

$$ 1 = \sum_{k=1}^{\infty} \varepsilon C_k = \sum_{k=1}^{\infty} \varepsilon C(\varepsilon k, \bar{\tau}) \approx \int_{0}^{\infty} C(x, \bar{\tau}) dx, \qquad (122) $$

which implies that $\int_0^\infty C(x, \bar{\tau}) dx$ is constant for all $\bar{\tau}$ and integration of (114) yields (115).
On the other hand, the first condition in (6) implies an additional normalisation condition, namely

$$1 = \varepsilon(V+W) + \varepsilon^2 \sum_{j=1}^{\infty} j C_j = \varepsilon(V+W) + \sum_{j=1}^{\infty} \varepsilon j C(\varepsilon j, \bar{\tau}) \varepsilon \\ \approx \varepsilon(V+W) + \int_{0}^{\infty} x C(x, \bar{\tau}) dx.$$

Using that V and W are close to 1 we then obtain the following normalisation condition as $\varepsilon \to 0$

$$\int_{0}^{\infty} x C(x, \bar{\tau}) dx = 1. \qquad (123)$$

The problem (114)–(116) yields the evolution of the cluster concentrations during Phase IV. The steady states of the system (114)–(115) have the form

$$C_s(x) = K \exp(-Kx), \quad x > 0,$$

where $K > 0$ is an arbitrary constant. Notice that $C_s(x)$ satisfies the normalisation condition (122) for any $K > 0$. In order to determine K we use the normalisation condition (123). We calculate

$$\int_{0}^{\infty} x C_s(x) dx = \frac{1}{K}.$$

Thus, (122) implies $K = 1$. Therefore the equilibrium distribution of clusters is given by

$$C_s(x) = \exp(-x), \quad x > 0. \qquad (124)$$

<figure><img src="image_30.png" /></figure>

# 7. Discussion

## 7.1. Other choices of initial concentrations

It is worth to note that choosing initial cluster distributions different from the ones considered in the previous subsections, i.e. the energy at initial time $E_0$ of order 1 and the characteristic length at initial time $L_0$ of order 1 (or more generally smaller than $1/\varepsilon$), the dynamics of the system differ from those described above in the four consecutive phases.

**Initial condition widely spread along the size axis.** Certainly, we can skip some of the phases starting with initial cluster contributions having, say $E_0 \approx 1$ and $L_0 \approx 1/\varepsilon$. This would result in an evolution without Phase I, starting directly in Phase II. If, in addition to $L_0 \approx 1/\varepsilon$ we have also $E_0 \approx \varepsilon^2$, we could have evolutions starting directly in Phase III, skipping the two previous phases.

**Initial condition concentrated on the small sizes.** If we take as starting point values of $(w^0, v^0)$ for which the concentrations $c_k$ tend to move towards smaller values of $k$, we can
obtain very large changes of the energy just in the early transient states. For instance, if we begin with $E_0$ close to 1 with $w^0 = v^0 < \varepsilon$, and concentrations $c_k^0$ concentrated in values of k of order one, we obtain a concentrations wave moving towards lower values of k. This results in a large increase of $c_1$ and as a consequence large changes of the Energy E. Due to this the initial, transient dynamics of the whole system can result in values of $(w, v)$ that differ significantly from a LV orbit. After the values of $(w, v)$ reach the line $w = v > \varepsilon$ we obtain a dynamics that can be described as indicated above, for suitable values of E.

**Initial condition concentrated in a dirac mass along the size axis.** Moreover, we can make choices of the initial cluster concentrations which differ more drastically of the previous phases, because they cannot be characterized in a meaningful manner by a single characteristic length $L_0$. This would be the case, for instance with distributions with the form

$$c_{j,0} = \varepsilon \delta_{j,R_{\varepsilon}} \quad \text{where} \quad R_{\varepsilon} = \frac{a}{\varepsilon}, \quad \text{for some } a \in (0, 1).$$

We will assume also that $v$ and $w$ are of order one. In this case we have approximately $E_0 \simeq 1 - a > 0$ for small $\varepsilon$. We will assume then that initially $v = w \simeq \frac{1-a}{2}$ by definiteness. On the other hand, it is not clear what should be the definition of $L_0$. Taking into account the set of values where the mass of the clusters is concentrated we should take $L_0 = R_{\varepsilon} = \frac{a}{\varepsilon}$. On the other hand, there is not dispersion in the concentration distributions and therefore it would be also reasonable to assume that $L_0 = 1$. Actually, the evolution of the cluster concentrations in this case differs from the one described in the previous sections. For these initial concentrations, we obtain oscillations of the concentrations $c_j$ in the space of cluster sizes $j$ in a manner similar to the one described in sections 2.3, 4.2, while at the same time diffusion in the space of clusters takes place (cf (36)). During the first oscillations, the values of $c_1$ are exponentially small in $\frac{1}{\varepsilon}$ and, due to this, the energy $E_n$ which characterizes each of the LV cycles remains almost constant. The diffusion in the space of cluster sizes increases slowly the width of the cluster distributions. The values of $c_1$ become significant (i.e. non-exponentially small), after $\frac{1}{\varepsilon}$ LV cycles. Actually, after this number of cycles, the distribution of clusters has a characteristic length of order $\frac{1}{\varepsilon}$ and the corresponding evolution becomes similar to the one described in Phase II, with the only difference that the cluster concentrations is not necessarily given by the Gaussian distribution $\psi(x)$ in (55). We then obtain an evolution similar to the one described in section 4.2, but where an additional evolution of the concentrations $\phi_n(x)$ by means of the iteration (39) must be included (cf figure 11). During this modified Phase II, the energy $E_n$ decreases until reaching values of order $\varepsilon$, when the system starts to evolve according to the mechanism described in Phase III, and eventually, the concentrations approach to the equilibrium as described in Phase IV.

## 7.2. Concluding remarks

In this paper, using asymptotic and numerical methods, we have described the behaviour of the solutions of the system of equations (7)-(9). This system of equations couples the dynamics of the classical LV oscillator with a generalized Becker-Döring model that describes the concentrations of clusters with each given length. The coupled system was introduced in [10]. This paper describes the way in which the chemical concentrations $v, w, \{c_k\}_{k \in \mathbb{N}}$ approach to
<figure><img src="image_31.png" /><figcaption>**Figure 11.** Numerical simulation of the evolution of the cluster size distribution and the energy taking as initial condition a Dirac mass located at the size $R_{\epsilon} = 25$ and $v = w \approx \frac{1-a}{2}$. The time evolution is indicated by the colours in figure (b). The figure (b) is truncated on the y-axis vertical axis at 0.001 and on the horizontal axis at 100. In the figure (a), the horizontal red lines illustrate the thresholds when the energy would be inferior to $\epsilon$.</figcaption></figure>

their equilibrium values if the parameter $\epsilon$ is small. This equilibrium is reached by means of a mechanism in which the concentrations of $v$, $w$ oscillate, with decreasing amplitude around a centre $(\epsilon, \epsilon)$, with the monomers spreading in the set of chemical concentrations $\{c_k\}_{k \in \mathbb{N}}$ in increasingly larger values of cluster sizes $k$, until reaching sizes of order $k \approx \frac{1}{\epsilon}$. In this paper we have derived a set of equations that describes the decrease of the amplitude of the oscillations in the space $(v, w)$ as well as a sequence of iterative mappings that describes the evolution of the chemical concentrations $\{c_k\}_{k \in \mathbb{N}}$. A characteristic feature of the mechanism that yields the damping of oscillations that we have obtained is the onset of some oscillations of the whole set of chemical concentrations $\{c_k\}_{k \in \mathbb{N}}$ in the space of cluster sizes $k$.

In the course of the analysis, we have raised a series of interesting asymptotic problems, e.g. the trend of the iterative system (64)–(65)–(66) towards its steady state (56)–(58), or the study of the nonlinear free-boundary problem (114)–(116), or yet a rigorous justification for the heuristic description of the ‘collision’ of the cluster concentration waves with regions where $k \approx 1$, in the fast-moving regions $t \in [T_n - \Delta t, T_n + \Delta t]$, see lemma 3.4. Their proof is let for future work [11].

In this paper we have assumed that the reaction coefficients are constant and we have modified the time scale to normalize them as one. It would be relevant to understand if the damping mechanism for the chemical oscillations obtained in this paper is also valid for more general choices of the chemical coefficients. Constant coefficients are well-adapted for polymerisation by one or the two ends of fibrils, for instance, whereas linear or affine coefficients could take into account secondary nucleation, i.e. lateral polymerisation, and other more complex laws could take into account varying geometries of the clusters [21, 25]. In the case considered in this paper, the model can be rewritten in an obvious manner as a perturbation of the LV model.
However, in the case of more general coefficients, it is not clear if the model can be rewritten as a perturbation of LV or some more involved oscillator.

## Data availability statement

The data that support the findings of this study are openly available at the following URL/DOI: https://github.com/mmezache/BiMono_Becker_Doring.

## Appendix A. Extra computations for the LV system

In this appendix, we prove the asymptotic expansions of lemmas 2.2–2.5. We discuss also various orders of magnitude, which are useful to get some insight in the oscillatory behaviour of the chemical concentrations v, w.

**Proof of lemma 2.2: from 0 to t₁.** We develop asymptotically (14) as

$$ \frac{dv}{dt} = -vw + \varepsilon v = -v(1 - v + O(\varepsilon \log(\varepsilon))). $$

Hence,

$$ \frac{1}{v(1-v)} \frac{dv}{dt} = \frac{d}{dt} \left( \log \left( \frac{v}{1-v} \right) \right) = \frac{d}{dt} \log \left( -1 + \frac{1}{1-v} \right) = -1 + O(\varepsilon \log(\varepsilon)) $$

and we deduce

$$ \log \left( -1 + \frac{1}{1-v} \right) = -t + O(\varepsilon \log(\varepsilon)t), \quad -1 + \frac{1}{1-v} = e^{-t(1+O(\varepsilon \log(\varepsilon)))}, $$

so finally

$$ v(t) = 1 - \frac{1}{1 + e^{-t}} + O(e^{-t} \varepsilon \log(\varepsilon)(t+1)). $$

We then compute that $v(t_1) = \varepsilon$ if $e^{-t_1} = \varepsilon + O(e^{-t_1 \varepsilon \log(\varepsilon)(t_1+1)})$, which implies

$$ t_1 \sim -\log(\varepsilon). $$

The equation for w may now be expanded:

$$ \frac{dw}{dt} = vw - \varepsilon w = w \left( 1 - \frac{1}{1 + e^{-t}} + O(\varepsilon \log(\varepsilon)) \right), $$

hence

$$ \frac{d}{dt} \log(w) = 1 - \frac{e^t}{e^t + 1} + O(\varepsilon \log(\varepsilon)) = \frac{d}{dt} (t - \log(e^t + 1)) + O(\varepsilon \log(\varepsilon)), $$

so by integrating we get

$$ \log w(t) = -\log 2 + t - \log(e^t + 1) + \log(2) + O(\varepsilon \log(\varepsilon)t) $$
and finally

$$w(t) = \frac{e'}{1+e'} (1+O(\varepsilon \log(\varepsilon) t)),$$

which provides us with $w(t_1) = 1 + O(\varepsilon(\log \varepsilon)^2)$, and gives (20). More precisely, we compute from the energy

$$1 = w(t_1) - \varepsilon - \varepsilon \log \left( \frac{w(t_1)}{\varepsilon} \right) \implies w(t_1) = 1 - \varepsilon \log(\varepsilon) + \varepsilon + o(\varepsilon).$$

**Proof of lemma 2.3: From $t_1$ to $t_2$.** During this second phase, $v$ goes from $\varepsilon$ to exponentially small values, whereas $w$ goes from approximately 1 to $\varepsilon$. We thus define an intermediate time $t_{1,2}$ by $v(t_{1,2}) := \varepsilon^{3/2}$, and the energy equality becomes

$$1 = w(t_{1,2}) - \varepsilon \log \left( \frac{w}{\sqrt{\varepsilon}} \right) - 2\varepsilon + O(\varepsilon^{3/2}) \implies w(t_{1,2}) = 1 - \frac{\varepsilon}{2} \log(\varepsilon) + 2\varepsilon + O(\varepsilon^{3/2}),$$

from which we deduce that on $(t_1, t_{1,2})$,

$$\frac{dv}{dt} = -v(w - \varepsilon) = -v(1 + O(\varepsilon \log(\varepsilon))) \implies v(t) = \varepsilon e^{-(t-t_1)} + O(\varepsilon \log(\varepsilon)(t-t_1)),$$

hence we have

$$\varepsilon^{3/2} = \varepsilon e^{-(t_{1,2}-t_1)+O(\varepsilon \log(\varepsilon)(t_{1,2}-t_1))} \implies t_{1,2} - t_1 \sim -\frac{1}{2} \log(\varepsilon)$$

and on $(t_{1,2}, t_2)$ we have $v = O(\varepsilon^{3/2}) \ll \varepsilon$ hence

$$\frac{dw}{dt} = -w\left(\varepsilon + O\left(\varepsilon^{3/2}\right)\right) \implies w(t) = \left(1 - \frac{\varepsilon}{2}\log(\varepsilon) + 2\varepsilon + O\left(\varepsilon^{3/2}\right)\right)e^{-\varepsilon(t-t_{1,2})},$$

we then use $w(t_2) = \varepsilon$ to obtain

$$\begin{aligned} \varepsilon &= \left(1 - \frac{\varepsilon}{2}\log(\varepsilon) + 2\varepsilon + O\left(\varepsilon^{3/2}\right)\right) e^{-\varepsilon(t_2-t_{1,2})} \\ \implies t_2 - t_{1,2} &\sim -\frac{1}{\varepsilon}\log(\varepsilon), \end{aligned} \qquad (125)$$

so that

$$t_2 - t_1 \sim -\frac{1}{\varepsilon} \log(\varepsilon).$$

If we need an intermediate time $t_2^*$ such that $t_2^* - t_{1,2} \sim \frac{1}{\sqrt{\varepsilon}}$, useful for the typical symetric time interval during which mass accumulates in $c_1$, we see that $w(t_2^*) \sim e^{-\sqrt{\varepsilon}} \sim 1$, justifying the approximations carried out for $\Delta t$ in lemma 3.4.

We can also compute explicitly

$$\frac{dv}{dt} = -v(w - \varepsilon) = -v\left(\left(1 - \frac{\varepsilon}{2}\log(\varepsilon) + 2\varepsilon + O\left(\varepsilon^{3/2}\right)\right)e^{-\varepsilon(t-t_{1,2})} - \varepsilon\right),$$
which implies

$$v(t) = \varepsilon^{3/2} \exp \left( - \left( \frac{1}{\varepsilon} - \frac{1}{2} \log(\varepsilon) + 2 + O(\varepsilon^{1/2}) \right) \left( 1 - e^{-\varepsilon(t-t_{1,2})} \right) + \varepsilon(t-t_{1,2}) \right).$$

Taken in $t = t_2$ and using (125) we get

$$\begin{aligned} v(t_2) &= \varepsilon^{3/2} \exp \left( -\frac{1}{\varepsilon} + \frac{1}{2} \log(\varepsilon) - 2 + O(\varepsilon^{1/2}) + \frac{1}{\varepsilon} \varepsilon - \log(\varepsilon) \right) \\ &= \varepsilon^{3/2} \exp \left( -\frac{1}{\varepsilon} - \frac{1}{2} \log(\varepsilon) - 1 + O(\varepsilon^{1/2}) \right) \sim \varepsilon e^{-1} e^{-\frac{1}{\varepsilon}}, \end{aligned}$$

formula which we can also obtain from the conservation of energy:

$$\begin{aligned} E &= v + w - 2\varepsilon - \varepsilon \log \left( \frac{vw}{\varepsilon^2} \right) \\ 1 &= v + \varepsilon - 2\varepsilon - \varepsilon \log \left( \frac{vw}{\varepsilon^2} \right) \\ 1 &= -\varepsilon - \varepsilon \log \left( \frac{v\varepsilon}{\varepsilon^2} \right) + O(\varepsilon^{3/2}) \\ v &\sim \varepsilon e^{-1} \exp \left( -\frac{1}{\varepsilon} \right) \end{aligned}$$

**Proof of lemma 2.4: from $t_2$ to $t_3$.** We divide the trajectory connecting $(v(t_2), w(t_2))$ with $(v(t_3), w(t_3))$ into two (symmetric) parts, by defining $t_{2,3}$ by

$$v(t_{2,3}) = w(t_{2,3}).$$

For $t \in [t_2, t_{2,3}]$ we have $v(t) \le w(t)$, and for $t \in [t_{2,3}, t_3]$ we have $v(t) \ge w(t)$. During the first interval we have that $v$ is much smaller than $\varepsilon$ and we approximate (14) as

$$\frac{dv}{dt} = -vw + v\varepsilon, \quad \frac{dw}{dt} = -\varepsilon w.$$

We then obtain the approximations

$$\begin{aligned} w &\sim \varepsilon \exp(-\varepsilon(t-t_2)) \\ v &\sim v(t_2) \exp(\varepsilon(t-t_2) - [1 - \exp(-\varepsilon(t-t_2))]) \end{aligned}$$

for $t \in [t_2, t_{2,3}]$. We can determine $t_{2,3}$ by means of the condition $v(t_{2,3}) = w(t_{2,3})$ applied to the energy conservation: since $v = w = O(e^{-1/\varepsilon})$, we compute

$$1 = 2v - 2\varepsilon - \varepsilon \log \left( \frac{v^2}{\varepsilon^2} \right) = -2\varepsilon \left( 1 + \log \left( \frac{v}{\varepsilon} \right) \right) + O \left( e^{-\frac{1}{\varepsilon}} \right).$$

Hence $v(t_{2,3}) = w(t_{2,3}) \sim \varepsilon e^{-1} e^{-\frac{1}{2\varepsilon}}$, and using $w(t_{2,3}) \sim \varepsilon e^{-\varepsilon(t_{2,3}-t_2)}$ we get

$$t_{2,3} - t_2 \sim \frac{1}{2\varepsilon^2}.$$
We can now compute $(t_3 - t_{2,3})$ using a same argument, due to the symmetry of the equation. Notice that the energy formula (11) as well as $v(t_3) = \varepsilon$ imply $w(t_3) \sim \varepsilon e^{-1} \exp(-\frac{1}{\varepsilon})$. In the interval $t \in [t_{2,3}, t_3]$ we can approximate (14) as

$$ \frac{dv}{dt} = \varepsilon v, \quad \frac{dw}{dt} = vw - \varepsilon w, $$

because $w \ll \varepsilon$ for $t \in [t_{2,3}, t_3]$. We then obtain the following approximation $(t_3 - t_{2,3}) \sim \frac{1}{2\varepsilon^2}$ for small $\varepsilon$. Then

$$ t_3 \sim t_2 + \frac{1}{\varepsilon^2} \quad \text{as } \varepsilon \to 0. $$

**Proof of lemma 2.5: From $t_3$ to $t_4$.**

For $t > t_3$ both $v$ and $w$ increase. This stage can be analyzed in a manner symmetric to the one from $t_1$ to $t_2$. During this phase $v$ becomes of order one. This happens in times of order $\frac{1}{\varepsilon}$ (up to a logarithmic correction). We can define $t_4$ by means of $w(t_4) = \varepsilon$. During this phase we can do the same computations as in $(t_1, t_2)$ so that we obtain:

$$ t_4 - t_3 \sim \frac{1}{\varepsilon} \log \left( \frac{1}{\varepsilon} \right). $$

**Proof of lemma 2.6: From $t_4$ to $t_5$.**

Finally, the trajectory returns to the line $\{w = v\}$. Using the conservation of energy formula (11) we can approximate the evolution of $(w, v)$ by means of the equation

$$ \frac{dw}{dt} = w(1-w), \quad v = 1-w $$

We have defined $t_5$ by means of $w(t_5) = v(t_5) > \varepsilon$. Using that $w(t_4) = \varepsilon$ we obtain the following approximation, symmetric to the interval $(0, t_1)$:

$$ t_5 - t_4 \sim \log \left( \frac{1}{\varepsilon} \right). $$

**Proof of lemma 2.7.** Using (17) and (18) we obtain

$$ \begin{aligned} D(E, \varepsilon) &= \frac{E}{2} \int_{0}^{T(E, \varepsilon)} (w+v) \left(Es; 1, \frac{\varepsilon}{E}\right) ds = \frac{1}{2} \int_{0}^{ET(E, \varepsilon)} (w+v) \left(s; 1, \frac{\varepsilon}{E}\right) ds \\ &= \frac{1}{2} \int_{0}^{T(1, \frac{\varepsilon}{E})} (w+v) \left(s; 1, \frac{\varepsilon}{E}\right) ds = D\left(1, \frac{\varepsilon}{E}\right). \end{aligned} $$

**Proof of lemma 2.8**

The equivalent $T(1, \varepsilon) \sim \frac{1}{\varepsilon^2}$ comes directly from lemmas 2.2 to 2.6. It is known that for the LV system, the mean value of $v$ and $w$ satisfies the equality

$$ \frac{1}{\varepsilon} \int_{0}^{T} v(s) ds = \frac{1}{\varepsilon} \int_{0}^{T} w(s) ds = \varepsilon, $$
hence their equivalent and the one for $D(1, \varepsilon)$. This proves (35). The maximum value of $Y$ is given by

$$Y_{\max}(t) = \int_{0}^{\frac{T}{2}} (w(s) - v(s)) ds = \int_{0}^{t_{2,3}} (w(s) - v(s)) ds$$

and with lemmas 2.2 to 2.6 we get

$$ 
\begin{aligned}
Y(t_1) &= \int_0^{t_1} (w(s) - v(s)) ds & \sim \int_0^{t_1} ds & \sim -\log(\varepsilon), \\
\int_{t_1}^{t_{1,2}} (w(s) - v(s)) ds & \sim t_{1,2} - t_1 \sim -\frac{1}{2} \log(\varepsilon), \\
\int_{t_{1,2}}^{t_{2,3}} (w(s) - v(s)) ds & \sim \int_{t_{1,2}}^{t_2} e^{-\varepsilon(s-t_{1,2})} ds \sim \frac{1}{\varepsilon}, \\
\int_{t_2}^{t_{2,3}} (w(s) - v(s)) ds & \le \int_{t_2}^{t_{2,3}} w(s) ds \sim \int_{t_2}^{t_{2,3}} \varepsilon e^{-\varepsilon(t-t_2)} ds \le 1.
\end{aligned}
$$

Finally, we see that the maximal contribution to the displacement of $Y$ is linked to the part $(t_{1,2}, t_2)$, and we have $Y_{\max} \sim \frac{1}{\varepsilon}$.

## Appendix B. Linearised stability around the positive steady state

In this appendix we study the linear stability of the steady states (10).

**Proposition B.1 (Linear stability around the positive steady state).** For $M=1$ and $\varepsilon \ll 1$, the unique positive steady state defined by (10) is locally asymptotically stable.

**Proof.** We provide here the main lines for the linear stability analysis, letting a fully rigorous and complete study for future work. We introduce new variables by means of

$$c_j = \bar{c}_j(1 + \varphi_j), \quad j \ge 1, \quad v = \bar{v}(1 + V), \quad w = \bar{w}(1 + W), \quad (126)$$

and recall that the positive steady state $(\bar{v}, \bar{w}, \bar{c}_j)$ is defined by (10), namely

$$ 
\begin{aligned}
\theta &:= 1 - \frac{1}{2} \left( \frac{1}{\bar{\varepsilon}} - \sqrt{\left(\frac{1}{\bar{\varepsilon}} - 2\right)^2 + 4} \right) \sim \varepsilon, \\
\bar{c}_1 &:= \varepsilon\theta \sim \varepsilon^2, \quad \bar{c}_i := (1 - \theta)^{i-1}\bar{c}_1, \quad \bar{v} := \varepsilon, \quad \bar{w} := \varepsilon - \bar{c}_1 = \varepsilon(1 - \theta).
\end{aligned}
$$

Assuming $|\varphi_j| + |V| + |W| \ll 1$, neglecting higher order terms, we get the following linearised system:

$$ \frac{dV}{dt} = -\varepsilon(1-\theta)W - \varepsilon\theta\varphi_1, \quad (127) $$

$$ \frac{dW}{dt} = \varepsilon V, \quad (128) $$

$$ \frac{d\varphi_1}{dt} = -\varepsilon(1-\theta)(W+\varphi_1-V-\varphi_2), \quad (129) $$

$$ \frac{d\varphi_j}{dt} = \varepsilon(\varphi_{j-1} - \varphi_j + W - V) - \varepsilon(1-\theta)(\varphi_j - \varphi_{j+1} + W - V). \quad (130) $$

Introducing the time variable $\tau = \varepsilon t$ we get

$$ \frac{dV}{d\tau} = -(1-\theta)W - \theta\varphi_1, \quad (131) $$
$$ \frac{dW}{d\tau} = V, \qquad (132) $$

$$ \frac{d\varphi_1}{d\tau} = -(1-\theta)(W+\varphi_1-V-\varphi_2), \qquad (133) $$

$$ \frac{d\varphi_j}{d\tau} = (\varphi_{j-1} - \varphi_j + W - V) - (1-\theta)(\varphi_j - \varphi_{j+1} + W - V). \qquad (134) $$

We look for solutions to the eigenproblem, hence of the form $(V(t), W(t), \varphi_j(t)) = e^{\lambda\tau}(V, W, \varphi_j)$ with $\lambda \in \mathbb{C}$. We get the system:

$$ \lambda V = -(1-\theta)W - \theta\varphi_1, \qquad (135) $$

$$ \lambda W = V, \qquad (136) $$

$$ \lambda\varphi_1 = -(1-\theta)(W + \varphi_1 - V - \varphi_2), \qquad (137) $$

$$ \lambda\varphi_j = (\varphi_{j-1} - \varphi_j + W - V) - (1-\theta)(\varphi_j - \varphi_{j+1} + W - V). \qquad (138) $$

We recall that $\theta \sim \varepsilon$ and, thus, $\theta \ll 1$ as soon as $\varepsilon \ll 1$, hence we treat the solutions to (135)–(138) in a perturbative manner with respect to $\theta$.

At the limit $\theta = 0$, we get

$$ \lambda V = -W, \qquad (139) $$

$$ \lambda W = V, \qquad (140) $$

$$ \lambda\varphi_1 = -W - \varphi_1 + V + \varphi_2, \qquad (141) $$

$$ \lambda\varphi_j = \varphi_{j-1} - 2\varphi_j + \varphi_{j+1}. \qquad (142) $$

We notice that we have two eigenvalues $\lambda_0 = \pm i$, with $V = \pm iW$, linked to an oscillatory regime of period $2\pi$ (coherent with the period found at the end of Phase III), and, in the case $V = W = 0$, a continuous spectrum similar to the one of the discrete heat equation, whose admissible generalised eigenvalues all have a negative real part except $\lambda = 0$ (associated to the constant sequence $\varphi_j = \varphi_1$ for all $j$).

**Case (a):** $\lambda = \pm i$, $V = \pm iW$. We only detail the computations for $\lambda = i$ since the case $\lambda = -i$ follows from taking the conjugate. For $W_0 = -i$, $V_0 = 1$ we obtain

$$ \varphi_2 = (i+1)\varphi_1 - (1+i), \\ \varphi_{j+1} = (i+2)\varphi_j - \varphi_{j-1}, \quad j \ge 2. $$

We recognize a linear recurrent sequence of order 2 for $j \ge 2$, whose characteristic equation is

$$ r^2 - (i+2)r + 1 = 0 \implies r_{\pm} = \frac{i+2 \pm \sqrt{(i+2)^2 - 4}}{2} = 1 \pm \frac{i}{2} \pm \sqrt{i - \frac{1}{4}}. $$

Since $|r_{-}| < 1$ and $|r_{+}| > 1$, the only admissible value is $r_{-}$, that we denote $r = r_{-}$ for simplicity. We finally have

$$ V_0 = 1, \quad W_0 = -i, \quad \varphi_{j,0} = \frac{1+i}{1-r} \mu^j, \quad r = 1 + \frac{i}{2} - \sqrt{i - \frac{1}{4}}. \qquad (143) $$
**Case (b):** $W = V = 0$. Then (141) implies $(\lambda + 1)\varphi_1 = \varphi_2$, and (142) implies for $j \ge 2$

$$ 
\begin{aligned}
\lambda\varphi_j &= \varphi_{j-1} - 2\varphi_j + \varphi_{j+1}, & j \ge 2, \\
\lambda\varphi_1 &= -(\varphi_1 - \varphi_2).
\end{aligned}
$$

We recognize the discrete heat equation. As previously, the characteristic equation is

$$ r^2 - (\lambda + 2)r + 1 = 0 \implies r_{\pm} = \frac{\lambda + 2 \pm \sqrt{(\lambda + 2)^2 - 4}}{2} = \frac{\lambda + 2 \pm \sqrt{\lambda(\lambda + 4)}}{2} $$

We see that for $\lambda = 0$ we have $r_{\pm} = 1$ linked to the constant sequence $\varphi_j = 1$.

If $\lambda \ne 0$, generalised eigenvectors are given by linear combinations of the two sequences $(r_{\pm}^{j-1})_{j \ge 1}$. Using the boundary condition at $j = 1$ and $r_{+} + r_{-} = \lambda + 2$, we obtain

$$ \varphi_j = r_{+}^{j-1} - \frac{1 - r_{-}}{1 - r_{+}} r_{-}^{j-1}, $$

so that to ensure $(\varphi_j) \in \ell^{\infty}(\mathbb{C})$ we need $|r_{\pm}| \le 1$. Since $r_{+}r_{-} = 1$, we define $r_{\pm} = e^{\pm i\beta}$ with $\beta \in [0, 2\pi)$, and compute

$$ \lambda = e^{i\beta} + e^{-i\beta} - 2 = 2(\cos(\beta) - 1) \in [-4, 0]. $$

To conclude: We have found that no admissible (generalised) eigenvalue has a positive real part, and the only eigenvalues with a nonnegative real part are $\lambda = 0$, $\lambda = i$ and $\lambda = -i$.

In the case $\theta \ll 1$, we consider a perturbation of the eigenvalue $\lambda_0 = i$ (the case $\lambda_0 = -i$ is similar) and its eigenvector given by (143) by writing the asymptotics

$$ \lambda = \lambda_0 + \theta\lambda_1, \quad V = V_0 + \theta V_1, \quad W = W_0 + \theta W_1, \quad \varphi_j = \varphi_{j,0} + \theta\varphi_{j,1}, $$

then keeping only the first order terms in $\theta$, the system (135)–(138) becomes

$$ 
\begin{aligned}
\lambda_0 V_1 + \lambda_1 V_0 &= -W_1 + W_0 - \varphi_{1,0}, & \implies iV_1 + \lambda_1 &= -W_1 - i - \varphi_{1,0} \\
\lambda_0 W_1 + \lambda_1 W_0 &= V_1, & \implies iW_1 - i\lambda_1 &= V_1 \\
\implies -W_1 + 2\lambda_1 &= -W_1 - i - \varphi_{1,0} & \implies \lambda_1 &= -\frac{1}{2} \left( i + \frac{1+i}{1-r} r \right),
\end{aligned}
$$

and we check numerically that $\mathrm{Re}(\lambda_1) < 0$, leading to damping oscillations.

We now consider the stability of the continuous spectrum, when $\lambda_0 = 2(\cos(\beta) - 1)$. Detailed computations are necessary for a fully rigorous justification, and are left for future work. In brief, we notice that as soon as $\theta \lesssim \beta$, we have $\mathrm{Re}(\lambda) = \lambda_0 + \theta\mathrm{Re}(\lambda_1) < 0$, so this part of the spectrum remains stable. The more delicate part is thus when $\beta \lesssim \theta \ll 1$. We can carry out an asymptotic expansion as before, and conclude that the continuous spectrum remains in the half-plane with nonpositive real parts.

To conclude, we note that the oscillations have a period of order $2\pi$ in the variable $\tau$ and the damping occurs in a time of order $1/\theta \sim 1/\epsilon$, which corresponds to a period of order $2\pi/\epsilon$ and a damping as $1/\epsilon^2$ in the original variable $t$: this is in line with the analysis carried out for the end of Phase III. $\square$
# Appendix C. Numerical computations

In this appendix, we give details on the numerical simulation used to illustrate the evolution of the system (2)–(4). The code is written in MATLAB and is published online (https://github.com/mmezache/BiMono_Becker_Doring). The numerical study has an illustrative purpose and is divided based on the framework considered for the model, either the size-continuous framework or the size-discrete framework.

## C.1. The continuous model

We recall the approximation of the model by a diffusion-convection equations (27)–(28) describing the evolution of clusters during Phases I and II:

$$ \partial_t c(j,t) + \tilde{V}(t) \partial_j c(j,t) = \frac{d(t)}{2} \partial_j^2 c(j,t), \quad j \in \Omega \subset \mathbb{R}_+, t \in [0, T) $$

where

$$ \tilde{V}(t) = w(t) - v(t), \quad d(t) = w(t) + v(t) $$

and where $v$ and $w$ are the solution of the LV system (14):

$$ \left\{ \begin{array}{l} \frac{d}{dt}v(t) = v(t)(\epsilon - w(t)), \quad v(0) = 0.6, \\ \\ \frac{d}{dt}w(t) = w(t)(v(t) - \epsilon), \quad w(0) = 0.6, \end{array} \right. $$

where $\epsilon = \int_{\Omega} c(j,t) dj$.

## Computation of the monomers

The system (14) depends on the parameters $\epsilon \ll 1$ which makes it stiff when one attempts to solve it numerically. The main reason for the instability of the numerical schemes is that the trajectories of $v$ and $w$, which form the following closed-curve

$$ E = v(t) + w(t) - 2\epsilon - \epsilon \log \left( \frac{v(t)w(t)}{\epsilon^2} \right), $$

pass through strictly positive values but very close to 0. The issue is then to find the balance between the accuracy of the numerical schemes and preserving the positivity properties (see for instance [6]). To illustrate this phenomenon, let $\epsilon = 10^{-2}$ and we define the function

$$ g : x \in \mathbb{R}_+^* \mapsto 2(x - \epsilon - \epsilon(\log(x) - \log(\epsilon))) - 1. $$

The function $g$ corresponds to the Hamiltonian $E$ (11) minus 1. Then, the root $x^* \in (0, \epsilon)$ such that $g(x^*) = 0$ gives insights on the minimal values taken by $v$ and $w$. We have

$$ \begin{aligned} g(x^*) &= 0 \\ -2\epsilon \log(x^*) &= 1 - 2\epsilon \log(\epsilon) + 2\epsilon - 2x^* \\ -2\epsilon \log(x^*) &\ge 1 - 2\epsilon \log(\epsilon) \\ x^* &\le \exp\left(\frac{1-2\epsilon \log(\epsilon)}{-2\epsilon}\right) \le 2 \times 10^{-24}. \end{aligned} $$
Since the solution of system (14) goes very close to 0, the numerical scheme needs to be capable of preserving the positivity. In order to ensure the positivity of the numerical solutions of system (14), the structure of the LV system that forces the solutions to remain positive can be used. We apply the following change of variable $x = \log(v)$ and $y = \log(w)$. Hence, we have

$$ \left\{ \begin{array}{l} \frac{d}{dt}x(t) = \epsilon - e^{y(t)}, \quad x(0) = \log(0.6), \\ \\ \frac{d}{dt}y(t) = e^{x(t)} - \epsilon, \quad y(0) = \log(0.6). \end{array} \right. \qquad (144) $$

The numerical solution of (144) is then computed with the high accuracy Runge-Kutta scheme at the 8th order.

### Numerical scheme for the diffusion-convection equation

We propose to use an implicit scheme to solve the diffusion-convection equation (27) with a centred difference operator for the convection term. Consider a constant step time discretisation of the interval $[0, T]$ $t_n = n\Delta t, n = 0, ..., N_1$, where $T = t_{N_1}$ and $\Delta t$ is the time step. Consider a constant step 'space' discretisation of the interval of sizes $[0, L]$ $j_k = k\Delta j, k = 0, ..., N_2$, where $L = j_{N_2}$ and $\Delta j$ is the space step. We look for an approximation $c_k^n$ of $c(k\Delta j, n\Delta t)$. The scheme is then the following:

$$ \frac{c_k^{n+1} - c_k^n}{\Delta t} + \tilde{V}^{n+1} \frac{c_{k+1}^{n+1} - c_{k-1}^{n+1}}{2\Delta j} - \frac{d^{n+1}}{2} \frac{c_{k+1}^{n+1} - 2c_k^{n+1} + c_{k-1}^{n+1}}{(\Delta j)^2} = 0, \quad (145) $$

for $n = 0, ..., N_1 - 1$ and $k = 1, ..., N_2 - 1$. The approximation of the integral terms are given by

$$ \epsilon^n = \Delta j \sum_{k=1}^{N_2} c_k^n, \qquad (146) $$

$$ M^n = \Delta j \sum_{k=1}^{N_2} j_k c_k^n. \qquad (147) $$

### Conservation of the total concentration of clusters

**Proposition C.1.** The numerical scheme defined in (145) is conservative for the quantity $\epsilon$ if the following equation holds for the boundary conditions:

$$ \begin{aligned} 0 = c_0^n & - (a^{n+1} + b^{n+1} + 1) c_0^{n+1} - (a^{n+1} - b^{n+1}) c_1^{n+1} + c_{N_2}^n \\ & + (a^{n+1} - b^{n+1} - 1) c_{N_2}^{n+1} + (a^{n+1} + b^{n+1}) c_{N_2-1}^{n+1}, \end{aligned} \qquad (148) $$

where

$$ a^{n+1} = \frac{\Delta t \tilde{V}^{n+1}}{2\Delta j} \quad \text{and} \quad b^{n+1} = \frac{\Delta t d^{n+1}}{2(\Delta j)^2}. $$
**Proof.** For the sake of clarity, we use the following notations $a^{n+1} = a$ and $b^{n+1} = b$. The approximation of the integral term is defined by (146), dividing by the space step, we have

$$
\begin{align*}
\sum_{k=0}^{N_2} c_k^n &= c_0^n + \sum_{k=1}^{N_2-1} c_k^n + c_{N_2}^n, \\
&= c_0^n + \sum_{k=1}^{N_2-1} \left( -(a+b) c_{k-1}^{n+1} + (1+2b) c_k^{n+1} + (a-b) c_{k+1}^{n+1} \right) + c_{N_2}^n, \\
&= c_0^n + c_{N_2}^{n+1} + \sum_{k=2}^{N_2-2} c_k^{n+1} - (a+b) (c_0^{n+1} + c_1^{n+1}) \\
&\quad + (1+2b) (c_1^{n+1} + c_{N_2-1}^{n+1}) + (a-b) (c_{N_2-1}^{n+1} + c_{N_2}^{n+1}), \\
&= c_0^n - (a+b+1) c_0^{n+1} - (a-b) c_1^{n+1} \\
&\quad + c_{N_2}^n + (a-b-1) c_{N_2}^{n+1} + (a+b) c_{N_2-1}^{n+1} + \sum_{k=0}^{N_2} c_k^{n+1}.
\end{align*}
$$

Hence, if the equation (148) holds then the total cluster concentration defined by (146) is the same at each time iteration. □

As a direct result of proposition C.1, the boundary conditions can be chosen as the following in order to keep the tridiagonal structure of the matrix for the numerical scheme (145) :

$$
\left\{
\begin{aligned}
c_0^n &= (1 + 2a^{n+1}) c_0^{n+1} + \Delta j (a^{n+1} - b^{n+1}) \frac{c_1^{n+1} - c_0^{n+1}}{\Delta j}, \\
c_{N_2}^n &= (1 - 2a^{n+1}) c_{N_2}^{n+1} + \Delta j (a^{n+1} + b^{n+1}) \frac{c_{N_2}^{n+1} - c_{N_2-1}^{n+1}}{\Delta j}.
\end{aligned}
\right.
$$

**Discussion about the total mass of the system (14)–(27)**

The system composed by equations (14) and (27) is a continuous approximation of the behaviour of the more complex polymerisation/depolymerisation chemical reactions with two monomers and the catalytic depolymerisation. As such, the system does not ensure the conservation of the total mass. For instance, let $\Omega = \mathbb{R}_+$, assume sufficient regularity conditions for $(c, v, w)$ and the correct boundary conditions such that $\epsilon$ remains constant. Then

$$
\begin{align*}
\frac{d}{dt} \left( v(t) + w(t) + \int_{\Omega} j c(j,t) dj \right) &= -\tilde{V}(t) \epsilon + \tilde{V}(t) \int_{\Omega} c(j,t) dj - \frac{d(t)}{2} \int_{\Omega} \partial_j c(j,t) dj, \\
&= \frac{d(t)}{2} c(0,t).
\end{align*}
$$

This implies that the solution of (27) has a homogeneous Dirichlet boundary conditions. Moreover, the conservation of the total concentration of clusters imposes a homogeneous Robin boundary conditions. Taking into account the two boundary conditions, the solution is locally constant and equal to 0 in a neighbourhood of $\partial\Omega = 0$. Hence, the numerical scheme proposed above is conservative only for the total concentration of clusters. The correction imposed by the mass conservation law in the equations for the monomers is neglected because of its small order of magnitude in the continuous approximation. The total mass increases at each time step.

Using the previous notations, consider $N_2 = \infty$, for the sake of simplicity. The discretization scheme (145), with the choice of numerical quadratures (146) and (147) imply that the total
mass is increasing by the following amount at each iteration

$$ \frac{d^{n+1} - \Delta j \tilde{V}^n + 1}{2} c_0^{n+1}. $$

One thing to note is that the purpose of the discretization of the continuous approximation is to illustrate the transient dynamics of the clusters over one cycle of the LV oscillations. Hence, the increase of the total mass does not significantly affect the interpretation that can be made of the graphs.

## Numerical simulation of the advection-diffusion equation

The numerical simulations of (27) and (14) are illustrated in figures 3 and 4. The parameters are $T = 10^4$, $\Delta t = 0.05$, $L = 250$ and $\Delta j = 0.5$. The initial conditions are $v(0) = w(0) = 0.6$ and

$$ c_0(x) = e^{-\frac{x^2}{2\sigma^2}} \mu \sqrt{\frac{2}{\pi\sigma}} $$

where $\sigma = 10$ and $\mu = 0.02$. Hence, $\epsilon^0 \approx 0.0212$ and $v^0 + w^0 + M^0 \approx 1.2503$.

### C.2. The discrete model

The numerical computation of the discrete model is using the 8th order Runge Kutta scheme. In order to ensure the nonnegativity of the concentrations of monomers, we apply the change of variables

$$ v \to e^v \quad \text{and} \quad w \to e^w. $$

The maximal size of the clusters is set to 500. The code can be found on the following link: https://github.com/mmezache/BiMono_Becker_Doring.

## References

[1] Armiento A, Moireau P, Martin D, Lepejova N'a, Doumic M and Rezaei H 2017 The mechanism of monomer transfer between two structurally distinct PrP oligomers *PLoS One* 12 e0180538

[2] Ball J M, Carr J and Penrose O 1986 The Becker-Döring cluster equations: basic properties and asymptotic behaviour of solutions *Commun. Math. Phys.* 104 657-92

[3] Beal D M, Tournus M, Marchante R, Purton T J, Smith D P, Tuite M F, Doumic M and Xue W-F 2020 The division of amyloid fibrils: systematic comparison of fibril fragmentation stability by linking theory with experiments *Science* 23 101512

[4] Becker R and Döring W 1935 Kinetische Behandlung der Keimbildung in übersättigten Dämpfen *Ann. Phys., Lpz.* 416 719-52

[5] Bender C M and Orszag S A 2013 *Advanced Mathematical Methods for Scientists and Engineers I: Asymptotic Methods and Perturbation Theory* (Springer)

[6] Blanes S, Iserles A and Macnamara S 2022 Positivity-preserving methods for ordinary differential equations *ESAIM: Math. Modelling Numer. Anal.* 56 1843-70

[7] Budzinskiy S S, Matveev S A and Krapivsky P L 2021 Hopf bifurcation in addition-shattering kinetics *Phys. Rev. E* 103 L040101

[8] Calvez V, Hawkins R J, Meunier N and Voituriez R 2012 Analysis of a nonlocal model for spontaneous cell polarization *SIAM J. Appl. Math.* 72 594-622

[9] Calvez V, Meunier N and Voituriez R 2010 A one-dimensional keller-segel equation with a drift issued from the boundary *C. R. Math.* 348 629-34

[10] Doumic M, Fellner K, Mezache M and Rezaei H 2019 A bi-monomeric, nonlinear Becker-Döring-type system to capture oscillatory aggregation kinetics in prion dynamics *J. Theor. Biol.* 480 241-61
[11] Doumic M, Fellner K, Mezache M, and Velázquez J J L 2024 Analysis of some discrete and non-local PDE models arising in the study of a model of prion polymerization accepted

[12] Fornara B, Igel A, Beringue V, Martin D, Sibille P, Pujo-Menjouet L and Rezaei H 2024 The dynamics of prion spreading is governed by the interplay between the non-linearities of tissue response and replication kinetics *iScience* **27** 111381

[13] Fortin J-Y and Choi M 2023 Stability condition of the steady oscillations in aggregation models with shattering process and self-fragmentation *J. Phys. A: Math. Theor.* **56** 385004

[14] Guckenheimer J and Holmes P 2013 *Nonlinear Oscillations, Dynamical Systems and Bifurcations of Vector Fields* vol 42 (Springer)

[15] Hawkins R J, Benichou O, Piel M and Voituriez R 2009 Rebuilding cytoskeleton roads: Active-transport-induced polarization of cells *Phys. Rev. E* **80** 040903

[16] Honoré S, Hubert F, Tournus M and White D 2019 A growth-fragmentation approach for modeling microtubule dynamic instability *Bull. Math. Biol.* **81** 722–58

[17] Jabin P-E and Niethammer B 2003 On the rate of convergence to equilibrium in the Becker-Döring equations *J. Differ. Equ.* **191** 51–543

[18] Kevorkian J and Cole J D 2013 *Perturbation Methods in Applied Mathematics* vol 34 (Springer)

[19] Niethammer B, Pego R L, Schlichting A and Velázquez J J L 2022 Oscillations in a Becker-Döring model with injection and depletion *SIAM J. Appl. Math.* **82** 1194–219

[20] Pego R L and Velázquez J J L 2020 Temporal oscillations in Becker-Döring equations with atomization *Nonlinearity* **33** 1812

[21] Prigent S, Ballesta A, Charles F, Lenuzza N, Gabriel P, Matar Tine L, Rezaei H and Doumic M 2012 An efficient kinetic model for assemblies of amyloid fibrils and its application to polyglutamine aggregation *PLoS One* **7** e43273

[22] Torrent J, Martin D, Noinville S, Yin Y, Doumic M, Moudjou M, Beringue V and Rezaei H 2019 Pressure reveals unique conformational features in prion protein fibril diversity *Sci. Rep.* **9** 1–11

[23] Vasseur A, Poupaud F, Collet J-F and Goudon T 2002 The Becker-Döring system and its Lifshitz-Slyozov limit *SIAM J. Appl. Math.* **62** 1488–500

[24] Velázquez J J L 1998 The Becker-Döring equations and the Lifshitz-Slyozov theory of coarsening *J. Stat. Phys.* **92** 195–236

[25] Xue W-F, Homans S W and Radford S E 2008 Systematic analysis of nucleation-dependent polymerization reveals new insights into the mechanism of amyloid self-assembly *PNAS* **105** 8926–31

[26] Xue W-F and Radford S E 2013 An imaging and systems modeling approach to fibril breakage enables prediction of amyloid behavior *Biophys. J.* **105** 2811–9