text:
# A simple asymmetric evolving random network

## 1 Motivations and results.

Evolving random graphs have recently attracted attention, see e.g. refs  and references therein. This interest is mainly motivated by concrete problems related to the structure of communication or biological networks. Experimental data are now available in many contexts .

In these examples, the asymmetry and the evolving nature of the networks are likely to be important ingredients for deciphering their statistical properties. It is however far from obvious to find solvable cases that would possibly account for some relevant features of, say, the regulating network of a genome. Although biology has strongly influenced our interest in evolving networks, the model we solve is not based on realistic biological facts but it nevertheless incorporates asymmetry and chronological order. Understanding such simple evolving graphs may help understanding biological networks, at least by comparison and opposition.

We were initially motivated by the study of the yeast genetic regulatory network presented in ref.. The authors studied in and out degree distributions and discovered a strong asymmetry: a single gene may participate to the regulation of many other genes – the law for out-degrees seems to be large –, but each genes is only regulated by a few other genes – the law for in-degrees seems to have finite moments. This is why we consider oriented evolving random graphs in the sequel. A biological interpretation for the asymmetry is that the few promoter-repressor sites for each gene bind only to specific proteins, but that along the genome many promoter-repressor sites are homologous. However, this does not predict the precise laws. An understanding of the same features from a purely probabilistic viewpoint would be desirable as well.

The recent experimental studies dealt with global statistical properties of evolving graphs, i.e. when the evolving network is observed at some fixed time with the ages of different vertices and edges not taken into account. There are simple experimental reasons for that : to keep track of the ages would in many cases dramatically reduce the statistics, and in other cases this information is even not available. Our second motivation is a better understanding of the local-in-time statistical properties of evolving networks. This helps dating or assigning likely ages to different structures of the networks. As we shall later see, the global analysis, which is like a time average, gives a distorted view of the real structure of the networks. We shall present a detailed analysis of local-in-time features in our model.

The model we study is the natural evolving cousin of the famous Erdös-Renyi random graphs . Starting from a single vertex at time $`1`$, a new vertex is created at each time step – so that at time $`t`$, the size of the system, i.e. the number of vertices, is $`t`$ –, and new oriented edges are created with specified probabilistic rules. A tunable parameter $`\alpha `$ ranging from $`0`$ to $`\mathrm{}`$ describes asymptotically the average number of incoming edges on a vertex. Precise definitions are given in the next section.

Our main results are the following :

From very simple rules, we see an asymmetry emerging. The global in and out degree distributions are different. We also compute the local profiles of in and out degree distributions, and comment on the differences.

We make a detailed global analysis for the structure and sizes of the connected components. We use generating function methods to write down a differential equation that implies recursion relations for the distribution of component sizes, see eqs.(11,14).

A salient global feature of the model is a percolation phase transition at a critical value of the average connectivity. Below this value, no single component contains a finite fraction of the sites in the thermodynamic limit, i.e. in the large $`t`$ limit. However, a slightly unusual situation occurs in that below the transition the system contains components whose sizes scale like a power of the total size of the graph, see eq.(26). Correspondingly, the probability distribution for component sizes has an algebraic queue, see eq.(23), and its number of finite moments jumps at specific values of the average connectivity. Above the transition, this probability distribution becomes defective, but its decrease is exponential, see eq.(43). The transition is continuous. Close to the threshold, the fraction of sites in the giant component – the percolation cluster – has an essential singularity, see eq.(41). We argue that this result is universal, with the meaning used in the study of critical phenomena. The essential singularity at the percolation threshold had already been observed numerically by  in a different model which we show to be in the same universality class as ours for the percolation transition, and computed analytically for another class of models in .

We then turn to the study of local-in-time profiles of connected components. Guided by a direct enumeration based on tree combinatorics, we show that they satisfy recursion relations, and we give the first few profiles (isolated vertices, pairs, triples) explicitly. The profile of the giant component is given by a differential equation, from which we extract the singularity in the far past and the critical singularity in the present – see eqs(54,55). In particular the giant component invades all the time slices of the graph above the transition. One strange feature of profiles, which would deserve a good explanation, is that in several instances the formal parameter involved in generating functions for global quantities is simply traded for the relative age to obtain interesting local-in-time observables, see eqs.(52,56).

We have compared our analytical results with numerical simulation whenever possible.

While polishing this paper, we became aware of , whose goals overlap partly with ours. When they can be compared, the results agree.

## 2 The model.

We construct evolving random graphs with the following rules:
(i) We consider a triangular array of independent random variables $`\mathrm{}_{i,j}`$, $`1i<j`$, where $`\mathrm{}_{i,j}`$ takes value $`1`$ with probability $`p_j[0,1]`$ and value $`0`$ with probability $`q_j1p_j`$.
(ii) We start from the graph made of single vertex at initial time $`t=1`$. At time $`t`$, $`t2`$, a vertex with label $`t`$ is added together with the directed edges $`[jt]`$ for which $`\mathrm{}_{j,t}=1`$. We shall often take the viewpoint that the (biased) coin tossings defining $`\mathrm{}_{i,t}`$ are done at time $`t`$.

We shall assume that $`p_t\alpha /t`$ at large time $`t`$, with $`\alpha `$ a parameter which we shall identify as half the average connectivity. This choice ensures the convergence of various distributions to stationary measures, most of them being independent of the precise values of the early probabilities.

By construction all edges arriving at a given vertex are simultaneously created at the instant of creation of this vertex. As a consequence, these graphs are not only oriented but chronologically oriented – this is unrealistic from the biological viewpoint.

## 3 Edge distributions.

In this section we give the incoming and outgoing edge distributions.
Let $`\mathrm{}_j^{}(t)`$ be the number of incoming edges at the vertex $`j`$, and $`\mathrm{}_j^+(t)`$ be the number of outgoing edges at this vertex at time $`t`$.
Let $`v_k^{}(t)`$ be the number of vertices with $`k`$ incoming edges, and $`v_k^+(t)`$ be the number of vertices with $`k`$ outgoing edges at time $`t`$.

We may look either at the edge distributions at a given vertex, or we may look at the edge distributions defined by gathering averaged histograms over whole graphs. The former are specified by their generating functions,

$$z^{\mathrm{}_j^\pm (t)}$$

where $`\mathrm{}`$ denotes expectation value. It may depend on the specified vertex labeled by $`j`$. The latter is defined by the generating functions,

$`𝒱_t^\pm (z){\displaystyle \frac{1}{t}}{\displaystyle \underset{0kt}{}}v_k^\pm (t)z^k={\displaystyle \frac{1}{t}}{\displaystyle \underset{1jt}{}}z^{\mathrm{}_j^\pm (t)}`$ (1)
We remark that this global histogram distribution is the average of the local-in-time quantity $`z^{\mathrm{}_j^\pm (t)}`$. Since at time $`t`$ the total number of vertices is $`t`$, $`𝒱_t^\pm (z)`$ is properly normalized, $`𝒱_t^\pm (1)=1`$, to define an averaged probability distribution function, independent of the vertices, for the incoming or outgoing edge variables $`\mathrm{}^\pm `$:

$$\underset{k}{}z^k\mathrm{Prob}_{(\mathrm{}^\pm =k)}𝒱_t^\pm (z)$$

Incoming vertices. The number of incoming edges $`\mathrm{}_j^{}(t)=_i\mathrm{}_{i,j}`$ at vertex $`jt`$ asymptotically possesses a Poisson distribution since

$$z^{\mathrm{}_j^{}(t)}=\underset{i<j}{}[q_j+zp_j]_j\mathrm{}\mathrm{exp}(\alpha (z1))$$

The convergence of this distribution justifies our choice of asymptotic probabilities $`p_j\alpha /j`$. Only the vertices whose ages $`j`$ scale with the age of the graph, i.e. with $`j/t=\sigma `$ fixed, $`0\sigma 1`$, give non trivial contributions at large time to the averaged histogram (1) and

$`𝒱_t^{}(z)_t\mathrm{}\mathrm{exp}(\alpha (z1))`$ (2)
This expression may also be retrieved by looking at the evolution equation of $`𝒱_t^{}(z)`$. Indeed, consider adding the new vertex at time $`t`$. Since the edges are oriented from older to younger vertices 111A vertex is older than another if it appeared before in the evolution, i.e. if it corresponds to a smaller value of $`\sigma `$., we have $`t𝒱_t^{}(z)=(t1)𝒱_{t1}^{}(z)+z^{\mathrm{}_t^{}(t)}`$ from the second definition in eq.(1). This is equivalent to

$`t𝒱_t^{}(z)=(t1)𝒱_{t1}^{}(z)+(q_t+zp_t)^{t1}`$
As $`(q_t+zp_t)^{t1}e^{\alpha (z1)}`$ at large time, the stationary limit is given by eq.(2). This yields a Poissonian distribution with probabilities

$`\mathrm{Prob}_{(\mathrm{}^{}=k)}=e^\alpha {\displaystyle \frac{\alpha ^k}{k!}}`$ (3)
Outgoing vertices. At a given vertex $`jt`$, with $`j/t=\sigma `$ fixed, the number of outgoing edges $`\mathrm{}_j^+(t)=_{it}\mathrm{}_{j,i}`$ at vertex $`j`$ also have a Poisson distribution at large time,

$$z^{\mathrm{}_{\sigma t}^+(t)}=\underset{\sigma t=j<it}{}[q_i+zp_i]_t\mathrm{}\mathrm{exp}(\alpha \mathrm{log}\sigma (z1))$$

but with a parameter $`\alpha \mathrm{log}(1/\sigma )`$ depending on the age of the vertex. Approximating at large time the sum over $`j`$ in eq.(1) by an integral over $`\sigma `$ gives the histogram distribution:

$`𝒱_t^+(z)_t\mathrm{}{\displaystyle _0^1}𝑑\sigma z^{\mathrm{}_{\sigma t}^+}(t)={\displaystyle \frac{1}{1+\alpha (1z)}}`$ (4)
As for incoming vertices, this formula follows from the evolution equation for $`𝒱_t^+(z)`$. Indeed, since the numbers of outgoing edges $`\mathrm{}_j^+(t)`$ from vertex $`j`$ at time $`t`$ and $`t1`$ differ by $`\mathrm{}_{j,t}`$ we have $`z^{\mathrm{}_j^+(t)}=z^{\mathrm{}_j^+(t1)}z^{\mathrm{}_{j,t}^+}`$. From definition (1) this gives

$`t𝒱_t^{}(z)=1+(t1)𝒱_{t1}^{}(z)(q_t+zp_t)`$
where the first term is the contribution of the newly added vertex at time $`t`$. The stationary limit is given by eq.(4). This is a geometric distribution, slightly larger than the Poisson distribution, with probabilities

$`\mathrm{Prob}_{(\mathrm{}^+=k)}={\displaystyle \frac{\alpha ^k}{(1+\alpha )^{k+1}}}`$ (5)
Mixed distribution. Let $`v_{k_+,k_{}}(t)`$ be the number of vertices with $`k_+`$ outgoing and $`k_{}`$ incoming edges at time $`t`$. As in eq.(1), the generating function for the mixed histogram distribution is defined by

