Ensemble nonequivalence and Bose-Einstein condensation in weighted networks

The two ensembles can be made conjugate to each other by choosing a specific value $E^{*}$ and simultaneously setting the total energy $E(\mathbf{W})$ of each realized configuration $\mathbf{W}$ in the microcanonical ensemble equal to $E(\mathbf{W})\equiv E^{*}$ and the inverse temperature $\beta$ in the canonical ensemble to the corresponding value $\beta\equiv\beta^{*}$ such that the resulting average value $\langle E\rangle_{\beta^{*}}$ of the fluctuating total energy under the canonical probability $P_{\textrm{can}}(\mathbf{W}|\beta^{*})$ equals $E^{*}$, i.e.

$$\langle E\rangle_{\beta^{*}}\equiv E^{*}.$$

(2)

For systems with finite size, the two conjugate ensembles are unavoidably different, because in the microcanonical ensemble the hardness of the constraint implies extra dependencies among the state of the microscopic constituents with respect to the canonical case. However, in the thermodynamic limit (i.e. when the number of units in the system goes to infinite) and under certain ‘natural’ circumstances, the two associated probabilistic descriptions are expected to become effectively equivalent (i.e. the canonical fluctuations and the microcanonical dependencies are both expected to play an asymptotically vanishing role) as a result of some form of the law of large numbers. This idea, which dates back to Gibbs himself Gibbs (1902) and has continued to attract a lot of interest until presently Campa et al. (2009); Touchette (2015); Squartini et al. (2015a); Squartini and Garlaschelli (2017a), goes under the name of ensemble equivalence (EE). When EE holds, one can treat the two ensembles as asymptotically interchangeable, and hence use any of them based on mathematical or computational convenience. For instance, analytical calculations are signficantly easier in the canonical ensemble, while numerical randomizations of an initial configuration can be carried out more naturally in the microcanonical ensemble.

Most statistical physics textbooks still convey the message that EE is expected to hold in general as a sort of principle at the basis of ensemble theory. The possible breakdown of EE, also known as ensemble nonequivalence (EN), is still not discussed systematically in the literature. However, several observations of EN have been documented over the past decades Ellis et al. (2000); Blume et al. (1971); Barré et al. (2001); Ellis et al. (2004); Lynden-Bell (1999); Chavanis (2003); d’Agostino et al. (2000); Barré and Gonçalves (2007); Radin and Sadun (2013); Ellis et al. (2002); Kastner (2010); Campa et al. (2009); Squartini et al. (2015a). These observations motivated various efforts aimed at elucidating both the possible physical mechanisms at the origin of EN and, in parallel, its proper mathematical definition(s).

Traditionally, the ‘natural’ circumstances generally invoked to ensure EE mainly concern the presence of (loosely speaking) ‘at most weak’ interactions between the constituents of the system. This condition is automatically realized when the system consists of independent units or units with short-range interactions and sufficiently high temperature (to stay away from possible low-temperature phases with broken symmetries, for which the canonical average value of the energy is no longer the typical value). Indeed, violations of EE have been documented in presence of long-range interactions (e.g. in gravitational systems) or phase transitions (e.g. in interacting spin systems) Ellis et al. (2000); Blume et al. (1971); Barré et al. (2001); Ellis et al. (2004); Lynden-Bell (1999); Chavanis (2003); d’Agostino et al. (2000); Barré and Gonçalves (2007); Radin and Sadun (2013); Campa et al. (2009); Touchette (2015).

However, recent research on complex systems encountered beyond the usual realm of physics has found an additional mechanism that can break EE, even in presence of weak (or no) interactions: namely, the presence of an extensive number of local constraints Squartini et al. (2015a); Squartini and Garlaschelli (2017a). This situation is frequently found in networks with constraints on the number of links (degree) of each node. More specifically, the binary configuration model Squartini and Garlaschelli (2017b) is a widely used null model of graphs with a given degree sequence, i.e. a given vector of node degrees. The model captures many properties found in real-world networks, because the local character of the degree constraint can accomodate the strong structural heterogeneity typically observed across nodes in real networks. Unlike the traditional thermodynamic example where the total energy (and possibly a small, finite number of additional macroscopic properties) is a global and unique constraint for the system, networks with given node degrees are systems with as many constraints as the number of fundamental units, i.e. where constraints are extensive in number and local in nature. This situation has been found to break the equivalence of the corresponding canonical and microcanonical ensembles, even without long-range interactions or phase transitions Squartini et al. (2015a); Garlaschelli et al. (2016); Roccaverde (2019). Notably, since systems with local constraints are generic models for virtually any heterogenous system, the new mechanism significantly widens the range of real-world cases for which EE may break down. This novel result deserves further research.

Besides the aforementioned advances in the study of the physical mechanisms at the origin of EN, significant progress has been made in the mathematical characterization of EN as well. Traditionally, the informal criterion Gibbs (1902) used to test whether two conjugate ensembles are equivalent is checking whether the relative fluctuations of the constraint in the canonical ensemble, e.g. the ratio $\sigma_{\beta^{*}}(E)/E^{*}$ of the canonical standard deviation to the average value of the energy, vanish in the thermodynamic limit. If this happens, then the canonical fluctuations of the total energy are negligible with respect to the total energy itself and, intuitively, the energy in the canonical ensemble can be though of as an effectively deterministic quantity, very much like in the conjugated microcanonical ensemble. Similarly, the extra dependencies among the microscopic constituents in the microcanonical ensemble are expected to play a smaller and smaller role, thus coming closer to the canonical case.

More recent approaches have considered different rigorous definitions of EE, which can be beautifully summarized Touchette (2015) as the following three notions: thermodynamic equivalence (convexity of the microcanonical entropy density), macrostate equivalence (equality of the expected values of macroscopic quantities under the two ensembles), and measure equivalence (vanishing of the relative entropy density between the microcanonical and canonical probability distributions). Under mild conditions, these notions have been shown to be equivalent Touchette (2015). In this paper, we use measure equivalence as it is more transparently related to the ensemble probabilities.

As a useful result, research on the relationship between statistical physics and combinatorics has revealed that the relative entropy between the microcanonical and the canonical probability distributions is, under certain conditions, asymptotically proportional to the logarithm of the determinant of the covariance matrix of the effective constraints in the canonical ensemble Squartini and Garlaschelli (2017a). The effective constraints are those that are neither redundant, i.e. trivially replicating other constraints, nor degenerate, i.e. deterministically restricting the canonical and microcanonical configurations in exactly the same way. For instance, formally imposing ‘two’ constraints where one is the total energy $E$ and the other one is twice the total energy $2E$ is clearly a redundant choice: the effective number of constraints is just one in this case. As another example, if in addition to the energy $E$ we impose its square value $E^{2}$, then for those values of the Lagrange multipliers such that $\langle E^{2}\rangle=\langle E\rangle^{2}$ the variance of the energy will be zero also in the canonical ensemble: $E$ will become degenerate and deterministically equal to its imposed value in both ensembles, so not contributing any difference between the two (by contrast, for parameter values such that $\langle E^{2}\rangle>\langle E\rangle^{2}$ there are no allowed configurations in the microcanonical ensemble because the hard values of $E^{2}$ and $E$ become conflicting, thereby breaking the equivalence with the canonical one). In general, if the problem is not ill-posed from the beginning, the number of effective constraints coincides with the number of enforced constraints. However, it may happen that some constraints become ineffective for certain degenerate values of the parameters. In any case, for a given parameter value the number of effective constraints coincides with the rank of the covariance matrix of all imposed constraints Squartini and Garlaschelli (2017a).

