{'fraction_non_alphanumeric': 0.042908801002249364, 'fraction_numerical': 0.012812847014606646, 'mean_word_length': 4.287069095288273, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Kinetic Theory of Random Graphs: from Paths to Cycles

Structural properties of evolving random graphs are investigated. Treating linking as a dynamic aggregation process, rate equations for the distribution of node to node distances (paths) and of cycles are formulated and solved analytically. At the gelation point, the typical length of paths and cycles, l, scales with the component size k as l ~ k^{1/2}. Dynamic and finite-size scaling laws for the behavior at and near the gelation point are obtained. Finite-size scaling laws are verified using numerical simulations.


I. INTRODUCTION
A random graph is a set of nodes that are randomly joined by links. When there are sufficiently many links, a connected component containing a finite fraction of all nodes, the so-called giant component, emerges. Random graphs, with varying flavors, arise naturally in statistical physics, chemical physics, combinatorics, probability theory, and computer science [1,2,3,4,5].
Several physical processes and algorithmic problems are essentially equivalent to random graphs. In gelation, monomers form polymers via chemical bonds until a giant polymer network, a "gel", emerges. Identifying monomers with nodes and chemical bonds with links shows that gelation is equivalent to the emergence a giant component [6,7,8]. A random graph is also the most natural mean-field model of percolation [9,10]. In computer science, satisfiability, in its simplest form, maps onto a random graph [11]. Additionally, random graphs are used to model social networks [12,13].
Random graphs have been analyzed largely using combinatorial and probabilistic methods [3,4,5]. An alternative statistical physics methodology is kinetic theory, or equivalently, the rate equation approach. The formation of connected components from disconnected nodes can be treated as a dynamic aggregation process [14,15,16,17]. This kinetic approach was used to derive primarily the size distribution of components [18,19,20].
Recently, we have shown that structural characteristics of random graphs can be analyzed using the rate equation approach [21]. In this study, we present a comprehensive treatment of paths and cycles in evolving random graphs. The rate equation approach is formulated by treating linking as a dynamic aggregation process. This approach allows an analytic calculation of the path length distribution. Since a cycle is formed when two connected nodes are linked, the path length distribution yields the cycle length distribution. More subtle statistical properties of cycles in random graphs can be calculated as well.
In particular, the probability that the system contains no cycles and the size distribution of the first, second, etc. cycles are obtained analytically.
We focus on the behavior near and at the phase transition point, namely, when the gel forms. We show that the path and the cycle length distribution approach selfsimilar distributions near the gelation transition. At the gelation point, these distributions develop algebraic tails.
The exact results obtained for an infinite system allow us to deduce scaling laws for finite systems. Using heuristic and extreme statistics arguments, the size of the giant component at the gelation point is obtained. This size scale characterizes the size distribution of components and it leads to a number of scaling laws for the typical path size and cycle size. Extensive numerical simulations validate these scaling laws for finite systems.
The rest of the paper is organized as follows. First, the evolving random graph process is introduced (Sec. II), and then the size distribution of all components is analyzed in Sec. III. Statistical properties of paths are derived in Sec. IV and then used to obtain statistical properties of all cycles (Sec. V) and of the first cycle (Sec. VI). We conclude in Sec. VII. Finally, in an appendix, some details of contour integration used in the body of the paper are presented.

II. EVOLVING RANDOM GRAPHS
A graph is a collection of nodes joined by links. In a random graph, links are placed randomly. Random graphs may be realized in a number of ways. The links may be generated instantaneously (static graph) or sequentially (evolving graph); additionally a given pair of nodes may be connected by at most a single link (simple graph) or by multiple links (multi-graph).
We consider the following version of the random graph model. Initially, there are N disconnected nodes. Then, a pair of nodes is selected at random and a link is placed between them (Fig. 1). This linking process continues ad infinitum and it creates an evolving random graph. The process is realized dynamically. Links are generated with a constant rate in time, set equal to (2N ) −1 without loss of generality. There are no restrictions associated with the identity of the two nodes. A pair of nodes may be selected multiple times, i.e., a multi-graph is created. Additionally, the two nodes need not be different, so selfconnections are allowed.
At time t, the total number of links is on average N t/2, the average number of links per node (the degree) is t, and the average number of self-connections per node is N −1 t/2. Therefore, whether or not self-connections are allowed is a secondary issue. Since the linking process is completely random, the degree distribution is Poissonian with a mean equal to t.

III. COMPONENTS
The evolving random graph model has several virtues that simplify the analysis. First, the linking process is completely random as there is no memory of previous links. Second, having at hand a continuous variable (time) allows us to use continuum methods, particularly the rate equation approach. This is best demonstrated by determination of the size distribution of connected components.
As linking proceeds, connected components form. When a link is placed between two distinct components, the two components join. For example, the latest link in Fig. 1 joins two components of size i = 2 and j = 4 into a component of size k = i + j = 6. Generally, there are i × j ways to join disconnected components. Hence, components undergo the following aggregation process Two components aggregate with a rate proportional to the product of their sizes.

A. Infinite Random Graph
Let c k (t) be the density of components containing k nodes at time t. In terms of N k (t), the total number of components with k nodes, then c k (t) = N k (t)/N . For finite random graphs, both N k (t) and c k (t) are random variables, but in the N → ∞ limit the density c k (t) becomes a deterministic quantity. It evolves according to the nonlinear rate equation (the explicit time dependence is dropped for simplicity) The initial condition is c k (0) = δ k,1 . The gain term accounts for components generated by joining two smaller components whose sizes sum up to k. The second term on the right-hand side of Eq. (2) represents loss due to linking of components of size k to other components. The corresponding gain and loss rates follow from the aggregation rule (1). The rate equations can be solved using a number of techniques. Throughout this investigation, we use a convenient method in which the time dependence is eliminated first. Solving the rate equations recursively yields c 1 = e −t , c 2 = 1 2 te −2t , c 3 = 1 2 t 2 e −3t , etc. These explicit results suggest that c k (t) = C k t k−1 e −kt . Substituting this form into (2), we find that the coefficients C k satisfy the recursion relation subject to C 1 = 1. This recursion is solved using the generating function approach. The form of the righthand side of Eq. (3) suggests to utilize the generating function of the sequence kC k rather than C k , i.e., G(z) = k k C k e kz . Multiplying Eq. (3) by k e kz and summing over all k, we find that the generating function satisfies the nonlinear ordinary differential equation Integrating this equation, z = ln G − G + A and using the asymptotics G → e z as z → −∞ fixes the constant A = 0. Thus, we arrive at an implicit solution for the generating function The coefficients C k can be extracted from (5) via the Lagrange inversion formula, or using contour integration as detailed in Appendix A. Substituting r = 1 in Eq. (A1) reproducing the well-known result for the size distribution [18,19] In the following, we shall often use the generating function for the size distribution c(z, t) = k k c k (t)e kz . This generating function is readily expressed via the auxiliary generating function G(z) = k k C k e kz : Let us consider the fraction of nodes in finite components, M 1 = k k c k (t). This quantity is merely the first moment of the size distribution (hence the notation). Equivalently M 1 = c(z = 0, t). From (7) we find M 1 = τ /t with τ = G(ln t − t). Using (5), we express τ through t: For t < 1, there is a single root τ = t, and all nodes reside in finite components, M 1 = 1. For t > 1 the physical root satisfies τ < t and only a fraction of the nodes resides in finite components, M 1 < 1. Thus, at time t = 1, the system undergoes a gelation transition with a finite fraction of the nodes contained in infinite components. We term this time the gelation time, t g = 1. In the late stages of the evolution t ≫ 1, one has τ ≃ te −t and M 1 ≃ c 1 = e −t , so the system consists of a single giant component and a small number of isolated nodes. The behavior at and near the transition point are of special interest. The critical behavior of the component size distribution is echoed by other quantities as will be shown below. Size distributions become algebraic near the critical point. Moreover, there is a self-similar behavior as a function of time (dynamical scaling) and as a function of the system size (finite-size scaling).

B. Finite Random Graphs
In the previous subsection, we applied kinetic theory to an infinite system. This approach can be extended to finite systems. Unfortunately, such treatments are very cumbersome [22,23]. Since the number of components is finite, the fluctuations are no longer negligible, and instead of a deterministic rate equation approach, a stochastic approach is needed. Here we follow an alternative path, employing the exact infinite system results in conjunction with scaling and extreme statistics arguments.
The characteristic size of components at the gelation point exhibits nontrivial dependence on the system size. This is conveniently seen via the cumulative size distribution. The size of the largest component in the system, k g , is estimated from the extreme statistics criterion, N k≥kg c k ∼ 1, to be The largest component in the system grows sub-linearly with the system size [3]. The time by which this component emerges approaches unity for large enough systems as follows from the diverging characteristic size scale The maximal component size (11) The scaling function has the following extremal behaviors The small-ξ behavior corresponds to sizes well below the characteristic size and thus reflects the infinite system behavior (9). The large-ξ behavior was obtained numerically with γ ∼ = 3. To appreciate the large-ξ asymptotic, let us estimate the probability that the system managed to generate the largest possible component of size N/2 at time t = 1. The lower bound for this probability can be established via a "greedy" evolution which assumes that after k linking events the graph is composed of a tree of size k + 1 and N − k − 1 disconnected nodes. Such evolution occurs with probability that scales as e −N . While this lower bound is not necessarily optimal, it suggests that the actual probability is exponentially small. The scaling variable To check the critical behavior in finite systems, we performed numerical simulations. In the simulations, N/2 links are placed randomly and sequentially among the N nodes as follows. A node is drawn randomly, and then another node is drawn randomly. Last, these two nodes are linked. Self-connections are therefore allowed. The simulations differ slightly from the above random graph model in that the number of links is not a stochastic variable. For large N , this simulation is faithful to the evolving random graph model because the number of links is self-averaging. The simulation results are consistent with the postulated finite-size scaling form (13). We note that the scaling function Ψ c (ξ) converges slowly as a function of N . The simulations reveal an interesting behavior of the finite-size scaling function. The function c k (N ) has a "shoulder" -a non-monotonic behavior compared with the pure algebraic behavior (9) characterizing infinite systems (Fig. 2). The properly normalized scaling function (2πξ 5 ) 1/2 Ψ c (ξ) is a non-monotonic function of ξ (Fig. 3). Obtaining the full functional form of the scaling function Ψ c (ξ) remains a challenge. A very similar shoulder has been observed for the degree distribution of finite random networks generated by preferential attachment [24,25,26,27].

IV. PATHS
Structural characteristics of components can be investigated in a similar fashion. By definition, every two nodes in a component are connected. In other words, there is a path consisting of adjacent links between two such nodes. We investigate statistical properties of paths in components. Characterization of paths yields useful information regarding the connectivity of components as well as internal structures such as cycles.
For every node in the graph, there are (generally) multiple paths that connect it with all other nodes in the respective component. With new links, new paths are formed. For every pair of paths of lengths n and m originating at two separate nodes, a new path is formed as follows n, m −→ n + m + 1.
In Fig. 1, linking two paths of respective lengths n = 1 and m = 2 generates a path of length n + m + 1 = 4. Thus, paths also undergo an aggregation process. However, this aggregation process is simpler than (1) because the aggregation rate is independent of the path length. Let q l (t) be the density of distinct paths containing l links at time t. By distinct we mean that the two paths connecting two nodes are counted separately. By definition, q 0 (t) = 1. The rest of the densities grow according to the rate equation for l > 0. The initial condition is q l (0) = δ l,0 . This rate equation reflects the uniform aggregation rate. Another notable feature is the lack of a loss term -once a path is created, it remains forever. Solving recursively gives q 1 = t, q 2 = t 2 , etc. By induction, the path length density is Indeed, this expression satisfies both the rate equation and the initial condition. The first quantity q 1 = t is consistent with the facts that the link density is equal to t/2 and that every link corresponds to two distinct paths of length one. The above path density represents an aggregate over all nodes and all components. Characterization of path statistics in a component of a given size is achieved via p l,k , the density of paths of length l in components of size k. Note the obvious length bounds 0 ≤ l ≤ k − 1 and the sum rule l p l,k = k 2 c k reflecting that there are k 2 distinct paths in a component of size k (every pair of nodes is connected). The density of the linkless paths is p 0,k = kc k , because kc k is the probability that a node belongs to a component of size k.
We have seen that components and paths form via the aggregation processes (1) and (15), respectively.
The joint distribution p l,k therefore undergoes a biaggregation process [28]. In the present case, where the first index corresponds to the path length and the second to the component size. The joint distribution evolves according to the rate equation The initial conditions are p l,k (0) = δ k,1 δ l,0 . The first term on the right-hand side of Eq. (19) describes newly formed paths due to linking. The last two terms correspond to paths that do not contain the newly placed link.
We now repeat the steps used to determine the size distribution. The time dependence is eliminated using the ansatz p l,k = P l,k t k−1 e −kt . The corresponding coefficients P l,k satisfy the recursion The generating function P l (z) = k P l,k e kz satisfies the recursion relation (1 − G) dP l dz = n+m=l−1 P n P m + P l for l > 0. Dividing this equation by (4) yields for l > 0. As noted above P 0,k = kC k , so P 0 (z) = G(z). Solving Eq. (21) recursively gives P 1 = G 2 , P 2 = G 3 , etc. In general, This solution can be validated directly. The time dependent generating function p l (z) = k p l,k e kz is therefore p l (z) = t −1 G l+1 (z +ln t−t). The total density of paths of length l, p l (z = 0) = t l , coincides with (17) prior to the gelation transition (t < 1) because all components are finite. However, the total number of paths is reduced, p l (z = 0) = t −1 τ l+1 , past the gelation time (t > 1). One may also obtain the bivariate generating function p(z, w) = l,k p l,k w l e kz . Using (22) one gets The total density of paths in finite components is of course g = l,k p l,k , so g ≡ p(z = 0, w = 1). Generally, g = τ t(1−τ ) ; for t < 1 the total density of paths is The coefficients are found via the contour integration P l,k = (2πi) −1 dy P l y −k−1 (see Appendix A). Substituting r = l + 1 in Eq. (A1) yields P l,k = (l + 1) k k−l−2 (k−l−1)! .
As a result, the density of paths of length l in components of size k is Comparing (24) and (6) we notice that the densities of the two shortest paths satisfy p 0,k = kc k and p 1,k = 2(k − 1)c k . The latter reflects that there are k − 1 links in a tree of size k and that with unit probability all components are trees (as discussed in the next section). Note also that the longest possible path, l = k − 1, corresponds to linear (chain-like) components. According to Eq. (24), the density of such paths is p k−1,k = t k−1 e −kt . This density decays exponentially with length, so these components are typically small, their length being of the order one.
The path length density can be simplified in the large k-limit by considering the properly normalized ratio of factorials Using the Stirling formula, in the limits k ≫ 1 and l ≫ 1, the path density becomes As was the case for the component size distribution, the path length density is self-similar in the vicinity of the gelation point, Thus, the characteristic path length diverges near the gelation point, l ∼ (1 − t) −1 . At the critical point, the path length density becomes It is evident that the typical path length scales as square root of the component size For finite systems, the scaling law for the typical path length (28) combined with the characteristic component size (11) leads to the following characteristic path length One can deduce several other scaling laws and finite-size scaling functions underlying the path density. For example, substituting the gelation time 1 − t g ∼ N −1/3 into the total number of paths g = (1 − t) −1 yields g ∼ N 1/3 .

V. CYCLES
Each component has a certain number of nodes and links. The complexity of a component is defined as the number of links minus the number of nodes. Components with complexity −1 are trees; components with complexity 0 and 1 are termed unicyclic and bicyclic correspondingly. Finite components are predominantly trees. We have seen that the overall number of links is proportional to N and that the overall the number of self-links is of the order unity. The overall numbers of trees and of unicyclic components mirror this behavior. Generally, the number of components of complexity R is proportional to N −R (this result is well-known, see e.g. [5,21] and especially [29]). Therefore, it suffices to characterize trees and unicyclic components only.
Each unicyclic component contains a single cycle. Cycles are an important characteristic of a graph [30,31]. In this section, we analyze cycles and unicyclic components using the rate equation approach. We first note that cycles in random graphs were also studied using various other approaches: Janson [32,33] employs probabilistic and combinatorial techniques; Marinari and Monasson [31] assign an Ising spin to each node and deduce certain properties of loops from the partition function of the Ising model; Burda et al [34] modify a random graph model to favor the creation of short cycles, and examine the model using a diagrammatic technique. A number of authors also studied cycles on information networks like the Internet (see [35] and references therein).

A. Infinite System
There is a significant difference between the distribution of trees and unicyclic components. In the thermodynamic limit, the number of trees is extensive and as a result, it is a deterministic, or a self-averaging quantity. The number of unicyclic components is not extensive, but rather of the order unity; as a result it is a random quantity with a nontrivial distribution even for infinite random graphs. In what follows, we study the average number of unicyclic components of a given size or cycle length.
The average number of cycles follows directly from the path length density. Quite simply, when the two extremal nodes in a path are linked, a cycle is born. Let the number of cycles of size l at time t be w l (t). It grows according to the rate equation The right-hand side equals the link creation rate 1/(2N ) times the total number of paths N q l−1 ; indeed, the total number of cycles of a given length is of the order one. The cycle length distribution is In particular, at the gelation point, the cycle length distribution is inversely proportional to the cycle length [5] w l = (2l) −1 .
This result can alternatively be obtained using combinatorics.
To characterize cycles in a given component size, we consider the joint distribution u l,k , the average number of unicyclic components of size k containing a cycle of length l with 1 ≤ l ≤ k. This joint distribution evolves according to the linear rate equation for l ≥ 1. Initially there are no cycles, and therefore u l,k (0) = 0. Eliminating the time dependence via the substitution u l,k = U l,k t k e −kt , the coefficients satisfy the recursion Using the generating function U l (z) = k e kz U l,k this recursion is recast into the differential equation Dividing by (4), we obtain Integrating this equation yields the generating function Consequently, the cycle length distribution (in finite components only) is p l = τ l 2l , in agreement with (31) prior to the gelation time (t < 1).
Additionally, the joint generating function defined as u(z, w) = l,k e kz w l u l,k is given by As for paths, statistics of cycles are directly coupled to statistics of components via the generating function G(z). The total number of unicyclic components of finite-size h = l,k u l,k is therefore Below the gelation point, h(t) = 1 2 ln 1 1−t , for t < 1. The total number of unicyclic components can alternatively be obtained by noting that (i) it satisfies the rate equation dh/dt = 1 2 k k 2 c k = 1 2 M 2 , and (ii) the second moment of the size distribution is M 2 = (1 − t) −1 for t < 1 as follows from (7).
The coefficients underlying the cycle distribution are found using contour integration. Indeed, writing U l,k = (2πi) −1 U l y −k−1 dy and substituting r = l in (A1) gives U l,k = 1 2 k k−l−1 (k−l)! [4]. The cycle length-size distribution is therefore The smallest cycle, l = 1, is a self-connection, and the average number of such cycles is u 1,k = t 2 kc k . The largest cycles are rings, l = k, and their total number is on average u k,k = 1 2k t k e −kt . The large-k behavior of the cycle length distribution is found following the same steps leading to (25) This distribution is self-similar in the vicinity of the gelation transition, u l, We see that the cycle length is characterized by the same scale as the path length, l ∼ (1 − t) −1 . At the gelation point, the distribution is Fixing the component size, the typical cycle length behaves as the typical path length, l ∼ k 1/2 . The size distribution of unicyclic components is found from the joint distribution v k = l u l,k . Using (39) we get [21] v k (t) = 1 2 This distribution can alternatively be derived from the linear rate equation This equation is obtained from (33) using the equality k 2 c k = l p l,k . It reflects that linking a pair of nodes in a component generates a unicyclic component. Integrating (42) over the cycle length, the critical size distribution of unicyclic components has an algebraic tail

B. Finite Systems
We turn now to finite systems, restricting our attention to the gelation point. The total number of unicyclic components is obtained by estimating h (N, t g ). Substituting (12) into (38) shows that the average number of unicyclic components (and hence, cycles) grows logarithmically with the system size ( Fig. 4) h(N ) ≃ 1 6 ln N. Comparing the path length distribution (27) and the cycle length distribution (42), we conclude that the characteristic cycle length and the characteristic path length obey the same scaling law, l ∼ N 1/3 . This implies that the cycle length distribution in a finite system of size N , w l (N ), obeys the finite-size scaling law Numerical simulations confirm this behavior (Fig. 5).
In the simulations, analysis of cycle statistics requires us to keep track of all links. Cycles are conveniently identified using the standard "shaving" algorithm. Dangling links, i.e., links involving a single-link node are removed from the system sequentially. The link removal procedure is carried until no dangling links remain. At this stage, the system contains no trees. Simple cycles are those components with an equal number of links and nodes.
The extremal behaviors of the finite-size scaling function are as follows The small-η behavior follows from (32). Shown is 2ηΨw(η) versus η obtained using systems with size N = 10 4 , 10 5 , and 10 6 . The data represents an average over 10 6 independent realizations.
links, the largest possible cycle has length l = N/2. Its likelihood w(n, 2n) is obtained using combinatorics There are 2n n ways to choose the nodes participating in the cycle and the next term is the number of ways to arrange them in a cycle. The corrective factor 2n accounts for rotation and reflection symmetries. The last term is the probability that each pair of consecutive nodes are linked. The large-n asymptotic behavior is Therefore, w(n, 2n) ∼ exp(−C N ). Substituting l ∼ N into the scaling form (47) leads to the super-exponential behavior Ψ w (η) ∼ exp(−Cη 3/2 ), see Fig. 6. Typically, cycles are of size N 1/3 . The average moments l(N ) = l l w l (N )/ l w l (N ) reflect this law. However, the algebraic divergence, w l ∼ l −1 , leads to a logarithmic correction as follows from (46)-(48): The behavior of the average cycle length is verified numerically (Fig. 7).
Finite-size scaling of other cycle statistics such as the joint distribution can be constructed following the same procedure. For example, the size distribution of unicyclic components should follow the scaling form The scaling function diverges Ψ v (ξ) ≃ (4ξ) −1 for ξ → 0.

VI. THE FIRST CYCLE
The above statistical analysis of cycles characterizes the average behavior but not necessarily the typical one because the number of cycles is a fluctuating quantity. There are numerous interesting features concerning cycles that are not captured by the average number of cycles. For instance, what is the probability that the system does not contain a cycle up to time t? It suffices to answer this question in the pre-gel regime as the giant component certainly contains cycles.
Let s 0 (t) be the (survival) probability that the system does not contain a cycle at time t. The cycle production rate is J = dh dt = 1 2(1−t) . The number of cycles is finite in the pre-gel regime, since cycles are independent of each other in the N → ∞ limit. This assertion (supported by numerical simulations, see Fig. 8) implies that the cycle production process is completely random. The cycle production rate characterizes the survival probability s 0 as follows The initial condition is s 0 (0) = 1. As a result, the survival probability is for t ≤ 1. The survival probability vanishes beyond the gelation point, s 0 (t) = 0 for t > 1. This reiterates that in the thermodynamic limit, a cycle is certain to form prior to the gelation transition [5].
Since the number of cycles produced is of the order of one in the pre-gel regime, one may expect that statistical properties of cycles strongly depend on their generation number or alternatively on their creation time. This is manifested by the first cycle. The quantity dt s 0 dw l dt is the probability that (i) the system contains no cycles at time t, (ii) a cycles is produced during the time interval (t, t + dt), and (iii) its length is l. Summing these probabilities gives the probability that the first cycle produced sometimes during the pre-gel regime has length l: Summing these quantities, we verify the normalization The length distribution of the first cycle can be expressed in terms of the beta function f l = 1 2 B(3/2, l) or alternatively The probability distribution f l has an algebraic tail, with C = √ π 4 for l ≫ 1. The tail exponent characterizing the distribution of the first cycle is larger compared with the exponent characterizing all cycles, reflecting the fact that the first cycle is created earlier.
Similarly, one can obtain additional properties of the first cycle. We mention the probability F k that the first unicyclic component has size k, This integral can be expressed in terms of the confluent hypergeometric function. Its asymptotic behavior can be readily found by noting that the integrand has a sharp maximum in the region 1 − t ∼ k −1/2 leading to I k ≃ 2 −1/4 Γ(3/4)k −3/4 e −k . Using this in conjunction with the Stirling's formula, the size distribution has the algebraic tail with C = 2 −7/4 π −1/2 Γ(3/4) for k ≫ 1. Under the assumption that cycle production is completely random, the number of cycles obeys Poisson statistics. The probability that there are n cycles, s n , then satisfies the straightforward generalization of Eq. (53), viz. dsn dt = J[s n−1 − s n ] with the initial condition s n (0) = δ n,0 . The solution is the Poisson distribution s n = h n n! e −h , see Fig. 8. Explicitly, the distribution reads (60) The cumulative distribution S n (t) = s 0 (t) + . . . + s n (t) is plotted in Fig. 9. The Poisson distribution (60) can also be used to calculate f n,l the size distribution of the nth cycle. We merely quote the large-l tail behavior Indeed, summation over the cycle generation reproduces the overall cycle distribution (32).
In finite systems, it is possible that no cycle are created at the gelation time. This probability decreases algebraically with the system size, as seen by substituting (12) into (54) This prediction agrees with simulations, see Fig. 10. In practice, this slow decay indicates that a relatively large system may contain no cycles after N/2 links are placed. Generally, the probability that there is a finite number of cycles increases with the number of cycles The length distribution of the first cycle is characterized by the same l ∼ N 1/3 size scale as does the overall cycle distribution. We focus on the behavior of the moments This behavior is obtained from the distribution (57) that should be integrated up to the appropriate cutoff, i.e., l n ∼ N 1/3 1 dl l n l −3/2 . As a result, the average size of the first cycle is much smaller than the characteristic cycle size l ∼ N 1/6 . Moments corresponding to the size of the first unicyclic component grow as follows as obtained from (59). Consequently, the average size of the first unicyclic component is smaller than the characteristic component size, k ∼ N 1/2 .

VII. CONCLUSIONS
In summary, we have extended the kinetic theory description of random graphs to structures such as paths and cycles. Modeling the linking process dynamically leads to an aggregation process for both components and paths. The density of paths in finite components is coupled to the component size distribution via nonlinear rate equations while the average number of cycles is coupled to the path density via linear rate equations. Both path and cycle length distributions are coupled to the component size distribution.
Generally, size distributions decay exponentially away from the gelation point, but at the gelation time, algebraic tails emerge. As the system approaches this critical point, the size distributions follow a self-similar behavior characterized by diverging size scales.
The kinetic theory approach is well-suited for treating infinite systems. The complementary behavior for finite systems can be obtained from heuristic scaling arguments. This approach yields scaling laws for the typical component size, path length, and cycle length at the gelation point. These scaling laws can be formalized using finite-size scaling forms, i.e., self-similarity as a function of the system size, rather than time. Obtaining the exact form of these scaling functions is a nice challenge in particular for the most fundamental quantity, the component size distribution that is characterized by a non-monotonic scaling function.
The kinetic theory approach seems artificial at first sight. Indeed, graphs are discrete in nature and therefore combinatorial approaches appear more natural. Yet, once the rate equations are formulated, the analysis is straightforward. Utilizing the continuous time variable allows us to employ powerful analysis tools. Moreover, some of the kinetic theory results are less cumbersome compared with the combinatorial results.
The same methodology can be expanded to analyze other features of random graphs. For example, correlations between the node degree and the cluster size can be analyzed using bi-aggregation rate equations [36]. It is quite possible that structural properties in other aggregation processes, for example, polymerization with a sum kernel [17], and in other variants of random graphs such as small-world networks [37] can be analyzed using kinetic theory.
One could try to utilize kinetic theory to probe the distribution of various families of subgraphs. We have limited ourselves to cycles since they, alongside with trees, do appear in random graphs while more interconnected families of subgraphs are very rear [29]. Yet in biological and technological networks certain interconnected families of subgraphs do appear. Such populated families of subgraphs, motifs, are believed to carry information processing functions [38,39]. It will be interesting to use kinetic theory to analyze motifs in special random graphs.
This research was supported by the DOE (W-7405-ENG-36).

APPENDIX A: CONTOUR INTEGRATION
Let A(z) = k A k e kz be the generating function of the coefficients A k . For the family of generating func-tions A(z) = G r (z) with G(z) satisfying Ge −G = e z , the coefficients A k can be obtained via contour integration in the complex y plane where y = e z as follows Since Ge −G = e z , it is convenient to perform the integration in the complex G plane. In writing the third line, we used dy dG = (1 − G)e −G .