$`𝒱_t(z_+,z_{}){\displaystyle \frac{1}{t}}{\displaystyle \underset{k_+,k_{}}{}}v_{k_+,k_{}}(t)z_+^{k_+}z_{}^k_{}={\displaystyle \frac{1}{t}}{\displaystyle \underset{1jt}{}}z_+^{\mathrm{}_j^+(t)}z_{}^{\mathrm{}_j^{}(t)}`$
By construction the outgoing and incoming edges variables $`\mathrm{}_j^\pm (t)`$ are statistically independent for $`j`$ fixed, so that the last expectation values factorize. As above we may derive an evolution equation by evaluating the contribution of the newly added vertex at time $`t`$. This yields:

$`t𝒱_t(z_+,z_{})=(q_t+z_{}p_t)^{t1}+(t1)𝒱_{t1}(z_+,z_{})(q_t+z_+p_t)`$
Its stationary limit is factorized:

$`𝒱(z_+,z_{})={\displaystyle \frac{e^{\alpha (z_{}1)}}{1+\alpha (1z_+)}}`$
Outgoing and incoming edges are statistically independent at large time.

## 4 Cluster distributions.

In this section, we present the main relations governing the probability distributions of connected components of the graphs. Two vertices belong to the same connected component if they can be joined by a path made of edges, without any reference to orientation. This definition ensures that the property of being in the same connected component is an equivalence relation, but does not ensure that two points in a connected component can be joined by an oriented path.

To partly avoid repetitions, the term cluster is used as a synonymous for connected component in the sequel.

Intuitively, the fact that the network is fragmented can be understood as follows : when a vertex $`t_0`$ is created, it has a finite probability to be isolated, and the probability that none of the vertices $`t_0+1,\mathrm{},t`$ connects to vertex $`t_0`$ is $`q_{t_0+1}\mathrm{}q_t`$ which scales as $`(t_0/t)^\alpha `$. This quantity remains finite as long as $`t_0/t`$ does. This argument shows that there are isolated vertices in the system. A small extension of the argument shows that there are also finite components and that young vertices are more likely to be in small components than old ones. This will be made more rigourous in the study of profiles, see section 8.

Let $`N_k(t)`$ be the number of connected components with $`k`$ vertices at time $`t`$ and let $`N_t(z)`$ be the generating function,

$$N_t(z)=\underset{k1}{}N_k(t)z^k$$

By definition, $`_kN_k(t)`$ is the number of components and $`_kN_k(t)k`$ the total number of vertices, $`_kN_k(t)k=t`$ at any finite time.

Let us write an evolution equation for $`N_t(z)`$. At time $`t+1`$, we add the vertex with label $`t+1`$ which may then be connected to $`n_k(t)`$ connected components of size $`k`$. This creates a new component of size $`1+_kn_k(t)k`$, but also removes $`n_k(t)`$ components of size $`k`$. Thus, at time $`t+1`$ we have:

$`N_k(t+1)=N_k(t)n_k(t)+\delta _{k;1+_pn_p(t)p}`$ (6)
with $`\delta _{j;k}`$ the Kronecker symbol. Alternatively,

$`N_{t+1}(z)=N_t(z){\displaystyle \underset{k1}{}}n_k(t)z^k+z^{1+_kn_k(t)k}`$ (7)
As is apparent from this formulation, the transition probability from a given $`N_t(z)`$ to a given $`N_{t+1}(z)`$ can be given in closed form. To be precise, the admissible $`N_t(z)`$’s (describing the accessible distributions of components at time $`t`$) are polynomials with integral non-negative coefficients, whose derivative at $`z=1`$ have value $`t`$. Now suppose $`N_t(z)`$ and $`N_{t+1}(z)`$ are admissible. If the difference $`N_{t+1}(z)N_t(z)`$ cannot be written as $`_{k1}n_k(t)z^k+z^{1+_kn_k(t)k}`$ for some set of nonnegative integers $`n_k(t)`$, the transition is forbidden. If it can, then the $`n_k(t)`$’s are uniquely defined and the transition probability is

$$\mathrm{Prob}(N_t(z)N_{t+1}(z))=\underset{k}{}\left(\genfrac{}{}{0pt}{}{N_k(t)}{n_k(t)}\right)q_{t+1}^{k(N_k(t)n_k(t))}(1q_{t+1}^k)^{n_k(t)}$$

The meaning of this equation is simple. At time $`t+1`$, the new vertex is added, and for each of the former $`t`$ points a (biased) coin is tossed to decide the value of the edge variables $`\mathrm{}_{j,t+1}`$. The tossings are independent with the same law, so the probability that the new point does not attach to a given component of size $`k`$ is $`q_{t+1}^k`$, and distinct components are independent. Hence for each $`k`$ one makes $`N_k(t)`$ independent Bernoulli trials with failure probability $`q_{t+1}^k`$, and the transition from $`N_t(z)N_{t+1}(z)`$ requires exactly $`n_k(t)`$ successes. This shows that the graph evolution is a (time inhomogeneous) Markov process on the space of components distributions, a fact that we shall use for the purpose of numerical simulations.

This explicit representation of the transition probability could be used to average equation (7). Alternatively, one can represent the number $`n_k(t)`$ of components of size $`k`$ which are connected to the new vertex in terms of the edge variables $`\mathrm{}_{i,j}`$ as

$$n_k(t)=\underset{[k]=1}{\overset{N_k(t)}{}}[1\underset{j[k]}{}(1\mathrm{}_{j,t+1})]$$

where $`[k]`$ runs over connected components of size $`k`$. Since the edge variables $`\mathrm{}_{j,t+1}`$ are statistically independent of the earlier edge variables, $`\mathrm{}_{j,k}`$ with $`kt`$, and therefore also independent of the $`N_k(t)`$’s, we have for any $`w`$,

$`w^{n_k(t)}=[q_{t+1}^k+(1q_{t+1}^k)w]^{N_k(t)}`$ (8)
In particular,

$`n_k(t)=(1q_{t+1}^k)N_k(t)`$ (9)
We can now take the average value of eq.(7) to get

$`N_{t+1}(z)`$ $`=`$ $`N_t(zq_{t+1})+z{\displaystyle \underset{k1}{}}[q_{t+1}^k+(1q_{t+1}^k)z^k]^{N_k(t)}`$ (10)
$`=`$ $`N_t(zq_{t+1})+zq_{t+1}^t{\displaystyle \underset{k1}{}}[1+(q_{t+1}^k1)z^k]^{N_k(t)}`$
Assume, but this will be justified later using the asymptotic behavior of the $`p_j`$’s, that at large time and for fixed component size $`N_k(t)/t`$ is self-averaging, meaning that

$$N_k(t)/tC_k+o(1),\mathrm{as}t\mathrm{}$$

with $`C_k`$ equals to its averaged value with probability one and with the remaining other $`o(1)`$ terms random. This in particular implies that for any finite size $`k`$ the number of connected components of size $`k`$ scales thermodynamically with the graph size. If both sides of eq.(10) are expended in powers of $`z`$, a given degree involves only components of bounded ($`t`$-independent) size. So, order by order in $`z`$, self-averaging applies and we conclude that taking $`1q_t\alpha /t`$, the following is an accurate approximation for the last term in eq.(10) at large times:

$$q_{t+1}^t\underset{k1}{}[1+(q_{t+1}^k1)z^k]^{N_k(t)}_t\mathrm{}\mathrm{exp}(\alpha +\alpha \underset{k1}{}kz^kC_k).$$

The averaged evolution equation (10) then gives a deterministic differential equation

$`0=C(z)\alpha z_zC(z)+z\mathrm{exp}(\alpha +\alpha z_zC(z))`$ (11)
for the generating function

$$C(z)\underset{k1}{}z^kC_k.$$

The function

$$z_zC(z)=\underset{k1}{}z^kP_k,P_kkC_k$$

has a direct probabilistic interpretation. Indeed, $`P_k`$ is the fraction of points in clusters of size $`k`$, or equivalently the probability that a randomly chosen vertex belongs to a connected component of size $`k`$.

By construction $`C(z)`$ is a Taylor series with positive coefficients. The series $`_kkC_k=_kP_k`$ is convergent because it counts the fraction of points in finite clusters which is $`1`$. As a consequence the radius of convergence of $`C(z)`$ is at least $`1`$: if we denote by $`R(\alpha )`$ the radius of convergence $`C(z)`$, we know that $`R(\alpha )1`$ for $`0\alpha <+\mathrm{}`$.

We show in Appendix A that analogous methods can be used to describe the connected components of the Erdös-Renyi random graph model.

Most of our subsequent analysis will be based on (11). It recursively determines the $`C_k`$’s. The first few are:

$`C_1`$ $`=`$ $`{\displaystyle \frac{e^\alpha }{(\alpha +1)}},C_2={\displaystyle \frac{\alpha e^{2\alpha }}{(\alpha +1)(2\alpha +1)}},`$ (12)
$`C_3`$ $`=`$ $`{\displaystyle \frac{\alpha ^2(6\alpha +5)e^{3\alpha }}{2(\alpha +1)^2(2\alpha +1)(3\alpha +1)}},\mathrm{}`$
More generally, the combinations $`C_ke^{k\alpha }`$ are rational functions in $`\alpha `$. They are more efficiently computed using a second order differential equation which follows from derivation of the logarithm of eq.(11):

$`z_z(C+\alpha z_zC)=(C+\alpha z_zC)(1+\alpha (z_z)^2C)`$ (13)
This leads to

$`(n1)(1+n\alpha )C_n={\displaystyle \underset{k+l=n}{}}k^2(1+\alpha l)C_kC_l`$ (14)
A general feature of the recursion relations is that the rational function $`C_ke^{k\alpha }`$ has no poles except possibly at the points $`1,1/2,\mathrm{},1/k`$.

Eq.(10) is not a closed formula for the probability distributions of the number of connected components. But as shown in Appendix B, one may improve slightly the above argument to obtain a closed system (see eq.58) which can be used to prove systematically that the variables $`C_k`$ are self-averaging, a task we perform for $`C_1`$ and $`C_2`$ in Appendix B. The proof becomes more and more tedious for large $`k`$. But we are confident that the self averaging property is general because of the close relationship with the Erdös-Renyi model presented in Appendix A.

## 5 Percolation transition.

Up to now, our analysis of connected components has always concentrated on finite components: the thermodynamic limit $`t\mathrm{}`$ has been taken while keeping $`k`$, the component size, arbitrary but fixed. We have argued that in this regime, the number of components of size $`k`$ is proportional to $`t`$. This has lead to a satisfactory description of $`C_k`$. Our arguments did not use any hypothesis on whether or not only finite components play a role in the thermodynamic limit. However, as we observed already, the sum $`_kP_k`$ measures exactly the fraction of the sites that are either isolated or in components of size $`2`$ or $`3`$ or $`\mathrm{}`$. To rephrase it more vaguely, $`_kP_k`$ counts the fraction of vertices that belong to finite connected components. If only clusters of finite size contribute to the thermodynamic limit, then $`_kP_k=1`$. Else, by standard physical arguments, a single giant component of size $`t(1_kP_k)`$ accounts for the deficit. The giant component is also called the percolating cluster. Its relative size, which we denote by $`P_{\mathrm{}}`$, is