Since nonequivalence in the measure sense corresponds to the (super)extensivity of the relative entropy, studying the asymptotic behaviour of the determinant of the (effective) covariance matrix is enough in order to assess ensemble nonequivalence. It is worth noticing that, if there is a single constraint (say, the total energy $E$), then the determinant of the covariance matrix coincides with the corresponding variance $\sigma^{2}_{\beta^{*}}(E)$ and the relative entropy grows asymptotically (under the necessary hypotheses) as $\ln\sigma^{2}_{\beta^{*}}(E)=2\ln\sigma_{\beta^{*}}(E)$. On the other hand, since $E$ is a global constraint, it is generally extensive in the number $n$ of units of the system. Therefore the vanishing of the relative fluctuations, i.e. $\sigma_{\beta^{*}}(E)=o(E)$, implies that the relative entropy is subextensive, i.e. the relative entropy density vanishes in the thermodynamic limit. This suggests that, in presence of a global constraint, the vanishing of the relative fluctuations implies ensemble equivalence.¹¹1Note that the converse is not necessarily true. However, observing ensemble equivalence and non-vanishing relative fluctuations simultaneously requires some rather uncommon circumstances: for instance, if $\sigma_{\beta^{*}}(E)$ grows like $E^{\alpha}$ with $\alpha\geq 1$ while the entropy still grows like $E$, then the relative entropy is still subextensive, while the relative fluctuations do not vanish. How this picture changes in presence of an extensive number of local constraints has not been investigated yet. In particular, whether the vanishing of relative fluctuations still implies ensemble equivalence remains an open question.

This paper connects to both lines of research described above (physical mechanisms and mathematical definitions for EN) and its aim is therefore twofold. On the one hand, we aim at investigating for the first (to the best of our knowledge) time the phenomenon of EN in a model system that combines the presence of an extensive number of local constraints with a phase transition. On the other hand, we aim at understanding whether the intuitive criterion of vanishing relative fluctuations of the constraints still ensures EE in this more general setting.

Concretely, we consider the weighted configuration model Squartini and Garlaschelli (2017b), namely a model of weighted (as opposed to binary) networks with given strength sequence, i.e. with given values of the strength (sum of the weights of incident links) of each node. The weighted character of the model allows for the emergence of a phase transition that is impossible to observe in the corresponding binary configuration model, namely Bose-Einstein Condensation (BEC) Park and Newman (2004); Garlaschelli and Loffredo (2009). For the sake of clarity, it is worth mentioning here that, although a form of BEC in networks was identified for the first time in growing binary graphs Bianconi and Barabási (2001), the notion we refer to here refers to static networks and as such can only occur in weighted networks Park and Newman (2004). Indeed, while the configuration model for weighted networks obeys Bose-Einstein statistics, the configuration model for binary networks obeys Fermi-Dirac statistics Park and Newman (2004); Garlaschelli and Loffredo (2009); Garlaschelli et al. (2013). BEC can arise in our model by appropriately tuning the strength sequence. In particular, we are going to show that we can make the strength sequence temperature-dependent and generate BEC by picking a sufficiently low temperature, below a certain critical value. The simplest such setting is one where the network has a ‘core-periphery’ structure, with BEC appearing in the core.

We find that, for all temperatures and irrespective of whether BEC emerges, the canonical and microcanonical ensembles are always nonequivalent as signalled by a nonvanishing relative entropy density. On the other hand, the relative fluctuations of all the constraints vanish when BEC is absent, while some of them do not vanish when BEC is present. This shows that the relative fluctuations cannot distinguish between equivalence and nonequivalence of the ensembles in this more general case where multiple constraints are present. In fact, what they do is detecting the presence of BEC. Therefore the traditional criterion for EE based on the vanishing of the relative fluctuations is no longer valid in presence of an extensive number of constraints, even when applied simultaneously to all constraints. These results enrich our understanding of the phenomenology of EN and shed more light on its relationship with both the extensivity of the constraints and the presence of phase transitions.

The remainder of this paper is organized as follows. In Sec. II we rigorously define the canonical and microcanonical ensembles of weighted networks with given strength sequence. In Sec. III we introduce two criteria for the (non)equivalence of the ensembles, one based on the relative entropy between the corresponding probability distributions (measure equivalence) and one based on the relative fluctuations of the constraints. In Sec. IV we study in detail a model defined by the simplest family of strength sequences, driven by a temperature parameter, such that we can observe both a BEC and a non-BEC phase. In Sec. V we offer our conclusions. Finally, the Appendix contains useful calculations needed to establish the scaling of the relative entropy in all the regimes considered.

Weighted networks are widely used to describe systems with a large number of components and heterogeneous patterns of interaction Newman et al. (2006). We represent a possible configuration of a weighted network with $n$ nodes as an $n\times n$ weighted adjacency matrix $\mathbf{W}$. Each entry of the matrix $w_{ij}$ ($1\leq i\leq n,1\leq j\leq n)$ denotes the weight of the link between node $i$ and node $j$, which is taken from the set $\mathbb{N}$ of natural numbers (including zero, which corresponds to the absence of a link between $i$ and $j$). In this work, we only consider undirected networks without self-loops, thus the weighted matrix $\mathbf{W}$ is a symmetric matrix ($w_{ij}=w_{ji}$ for all $i,j$) and its diagonal elements are zero ($w_{ii}=0$ for all $i$). The number of independent entries of each such matrix is therefore ${n\choose 2}=n(n-1)/2$.

An ensemble of weighted networks on $n$ nodes is the discrete (infinite) set $\mathcal{W}_{n}=\mathbb{N}^{n(n-1)/2}$ of all available configurations for the matrix $\mathbf{W}$ and a probability distribution $P(\mathbf{W})$ over $\mathcal{W}_{n}$ that is specified by a given vector $\vec{C}(\mathbf{W})$ of constraints, which can be enforced either as a soft constraint (canonical ensemble) or as a hard constraint (microcanonical ensemble) Squartini and Garlaschelli (2017b). So the matrix $\mathbf{W}$ is a possible outcome of a random variable. We will consider the weighted configuration model, for which the constraints are the strengths of all nodes, i.e. the strength sequence $\vec{C}(\mathbf{W})=\vec{s}(\mathbf{W})$, where the strength $s_{i}(\mathbf{W})$ of node $i$ is a local sum of all the link weights that connect $i$ to its neighbours in the particular network $\mathbf{W}$:

$$s_{i}(\mathbf{W})=\sum_{j=1}^{n}w_{ij},\qquad i=1,n.$$

(3)

Clearly, the number of scalar constraints is $n$, which coincides with the number of nodes, so this model is a perfect example of a system subject to an extensive number of local constraints. In fact, it is the weighted counterpart of the binary configuration model, where EN driven by local constraints was observed for the first time Squartini et al. (2015a).

A crucial consequence of the presence of local constraints in weighted networks is illustrated in Fig.1, where we show an example of different networks with the same strength sequence. For a given choice of $\vec{s}^{*}$, in different realizations of the network each node $i$ can have different neighbours and different distributions of weights on the links that connect it to those neighbours. In particular, the strength $s_{i}^{*}$ can be more or less concentrated on specific neighbours (a property that is usually quantified by the so-called disparity Serrano et al. (2006)). However, more homogeneous choices of $\vec{s}^{*}$ unavoidably result in less concentrated link weights, while more heterogeneous choices of $\vec{s}^{*}$ impose more concentrated link weights. This fact will allow us to consider (in Sec. IV) different structural regimes ranging between two extreme limits: a constant (infinite-temperature) strength sequence implying that on average each node is connected to its neighbours in an equally strong way, and a ‘step-like’ (zero-temperature) strength sequence implying an extreme concentration of link weights among a small subset of the $n$ nodes. In between these two limits, a certain critical temperature separates a ‘non-condensed’ (high-temperature) phase from a ‘condensed’ (low-temperature) phase featuring the properties of BEC.