$`P_{\mathrm{}}1{\displaystyle \underset{k1}{}}P_k=1_zC(1)`$ (15)
To discriminate between the two situations, $`P_{\mathrm{}}`$ may be computed numerically by evaluating a large number of coefficients $`P_k`$ using the recursion relation (14). As we shall see later, the convergence of the series is slow for the whole relevant range of $`\alpha `$. The result of such a partial summation $`_{kk_{max}}P_k`$ is plotted in Fig.1 for $`k_{max}=2^{11},\mathrm{},2^{15}`$. It reveals a phase transition at a value $`\alpha _c`$ between $`.24`$ and $`.29`$, going from a regime where finite components contain all vertices to a regime where they do not. Below $`.24`$ and above $`.29`$, the plots corresponding to different values of $`k_{max}`$ are hard to distinguish, but in the transition region, large values of $`k`$ make substantial contributions to $`P_k`$. The transition is also manifest on an analogous study of $`kP_k`$. We shall show later that $`1/4`$ is the exact threshold in this model and that $`kP_k`$ is discontinuous at the transition.

The growth of $`P_{\mathrm{}}`$ just above the threshold seems to start with many vanishing derivatives, in strong contrast with what happens in the Erdös-Renyi random graph, for which the growth of the giant component is linear close to the transition. This can be related to the following observations:

– As we have recalled in Appendix A, in the Erdös-Renyi model the components of size $`k`$ occupy a fraction $`\frac{k^{k1}}{k!}\alpha ^{k1}e^{k\alpha }`$. As a function of $`\alpha `$, this fraction has a single maximum at $`\alpha =\frac{k1}{k}`$. These values accumulate at $`\alpha =1^{}`$, the well-known transition point for the standard random graph. Then for $`\alpha 1`$, the fraction of sites occupied by components of size $`k`$ decreases with a finite slope for all $`k`$’s, and so does the sum, so that the growth of the giant component is linear close to the transition.

– In the model studied in this paper, the behavior of $`P_k(\alpha )`$ as a function of $`\alpha `$ for generic $`k`$ is not so easy to get at. However, a simple numerical analysis leads to the following picture: $`P_1`$ is a decreasing function of $`\alpha `$, but for $`k>1`$, $`P_k`$ has a single maximum, at say $`\alpha _k`$. This sequence starts with $`\alpha _2.241`$, $`\alpha _3.311`$ $`\alpha _4.341`$, is maximum for $`k=12`$ with $`\alpha _{12}.375`$ and then decreases very slowly ($`\alpha _{100}.338`$, $`\alpha _{1000}.301`$, $`\alpha _{10000}.282`$), apparently getting closer and closer to $`1/4`$. At the transition, most $`P_k(\alpha )`$’s are still increasing, and close enough to the transition a finite but large number of them is still increasing. So subtle compensation mechanisms can take place, leading possibly to the vanishing of (infinitely) many derivatives of $`P_{\mathrm{}}`$ at $`\alpha =1/4`$.

This is confirmed in the following subsections, which are also devoted to a more precise description of the distributions of finite and infinite clusters. To summarize:

$`P_{\mathrm{}}`$ $`=`$ $`0,\mathrm{if}\alpha 1/4,`$
$`P_{\mathrm{}}`$ $`>`$ $`0,\mathrm{if}\alpha >1/4.`$
## 6 Behavior below the transition, $`\alpha <\alpha _c`$.

We turn to the examination of the consequences of eq.(11), the equation that determines the generating function for the number of clusters of given finite size.

We know that $`R(\alpha )`$, the radius of convergence of $`C(z)`$, is at least $`1`$. To analyze the behavior of $`C(z)`$ around $`z=1`$, we define

$`F(\tau )\alpha \tau \alpha _\tau Y(\tau ),Y(\tau )=C(e^\tau ).`$ (16)
From eq.(11) $`F`$ satisfies

$`\alpha (1e^F)_\tau F=F\tau `$ (17)
Below the transition, there is no percolating cluster so that $`_\tau Y(0)=_kkC_k=1`$ or alternatively,

$`{\displaystyle \underset{k1}{}}P_k=1,\alpha <1/4`$ (18)
This makes clear that the normalized positive numbers $`P_k`$ are the probabilities for a vertex to be in a connected component of size $`k`$. For $`F(\tau )`$ this translates into the boundary condition

$$F(0)=0,\alpha <1/4$$

### 6.1 Scaling laws.

We first look for a formal solution $`F_{tr}(\tau )`$ of eq.(17) in the form of a Taylor series in $`\tau `$:

$`F_{tr}(\tau )=\tau \alpha {\displaystyle \underset{m1}{}}{\displaystyle \frac{\mu _m}{m!}}\tau ^m`$ (19)
The Taylor coefficients are the moments of the measure $`P_k`$: $`\mu _m=_kk^mP_k`$. For example $`\mu _1`$ is the average proportion of vertices per component. As eq.(13) the differential equation eq.(11) may be turned into a second order differential equation for $`Y(\tau )`$,

$`{}_{\tau }{}^{}Y+\alpha ^2{}_{\tau }{}^{}Y=(Y+\alpha {}_{\tau }{}^{}Y)(1+\alpha {}_{\tau }{}^{}{}_{}{}^{2}Y)`$ (20)
Eq.(11) with $`{}_{\tau }{}^{}Y(0)=1`$ gives $`Y(0)=1\alpha `$. Eq.(20) then allows us to recursively compute the $`\mu _m`$’s:

$`\mu _1`$ $`=`$ $`{\displaystyle \frac{4}{(1+b)^2}}={\displaystyle \frac{12\alpha \sqrt{14\alpha }}{2\alpha ^2}}`$ (21)
$`{\displaystyle \frac{1}{2}}(n1(n+1)b)\mu _n`$ $`=`$ $`(1+\alpha \mu _1)\mu _{n1}\alpha {\displaystyle \underset{\genfrac{}{}{0pt}{}{k+l=n+1}{k,l2}}{}}\left({\displaystyle \genfrac{}{}{0pt}{}{n}{k}}\right)\mu _l(\mu _{k1}+\alpha \mu _k)`$
with

$`4\alpha =1b^2,b=\sqrt{14\alpha }`$
The square root singularity in the expression of $`\mu _1`$ indicates that the initial boundary condition (18) becomes pathological at $`\alpha =1/4`$ and thus signals the percolation transition. The recursion relation for the higher Taylor coefficients shows that $`\mu _n`$ possesses a pole at $`b=(n+1)/(n1)`$. It actually changes sign, from positive to negative, across the pole. As a function of $`n`$, the $`\mu _n`$’s have a simple pole at $`n_{}(b)`$,

$`\mu _n{\displaystyle \frac{\mathrm{reg}.}{nn_{}(b)}},n_{}(b)={\displaystyle \frac{1+b}{1b}},n2`$ (22)
where the numerator does not vanish at $`(n1)=b(n+1)`$. Since $`\mu _n=_kk^nP_k`$ at least for $`n<n_{}(b)`$, this simple pole codes for the asymptotic behavior of $`P_k`$ at large $`k`$. Indeed, recall that if $`P_kk^\nu `$ as $`k\mathrm{}`$ then $`_kk^xP_k`$ diverges as $`x\nu 1`$ from below with a simple pole at this value. Thus assuming that we may extend the simple pole (22) to non-integers value of $`n`$ we learn that

$`P_k_k\mathrm{}\mathrm{const}.k^{\nu (\alpha )}`$ (23)
with

$`\nu (\alpha )=1+n_{}(b)={\displaystyle \frac{1+\sqrt{14\alpha }}{2\alpha }}`$ (24)
This means in particular that the probability distribution $`P_k`$ (which describes the proportion of the system occupied by clusters of size $`k`$) is large and only the moments $`\mu _n`$ for $`n<n_{}(b)`$, exist.

The fact that the Taylor coefficients $`\mu _n`$, computed from the recursion relation (21), cannot consistently be interpreted as moments of a probability distribution for $`n>n_{}(b)`$ indicates that the expansion (19) is only up to $`o(\tau ^n_{})`$ terms. Indeed, the differential equation (17) is compatible with an expansion of $`F(\tau )`$ for $`\tau <0`$, $`z=e^\tau <1`$, of the form:

$`F(\tau )=\tau \alpha {\displaystyle \underset{q,p0}{}}{\displaystyle \frac{()^q}{q!}}y_{p,q}(\tau )^{q+pn_{}(b)},`$
with $`y_{0,0}=0`$ and $`y_{0,m}=\mu _m`$. As a function of the complex variable $`z`$, $`_zC(z)`$ has thus a branch cut starting at $`z=1`$:

$`_zC(z)\mathrm{const}.(1z)^{n_{}(b)}+\mathrm{},\mathrm{around}z=1`$ (25)
For $`n_{}(b)`$ an integer, this formula should becomes $`_zC(z)\mathrm{const}.(1z)^{n_{}(b)}\mathrm{log}(1z)+\mathrm{}`$. The cut implies the asymptotic behavior (23) for $`P_k=_0\frac{dz}{2i\pi }_zC(z)z^k`$.

At the transition but from below, i.e. $`\alpha =1/4^{}`$, there are logarithmic correction to the scaling behavior (23), and the branch cut equation is $`_zC(z)\mathrm{const}.(1z)/\mathrm{log}(1z)+\mathrm{}`$. This ensures that the first moment is finite, and its value, computed from (21), is $`\mu _1|_{\alpha =1/4^{}}=4`$

### 6.2 Scaling domain, $`\alpha <\alpha _c`$.

The scaling law (23) may be linked to the typical growth rate of large clusters in the system. For concreteness, consider the component of vertex $`t^{}`$ for any given $`t^{}`$. For very large $`t`$, the number of arrows emerging from vertex $`t^{}`$ grows like $`\alpha \mathrm{log}t`$. Then we infer that the size of the genealogical tree of $`t^{}`$ will grow like $`t^\alpha `$ (under the hypothesis that the genealogical tree is indeed tree-like, a reasonable assumption for small $`\alpha `$). This counting of descendants gives a crude lower bound for the size of the connected component of $`t^{}`$. Hence we expect that the system contains components whose sizes grow like a power of $`t`$. To estimate this power, we argue as follows.

Consider a given large cluster of size $`k(t)t`$ at time $`t`$:

i) $`k(t+1)k(t)`$ is $`0`$ with probability $`q_{t+1}^{k(t)}`$ and $`1+_ppn_p(t)`$ with probability $`1q_{t+1}^{k(t)}\alpha k(t)/t`$ times the probability that vertex $`t+1`$ connects to $`n_p(t)`$ clusters of size $`p`$ apart from the large cluster;

ii) removing the given large cluster does not change the thermodynamical properties of the graph, so the probability that vertex $`t+1`$ connects to $`n_p(t)`$ clusters of size $`p`$ apart from the large cluster is simply the probability that vertex $`t+1`$ connects to $`n_p(t)`$ clusters of size $`p`$. Hence for large $`t`$, from eq.(9), $`_ppn_p(t)\alpha _kkP_k=\alpha \mu _1`$;

iii) suppose we add $`\delta t`$ new vertices with $`\delta tt`$ but $`k(t)\delta tt`$. Between time $`t`$ and time $`t+\delta t`$, many new clusters have be connected to the given large components so $`1k(t+\delta t)k(t)k(t)`$, but this has not changed the thermodynamical properties of the graph. Hence we can average the equation in i) to get a deterministic equation

$$k(t+\delta t)k(t)\alpha k(t)(1+\alpha \mu _1)\delta t/t.$$

This leads to

$$k(t)t^{1/\nu }\text{ with }1/\nu =\alpha (1+\alpha \mu _1)=2\alpha /(1+b).$$
(26)
We find that the growth rate of large clusters is universal. As expected, their growth exponent is larger than $`\alpha `$, the genealogical tree growth, because it takes not only descendants into account but the whole component. The difference is maximum at the transition, where $`\nu =2`$.

The fact that the same exponent, $`\nu `$, governs the asymptotic behavior of $`P_k`$ at large $`k`$ and the size of large clusters for large $`t`$ can be understood directly as follows. First consider one realization of the random graph for a given $`t`$ and suppose that there is a single component of maximal size, say $`K`$. Observe that $`_{kl}kN_k(t)`$ is by definition the number of points in components of size larger than $`l`$. This number is strictly larger than $`l`$ if $`l<K`$, but vanishes if $`l>K`$. So $`K`$ is characterized by the identity $`_{kK}kN_k(t)=K.`$

From this we infer by taking the average that for large $`t`$ the relation

$$\underset{kk(t)}{}kN_k(t)k(t)$$

gives a sensible characterization for $`k(t)`$, the order of magnitude of the size of large components in the graph. We write $`_{kk(t)}kN_k(t)=t_{k<k(t)}kN_k(t)`$ and use that for $`\alpha <1/4`$, $`_kP_k=1`$ to write $`_{kk(t)}tP_k+_{k<k(t)}tP_kkN_k(t)k(t)`$. For large $`k(t)`$, the asymptotics of the first sum is $`_{kk(t)}tP_ktk(t)^{\nu +1}`$. The second sum is made of finite size corrections. If we assume that these that these are not too large, we conclude that $`tk(t)^{\nu +1}k(t)`$, i.e. that $`k(t)t^{1/\nu }`$.

In fact, experience from finite size scaling suggests that the two sums give contributions of the same order of magnitude. The idea is that $`tP_k`$ is the main contribution to $`kN_k(t)`$ not only when $`tP_k`$ is of order $`t`$, but even when simply $`tP_k1`$ so that self averaging remains valid. This means that $`tP_kkN_k(t)`$ and $`kN_k(t)`$ become of the same order of magnitude only when $`k`$ is so large that $`tP_k1`$. So again we see that $`k(t)`$ is characterized by $`tP_{k(t)}1`$, and we conjecture that $`tk^\nu `$ is a scaling variable.

## 7 Behavior above the transition, $`\alpha >\alpha _c`$.

Above the transition there is an giant component. Let $`P_{\mathrm{}}`$ be its relative size. By definition, $`{}_{\tau }{}^{}Y(0)=1P_{\mathrm{}}`$ or alternatively

$`{\displaystyle \underset{k1}{}}P_k+P_{\mathrm{}}=1,\alpha >1/4`$ (27)
This makes clear that the $`P_k`$’s and $`P_{\mathrm{}}`$ define the probability distribution of vertices among the clusters of different sizes, with $`P_{\mathrm{}}`$ the probability for a vertex to be in the percolating cluster.

The size of the giant component $`k_{\mathrm{}}`$ increases linearly with time : $`k_{\mathrm{}}tP_{\mathrm{}}`$ for $`t`$ large. The slope may be evaluated as follows. Imagine adding a new vertex at time $`t+1`$. It is connected to the percolating cluster with probability $`1q_{t+1}^k_{\mathrm{}}`$ so that

$$k_{\mathrm{}}(t+1)(1q_{t+1}^k_{\mathrm{}})(k_{\mathrm{}}(t)+1+\underset{p}{}pn_p(t))+q_{t+1}^k_{\mathrm{}}k_{\mathrm{}}(t)$$

where, as in eq.(6), $`n_p(t)`$ are the numbers of components of size $`p`$ connected to the new vertex. This can be rearranged as $`k_{\mathrm{}}(t+1)k_{\mathrm{}}(t)(1q_{t+1}^{k_{\mathrm{}}(t)})(1+_ppn_p(t))`$. The quantity $`P_{\mathrm{}}`$ is self averaging, so using $`q_t1\alpha /t`$ and $`k_{\mathrm{}}tP_{\mathrm{}}`$ at large time we infer $`P_{\mathrm{}}=(1e^{\alpha P_{\mathrm{}}})(1+_ppn_p(t))`$ which by use of eq.(9) leads to

$$P_{\mathrm{}}=(1e^{\alpha P_{\mathrm{}}})(1+\alpha \mu _1)$$
(28)
As shown below the above equation is exact, but it does not determine $`P_{\mathrm{}}`$ as this requires knowing the moment $`\mu _1`$.

Above the transition, $`F(\tau )`$ still satisfies the differential equation (17) but with a different boundary condition:

$`F(0)=\alpha P_{\mathrm{}},\alpha >1/4`$ (29)
with $`P_{\mathrm{}}`$ vanishing at the transition i.e. as $`\alpha 1/4^+`$. This modifies the behavior of its moments. So let us expand $`F(\tau )`$ in Taylor series around the origin:

$`F(\tau )=\alpha P_{\mathrm{}}\tau \alpha {\displaystyle \underset{m1}{}}{\displaystyle \frac{\mu _m}{m!}}\tau ^m`$
with $`\mu _m=_kk^mP_k`$. Eq.(17) implies eq.(28) which is thus exact. The second order differential equation (20) does not fix $`P_{\mathrm{}}`$ but determines recursively the $`\mu _m`$’s which depend parametrically on $`P_{\mathrm{}}`$.

It remains to decipher what the behavior of the size of the giant component is, at least close to the transition. This will follow from an analysis of the behavior of $`F(\tau )`$ close to its singularities. As we are going to show, $`F(\tau )`$ possesses square root branch cut at a point $`\tau _c>0`$ but exponentially close to the origin. This means that above the transition the singularity of the function $`_zC(z)`$ is at $`z_c=e^{\tau _c}>1`$, with $`z_c`$ depending on $`\alpha `$, a behavior which has to be compared with the cut at $`z=1`$ below the transition, eq.(25). More precisely, we show that

$$\mathrm{log}P_{\mathrm{}}\pi /\beta ,\tau _cP_{\mathrm{}}/8e\text{for }\alpha 1/4^+$$
(30)
with $`\beta `$ defined by

$$4\alpha \beta ^2+1,\beta \sqrt{4\alpha 1}.$$

### 7.1 Preliminaries.

To prove eq.(30) we shall look at the behavior of $`F`$ in the neighborhood of three different points:

– at the origin $`\tau =0`$ where $`F(\tau )`$ takes the boundary value (29),

– at the branch point $`\tau _c>0`$ with $`F(\tau _c)=0`$, $`{}_{\tau }{}^{}F(\tau _c)=\mathrm{}`$ and,

– around the point $`\tau _d<0`$ specified by the condition $`F(\tau _d)+2\tau _d=0`$.

Let us first show $`\tau _c`$ and $`\tau _d`$ exist. The function we are interested in is defined by $`F=\tau +\alpha (1_{k1}P_ke^{k\tau })`$, and satisfies the differential equation eq.(17). So $`_\tau F<0`$ and $`F`$ decreases from $`+\mathrm{}`$ at $`\tau =\mathrm{}`$ to $`0`$ at a point $`\tau _c`$ which is non-negative since $`F(0)0`$. $`F`$ also satisfies the obvious inequalities $`0<F+\tau <\alpha `$ for $`\tau <0`$ and $`{}_{\tau }{}^{}{}_{}{}^{2}F<0`$. So there is a single point $`\tau _d`$, with $`\alpha <\tau _d<0`$, at which $`F(\tau _d)+2\tau _d=0`$ and $`F+2\tau `$ has the sign of $`\tau \tau _d`$. Both $`\tau _c`$ and $`\tau _d`$ vanish at the transition, $`\beta 0`$.

It turns out that most of these properties could be proved using only the differential equation. More precisely, take any solution $`E`$ of (17) with an initial condition at $`\tau _i<0`$ such that $`E(\tau _i)+\tau _i>0`$. The function $`E`$ can be extended on the right as a positive function on a maximal interval $`[\tau _i,\tau _f[`$. Then $`E`$ is strictly decreasing for $`\tau >\tau _i`$, $`lim_{\tau _f^{}}E=0`$, $`\tau _f0`$, and $`\tau _f>0`$ if $`\alpha >1/4`$. This discussion emphasizes the intrinsic role played by $`\alpha =1/4`$ in this problem.

Before proving eq.(30) let us get intuition from a simpler, but more universal, version of the differential equation (17). If we look at a region where $`F`$ is small – and we know that such regions exist generically – we can estimate $`1e^F`$ by $`Ff`$ with $`f`$ approximating $`F`$. This leads to the simplified, but still non-linear, differential equation

$$\alpha f{}_{\tau }{}^{}f+\tau +f=0.$$
(31)
It turns out that this simpler equation can be solved in closed form. Before doing that, let us further assume that $`\tau `$ is close to $`\tau _c`$ so that $`f`$ can be neglected compared to $`\tau `$. Then we get $`\alpha f{}_{\tau }{}^{}f\tau `$, with solution

$$f\sqrt{(\tau _c^2\tau ^2)/\alpha },$$

leading to the announced square root singularity. But this does not give informations on the location of the branch point.

Consider now the following functional of an arbitrary function $`f`$ of $`\tau `$

$$\frac{1}{2}\mathrm{log}(\alpha f^2+\tau f+\tau ^2)\frac{1}{\beta }\mathrm{arctan}\left(\frac{\beta f}{f+2\tau }\right).$$
(32)
To fix conventions, we specify the function $`\mathrm{arctan}`$ by demanding that it is continuous and takes value in $`]\pi /2,\pi /2[`$. The total derivative of this functional with respect to $`\tau `$ is

$$\frac{\alpha f{}_{\tau }{}^{}f+\tau +f}{\alpha f^2+\tau f+\tau ^2}.$$

It vanishes if $`f`$ is a solution of eq.(31), so we have indeed solved in closed form equation (17) for small $`F`$. As this is the domain we are interested in, it is tempting to argue that $`F`$ and $`f`$ should exhibit the same singular behaviour. This turns out to be true, but there are some subtleties because the limit $`\alpha 1/4^+`$ is singular: the size of the domain for which $`f`$ is a uniform approximation to $`F`$ shrinks to $`0`$.

### 7.2 Rigourous estimates

Our strategy is to use the invariant (32) for the approximate equation (31) to derive exact inequalities for $`F`$.

Let us observe that the functional (32) is singular at points where $`f+2\tau =0`$ where it has a jump of amplitude $`\pm \pi /\beta `$. With this in mind, we define the functional