We first discuss how to implement the strength sequence constraint mathematically in the canonical ensemble (soft constraint) Squartini and Garlaschelli (2017b). Recall that in the traditional canonical ensemble the inverse temperature $\beta^{*}$ is the only (scalar) parameter of the canonical probability distribution $P_{\textrm{can}}(\mathbf{W}|{\beta}^{*})$, conjugate to a certain (scalar) total energy $E^{*}$ in the corresponding microcanonical ensemble. Explicitly, $P_{\textrm{can}}(\mathbf{W}|{\beta}^{*})$ is the Boltzmann distribution $P_{\textrm{can}}(\mathbf{W}|{\beta}^{*})=e^{-\beta^{*}E(\mathbf{W})}/Z(\beta^{*})$ with inverse temperature $\beta^{*}=1/kT^{*}$. By contrast, in our setting the canonical distribution $P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})$ has to depend on an $n$-dimensional vector $\vec{\beta}^{*}$ of parameters, conjugate to the $n$-dimensional constraint $\vec{s}^{*}$ which, in turn, defines the microcanonical ensemble. This distribution is found by maximizing the Gibbs-Shannon entropy functional defined in Eq. (1) under the soft constraint

$$\langle\vec{s}\rangle_{\vec{\beta}^{*}}=\sum_{\mathbf{W}\in\mathcal{W}_{n}}P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})\vec{s}(\mathbf{W})\equiv\vec{s}^{*}$$

(4)

which generalizes the conjugacy condition in Eq. (2). The solution to the maximization problem sees $\vec{\beta}^{*}$ play the role of a vector of Lagrange multipliers coupled to the strength sequence $\vec{s}^{*}$, and is given by Squartini and Garlaschelli (2017b)

$$P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})=\frac{e^{-H(\mathbf{W},\vec{\beta}^{*})}}{Z(\vec{\beta}^{*})},$$

(5)

where the (network) Hamiltonian $H(\mathbf{W},\vec{\beta})=\vec{\beta}\cdot\vec{s}(\mathbf{W})$ is the linear combination of the node strengths, the partition function $Z(\vec{\beta})=\sum_{\mathbf{W}\in\mathcal{W}_{n}}e^{-H(\mathbf{W},\vec{\beta})}$ is a normalization constant, and $\vec{\beta}^{*}$ is the unique parameter value realizing Eq. (4). Note that Eq. (5) has still the form of the Boltzmann distribution, with the important caution that the inverse temperature has been reabsorbed into the Hamiltonian. Therefore, to keep the parallel with the traditional physical situation, in our setting the Hamiltonian should be thought of as the inverse temperature times the energy, and $\vec{\beta}^{*}$ as the inverse temperature times a vector of ‘fields’, each coupled to a different constraint. Clearly, since the probability in Eq. (5) must be dimensionless, the product $\vec{\beta}\cdot\vec{s}(\mathbf{W})$ must be dimensionless as well. In Sec. IV, we will notice that the Hamiltonian can be further reinterpreted as also incorporating a ‘chemical potential’ governing the expected weight of the links in the network Garlaschelli et al. (2013).

Notably, $P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})$ depends on $\mathbf{W}$ only through $\vec{s}(\mathbf{W})$. In particular, it gives the same value $P_{\textrm{can}}(\mathbf{W}^{*}|\vec{\beta}^{*})$ to any network $\mathbf{W}^{*}$ such that $\vec{s}(\mathbf{W}^{*})=\vec{s}^{*}$. Explicitly, given the definition of node strength in Eq. (3), the network Hamiltonian can be written for a generic value of $\vec{\beta}$ as

$$H(\mathbf{W},\vec{\beta})=\sum_{i=1}^{n}\sum_{i<j}(\beta_{i}+\beta_{j})w_{ij}.$$

(6)

The partition function can be easily shown Squartini and Garlaschelli (2017b) to be

$$Z(\vec{\beta})=\prod_{i=1}^{n}\prod_{i<j}\frac{1}{1-e^{-(\beta_{i}+\beta_{j})}}$$

(7)

provided that $\beta_{i}+\beta_{j}>0$ for all $i,j$ (otherwise, the model admits no solution). The canonical probability distribution therefore factorizes over pairs of nodes as

$$P_{\textrm{can}}(\mathbf{W}|\vec{\beta})=\prod_{i=1}^{n}\prod_{i<j}q_{ij}(w_{ij}|\vec{\beta}),$$

(8)

where

$$q_{ij}(w|\vec{\beta})=\frac{e^{-(\beta_{i}+\beta_{j})w}}{[1-e^{-(\beta_{i}+\beta_{j})}]^{-1}}$$

(9)

is the probability that the weight of the link between nodes $i$ and $j$ takes the particular value $w$. Therefore different pairs of nodes are statistically independent in the canonical ensemble (while they are not in the microcanonical one).

Note that $q_{ij}(w|\vec{\beta})$ is a geometric distribution Squartini et al. (2015b); Squartini and Garlaschelli (2017b) with expected value

	$\displaystyle\langle w_{ij}\rangle_{\vec{\beta}}$	$\displaystyle=$	$\displaystyle\sum_{w\in\mathbb{N}}w\,q_{ij}(w|\vec{\beta})$		(10)
		$\displaystyle=$	$\displaystyle\sum_{w\in\mathbb{N}}w\,\frac{e^{-(\beta_{i}+\beta_{j})w}}{[1-e^{-(\beta_{i}+\beta_{j})}]^{-1}}$
		$\displaystyle=$	$\displaystyle\frac{e^{-(\beta_{i}+\beta_{j})}}{1-e^{-(\beta_{i}+\beta_{j})}}$

(representing the expected weight of the link connecting nodes $i$ and $j$) and variance

	$\displaystyle\textrm{Var}_{\vec{\beta}}(w_{ij})$	$\displaystyle=$	$\displaystyle\frac{e^{-(\beta_{i}+\beta_{j})}}{[1-e^{-(\beta_{i}+\beta_{j})}]^{2}}$		(11)
		$\displaystyle=$	$\displaystyle\langle w_{ij}\rangle_{\vec{\beta}}(1+\langle w_{ij}\rangle_{\vec{\beta}}).$

As we will discuss in detail is Sec. IV, Eq. (10) has the form of Bose-Einstein statistics, where $\langle w_{ij}\rangle_{\vec{\beta}}$ plays the role of an expected occupation number for the state labeled by nodes $i$ and $j$. In an appropriate ‘low-temperature’ regime, BEC can emerge into the model through the divergence of the occupation number $\langle w_{ij}\rangle_{\vec{\beta}}$ for one (possibly degenerate) ‘ground state’ (corresponding to $\beta_{i}+\beta_{j}\to 0^{+}$), while the occupation number for the other states remains finite Park and Newman (2004); Squartini and Garlaschelli (2017b); Garlaschelli and Loffredo (2009). This discussion requires a series of considerations that we leave for later. For the moment, we notice that Eq. (10) allows us to determine the special value ${\vec{\beta}}^{*}$ corresponding to the given strength sequence ${\vec{s}}^{*}$. Summing over all nodes $j\neq i$, the average value of the strength of node $i$ is

$$\langle s_{i}\rangle_{\vec{\beta}}=\sum_{j\neq i}\langle w_{ij}\rangle_{\vec{\beta}}=\sum_{j\neq i}\frac{e^{-(\beta_{i}+\beta_{j})}}{1-e^{-(\beta_{i}+\beta_{j})}},$$