$`I(\alpha ^{},F(\tau ),\tau )`$ $``$ $`{\displaystyle \frac{1}{2}}\mathrm{log}(\alpha ^{}F^2+\tau F+\tau ^2)`$ (37)
$`\{\begin{array}{ccc}\frac{1}{\beta ^{}}\mathrm{arctan}\left(\frac{\beta ^{}F}{F+2\tau }\right)+\frac{\pi }{\beta ^{}}\hfill & \mathrm{if}& \hfill \tau <\tau _d<0\\ \frac{\pi }{2\beta ^{}}\hfill & \mathrm{if}& \hfill \tau =\tau _d\\ \frac{1}{\beta ^{}}\mathrm{arctan}\left(\frac{\beta ^{}F}{F+2\tau }\right)\hfill & \mathrm{if}& \hfill \tau _d<\tau <\tau _c\end{array}`$
In this definition, $`\alpha ^{}=(1+\beta ^2)/4>1/4`$ is a priori independent of $`\alpha `$, the value for which $`F(\tau )`$ is considered. The functional $`I(\alpha ^{},F(\tau ),\tau )`$ is a smooth function of $`\tau `$ on $`]\mathrm{},\tau _c[`$. Using the differential equation for $`F`$, its total derivative with respect to $`\tau `$ is found to be

$$\frac{dI}{d\tau }=\frac{\tau +F}{\alpha ^{}F^2+\tau F+\tau ^2}\left(1\frac{\alpha ^{}}{\alpha }\frac{F}{1e^F}\right).$$
(38)
If $`\alpha ^{}=\alpha `$, the right-hand side is always negative, so $`I`$ is decreasing, and comparing its values at $`\tau <\tau _d`$, at $`\tau _d`$, at $`0`$ and at $`\tau _c`$ we find

$$I(\alpha ,F(\tau ),\tau )>\mathrm{log}(\beta |\tau _d|e^{\pi /2\beta })>\mathrm{log}(\alpha ^{3/2}P_{\mathrm{}})\frac{1}{\beta }\mathrm{arctan}\beta >\mathrm{log}\tau _c.$$
(39)
We know that $`\tau _d>\alpha `$, so taking a fixed $`\tau <1/4`$, we can take the limit $`\alpha 1/4^+`$. At point $`\tau `$, $`F`$ is analytic in $`\alpha `$, so $`I(\alpha ,F(\tau ),\tau )+\pi /\beta `$ has a finite limit, and we get exponentially small upper-bounds for $`|\tau _d|,P_{\mathrm{}}`$ and $`\tau _c`$. For instance, $`|\tau _d|\beta e^{\pi /2\beta }`$ is bounded above when $`\alpha 1/4^+`$.

For $`\alpha ^{}<\alpha `$, $`\frac{dI}{d\tau }`$ changes sign at a point $`\tau ^{}`$ solution of $`1\frac{\alpha ^{}}{\alpha }\frac{F}{1e^F}=0`$. This point is unique as $`F`$ is decreasing and goes from $`+\mathrm{}`$ at $`\mathrm{}`$ to $`0`$ at $`\tau _c`$. Let us choose $`\alpha ^{}`$ such that $`\tau ^{}\tau _d`$ – and such choices exist. Then $`\frac{dI}{d\tau }`$ goes from negative to positive if $`\tau `$ increases so that $`I`$ increases for $`\tau `$ varying from $`\tau ^{}`$ to $`\tau _c`$. This leads to lower bounds

$$I(\alpha ^{},F(\tau ^{}),\tau ^{})<\mathrm{log}(\beta ^{}|\tau _d|e^{\pi /2\beta ^{}})<\mathrm{log}(\alpha _{}^{}{}_{}{}^{1/2}\alpha P_{\mathrm{}})\frac{1}{\beta ^{}}\mathrm{arctan}\beta ^{}<\mathrm{log}\tau _c.$$
(40)
Take $`\alpha ^{}/\alpha =(e^{2\tau _d}1)/2\tau _d`$, in such a way that $`\frac{dI}{d\tau }`$ vanishes exactly at $`\tau _d`$. Then $`\alpha ^{}=\alpha (1+\tau _d+\mathrm{})`$, and $`\beta ^{}=\beta (1+\tau _d/2\beta ^2+\mathrm{})`$ when $`\alpha 1/4^+`$. Comparing the lower bounds with the upper bounds obtained above, we see that

$$\beta |\tau _d|e^{\pi /2\beta }P_{\mathrm{}}/8e\tau _c\mathrm{for}\alpha 1/4^+.$$

To get a lower bound for $`|\tau _d|`$ we take $`\alpha ^{}`$ such that $`\frac{dI}{d\tau }`$ vanishes at $`\tau ^{}<\tau _d`$. Then $`\frac{1}{\beta ^{}}\mathrm{arctan}\left(\frac{\beta ^{}F}{F+2\tau }\right)>0`$ and $`\alpha ^{}F^2+\tau F+\tau ^2=\beta _{}^{}{}_{}{}^{2}F^2/4+(\tau +F/2)^2\beta _{}^{}{}_{}{}^{2}F^2/4`$ so we have a crude bound $`I(\alpha ^{},F(\tau ^{}),\tau ^{})+\pi /\beta ^{}>\mathrm{log}\beta ^{}F(\tau ^{})/2`$. Taking for instance $`\beta ^{}=\beta \beta ^2`$, so that $`\beta ^{}/\beta 1`$ and $`F(\tau ^{})0`$ as $`\beta 0^+`$. Reporting in $`\frac{dI}{d\tau }=0`$ yields $`F4\beta ^3`$. Thus, from (40), $`|\tau _d|\beta ^3e^{\pi /2\beta }`$ is bounded below when $`\alpha 1/4^+`$.

The above upper and lower bounds then imply

$$2\mathrm{log}\tau _d\mathrm{log}\tau _c\mathrm{log}P_{\mathrm{}}\pi /\beta \text{as}\beta 0^+.$$

We restate the physically most important result: up to an algebraic prefactor, the size of the percolation cluster close to percolation is

$$P_{\mathrm{}}e^{\pi /\sqrt{4\alpha 1}}.$$
(41)
In fact, we have obtained a better estimate. We expect that $`\tau _d\beta ^{1\gamma }e^{\pi /2\beta }`$ as a finite limit for $`\beta 0^+`$ for a certain $`\gamma `$, or equivalently

$$P_{\mathrm{}}\mathrm{const}\beta ^\gamma \mathrm{e}^{\pi /2\beta },\beta 0^+.$$
(42)
We have proved that if $`\gamma `$ exists, $`0\gamma 4`$. The asymptotic behavior of the probabilities $`P_k`$ is not the same below and above the transition as $`C(z)`$ does not have the same analytical properties on the two sides of the transition.

Since above the transition the branch point is located at $`z_c=\mathrm{exp}\tau _c>1`$, the $`P_k`$’s now decrease exponentially. More precisely, $`_zC(z)`$ possesses a square root branch point at $`z_c`$,

$$_zC(z)=\mathrm{const}.\sqrt{z_cz}+\mathrm{}$$

so that $`P_k=\frac{dz}{2i\pi }z^k_zC(z)`$ behave as

$`P_k_k\mathrm{}\mathrm{const}.k^{3/2}z_c^k\mathrm{const}.k^{3/2}e^{k\tau _c},`$ (43)
to be compared with eq.(23).

### 7.3 Universality.

We give now a universality argument suggesting that $`\gamma =0`$ and this is confirmed by solving numerically the differential equation (17).

Our invocation of universality rests on eq.(43), in which $`\tau _c`$, or equivalently (in the vicinity of the critical point) $`P_{\mathrm{}}/8e`$, controls the exponential decrease of the $`P_k`$’s and plays the role of a mass gap. This mass gap is exponentially small close to the transition, and by analogy we may argue that increasing $`\beta `$ is a marginally relevant perturbation of the percolating critical point. Introducing the Wilson-Callan-Symanzik beta function $`B(\beta )`$, we expect a relation

$$\tau _cP_{\mathrm{}}/8e\mathrm{exp}^\beta \frac{d\beta ^{}}{B(\beta ^{})}.$$
(44)
Comparison with our formula gives $`B(\beta )=\beta ^2/\pi +B_3\beta ^3+\mathrm{}`$. It is known from field theory that the coefficient $`B_3`$, which dictates the exponent of the algebraic prefactor, is universal.

Such universal features are controlled by the continuum limit. In our framework, the continuum region is reached at small $`\tau `$, and small $`F`$, and we expect that the continuum limit is governed by (31). Note that both $`\tau _c`$ and $`|\tau _d|`$ are exponentially small close to the transition, and $`F(\tau )`$ remains small for $`\tau `$ between these points. This a posteriori justifies looking at the approximate equation (31). Consequently, as far as universal quantities are concerned, the inequalities in (39) can be replaced by equalities. This leads to $`\gamma =0`$, or $`B_3=0`$ i.e.

$$P_{\mathrm{}}\mathrm{const}\mathrm{e}^{\pi /2\beta }$$

as announced.

As universality could suggest, eq(31) turns out to describe the continuum limit for a larger class of evolving networks than just the specific one we are studying. This is the case for the model studied in . We refer to the original paper for the definitions. It suffices to say that eq.(17) is replaced by

$$2\delta S_\tau S=S(e^\tau 1).$$

In this equation, $`1S`$ is a generating function, the coefficient of $`e^{k\tau }`$ giving the fraction of points in components of finite $`k`$, so $`S(0)`$ is the fraction of sites occupied by the giant component. The parameter $`\delta `$ is the average number of edges created at each time step, so this is is precisely the equivalent of our $`\alpha `$. To study $`S`$ close to $`\tau =0`$, the approximation is $`e^\tau 1\tau `$, and we retrieve (31), with $`\alpha `$ replaced by $`2\delta `$. The percolation threshold is $`\delta =1/8`$ (as expected, percolation thresholds are not universal), and the size of the infinite component behaves like $`e^{\pi /2\sqrt{2}(\delta 1/8)^{1/2}}`$. The prefactor $`\pi /2\sqrt{2}1.111`$ compares quite well with the numerical value $`1.132\pm 0.008`$ obtained in the original paper. This fact was already noted in .

### 7.4 Scaling regime, $`\alpha >\alpha _c`$.

We have accumulated evidence that the percolation phase transition is very similar to the Kosterlitz-Thouless phase transition in the XY model. For instance, the probabilities $`P_k`$ follow the scaling laws (23) below the transition with scaling exponents varying continuously with the parameter $`\alpha `$ while they decrease exponentially above the transition.

The discussion of universality in the previous section suggests to look for scaling functions describing the neighborhood of the critical point. For instance, from eq.(20) it follows that the moments $`\mu _m`$, $`m2`$, diverge as $`P_{\mathrm{}}0`$ with

$`\mu _m_{\alpha 1/4^+}{\displaystyle \frac{g_m}{P_{\mathrm{}}^{m1}}},g_1=12,g_2=64,\mathrm{}`$
In particular, $`\mu _1`$ is discontinuous at the transition, jumping from $`4`$ to $`12`$. These scaling relations imply that the function

$$G(x)P_{\mathrm{}}^1F(xP_{\mathrm{}})$$

has a finite limit, denoted by $`G_c(x)`$, as $`\alpha 1/4^+`$ for any fixed $`x`$,

$`G_c(x)=\alpha x\alpha {\displaystyle \underset{m1}{}}{\displaystyle \frac{g_m}{m!}}x^m`$ (45)
The differential equation (17) at $`\alpha =1/4^+`$ then reduces to $`G_c_xG_c+4(G_c+x)=0`$. With the boundary condition $`G_c(0)=1/4`$, this is integrated as:

$`\mathrm{log}\left(4(G_c+2x)\right)+{\displaystyle \frac{2x}{G_c+2x}}=0`$ (46)
The left-hand side is very reminiscent of the limit $`\alpha 1/4`$ of eq.(32). As a consequence $`G_c(x)`$ has a square root branch point at $`x_c=1/8e`$ at which $`G_c`$ vanishes and $`g_m/m!m^{3/2}(8e)^m`$ for $`m`$ large. Actually, a more precise computation based on the Lagrange formula presented in Appendix C yields the exact value of the scaling coefficients $`g_m`$:

$`g_{m+1}=m^m\mathrm{\hspace{0.17em}8}^{m+1},m1.`$ (47)
## 8 Chronological profiles

One can make a direct counting of the average number of connected components in the random graph at time $`t`$ which are copies of a given finite labeled graph. This will allow us to retrieve the results of section 4. But this also leads to a detailed local-in-time description of the connected components which illustrates the consequences of the chronological memory of our model.

### 8.1 Tree distributions.

Let $`G`$ be a labeled graph with vertices $`1,2,\mathrm{},k`$. We let $`m_j`$ be the number of edges connecting vertex $`j`$ to a vertex with smaller label, so that $`m=m_1+\mathrm{}+m_k`$ is the number of edges of $`G`$. We look for the average number of increasing maps $`v`$ from $`[1,\mathrm{},k]`$ to $`[1,\mathrm{},t]`$ such that the vertices $`v_1,\mathrm{},v_k`$ span a connected component of the random graph isomorphic to $`G`$. This number is the average number of connected components isomorphic to $`G`$ in the random graph.

By the rules of construction of the random graph, the probability that the vertices $`v_1,\mathrm{},v_k`$ with $`1v_1<\mathrm{}<v_k\mathrm{}t`$ span a connected component of the random graph isomorphic to $`G`$ is

$$\underset{i=1}{\overset{k}{}}\left((1q_{v_i})^{m_i}q_{v_i}^{v_i1m_i}\underset{v_i<w_i<v_{i+1}}{}q_{w_i}^i\right)$$

with the convention $`v_{k+1}t+1`$.

The average we look for is obtained by summing this expression over the $`v_i`$’s. For large $`t`$, using the asymptotic behavior $`p_j\alpha /j`$ for large $`j`$, the sum can be formally reinterpreted as a Riemann sum, leading to a contribution

$$t^{km}e^{k\alpha }\alpha ^m_{0\sigma _1\mathrm{}\sigma _k1}𝑑\sigma _1\sigma _1^{\alpha m_1}\mathrm{}𝑑\sigma _k\sigma _k^{\alpha m_k}.$$

In deriving this formula, we have not treated carefully the contribution of small values of the $`v_i`$’s. This is reflected in the fact that the integral can be divergent if some $`m_i`$’s are too large. However, the prefactor $`t^{km}`$ is the sign that in this case the sum over the $`v_i`$’s is nevertheless negligible in the large $`t`$ limit, as a more careful treatment would show.

The most salient feature of this formula is that only connected graphs with $`k=m+1`$ give a contribution proportional to $`t`$, i.e. contribute to $`C_k`$. Since $`k`$ is the number of vertices and $`m`$ the number of links, this relation characterizes trees as follows from the Euler formula: only trees contribute to the thermodynamic limit. A given tree on $`k`$ vertices with incoming degrees $`m_i`$, $`i=1,\mathrm{},k`$ gives a contribution

$$\frac{e^{k\alpha }\alpha ^{k1}}{(\alpha +1m_1)(2\alpha +2m_1m_2)\mathrm{}(k\alpha +km_1\mathrm{}m_k)}$$
(48)
to $`C_k`$. Observe that for a tree, $`m_1+\mathrm{}+m_ii1`$ for $`i=1,\mathrm{},k`$ so that all integrals are well-defined and finite for real non-negative $`\alpha `$. It is amusing to note also that the contribution of a single tree of size $`k`$ can contain poles in $`\alpha `$ at values that are not in the list $`1,1/2,\mathrm{},1/k`$. These poles have to cancel between different trees in the sum over trees of size $`k`$, because we know that they are absent in $`C_k`$. But we have no simple explanation for this cancellation.

This explicit formula makes it easy to show that $`C(z)`$ is (complex) analytic in $`\alpha `$ in a neighborhood of $`]0,+\mathrm{}[`$ for every $`z`$ such that $`|z|<1`$. Indeed, we know that for $`\alpha ]0,+\mathrm{}[`$, and $`|z|<1`$ the series $`_kP_k(\alpha )z^k`$ is absolutely convergent. But $`P_k(\alpha )`$ is a sum of non-negative contributions, each tree giving a contribution of the form (48). Now suppose that $`\alpha `$ is complex with positive real part. For each tree contribution and for fixed $`\mathrm{}\alpha `$, the modulus of

$`{\displaystyle \frac{e^{k\alpha }}{(\alpha +1m_1)(2\alpha +2m_1m_2)\mathrm{}(k\alpha +km_1\mathrm{}m_k)}}`$
is maximal when $`\mathrm{}\alpha =0`$ and then the expression is real and positive. This is because the statement is true for every factor. Taking the sum over trees we infer that $`P_k(\alpha )P_k(\mathrm{}\alpha )\left(\frac{|\alpha |}{\mathrm{}\alpha }\right)^{k1}`$. So the series $`_kP_k(\alpha )z^k`$ is absolutely convergent if $`|z\alpha |<\mathrm{}\alpha `$. This equation defines, for fixed $`|z|<1`$, a neighborhood of $`\alpha ]0,+\mathrm{}[`$ in which we have an absolutely convergent sum of analytic functions analytic of $`\alpha `$. Hence the sum is analytic in $`\alpha `$ as claimed.

To resum more explicitly the contribution of all trees of a given size, we need the generating function for labeled trees with given incoming degrees. Suppose more generally that we give a weight $`x_i`$ for each edge leaving vertex $`i`$ (i.e. connecting $`i`$ to a $`j>i`$) and a weight $`y_i`$ for each edge entering vertex $`i`$ (i.e. connecting $`i`$ to a $`j<i`$). The generating function $`𝒯`$ for weighted trees on $`n`$ vertices factorizes nicely as

$`𝒯`$ $`=`$ $`x_1(y_2x_1+y_2x_2+y_3x_2+\mathrm{}+y_kx_2)`$
$`(y_3x_1+y_3x_2+y_3x_3+y_4x_3+\mathrm{}+y_nx_3)\mathrm{}`$
$`(y_{n1}x_1+\mathrm{}+y_{n1}x_{n1}+y_nx_{n1})y_n.`$
This generalization of the famous Caley tree formula222Which states that there are $`n^{n2}`$ labeled trees on $`n`$ vertices. implies it immediately. It seems to be little known, although it is implicit in the mathematical literature . Gilles Schaeffer  provided us with a clean proof using a refined version of one of the standard proofs of the Caley tree formula, putting trees on $`k`$ vertices in one to one correspondence with applications from $`[1,\mathrm{},k]`$ to $`[1,\mathrm{},k]`$ fixing $`1`$ and $`k`$, see e.g. .

This formula can be specialized to $`x_i=1`$ and $`y_i=1/\sigma _i`$ for $`i=1,\mathrm{},k`$ to give

$`C_k`$ $`=`$ $`e^{k\alpha }\alpha ^{k1}{\displaystyle _{0\sigma _1\mathrm{}\sigma _k1}}𝑑\sigma _1\mathrm{}𝑑\sigma _k(\sigma _1\mathrm{}\sigma _k)^\alpha `$
$`({\displaystyle \frac{2}{\sigma _2}}+{\displaystyle \frac{1}{\sigma _3}}+\mathrm{}+{\displaystyle \frac{1}{\sigma _k}})({\displaystyle \frac{3}{\sigma _3}}+{\displaystyle \frac{1}{\sigma _4}}+\mathrm{}+{\displaystyle \frac{1}{\sigma _k}})\mathrm{}({\displaystyle \frac{k1}{\sigma _{k1}}}+{\displaystyle \frac{1}{\sigma _k}}){\displaystyle \frac{1}{\sigma _k}}.`$
If one integrates only over a subset of the $`\sigma `$’s, one gets marginal distributions. For instance, for $`k=1`$, we get that in the thermodynamic limit the fraction of sites with age close to $`t\sigma `$ that are isolated is $`e^\alpha \sigma ^\alpha `$. For $`k=2`$, if we integrate over $`\sigma _2`$, we get that the fraction of sites with age close to $`t\sigma `$ that are the older vertex of a tree on two vertices is $`e^{2\alpha }(\sigma ^\alpha \sigma ^{2\alpha })`$ while if we integrate over $`\sigma _1`$, we get that the fraction of sites with age close to $`t\sigma `$ that are the younger vertex of a tree on two vertices is $`e^{2\alpha }\sigma ^{2\alpha }\frac{\alpha }{\alpha +1}`$. The sum, $`e^{2\alpha }(\sigma ^\alpha \sigma ^{2\alpha }\frac{1}{\alpha +1})`$, gives the probability that a site with age close to $`t\sigma `$ belongs to a tree on $`2`$ vertices.

Our explicit representation in terms of trees shows that in this model, and at least for questions concerning connected components, the thermodynamic limit applies not only to the full system, but also to slices of fixed relative age $`\sigma =t^{}/t`$. In the next section we shall study $`\sigma `$-dependent profiles.

For small $`k`$’s, we have done all the integrals and checked the agreement with the value of $`C_k`$ obtained by the recursion relation. But a general proof valid for all $`k`$’s is lacking.

### 8.2 Dating finite components.

To illustrate the evolving nature of our model, we now determine the local-in-time distribution of the cluster sizes. This means determining for any given age interval what is the proportion of vertices of these ages which belong to clusters of given size.

Define $`p_k(t,t^{})`$ to be the probability that vertex $`t^{}`$ belongs to a component of size $`k`$ at time $`t`$. Guided by the previous tree representation, we infer that in the thermodynamical limit $`_{t^{}[t\sigma ,t(\sigma +d\sigma )]}p_k(t,t^{})t\rho _k(\sigma )d\sigma `$, with $`\rho _k(\sigma )`$ a deterministic function. By construction, $`_0^1\rho _k(\sigma )𝑑\sigma =kC_k`$, the total fraction of points that belong to components of size $`k`$.

The reasoning leading to (11) can be generalized: one writes down a recursion relation for $`p_k(t+1,t^{})`$ and then takes the average, a step justified by the explicit tree representation. The event that vertex $`t^{}`$ belongs to a component of size $`k`$ at time $`t+1`$ is the exclusive union of several events.

i) Vertex $`t^{}`$ belonged to a component of size $`k`$ at time $`t`$ and this component is not linked to the new vertex $`t+1`$. This has probability $`p_k(t,t^{})q_{t+1}^k`$.

ii) Vertex $`t^{}`$ belonged to a component of size $`l<k`$ at time $`t`$, this component is linked to the new vertex $`t+1`$, and together with the other components linked to $`t+1`$ (say $`n_m(t)`$ components of size $`m`$), it builds a component of size $`k=l+1+_mmn_m(t)`$. This has probability

$$p_l(t,t^{})(1q_{t+1}^l)\underset{m}{}\left(\genfrac{}{}{0pt}{}{N_m(t)\delta _{m,l}}{n_m(t)}\right)q_{t+1}^{m(N_m(t)\delta _{m,l}n_m(t))}(1q_{t+1}^m)^{n_m(t)}$$

To perform the explicit sum over $`l`$ and the $`n_m(t)`$’s, we introduce again generating functions and set $`p_{t,t^{}}(z)=_lp_l(t,t^{})z^l`$. This leads to

$`p_{t+1,t^{}}(z)`$ $`=`$ $`p_{t,t^{}}(zq_{t+1})`$
$`+\left({\displaystyle \underset{l}{}}p_l(t,t^{}){\displaystyle \frac{z^l(1q_{t+1}^l)}{q_{t+1}^l+z^l(1q_{t+1}^l)}}\right){\displaystyle \underset{m}{}}(q_{t+1}^m+z^m(1q_{t+1}^m))^{N_m(t)}`$
In the large $`t`$ limit this complicated formula simplifies if we use again the hypothesis of self-averaging and asymptotic independence. Defining

$`\rho (\sigma ,z)={\displaystyle \underset{k}{}}\rho _k(\sigma )z^k,`$ (50)
this leads to

$$\sigma _\sigma \rho =\alpha (1ze^{\alpha +\alpha z_zC})z_z\rho .$$

Together with the sum rule relating $`\rho `$ to $`C`$, this fixes completely the profiles $`\rho _k(\sigma )`$. A relation between $`\rho (1,z)`$ and $`C(z)`$ is obtained by integrating (8.2) for $`\sigma `$ between $`0`$ and $`1`$. Using the defining equation for $`C(z)`$, eq.(11), this leads to $`\rho (1,z)=(1+\alpha z_z)C(z)`$. Thus, we can summarize our knowledge on component time profiles with the four relations :

$`\sigma _\sigma \rho (\sigma ,z)`$ $`=`$ $`\alpha (1\rho (1,z))z_z\rho (\sigma ,z)`$
$`{\displaystyle _0^1}𝑑\sigma \rho (\sigma ,z)`$ $`=`$ $`z_zC(z)`$ (51)
$`\rho (1,z)`$ $`=`$ $`ze^{\alpha +\alpha _0^1\rho (\sigma ,z)𝑑\sigma }`$
$`\rho (1,z)`$ $`=`$ $`(1+\alpha z_z)C(z).`$
Expansion in powers of $`z`$ leads to recursion relations for the $`\rho _k(\sigma )`$’s. They can be shown to be alternating polynomials of degree $`k`$ in $`x=\sigma ^\alpha `$ with vanishing constant coefficient. The first few polynomials are

$`\rho _1`$ $`=`$ $`e^\alpha x`$
$`\rho _2`$ $`=`$ $`e^{2\alpha }\left(x{\displaystyle \frac{x^2}{\alpha +1}}\right)`$
$`\rho _3`$ $`=`$ $`e^{3\alpha }\left({\displaystyle \frac{3\alpha +2}{2(\alpha +1)}}x{\displaystyle \frac{2}{\alpha +1}}x^2+{\displaystyle \frac{3\alpha +2}{2(\alpha +1)^2(2\alpha +1)}}x^3\right)`$
Again, we have checked that for small $`k`$’s the values of $`\rho _k`$ as computed from iterated tree integrals and from the generating function coincide, but we have no general proof.

As an application, let us look at the profile of vertices of relative age close to $`\sigma `$ which are the youngest in components of size $`k`$. These are vertices that created, when they appeared, a component of size $`k`$ – this happens with probability $`\rho _k(1)=(1+k\alpha )C_k`$ – which was then left untouched for the rest of the evolution – this happens with probability $`\sigma ^{k\alpha }`$. So the distribution of vertices that are the youngest in their connected (finite) component is

$$\underset{k}{}(1+k\alpha )C_k\sigma ^{k\alpha }=(1+\alpha z_z)C(\sigma ^\alpha ).$$
(52)
In this expression, the formal parameter $`z`$ has been replaced by $`\sigma ^\alpha `$.

The giant component, when it exists, also has a well-defined profile to which we turn now.

### 8.3 Dating the percolating cluster.

We now determine the profile of the giant component.

So let $`\rho _{\mathrm{}}(\sigma )td\sigma `$ be the fraction of vertices whose ages are between $`\sigma t`$ and $`(\sigma +d\sigma )t`$ which belong to the giant component. By definition $`_0^1𝑑\sigma \rho _{\mathrm{}}(\sigma )=P_{\mathrm{}}`$. We know from the tree representation that for finite clusters the thermodynamic limit applies not only on the full graph, but also on time slices. The giant component is the complement of finite components, so we expect that the density $`\rho _{\mathrm{}}(\sigma )`$ is self-averaging as is the size of the percolating cluster. To derive an equation fixing this density we look for the probability for a site of age $`j=\sigma t`$ not to be in the percolating cluster. On the one hand, by definition of the density this probability is

$$\mathrm{Prob}(j=\sigma t[k_{\mathrm{}}])=1\rho _{\mathrm{}}(\sigma )$$

On the other hand, this probability may be evaluated by demanding the vertex $`j`$ not to be connected to the older and younger vertices of the percolating cluster:

$$\mathrm{Prob}(j=\sigma t[k_{\mathrm{}}])=\underset{k<j,k[k_{\mathrm{}}]}{}(1\frac{\alpha }{j})\underset{k>j,k[k_{\mathrm{}}]}{}(1\frac{\alpha }{k})$$

At large time, the first above product converges to $`\mathrm{exp}(\alpha _0^\sigma 𝑑\zeta \frac{\rho _{\mathrm{}}(\zeta )}{\sigma })`$ and the second to $`\mathrm{exp}(\alpha _\sigma ^1𝑑\zeta \frac{\rho _{\mathrm{}}(\zeta )}{\zeta })`$. So we get the relation:

$`1\rho _{\mathrm{}}(\sigma )`$ $`=`$ $`\mathrm{exp}\left(\alpha {\displaystyle _0^\sigma }𝑑\zeta {\displaystyle \frac{\rho _{\mathrm{}}(\zeta )}{\sigma }}\alpha {\displaystyle _\sigma ^1}𝑑\zeta {\displaystyle \frac{\rho _{\mathrm{}}(\zeta )}{\zeta }}\right)`$ (53)
$`=`$ $`\mathrm{exp}\left(\alpha {\displaystyle _0^1}𝑑\zeta \rho _{\mathrm{}}(\zeta )\mathrm{min}({\displaystyle \frac{1}{\sigma }},{\displaystyle \frac{1}{\zeta }})\right)`$
As $`\rho _{\mathrm{}}`$ is positive, $`_0^1𝑑\zeta \rho _{\mathrm{}}(\zeta )\mathrm{min}(\frac{1}{\sigma },\frac{1}{\zeta })`$ is an increasing function of $`\sigma `$. So $`\rho _{\mathrm{}}(\sigma )`$ is decreasing and has a right limit at $`0`$. Unless $`\rho _{\mathrm{}}0`$ (i.e. $`\alpha 1/4`$), the integral has a logarithmic divergence, and henceforth $`\rho _{\mathrm{}}(0)=1`$. More precisely,

$$\rho _{\mathrm{}}(\sigma )=1\sigma ^\alpha \mathrm{exp}\left(\alpha _0^1𝑑\zeta \frac{1\rho _{\mathrm{}}(\zeta )}{\zeta }\right)+\mathrm{}$$
(54)
This means that the early vertices belong to the giant component with probability $`1`$.

On the other hand, by definition $`_0^1𝑑\zeta \rho _{\mathrm{}}(\zeta )=P_{\mathrm{}}`$ so taking $`\sigma =1`$ in (53) we get that

$$\rho _{\mathrm{}}(1)=1e^{\alpha P_{\mathrm{}}}.$$
(55)
This means that the late vertices always belong to the giant component with a non vanishing probability – although this probability is exponentially small close to the threshold. Hence the giant component invades all time slices above the threshold, and the term percolation transition is appropriate even with this unusual interpretation. Our results are illustrated on Fig.2.

To conclude this discussion, let us observe that $`\rho _{\mathrm{}}`$ can be expressed in terms of the profiles of finite components that we studied before:

$$\rho _{\mathrm{}}(\sigma )=1\underset{k}{}\rho _k(\sigma )=1\rho (\sigma ,z=1)$$

with $`\rho (\sigma ,z)`$ defined in eq.(50). We end up with yet another equation constraining $`\rho (\sigma ,z)`$,

$$\rho (\sigma ,z=1)=\sigma ^\alpha \mathrm{exp}\left(\alpha +_0^1𝑑\zeta \rho _{\mathrm{}}(\zeta )\mathrm{min}(\frac{1}{\sigma },\frac{1}{\zeta })\right).$$

If we expand $`\rho (\sigma ,z=1)`$ in powers of $`\sigma ^{n\alpha }`$ as $`\rho (\sigma ,z=1)=_{k1}D_n(1+n\alpha )\sigma ^{n\alpha }`$ and define $`D(z)=_{k1}\frac{D_n}{n}z^n`$, we can check that

$$z_zD(z)+\alpha (z_z)^2D(z)=z\mathrm{exp}(\alpha D(z)+D(1)+z_zD(1))$$
(56)
In particular, by differentiation this implies that

$$(n1)(1+n\alpha )D_n=\underset{k+l=n}{}k^2(1+\alpha l)D_kD_l$$
(57)
The resemblance with equations (11) and (14) for the function $`C(z)`$ is rather striking, but we have not been able to use this more deeply. One difference is that (11) determines $`C_1`$ from scratch, whereas (56) does it in a two step process. If $`d_n`$ is the solution of this equation satisfying $`d_1=1`$, then $`D_n=d_nD_1^n`$, and then $`D_1`$ has to satisfy

$$(1+\alpha )D_1=\mathrm{exp}(\alpha +\underset{n1}{}\frac{1+n\alpha }{n}d_nD_1^n)$$

Another difference is that the sequence $`D_n`$ alternates in sign. Again, the formal parameter $`z`$ receives a simple physical interpretation $`z=\sigma ^\alpha `$ where $`\sigma `$ labels the relative date of birth of vertices.

## 9 Conclusions

In this study, we have described detailed global and local-in-time features of evolving random graphs with uniform attachment rules.

Concerning global properties, we have shown that the model has a percolation phase transition at $`\alpha =1/4`$. Below the transition, the system contains clusters whose sizes scale like $`t^{(1\sqrt{14\alpha })/2}`$. Above the transition, a single component, the giant component, grows steadily with time. We have shown that, close to the threshold, the fraction of sites in the giant component has an essential singularity and behaves as $`e^{\pi /\sqrt{4\alpha 1}}`$. The behaviors below and above the transitions are strongly reminiscent of the two dimensional $`XY`$ model, so that our model can be interpreted as some algorithmic equivalent of it, but it seems unlikely that a direct connexion exists. This analogy and further scaling properties we present call for an alternative renormalization group approach to the transition.