(12)

whence we can reformulate Eq. (4) as

$$\sum_{j\neq i}\langle w_{ij}\rangle_{\vec{\beta}^{*}}=\sum_{j\neq i}\frac{e^{-(\beta^{*}_{i}+\beta^{*}_{j})}}{1-e^{-(\beta^{*}_{i}+\beta^{*}_{j})}}\equiv s^{*}_{i}\quad i=1,n$$

(13)

which fixes the unique parameter value $\vec{\beta}^{*}$. Notably, this value is also the one that maximizes the (log-)likelihood Garlaschelli and Loffredo (2008); Squartini and Garlaschelli (2017b), i.e.

$$\vec{\beta}^{*}=\textrm{argmax}_{\vec{\beta}}\ln P_{\textrm{can}}(\mathbf{W}^{*}|\vec{\beta}),$$

(14)

where, again, $\mathbf{W}^{*}$ is any configuration such that $\vec{s}(\mathbf{W}^{*})=\vec{s}^{*}$. In general, it is not possible to write $\vec{\beta}^{*}$ explicitly as a function of $\vec{s}^{*}$. However, Eq. (13) or equivalently Eq. (14) can be efficiently solved numerically Squartini et al. (2015b); Squartini and Garlaschelli (2017b) using various algorithms that have been coded for this purpose max ; meh . In any case, a general property (that we will use later) is that for any two nodes $i$ and $j$ with the same expected strength ($s^{*}_{i}=s^{*}_{j}$), the corresponding parameters $\beta^{*}_{i}$ and $\beta^{*}_{j}$ obey the same equation in (13) and are therefore equal. In other words, $s^{*}_{i}=s^{*}_{j}$ implies $\beta^{*}_{i}=\beta^{*}_{j}$.

Once $\vec{\beta}^{*}$ is calculated, we can plug it back into Eq. (8) to finally obtain the probability distribution

$$P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})=\prod_{i=1}^{n}\prod_{i<j}\frac{e^{-(\beta^{*}_{i}+\beta^{*}_{j})w_{ij}}}{[1-e^{-(\beta^{*}_{i}+\beta^{*}_{j})}]^{-1}}$$

(15)

that characterizes the canonical ensemble entirely. For practical purposes, $\vec{\beta}^{*}$ can be inserted into Eq. (9) to obtain the link weight probability $q_{ij}(w|\vec{\beta}^{*})$, from which several expected network properties can be calculated very directly. For instance, besides the expected link weight $\langle w_{ij}\rangle_{\vec{\beta}^{*}}$, we can calculate the probability that nodes $i$ and $j$ are connected by a link, irrepective of the weight of the latter, as follows:

	$\displaystyle p^{*}_{ij}$	$\displaystyle\equiv$	$\displaystyle\langle\Theta(w_{ij})\rangle_{\vec{\beta}^{*}}$		(16)
		$\displaystyle=$	$\displaystyle\sum_{w=1}^{\infty}w\,q_{ij}(w|\vec{\beta}^{*})$
		$\displaystyle=$	$\displaystyle 1-q_{ij}(0|\vec{\beta}^{*})$
		$\displaystyle=$	$\displaystyle e^{-(\beta^{*}_{i}+\beta^{*}_{j})}$
		$\displaystyle=$	$\displaystyle\frac{\langle w_{ij}\rangle_{\vec{\beta}^{*}}}{1+\langle w_{ij}\rangle_{\vec{\beta}^{*}}},$

where $\Theta(x)$ denotes the Heaviside step function, defined as $\Theta(x)=1$ if $x>0$ and $\Theta(x)=0$ if $x\leq 0$. Note that, if $i$ and $j$ belong to the condensed state where the expected link weight $\langle w_{ij}\rangle_{\vec{\beta}^{*}}$ diverges ($\beta^{*}_{i}+\beta^{*}_{j}\to 0^{+}$), then they become deterministically connected, i.e. $p^{*}_{ij}\to 1^{-}$. By contrast, non-condensed states have $\langle w_{ij}\rangle_{\vec{\beta}^{*}}<\infty$, $\beta^{*}_{i}+\beta^{*}_{j}>0$, and $p^{*}_{ij}<1$. The analogy with BEC will be discussed in much more detail in Sec. IV.

Besides the structural properties, one of the key quantities that we will need in order to determine EE (or the lack thereof) is the resulting canonical entropy $S^{*}_{\textrm{can}}$, obtained by inserting Eq. (15) into Eq. (1):

	$\displaystyle S^{*}_{\textrm{can}}$	$\displaystyle\equiv$	$\displaystyle S[P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})]$		(17)
		$\displaystyle=$	$\displaystyle\langle H(\mathbf{W},\vec{\beta}^{*})\rangle+\ln Z(\vec{\beta}^{*})$
		$\displaystyle=$	$\displaystyle\vec{\beta}^{*}\cdot\langle\vec{s}(\mathbf{W})\rangle+\ln Z(\vec{\beta}^{*})$
		$\displaystyle=$	$\displaystyle\vec{\beta}^{*}\cdot\vec{s}^{*}+\ln Z(\vec{\beta}^{*})$
		$\displaystyle=$	$\displaystyle-\ln P_{\textrm{can}}(\mathbf{W}^{*}|\vec{\beta}^{*}).$

Note that the calculation of the canonical entropy $S^{*}_{\textrm{can}}$ of the entire weighted network ensemble only requires the knowledge of the probability of one generic network $\mathbf{W}^{*}$ with strength sequence $\vec{s}^{*}$, which is in turn directly calculated through Eq. (15).

We now come to the microcanonical ensemble. Its governing probability distribution $P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})$ can be obtained by maximizing the Gibbs-Shannon entropy functional in Eq. (1) under the hard constraint

$$\vec{s}(\mathbf{W})=\vec{s}^{*}$$

(18)

that applies to each network $\mathbf{W}$ realized (with positive probability) in the set $\mathcal{W}_{n}$. The solution is obviously the uniform probability distribution

$$P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})=\left\{\begin{array}[]{cc}\Omega_{\vec{s}^{*}}^{-1}&\vec{s}(\mathbf{W})=\vec{s}^{*}\\ 0&\vec{s}(\mathbf{W})\neq\vec{s}^{*}\end{array}\right.$$

(19)

where $\Omega_{\vec{s}^{*}}$ is the number of networks for which the hard constraint in Eq. (18) is realized. An implicit assumption throughout this paper is that the particular strength sequence $\vec{s}^{*}$ is graphic, i.e. it can be realized by at least one network, so that $\Omega_{\vec{s}^{*}}>0$. In this case as well, the (microcanonical) entropy is obtained by inserting Eq. (19) into Eq. (1):

	$\displaystyle S^{*}_{\textrm{mic}}$	$\displaystyle\equiv$	$\displaystyle S[P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})]$		(20)
		$\displaystyle=$	$\displaystyle\ln\Omega_{\vec{s}^{*}}$
		$\displaystyle=$	$\displaystyle-\ln P_{\textrm{mic}}(\mathbf{W}^{*}|\vec{s}^{*}),$

which is also known as Boltzmann entropy. Note that $\mathbf{W}^{*}$ has the same meaning as in Eq. (17), therefore both the canonical and microcanonical entropies are equal to minus the log of the corresponding probability, evaluated in any state $\mathbf{W}^{*}$ realizing the hard constraint in Eq. (18).