By describing local-in-time profiles, we have shown that they offer a more accurate vision of the specificities of evolving graphs. It would be desirable to generalize this approach by answering the following question : assume that a procedure to assign ages to the vertices of some evolving graphs has been given, what informations on the microscopic evolution rules of these graphs can be decoded from the knowledge of local-in-time statistics ?

Aknowledgements: We thank Gilles Schaeffer for his clarifying help in tree generating functions and Sergei Dorogovtsev for his comments and interest.

## 10 Appendix A.

We show how our arguments of section 4 can be modified to describe the Erdös-Renyi random graph model. In this model, one starts with $`n`$ points, and any two points are connected by an edge with probability $`\alpha /n`$ (so in this model all the points are equivalent). Then a limit $`n\mathrm{}`$ is taken. This famous model describes a static graph, but it can also be rephrased as an evolving graph in the following way : set $`t=n`$ and suppose that points are added one by one, from $`1`$ to $`t`$, each new point connecting to any previous one with probability $`\alpha /t`$. From this point of view, it can be seen that looking only at the first $`t^{}`$ vertices ($`t^{}t`$) amount to look at an Erdös-Renyi random graph with a modified connectivity parameter $`\alpha ^{}=\alpha t^{}/t`$. To get a recursion relation, we start from an Erdös-Renyi random graph of size $`t`$ with connectivity parameter $`\alpha `$. We add vertex $`t+1`$ and connect any older vertex to it with probability $`q_{t+1}=\alpha /t`$, so that the effective connectivity parameter for the graph on $`t+1`$ vertices is $`\alpha (t+1)/t`$. Then the derivation proceeds as before, with the little proviso that in eq.(10), on the left-hand side $`N_{t+1}(z)`$ has an effective connectivity parameter $`\alpha (t+1)/t`$ instead of $`\alpha `$. So when we take the thermodynamic limit, nothing changes on the right-hand side of (10), but an additional term $`\alpha _\alpha `$ contributes to the left-hand side. If $`c`$ denotes the analog of $`C`$ but for the Erdös-Renyi random graph, then

$$\alpha _\alpha c=c\alpha z_zc+z\mathrm{exp}(\alpha +\alpha z_zc).$$

This little modification in the equation has drastic consequences. The new equation has a single solution regular at $`\alpha =0`$, namely

$$c=\underset{k1}{}\frac{k^{k2}}{k!}\alpha ^{k1}e^{k\alpha }z^k$$

which is the well-known result. We conclude that our self-averaging hypothesis is valid for the Erdös-Renyi model. This makes it more plausible that it works for our original model as well, a fact also confirmed by numerical simulations.

## 11 Appendix B.

Our goal is to prove the following formula and to investigate a few of its consequences. We claim that

$`{\displaystyle \underset{k1}{}}w_k^{N_k(t+1)}={\displaystyle _0}{\displaystyle \frac{d\xi }{2i\pi }}{\displaystyle \underset{j1}{}}{\displaystyle \frac{w_j}{\xi ^j}}{\displaystyle \underset{m1}{}}[q_{t+1}^mw_m+(1q_{t+1}^m)\xi ^m]^{N_m(t)}`$ (58)
where the contour integral is around the origin. It is a Fokker-Planck equation for the Markov process formed by the $`N_k(t)`$’s, which could be used to prove systematically that the variables $`C_k`$ are self-averaging, a task we perform for $`C_1`$ and $`C_2`$ at the end of this appendix.

We start from eq.(6) which may be rewritten as

$`{\displaystyle \underset{k}{}}w_k^{N_k(t+1)}={\displaystyle \underset{j}{}}w_j^{N_j(t)n_j(t)}w_j^{\delta _{j;1+_ppn_p(t)}}`$
To compute the last term we insert the tautological identity

$$1=\underset{k1}{}\delta _{k;1+_ppn_p(t)}$$

in the r.h.s. to get

$`{\displaystyle \underset{k}{}}w_k^{N_k(t+1)}={\displaystyle \underset{k1}{}}w_k{\displaystyle \underset{j}{}}w_j^{N_j(t)n_j(t)}\delta _{k;1+_ppn_p(t)}`$
To compute the r.h.s. expectation value we use a contour integral representation of the Kronecker symbol

$$\delta _{k;1+_ppn_p(t)}=_0\frac{d\xi }{2i\pi }\xi ^{k+_ppn_p(t)}$$

This yields

$`{\displaystyle \underset{k}{}}w_k^{N_k(t+1)}={\displaystyle \underset{k1}{}}{\displaystyle _0}{\displaystyle \frac{d\xi }{2i\pi }}\left({\displaystyle \frac{w_k}{\xi ^k}}\right){\displaystyle \underset{p}{}}(\xi ^p/w_p)^{n_p(t)}w_p^{N_p(t)}`$
The r.h.s. can now be computed using eq.(8) and gives eq.(58).

This computation has a rather simple combinatorial reinterpretation. We set

$$F_t(w_1,w_2,\mathrm{})=\underset{k}{}w_k^{N_k(t)},$$

and observe that going from time $`t`$ to $`t+1`$, we add vertex $`t+1`$ and edges from the rest of the graph to $`t+1`$. Suppose that the component of vertex $`t+1`$ has size say $`k`$. This component was build by “eating” some components of the graph at time $`t`$. A component of size $`j`$ is swallowed with probability $`1q_{t+1}^j`$ and survives with probability $`q_{t+1}^j`$. On the other hand, if one expands

$$F_t(w_1q_{t+1}+\xi (1q_{t+1}),w_2q_{t+1}^2+\xi ^2(1q_{t+1}^2),\mathrm{})$$

in powers of $`\xi `$, the term of degree $`k1`$ enumerates all the possibilities to build the component of vertex $`t+1`$ with the correct probability. Defining $`\widehat{w}_l(\xi )w_lq_{t+1}^l+\xi ^l(1q_{t+1}^l)`$, summation over $`k`$ gives

$$F_{t+1}(w_1,w_2,\mathrm{})=\underset{k}{}w_k_0\frac{d\xi }{2i\pi }\xi ^kF_t(\widehat{w}_1(\xi ),\widehat{w}_2(\xi ),\mathrm{}),$$

as obtained previously.

The Fokker-Planck eq.(58) may also be formulated as a difference equation, similar to a discrete Schrodinger equation. Let us for instance specify it for $`w_k=_j\zeta _j^k`$ with $`\zeta _j`$, $`j=0,\mathrm{},M`$, a set of complex numbers and let

$`\mathrm{\Psi }_t(\zeta _0,\mathrm{},\zeta _M){\displaystyle \underset{k}{}}[\zeta _0^k+\mathrm{}+\zeta _M^k]^{N_k(t)}`$
This parametrization is similar to that used in matrix theory where $`w_k`$ may be thought of as the trace of the $`k^{\mathrm{th}}`$ power of matrix whose eigenvalues are the $`\zeta _j`$’s. The contour integral in eq.(58) can then be explicitly evaluated by deforming the integration contour to pick the simple pole contributions located at the points $`\xi =\zeta _j`$. This gives:

$`\mathrm{\Psi }_{t+1}(\zeta _0,\mathrm{},\zeta _M)=q_{t+1}^t{\displaystyle \underset{j}{}}\zeta _j\mathrm{\Psi }_t(\mathrm{},\zeta _j/q_{t+1},\mathrm{})`$ (59)
Eq.(58) or (59) may be used to prove that the numbers $`C_k`$ are self-averaging. Let us choose for example two parameters $`\zeta _0=1`$ and $`\zeta _1=z/t`$. Then, only the clusters of size $`1`$ give a non trivial contribution to $`\mathrm{\Psi }_t`$ at large time so that

$$\mathrm{\Psi }_t^{(1)}(z)\mathrm{\Psi }_t(1,z/t)=\underset{k}{}(1+z^k/t^k)^{tC_k}_t\mathrm{}e^{zC_1}$$

At large time, the difference equation (59) then reduces to a differential equation for $`\mathrm{\Psi }_t^{(1)}(z)`$,

$$(1+\alpha )_z\mathrm{\Psi }_t^{(1)}(z)=e^\alpha \mathrm{\Psi }_t^{(1)}(z)$$

It implies that $`\mathrm{log}\mathrm{\Psi }_t^{(1)}(z)`$ is linear in $`z`$ at large time which means that $`C_1`$ is self-averaging and stationary at large time. More precisely, integrating the above equation gives:

$$\mathrm{\Psi }_t^{(1)}(z)=\mathrm{exp}(zC_1(\alpha )),C_1(\alpha )=\frac{e^\alpha }{\alpha +1}$$

in agreement with eq.(12).

Similarly, to prove that $`C_2`$ possesses a finite self-averaging limit as $`t\mathrm{}`$ we choose three parameters $`\zeta _0=1`$, $`\zeta _1=\zeta _2=\sqrt{z/t}`$ so that only clusters of size two survive in $`\mathrm{\Psi }_t`$ at large time,

$$\mathrm{\Psi }_t^{(2)}(z)\mathrm{\Psi }_t(1,\sqrt{z/t},\sqrt{z/t})_t\mathrm{}e^{2zC_2}$$

Eq.(59) then gives a first order differential equation for $`\mathrm{\Psi }_t^{(2)}(z)`$ which implies that

$$\mathrm{\Psi }_t^{(2)}(z)=\mathrm{exp}(2zC_2(\alpha )),C_2(\alpha )=\frac{e^{2\alpha }\alpha }{(2\alpha +1)(\alpha +1)}$$

in agreement with eq.(12). Although we do not have a global argument this proof may clearly be extended to recursively prove self-averageness of any $`C_k`$ by choosing the parameters $`\zeta _j=\omega ^j`$ with $`\omega ^k=z/t`$.

## 12 Appendix C.

Here we present the proof of eq.(47). Let us first recall Lagrange formula. Consider a variable $`X`$ defined by the implicit relation $`f(X)=y`$ for some given analytic function $`f`$. The solution of this equation is supposed to be unique so that $`X`$ is function of $`y`$. Given another analytic function $`g(w)`$ we look for the Taylor series in $`y`$ of $`g(X(y))`$. This composed function may be presented as a contour integral:

$$g(X(y))=\frac{dz}{2i\pi }\frac{f^{}(z)}{f(z)y}g(z)$$

Expanding the integrated rational function in Taylor series in $`y`$ gives Lagrange formula:

$$g(X(y))=\underset{n0}{}y^n\frac{dz}{2i\pi }\frac{f^{}(z)}{f(z)^{n+1}}g(z)$$

Let $`H(x)4(G_c(x)+2x)`$. Eq.(46) translates into $`H(x)\mathrm{log}H(x)=8x`$ or $`W(y)e^{W(y)}=y`$ for $`y=8x`$ and $`H(x)=e^{W(y)}`$. We now apply Lagrange formula with $`f(X)=Xe^X`$ and $`g(w)=e^w`$. This gives

$`H(y)`$ $`=`$ $`{\displaystyle \underset{n0}{}}y^n{\displaystyle \frac{dz}{2i\pi }\frac{z+1}{z^{n+1}}e^{(1n)z}}`$
$`=`$ $`1+{\displaystyle \underset{n1}{}}{\displaystyle \frac{(1n)^{n1}}{n!}}y^n`$
With $`H(x)=4(G_c(x)+2x)`$ and $`y=8x`$, this proves eq.(47).