Note that, although the derivation of $P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})$ is formally much more direct than that of $P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})$ in the conjugate canonical ensemble, its explicit calculation is more challenging, as it requires the combinatorial enumeration of all the $\Omega_{\vec{s}^{*}}$ weighted networks with strength sequence $\vec{s}^{*}$ (as a side remark, it is precisely the local nature of $\vec{s}^{*}$ that makes the calculation of $\Omega_{\vec{s}^{*}}$ daunting). Here, we will employ a recently proposed saddle-point asymptotic formula, for a generic discrete system under a $K$-dimensional vector $\vec{C}^{*}$ of effective (see Sec. I.3) constraints, for the number $\Omega_{\vec{C}^{*}}$ of microcanonical configurations Squartini and Garlaschelli (2017a). The formula uses only conjugate canonical quantities, namely the canonical entropy $S^{*}_{\textrm{can}}$ and the $K\times K$ covariance matrix $\mathbf{\Sigma}^{*}$ among the $K$ constraints in the canonical ensemble, and reads Squartini and Garlaschelli (2017a)

$$\Omega_{\vec{C}^{*}}={\frac{e^{S^{*}_{\textrm{can}}}}{\sqrt{\det(2\pi\mathbf{\Sigma^{*}})}}\prod_{k=1}^{K}[1+O(1/\lambda^{*}_{k})]},$$

(21)

where $\{\lambda^{*}_{k}\}_{k=1}^{K}$ are the eigenvalues of $\mathbf{\Sigma}^{*}$. The symbol $O(x)$ indicates a quantity with a finite limit when divided by $x$ as $n\to\infty$, i.e. $O(x)$ is asymptotically of the same order as $x$. Note that, since covariance matrices are positive semidefinite, $\lambda^{*}_{k}\geq 0$ for all $k$. Moreover, since the constraints are assumed to be non-redundant, then $\lambda^{*}_{k}>0$ for all $k$ Squartini and Garlaschelli (2017a) (if some of the constraints were redundant, there would be certain zero eigenvalues rendering the above equation inapplicable; that is why the formula should be applied to a maximal set of $K$ non-redundant constraints). Finally, if these eigenvalues grow sufficiently fast as $n\to\infty$, then the product on the right hand side becomes irrelevant, in which case the knowledge of $S^{*}_{\textrm{can}}$ and $\det(2\pi\mathbf{\Sigma^{*}})$ is enough in order to characterize the asymptotics of $\Omega_{\vec{C}^{*}}$.

In our setting where $\vec{C}^{*}=\vec{s}^{*}$ and $K=n$ (node strengths are all mutually independent as it is not possible to guess any individual node strength from the knowledge of the other $n-1$ ones), we calculate the entries of $\mathbf{\Sigma^{*}}$ as

	$\displaystyle{\Sigma}^{*}_{ij}$	$\displaystyle\equiv$	$\displaystyle\textrm{Cov}_{\vec{\beta}^{*}}(s_{i},s_{j})$		(22)
		$\displaystyle=$	$\displaystyle\langle s_{i}s_{j}\rangle_{\vec{\beta}^{*}}-\langle s_{i}\rangle_{\vec{\beta}^{*}}\langle s_{j}\rangle_{\vec{\beta}^{*}}$
		$\displaystyle=$	$\displaystyle\left.\frac{\partial^{2}\ln Z(\vec{\beta})}{\partial{\beta_{i}}\partial{\beta_{j}}}\right|_{\vec{\beta}=\vec{\beta}^{*}}$

where $Z(\vec{\beta})$ is given by Eq. (7). An explicit calculation gives

	$\displaystyle{\Sigma}^{*}_{ii}$	$\displaystyle=$	$\displaystyle\textrm{Var}_{\vec{\beta}^{*}}(s_{i})$		(23)
		$\displaystyle=$	$\displaystyle\sum_{j\neq i}\frac{e^{-(\beta^{*}_{i}+\beta^{*}_{j})}}{[1-e^{-(\beta^{*}_{i}+\beta^{*}_{j})}]^{2}}$
		$\displaystyle=$	$\displaystyle\sum_{j\neq i}\langle w_{ij}\rangle_{\vec{\beta}^{*}}(1+\langle w_{ij}\rangle_{\vec{\beta}^{*}})$

for the diagonal entries (i.e. the variances of the constraints) and

	$\displaystyle{\Sigma}^{*}_{ij}$	$\displaystyle=$	$\displaystyle\textrm{Cov}_{\vec{\beta}^{*}}(s_{i},s_{j})$		(24)
		$\displaystyle=$	$\displaystyle\textrm{Var}_{\vec{\beta}^{*}}(w_{ij})$
		$\displaystyle=$	$\displaystyle\frac{e^{-(\beta^{*}_{i}+\beta^{*}_{j})}}{[1-e^{-(\beta^{*}_{i}+\beta^{*}_{j})}]^{2}}$
		$\displaystyle=$	$\displaystyle\langle w_{ij}\rangle_{\vec{\beta}^{*}}(1+\langle w_{ij}\rangle_{\vec{\beta}^{*}})\quad(i\neq j)$

for the off-diagonal entries (i.e. the covariances between distinct constraints).

We finally obtain

	$\displaystyle S^{*}_{\textrm{mic}}$	$\displaystyle=$	$\displaystyle\ln\Omega_{\vec{s}^{*}}$
		$\displaystyle=$	$\displaystyle S^{*}_{\textrm{can}}-\ln\sqrt{\det(2\pi\mathbf{\Sigma}^{*})}+\sum_{k=1}^{n}\ln[1+O(1/\lambda^{*}_{k})]$

where we have used $\sqrt{\det(2\pi\mathbf{\Sigma}^{*})}=\prod_{k=1}^{n}\sqrt{2\pi\lambda^{*}_{k}}$. Note that the eigenvalues $\{\lambda^{*}_{k}\}_{k=1}^{n}$ are positive, the $n$ node strengths being linearly independent constraints Squartini and Garlaschelli (2017a). In principle, in order to compute (the leading term of) Eq. (II.3) explicitly, we need to specify a value for $\vec{s}^{*}$, calculate the resulting matrix $\mathbf{\Sigma^{*}}$, and the eigenvalues of the latter. However, in Sec. III we show that the knowledge of the diagonal elements of the covariance matrix is enough for our purposes. This result is then used in Sec. IV when we consider specific choices of $\vec{s}^{*}$.

As we have anticipated, EE can be stated mathematically using three different notions, namely thermodynamic, macrostate and measure equivalence Touchette (2015). These definitions turn out to be, under mild assumptions, essentially equivalent Touchette (2015). Here, we use the definition in the measure sense, which has been recently formulated explicitly for binary network ensembles Squartini et al. (2015a); Garlaschelli et al. (2016); Roccaverde (2019); Squartini and Garlaschelli (2017a) and is based on the vanishing of a suitable relative entropy density between the microcanonical and canonical probability distributions. Our calculations generalize those results to the case of weighted networks, for which measure equivalence has not been studied yet.

In general, the relative entropy (or Kullback-Leibler divergence) between two distributions $P$ and $Q$, both having support over a discrete set $\mathcal{W}_{n}=\{\mathbf{W}\}$ of configurations in analogy with Eq. (1), is defined as

$$D[P||Q]\equiv\sum_{\mathbf{W}\in\mathcal{W}_{n}}P(\mathbf{W})\ln\frac{P(\mathbf{W})}{Q(\mathbf{W})}$$

(26)

and quantifies ‘how far’ the distribution $P$ is from the reference distribution $Q$ Kullback and Leibler (1951). When $P$ and $Q$ represent the microcanonical and canonical distributions respectively, it can be shown Squartini et al. (2015a); Squartini and Garlaschelli (2017a) that $D[P||Q]$ reduces to the difference between the canonical and microcanonical entropy, which can be both estimated on a single configuration realizing the hard constraint (as we have indeed shown in the previous section for our weighted network model). Moreover, its calculation asymptotically requires only the canonical covariance matrix $\mathbf{\Sigma^{*}}$ between the constraints.

We are now going to see how these general results apply to our specific case. A visual illustration of the idea behind measure equivalence for our ensembles of weighted networks with given strength sequence is provided in Fig. 2. Following Eq. (26), the relative entropy between $P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})$ and $P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})$ is defined as

	$\displaystyle D^{*}$	$\displaystyle\equiv$	$\displaystyle D[P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})||P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})]$		(27)
		$\displaystyle=$	$\displaystyle\sum_{\mathbf{W}\in\mathcal{W}_{n}}P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})\ln\frac{P_{\textrm{mic}}(\mathbf{W}|\vec{s}^{*})}{P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})}.$

By inserting Eq. (19) into Eq. (27) and using the fact that $P_{\textrm{can}}(\mathbf{W}|\vec{\beta}^{*})$ has the same value for any network $\mathbf{W}^{*}$ matching the hard constraint $\vec{s}(\mathbf{W}^{*})=\vec{s}^{*}$, we confirm that $D^{*}$ can be estimated pointwise on $\mathbf{W}^{*}$ as

$$D^{*}=\ln P_{\textrm{mic}}(\mathbf{W}^{*}|\vec{s}^{*})-\ln P_{\textrm{can}}(\mathbf{W}^{*}|\vec{\beta}^{*}).$$

(28)

Moreover, using Eqs. (17) and (20), we also confirm that it reduces to the entropy difference

$$D^{*}=S^{*}_{\textrm{can}}-S^{*}_{\textrm{mic}}.$$

(29)

Now, Eq. (II.3) immediately allows us to obtain

$$D^{*}=\ln\sqrt{\det(2\pi\mathbf{\Sigma}^{*})}+\sum_{k=1}^{n}\ln[1+O(1/\lambda^{*}_{k})]$$

(30)

which depends only on the eigenvalues of the canonical covariance matrix $\mathbf{\Sigma^{*}}$, whose diagonal and off-diagonal entries are given in Eqs. (23) and (24) respectively.

The definition of measure equivalence is the vanishing of the relative entropy density, i.e. of the ratio $D^{*}/n$, in the thermodynamic limit Touchette (2015). Explicitly, EE in the measure sense corresponds to

$$d^{*}\equiv\lim_{n\rightarrow\infty}\frac{D^{*}}{n}=0$$

(31)

or equivalently Squartini and Garlaschelli (2017a)

$$D^{*}=o(n),$$

(32)

where $o(n)$ denotes a quantity that, if divided by $n$, vanishes as $n\to\infty$. Equation (29) allows us to understand the above definition of macrostate EE as follows. The microcanonical entropy $S^{*}_{\textrm{mic}}$ is the logarithm of the number of accessible configurations under hard constraints, while the canonical entropy $S^{*}_{\textrm{can}}$ is the logarithm of the corresponding ‘effective’ number of configurations under soft constraints. EE requires that, as $n$ increases, the typical configurations of the system under the two ensembles become the same, i.e. that the (effective) numbers of configurations in the two ensembles become closer to each other. This cannot happen if, as we keep adding one more unit to the system, the difference between the two entropies (i.e. the relative entropy $D^{*}$) keeps increasing by an arbitrary amount. The criterion in Eq. (31) establishes that if the entropy difference per node, i.e. $D^{*}/n$, does not vanish as $n$ diverges, then the two ensembles cannot be equivalent.

Equation (32) implies that, in order to assess whether the system is under EE, we do not need the exact value of $D^{*}$, but only its leading order with respect to $n$. Then from Eq. (30) we see that, since the term $1+O(1/\lambda^{*}_{k})$ is at most of the same order as $\sqrt{2\pi\lambda^{*}_{k}}$, the presence of $O(1/\lambda^{*}_{k})$ does not affect ensemble (non)equivalence:

$$D^{*}=O\left(\sum_{k=1}^{n}\ln\lambda^{*}_{k}\right)=O(\ln\det\mathbf{\Sigma}^{*}).$$

(33)

So, in order to check whether Eq. (32) holds, it is ultimately enough to check whether

$$\ln\det\mathbf{\Sigma}^{*}=o(n).$$

(34)

On the other hand, since our hypothesis of non-redundant constraints implies $\lambda^{*}_{k}>0$ for all $k$ (as discussed in Sec. II.3), we see that the contribution of the term $\sum_{k=1}^{n}\ln[1+O(1/\lambda^{*}_{k})]$ to the relative entropy in Eq. (30) is at most $O(n)$. So if $\ln\sqrt{\det(2\pi\mathbf{\Sigma}^{*})}$ grows faster than $O(n)$ then Eq. (30) reduces to the stronger result

$$D^{*}\approx\ln\sqrt{\det(2\pi\mathbf{\Sigma}^{*})}\approx\frac{1}{2}\ln\det\mathbf{\Sigma}^{*},$$

(35)

i.e. the leading term of the relative entropy (not only its leading order) can be calculated exactly from $\ln\det\mathbf{\Sigma}^{*}$ (throughout this paper, the symbol “$x\approx y$” indicates that the leading term of $x$ and $y$ is asymptotically the same, i.e. the two quantities differ by a quantity that vanishes if divided by either $x$ or $y$ as $n\to\infty$). Depending on the regimes considered later on in the paper, different techniques for calculating (the leading order of) the determinant of the covariance matrix $\mathbf{\Sigma}^{*}$ can be used. We will discuss these techniques when needed, and refer to the Appendix for explicit calculations. We will show that, except in a certain zero-temperature limit, the conditions ensuring Eq. (35) are met and the leading order of the relative entropy can be calculated exactly.

We now consider the relative fluctuations of the constraints, whose behaviour in the thermodynamic limit is, historically, the traditional criterion used to check whether statistical ensembles are equivalent Gibbs (1902). In the standard situation, where the canonical and microcanonical ensembles are defined through a single scalar constraint on the total energy $E$, the relative fluctuations are captured by a single scalar quantity $r^{*}\equiv\sigma_{\beta^{*}}(E)/{E^{*}}$ representing the ratio of the canonical standard deviation of the energy to the mean energy itself. EE is then associated with the vanishing of $r^{*}$ as $n\to\infty$. In statistics, $r^{*}$ is called the coefficient of variation of the random variable $E$ with expected value $E^{*}$ and variance $\sigma^{2}_{\beta^{*}}(E)$

In the case of networks with local constraints, there are $n$ coefficients of variation to consider. They have been calculated for both the binary and the weighted versions of the configuration model Squartini et al. (2015b). In extreme summary, those results show that, in the binary case (where there is a constraint on the degree $k^{*}_{i}$ for each node $i=1,n$), an upper bound $1/\sqrt{k^{*}_{i}}$ for the relative fluctuation $r_{i}^{*}=\sigma_{\beta^{*}}(k_{i})/{k^{*}_{i}}$ can been established. By contrast, in the weighted case (corresponding to the model considered here with a constraint on the strength $s^{*}_{i}$ for each node $i=1,n$), the value $1/\sqrt{s^{*}_{i}}$ becomes a lower bound for the relative fluctuation $r_{i}^{*}=\sigma_{\beta^{*}}(s_{i})/{s^{*}_{i}}$ Squartini et al. (2015b). Those results have two implications.

First, in the binary case the only regime for which general conclusions can be drawn about the vanishing of the relative fluctuations is the so-called dense regime where the expected degree of all nodes diverges, hence $r_{i}^{*}\to 0$. In the opposite sparse regime where the average degree of all nodes is finite, we only have a finite upper bound for the relative fluctuations, but their actual value depends on the specific network. In general, however, the decreasing behaviour of the upper bound for the relative fluctuations in binary networks with increasing degrees is opposite to that of the relative entropy density, which increases as the expected degree increases Squartini et al. (2015a); Garlaschelli et al. (2016); Roccaverde (2019).

Second, in the weighted case we have a somewhat opposite situation where we can only conclude that, in the sparse regime where the expected strengths are finite (apart from possible hubs), the relative fluctuations do not vanish. By contrast, in the dense case where the expected strengths diverge, the relative fluctuations can in principle pick any value. Confusingly, even if in the weighted case the lower bound for the relative fluctuations of the strengths goes to zero for nodes with diverging strength $s_{i}^{*}$, previous results seem to indicate that the realized value of $r_{i}^{*}$ in networks with heterogeneous strength sequence actually increases for higher $s_{i}^{*}$ Squartini et al. (2015b). This behaviour suggests that weighted networks, due to the many possible ways in which weight can accumulate on links, can behave very differently from binary networks. This observation requires further research and strengthens our motivation for studying the relative fluctuations in weighted networks in different scenarios (ranging from homogeneous to heterogenous concentrations of link weights), in conjunction with the relative entropy and, in general, ensemble (non)equivalence.

Using Eq. (23), we can immediately calculate the standard deviation of each constraint $s_{i}$ around its expected value $s^{*}_{i}$ as

	$\displaystyle\sigma_{\vec{\beta}^{*}}(s_{i})$	$\displaystyle=$	$\displaystyle\sqrt{\textrm{Var}_{\vec{\beta}^{*}}(s_{i})}$		(36)
		$\displaystyle=$	$\displaystyle\sqrt{\sum_{i\neq j}\frac{e^{-(\beta^{*}_{i}+\beta^{*}_{j})}}{[1-e^{-(\beta^{*}_{i}+\beta^{*}_{j})}]^{2}}}$
		$\displaystyle=$	$\displaystyle\sqrt{\sum_{i\neq j}\langle w_{ij}\rangle_{\vec{\beta}^{*}}(1+\langle w_{ij}\rangle_{\vec{\beta}^{*}})}$
		$\displaystyle=$	$\displaystyle\sqrt{s^{*}_{i}+\sum_{i\neq j}\langle w_{ij}\rangle^{2}_{\vec{\beta}^{*}}}.$

The relative fluctuation of the strength is therefore

$$r_{i}^{*}=\frac{\sigma_{\vec{\beta}^{*}}(s_{i})}{s^{*}_{i}}=\sqrt{\frac{1}{s^{*}_{i}}+\frac{\sum_{i\neq j}\langle w_{ij}\rangle_{\vec{\beta}^{*}}^{2}}{(\sum_{i\neq j}\langle w_{ij}\rangle_{\vec{\beta}^{*}})^{2}}},$$

(37)

Since $\langle w_{ij}\rangle_{\vec{\beta}^{*}}\geq 0$ for all $i,j$, we have $(\sum_{i\neq j}\langle w_{ij}\rangle_{\vec{\beta}^{*}})^{2}\geq\sum_{i\neq j}(\langle w_{ij}\rangle_{\vec{\beta}^{*}})^{2}$, showing that $1/\sqrt{s_{i}^{*}}$ is indeed a lower bound for $r_{i}^{*}$. When studying the asymptotic behaviour of the relative fluctuations in the thermodynamic limit, we will be interested in whether the limit

$\displaystyle\rho_{i}^{*}\equiv\lim_{n\to\infty}r_{i}^{*}$

(38)

is zero (vanishing relative fluctuations) or positive (nonvanishing relative fluctuations) for each node $i=1,n$.

As we have already recalled, weighted networks with a constraint on the strength sequence obey Bose-Einstein statistics, as opposed to binary networks that obey Fermi-Dirac statistics Park and Newman (2004); Garlaschelli and Loffredo (2009); Garlaschelli et al. (2013). Indeed, inserting Eq. (6) into Eq. (5) we get the probability of a configuration for a gas of free particles in the grandcanonical ensemble²²2In the grandcanonical ensemble, both the total energy and the total number of particles are treated as soft constraints. With respect to the canonical ensemble, the appearance of the number of particles as an additional soft constraint requires the introduction of an extra Lagrange multiplier, the chemical potential. Interestingly, in the context of BEC a fourth (so-called ‘Maxwell’s Demon’) ensemble has also been introduced where the total number of particles is soft while the total energy is hard Navez et al. (1997)., where the pair $i,j$ labels an energy state, the weight $w_{ij}$ is the number of particles in that state (occupation number), and the sum $\beta_{i}+\beta_{j}$ can be interpreted as

$$\beta_{i}+\beta_{j}=\frac{\epsilon_{ij}-\mu(T)}{kT}.$$

(39)

In the latter expression, $\epsilon_{ij}$ represents the energy of the state, $1/kT$ is the inverse temperature, and $\mu(T)$ is the chemical potential (required to fix the same expected overall number of particles for all values of $T$) Garlaschelli et al. (2013). Indeed, as we discussed in Sec. II.2, in our setting the energy and temperature (and in this case, the chemical potential as well) are all reabsorbed into $\vec{\beta}$. Therefore we can interpret the link weight $w_{ij}$ as the number of ‘elementary particles’ of weight, i.e. the number of quanta of unit weight, populating the link between nodes $i$ and $j$ Park and Newman (2004); Garlaschelli and Loffredo (2009). The total number of such particles in the system is the total weight $W$ of all links in the network:

$$W(\mathbf{W})=\sum_{i=1}^{n}\sum_{j<i}w_{ij}=\frac{1}{2}\sum_{i=1}^{n}s_{i}(\mathbf{W}).$$

(40)

The ‘weighted’ property $w_{ij}\in\mathbb{N}$, which leads to Eqs. (7), (8) and (9), corresponds to the possibility that the same state (pair of nodes) is occupied by indefinitely many particles (subject to the average number dictated by the chemical protential), which is a property of bosons. By contrast, in binary networks one has to impose $w_{ij}\in(0,1)$, which is a property of fermions Park and Newman (2004). An extensive treatment of the role of chemical potential and temperature in binary networks can be found in Garlaschelli et al. (2013). Here, to properly interpret what the weighted model is doing, we should give a series of clarifications.

First, we should make a clear distrinction beween the $n$ ‘units’ of our system (i.e. the nodes of the network) and the $W$ ‘particles’ of weight that, as a formal analogy, can be interpreted as populating the links of the network. The former are the real constituents of our physical system, while the latter are a mathematical abstraction used to represent the nature of the interactions (links) between such constituents. If we imagine doubling the size of our network, we should imagine doubling the number $n$ of nodes: indeed, we can imagine the network ‘growing in size’ by adding one single node at a time, but we cannot imagine adding one single pair of nodes at a time, without actually adding $n$ new pairs. One should also not be tempted to regard node pairs as the fundamental units by the fact that, mathematically, the $n(n-1)/2$ variables $\{w_{ij}\}$ involving different pairs of nodes are independent random variables: actually, this only occurs in the canonical ensemble and would in any case not be true for more general choices of the constraints. Moreover, even the $n(n-1)/2$ independent node pairs in the canonical ensemble cannot be assigned independent values of the parameters, since there are only $n$ parameters corresponding to the Lagrange multipliers attached to each node. Explicit (and strong) consequences of this fact will be illustrated precisely in the context of BEC. Therefore the physical size of our system is $n$, and this is why in Eq. (31) we defined the relative entropy density as the relative entropy divided by $n$ in the first place. How the total weight $W$ varies with the system size $n$ depends on a specific property, i.e. on how we make the entries of $\vec{s}^{*}$ (and the resulting value of $W^{*}=\sum_{i=1}^{n}s^{*}_{i}/2$) scale with $n$. For instance, we may choose to be in the sparse regime where $\vec{s}^{*}$ remains finite as $n\to\infty$, or in the dense regime where $\vec{s}^{*}$ diverges as $n\to\infty$. As we show below, the latter is the relevant case for BEC to emerge.

Second, we stress that, irrespective of the above, we always consider a hard number $n$ of nodes, and this is why we compare (only) the canonical (soft value of $\vec{s}$) and microcanonical (hard value of $\vec{s}$) ensembles of networks, both for fixed $n$ (which sets the dimension of $\vec{s}$). We do not consider the grandcanonical ensemble of network configurations where $n$ is soft. The grandcanonical ensemble introduced in the aforementioned analogy with systems of bosons is a different one; it may be denoted as an ensemble of weight quanta in a network with fixed $n$ and originates from the fact that the Hamiltonian in Eq. (6), and consequently the total link weight (not $n$) in Eq. (40), is a fluctuating quantity in the canonical ensemble of network configurations. The fluctuations in the (rescaled) energy $H$ (canonical ensemble of network configurations) are seen as fluctuations in the particle number $W$ (grandcanonical ensemble of weight quanta) in the Bose-Einstein analogy. Fluctuations in the particle number have been the subject of many studies in the literature on BEC Navez et al. (1997); Holthaus et al. (1998); Mullin and Fernandez (2003); Chatterjee and Diaconis (2014); Tarasov et al. (2015); Crisanti et al. (2019). Note that, in both canonical and microcanonical ensembles, the individual link weights $\{w_{ij}\}$ are fluctuating quantities, despite the fact that the total link weight $W^{*}$ is constant in the microcanonical ensemble. Therefore the numbers of ‘weight particles’ of individual links fluctuate in both ensembles.

Third, while we necessarily discuss the (non)equivalence of the canonical and microcanonical ensembles in the thermodynamic limit $n\to\infty$, the total weight $W^{*}$ can (and, across the canonical ensemble, will in any case) be arbitrarily large even for finite $n$. Indeed, the phase transition that we are about to discuss (namely, BEC) does not per se require the limit $n\to\infty$, while it definitely requires the limit $W^{*}\to\infty$. Abstractly, these two limits (and the associated phenomena of EN and BEC respectively) may appear as mathematically unrelated. However, in practice they are physically related once the scaling of $\vec{s}^{*}$ with $n$ is specified. In particular, we are going to show that, in order to observe BEC, we need be in a dense regime where $W^{*}=O(n^{2})$. This ensures that, when taking the thermodynamic limit $n\to\infty$ in order to study ensemble (non)equivalence, we are automatically implying $W^{*}\to\infty$ so that we can check for BEC at the same time.

Last, we recall that $\vec{\beta}\cdot\vec{s}(\mathbf{W})$ has to be dimensionless in order to ensure that the probability is a number. Therefore, since $w_{ij}$ is dimensionless, so are $s_{i}^{*}$ and $\beta_{i}^{*}$. In turn, this implies that both sides of Eq. (39) must be dimensionless. On the other hand, when modelling a real system, the ‘energy’ $\epsilon_{ij}$ may represent any physical ‘cost’ associated to the link between nodes $i$ and $j$ (more precisely, it represents the cost of reinforcing $w_{ij}$ by a unit of weight) and may therefore carry its own unit of measure (e.g. it may depend on some distance between nodes $i$ and $j$). Necessarily, the chemical potential $\mu(T)$ carries the same units as the energy. As for the ‘temperature’ $T$, it may be chosen to be dimensionless as it merely represents a control parameter (this is the choice that we will make later); alternatively, it may carry the same units of the energy if it is useful that temperature and energy live on the same scale. Irrespective of this choice, in our setting the ‘Boltzmann constant’ $k$ is simply a constant that takes care of all dimensional units of measure and makes the ratio on the right hand side of Eq. (39) dimensionless.

With the above clarifications, we can finally go back to our model. In the traditional physical situation, in the canonical ensemble the energy $\epsilon_{ij}$ of each state $i,j$ is a constant and the temperature $T$ can be varied. Clearly, $\epsilon_{ij}$ is independent of $T$, while the chemical potential $\mu(T)$ is chosen, as a function of temperature, in order to realize the correct ($T$-independent) expected total number $\langle W\rangle^{*}\equiv W^{*}$ of particles for all values of $T$. In this ‘direct problem’, every state will therefore have an expected occupation number governed by $\epsilon_{ij}$, $T$ and $\mu(T)$. In our ‘inverse’ setting, $T^{*}$ and $\mu^{*}$ are instead reabsorbed into $\vec{\beta}^{*}$, which in turn is induced by the chosen value of the strength sequence (rather than the other way around). We should therefore regard $\vec{s}^{*}=\vec{s}^{*}(T)$ and $\vec{\beta}^{*}=\vec{\beta}^{*}(T)$ as $T$-dependent, while $W^{*}$ remains $T$-independent. This means that the chemical potential $\mu^{*}(T)$ should be such that

$$\sum_{i=1}^{n}{s}_{i}^{*}(T)=2W^{*}\quad\forall T\geq 0.$$

(41)

BEC emerges when, below a certain critical temperature $T_{c}$, the occupation number of the state with minimum energy $\epsilon_{\textrm{min}}=\textrm{min}_{i,j}\{\epsilon_{ij}\}$ (ground state), or of a finite number of states with lowest energy, becomes so large that it reaches a finite fraction of the total number $W^{*}$ of particles. Clearly, this requires the existence of at least two different energy levels (the ground state and at least one ‘excited’ state). Therefore the simplest way to obtain BEC in our model is by considering a strength sequence of the form

$$s^{*}_{i}(T)=\left\{\begin{array}[]{ll}s^{*}_{+}(T)&i=1,n_{+}\\ s^{*}_{-}(T)&i=n_{+}+1,n\end{array}\right.\quad s^{*}_{+}(T)\geq s^{*}_{-}(T),$$

(42)

i.e. by partitioning the $n$ nodes into two classes, which we call core and periphery: the core has a finite number

$$n_{+}=O(1)$$

(43)

of nodes, each having a ‘large’ strength $s^{*}_{+}(T)$, while the periphery has an extensive number

$$n_{-}=n-n_{+}=O(n)$$

(44)

of nodes, each having a ‘small’ strength $s^{*}_{-}(T)$. What we mean precisely by ‘small’ and ‘large’ will be clarified below. For the moment, we notice that the BEC phase ($T\!<\!T_{c}$) corresponds to picking a ‘condensed’ value of $\vec{s}^{*}(T\!<\!T_{c})$ such that, in the thermodynamic limit, the core takes up a finite fraction of the total weight $W^{*}$ of all links in the network, despite having a finite size. In particular, in the zero-temperature limit all the total weight $W^{*}$ is in the core. By contrast, the non-condensed phase $T\!>\!T_{c}$ is one where $\vec{s}^{*}(T\!>\!T_{c})$ is such that no individual link receives a finite fraction of $W^{*}$. In particular, the infinite-temperature limit should be such that the energy difference between ground and excited states becomes ineffective, i.e. $s^{*}_{+}(T\!\to\!\infty)=s^{*}_{-}(T\!\to\!\infty)$. The different phases can be efficiently monitored by introducing a temperature-dependent order parameter $Q^{*}(T)$, as we show below.