Bakry–Émery Curvature Sharpness and Curvature Flow in Finite Weighted Graphs: Implementation

Abstract

In this paper, we discuss the implementation of a curvature flow on weighted graphs based on the Bakry–Émery calculus. This flow can be adapted to preserve the Markovian property and its limits as time goes to infinity turn out to be curvature sharp weighted graphs. After reviewing some of the main results of the corresponding paper concerned with the theoretical aspects, we present various examples (random graphs, paths, cycles, complete graphs, wedge sums and Cartesian products of complete graphs, and hypercubes) and exhibit various properties of this flow. One particular aspect of our investigations is asymptotic stability and instability of curvature flow equilibria. The paper ends with a description of the Python functions and routines freely available in an ancillary file on arXiv or via github. We hope that the explanations of the Python implementation via examples will help users to carry out their own curvature flow experiments.

1. Introduction

This paper is concerned with computational aspects of a curvature flow on weighted graphs based on the Bakry–Émery calculus. This curvature flow was originally introduced in our paper [1]. The current paper focuses on properties for various concrete examples of graphs as well as stability investigations of the flow limits. Many observations in this paper are gained with the help of our curvature flow software, which is freely available as an ancillary file of [2] or via github at https://github.com/georgestagg/graph-curvature-server (accessed on 6 June 2023). A summary of our findings is presented in Section 4.

A weighted graph in this paper is a finite simple mixed combinatorial graph $G=(V,E)$ with vertex set V and edge set $E=E^1\cup E^2$ of one- and two-sided edges, together with a weighting scheme of transition rates $p_{xy}\ge 0$ for $x,y\in V$ (which can be represented by a generally non-symmetric matrix P after an enumeration of the vertices). Transition rates $p_{xy}$ can only be positive if $x=y$ or if there is a one- or two-sided edge from x to y. One-sided edges are denoted by ordered pairs $(x,y)\in E^1\subset V^2$, and two-sided edges are denoted by sets $\{x,y\}\in E^2$. One- or two-sided edges $(x,y)\in E^1$ or $\{x,y\}\in E^2$ with vanishing transition rates $p_{xy}=0$ are called degenerate and a weighted graph $(G,P)$ is called non-degenerate if it does not have degenerate edges. A weighted graph $(G,P)$ is called Markovian, if we have ${\sum }_{y\in V}{p}_{xy}=1$ for all $x\in V$. In this case, P is a stochastic matrix and we can view the transition rates $p_{xy}$ as transition probabilities of a lazy random walk (with non-zero laziness if there exists a vertex $x\in V$ with $p_{xx}>0$). Our curvature flow does not affect the underlying combinatorial graph G, but it changes the weighting scheme. We will focus on a version of our flow that preserves the Markovian property. In other words, starting with an initial Markovian weighted graph $(G,P_0)$, this flow will provide a family of Markovian weighting schemes ${\left\{P\left(t\right)\right\}}_{t\ge 0}$ with $P\left(0\right)=P_0$, depending on the continuous time parameter t.

Before we introduce our curvature flow, we need to briefly discuss the relevant Bakry–Émery curvature background. This curvature notion is based on the weighted Laplacian $\Delta ={\Delta }_P$, acting on functions $f:V\to \mathbb{R}$ as follows:

$${\Delta }_{P}f\left(x\right)=\sum _{y\in V}{p}_{xy}(f\left(y\right)-f\left(x\right)).$$

The Laplacian gives rise to the following symmetric bilinear “carré du champ operators” $\mathsf{\Gamma }$ and ${\mathsf{\Gamma }}_2$:

$$\begin{array}{ccc}\hfill 2\mathsf{\Gamma }(f,g)& =& \Delta \left(fg\right)-f\Delta g-g\Delta f,\hfill \\ \hfill 2{\mathsf{\Gamma }}_2(f,g)& =& \Delta \mathsf{\Gamma }(f,g)-\mathsf{\Gamma }(f,\Delta g)-\mathsf{\Gamma }(g,\Delta f).\hfill \end{array}$$

1.1. Bakry–Émery Curvature and Curvature Sharpness

Bakry–Émery curvature depends on a dimension parameter N and is well-defined at every vertex $x\in V$ that is not isolated, that is, there exists another vertex $y\in V$ with $p_{xy}>0$. For isolated vertices, there is some ambiguity on how to define their curvature, and we decided to assign to such a vertex the curvature value 0. (Another natural choice of curvature for an isolated vertex $x\in V$ would be $K_N\left(x\right)=\infty$ for all $N\in (0,\infty ]$. An argument for that choice is that an isolated vertex can be viewed as a discrete analogue of a limit of round spheres with radii shrinking to 0, whose curvatures would diverge to infinity.) The definition of Bakry–Émery curvature reads as follows:

Definition 1

(Bakry–Émery curvature). The Bakry–Émery curvature of a non-isolated vertex $x\in V$ for a fixed dimension $N\in (0,\infty ]$ is the supremum of all values $K\in \mathbb{R}$, satisfying the curvature-dimension inequality

$${\mathsf{\Gamma }}_2\left(f\right)\left(x\right)\ge \frac{1}{N}{(\Delta f\left(x\right))}^2+K\mathsf{\Gamma }\left(f\right)\left(x\right)$$

(1)

for all functions $f:V\to \mathbb{R}$. We use the simplified notation $\mathsf{\Gamma }\left(f\right)=\mathsf{\Gamma }(f,f)$ and ${\mathsf{\Gamma }}_2\left(f\right)={\mathsf{\Gamma }}_2(f,f)$. We denote the curvature at $x\in V$ by $K_N\left(x\right)=K_{P,N}\left(x\right)$. If $x\in V$ is isolated, that is, we have $p_{xy}=0$ for all $y\in V\backslash \left\{x\right\}$, we set $K_N\left(x\right)=K_{P,N}\left(x\right)=0$ for all $N\in (0,\infty ]$.

This curvature notion is motivated by Bochner’s identity (see, e.g., [3] (Prop. 4.15)), a fundamental pointwise formula in the smooth setting of n-dimensional Riemannian manifolds involving gradients, Laplacians, Hessians, and Ricci curvature. It was introduced for the smooth setting in [4]. The curvature was then reintroduced several times in the setting of graphs, see [5,6,7]. For further research about Bakry–Émery curvature on finite graphs, see, e.g., [8,9,10,11,12,13,14,15,16,17,18,19].

By the definition, the inequality

$${\mathsf{\Gamma }}_2\left(f\right)\left(x\right)\ge \frac{1}{N}{(\Delta f\left(x\right))}^2+K_N\left(x\right)\mathsf{\Gamma }\left(f\right)\left(x\right)$$

(2)

holds for every function f, and therefore, also holds for the combinatorial distance function $d(x,·)$. Here, $d(x,y)$ is the length of a shortest directed path from x to y (if there is no such path, we set $d(x,y)=\infty$). If we have equality at x in (2) for this particular function $f=d(x,·)$, we say that the vertex $x\in V$ is N-curvature sharp. Curvature sharpness will be particularly important in our considerations. Curvature sharpness was originally introduced in [20] (Definition 1.4). The (equivalent) definition given in this paper is inspired by [21] (Proof of Theorem 1.2). For more details about relations between different curvature sharpness definitions see Section 3 of [1].

Definition 2

(Curvature sharpness). Let $(G,P)$ be a weighted graph and $N\in (0,\infty ]$. A vertex $x\in V$ is called N-curvature sharp if we have

$${\mathsf{\Gamma }}_2\left(f\right)\left(x\right)=\frac{1}{N}{(\Delta f\left(x\right))}^2+K_N\left(x\right)\mathsf{\Gamma }\left(f\right)\left(x\right)$$

(3)

for the distance function $f=d(x,·)$. Moreover, a vertex $x\in V$ is curvature sharp if it is curvature sharp for some dimension $N\in (0,\infty ]$. A weighted graph $(G,P)$ is called curvature sharp, if every vertex of G is curvature sharp.

Note that each function $f:V\to \mathbb{R}$ with $\mathsf{\Gamma }\left(f\right)\left(x\right)\ne 0$ gives rise to an upper curvature bound $K_{P,N}^f\left(x\right)$ via the inequality (2). Namely, we have

$$K_N\left(x\right)\le K_{P,N}^f\left(x\right):=\frac{1}{\mathsf{\Gamma }\left(f\right)\left(x\right)}\left({\mathsf{\Gamma }}_2\left(f\right)\left(x\right)-\frac{1}{N}{(\Delta f\left(x\right))}^2\right).$$

(4)

A vertex $x\in V$ is therefore N-curvature sharp if its Bakry–Émery curvature $K_N\left(x\right)$ agrees with the specific upper curvature bound $K_{P,N}^{d(x,·)}\left(x\right)$. We would also like to mention the following monotonicity property of curvature sharpness: If $x\in V$ is N-curvature sharp that this vertex is also curvature sharp for any dimension $\le N$ (see [1] (Prop. 3.1)).

In the next subsection, we present an important reformulation of Bakry–Émery curvature at $x\in V$ using a specific matrix $Q\left(x\right)$, which will be important in the definition of the curvature flow.

1.2. Reformulation of Curvature via a Schur Complement

The combinatorial distance function allows us to define distance spheres and distance balls,

$$\begin{array}{ccc}\hfill S_r\left(x\right)& =& \{z\in V:d(x,z)=r\},\hfill \\ \hfill B_r\left(x\right)& =& \{z\in V:d(x,z)\le r\}.\hfill \end{array}$$

Let $x\in V$ be a non-isolated vertex. It turns out that the Bakry–Émery curvature $K_N\left(x\right)$ is determined locally, that is, can be derived solely from the information about the 2-ball

$$B_2\left(x\right)=\left\{x\right\}\cup S_1\left(x\right)\cup S_2\left(x\right).$$

More precisely, denoting $S_1\left(x\right)=\{y_1,\dots ,y_m\}$ and $S_2\left(x\right)=\{z_1,\dots ,z_n\}$, there exist a column vector $\Delta \left(x\right)$ and a symmetric matrix $\mathsf{\Gamma }\left(x\right)$ of size m and a symmetric matrix ${\mathsf{\Gamma }}_2\left(x\right)$ of size $m+n$ such that, for functions $f,g:V\to \mathbb{R}$ with $f\left(x\right)=g\left(x\right)=0$,

$$\begin{array}{ccc}\hfill \Delta f\left(x\right)& =& \Delta {\left(x\right)}^{\top }{\overrightarrow{f}}_m,\hfill \\ \hfill \mathsf{\Gamma }(f,g)\left(x\right)& =& {\overrightarrow{f}}_m^{\top }\mathsf{\Gamma }\left(x\right){\overrightarrow{g}}_m,\hfill \\ \hfill {\mathsf{\Gamma }}_2(f,g)\left(x\right)& =& {\overrightarrow{f}}_{m+n}^{\top }{\mathsf{\Gamma }}_2\left(x\right){\overrightarrow{g}}_{m+n},\hfill \end{array}$$

where ${\overrightarrow{f}}_m={(f\left(y_1\right),\dots ,f\left(y_m\right))}^{\top }$ and

$${\overrightarrow{f}}_{m+n}={(f\left(y_1\right),\dots ,f\left(y_m\right),f\left(z_1\right),\dots ,f\left(z_n\right))}^{\top },$$

and ${\overrightarrow{g}}_n,{\overrightarrow{g}}_{m+n}$, accordingly. Using the $(n+m)$-block decomposition

$${\mathsf{\Gamma }}_2\left(x\right)=\left(\begin{array}{cc}{\mathsf{\Gamma }}_2{\left(x\right)}_{S_1}& {\mathsf{\Gamma }}_2{\left(x\right)}_{S_1,S_2}\\ {\mathsf{\Gamma }}_2{\left(x\right)}_{S_2,S_1}& {\mathsf{\Gamma }}_2{\left(x\right)}_{S_2}\end{array}\right)$$

and employing the Schur complement

$$Q\left(x\right)={\mathsf{\Gamma }}_2{\left(x\right)}_{S_1}-{\mathsf{\Gamma }}_2{\left(x\right)}_{S_1,S_2}{\mathsf{\Gamma }}_2{\left(x\right)}_{S_2}^{†}{\mathsf{\Gamma }}_2{\left(x\right)}_{S_2,S_1}$$

for matrices ${\mathsf{\Gamma }}_2\left(x\right)$ with positive semidefinite ${\mathsf{\Gamma }}_2{\left(x\right)}_{S_2}$-blocks with $A^{†}$ the pseudoinverse of A (for a given matrix $A\in {\mathbb{R}}^{N\times M}$, its pseudoinverse $A^{†}\in {\mathbb{R}}^{M\times N}$ is defined by the following conditions: $A{A}^{†}A=A$, $A^{†}A{A}^{†}=A^{†}$, and $A{A}^{†}\in {\mathbb{R}}^{N\times N}$ and $A{A}^{†}\in {\mathbb{R}}^{M\times M}$ are both symmetric matrices), we can reformulate Bakry–Émery curvature at $x\in V$ for dimension $N\in (0,\infty ]$ as follows (see Section 2.4 in [1]):

Let $x\in V$ be a non-isolated vertex. $K_N\left(x\right)$ is then the maximum of all $K\in \mathbb{R}$ such that

$$Q\left(x\right)-\frac{1}{N}\Delta \left(x\right)\Delta {\left(x\right)}^{\top }-K\mathsf{\Gamma }\left(x\right)⪰0,$$

where $A⪰B$ means that $A-B$ is positive semidefinite.

This curvature translation was motivated originally by the aim to reformulate the computation of Bakry–Émery curvature as an eigenvalue problem (see [22,23,24]).

The symmetric matrix $Q\left(x\right)$ of a non-isolated vertex $x\in V$ is of size m and is—in the non-degenerate case—closely related to another symmetric matrix $A_{\infty }\left(x\right)$, which, in turn, can be viewed as a discrete counterpart of the Ricci curvature tensor at a point $x\in M$ of a Riemannian manifold $(M,g)$ (see formula (1.2) and Section 7 in [24]). In the case of a Markovian weighted graph, curvature sharpness at a vertex $x\in V$ can also be alternatively expressed with the help of the matrix $Q\left(x\right)$ as follows.

Theorem 1

(see Theorem 1.3 in [1]). Let $(G,P)$ be a Markovian weighted graph and $x\in V$ be a non-isolated vertex with $S_1\left(x\right)=\{y_1,\dots ,y_m\}$. Then, the following statements are equivalent:

(1) 

x is curvature sharp;

(2) 

x is curvature sharp for dimension $N=2$;

(3) 

We have

$$Q\left(x\right){\mathbf{1}}_m=\frac{1}{2}{K}_{P,\infty }^{d(x,·)}\left(x\right){\mathbf{p}}_x,$$

(5)

where ${\mathbf{p}}_x={(p_{x{y}_1},\dots ,p_{x{y}_m})}^{\top }$ and ${\mathbf{1}}_m$ is the all-one column vector of size m.

Note that the term $K_{\infty }^{d(x,·)}\left(x\right)$ in (5) is the upper curvature bound introduced in (4) (in the special case $N=\infty$), that is,

$$K_{P,\infty }^{d(x,·)}\left(x\right)=\frac{{\mathsf{\Gamma }}_2\left(d(x,·)\right)\left(x\right)}{\mathsf{\Gamma }\left(d\right(x,·\left)\right)\left(x\right)}.$$

1.3. Curvature Flow

Let $(G,P_0)$ be a fixed initial weighted graph with $N=\left|V\right|$. For every non-isolated vertex $x\in V$, the size of the corresponding symmetric matrix $Q\left(x\right)$ agrees with the degree of the vertex x, that is, the number of one- and two-sided edges emanating from x, and the entries of $Q\left(x\right)$ are determined by the transition rates of edges of the 2-ball $B_2\left(x\right)$. Our curvature flow associates to this initial data a smooth matrix-valued function $P:[0,\infty )\to {\mathbb{R}}^{N\times N}$ with $P\left(0\right)=P_0$. The corresponding symmetric Q-matrices at time $t\in [0,\infty )$ depend on the weighting schemes $P\left(t\right)$, and we denote them henceforth by $Q_x\left(t\right)$ for all $x\in V$. Our curvature flow is now defined as follows.

Definition 3

(Curvature flow). Let $(G,P_0)$ be a finite weighted graph. The associated curvature flow is given by the following differential equations for all non-isolated vertices $x\in V$ and all $t\ge 0$:

$$\begin{array}{ccc}\hfill p_{xx}^{\prime }\left(t\right)& =& 0,\hfill \end{array}$$

(6)

$$\begin{array}{ccc}\hfill {\mathbf{p}}_x^{\prime }\left(t\right)& =& -4Q_x\left(t\right){\mathbf{1}}_m+2C_x\left(t\right){\mathbf{p}}_x\left(t\right),\hfill \end{array}$$

(7)

where $S_1\left(x\right)=\{y_1,\dots ,y_m\}$ and

$${\mathbf{p}}_x\left(t\right)={(p_{x{y}_1}\left(t\right),\dots ,p_{x{y}_m}\left(t\right))}^{\top }.$$

In the case of an isolated vertex $x\in V$, its curvature flow is given by the simple equation $p_{xx}^{\prime }\left(t\right)=0$, that is $p_{xx}\left(t\right)$ is a constant function in t.

We note that the curvature flow Equation (6) guarantees that the diagonal entries of the weighting scheme do not change. The functions $C_x\left(t\right)$ in the curvature flow Equation (7) play the role of a normalization since, for the choice $C_x\equiv 0$, various transition rates will be unbounded as the time parameter $t\ge 0$ increases. Note that in the smooth case of closed Riemannian manifolds $(M,g_0)$, a suitable normalization leads to volume perseverance of $(M,g_t)$ under the Ricci curvature flow. Our aim is to preserve the Markovian property, and it was shown in our first paper that the curvature flow preserves this property if we choose the normalization functions

$$C_x\left(t\right)=K_{P\left(t\right),\infty }^{d(x,·)}\left(x\right).$$

(8)

Let us give the explicit formulas for the curvature flow Equation (7) for this particular choice of $C_x\left(t\right)$, where $y,y^{\prime },y^{″}$ always represent vertices in $S_1\left(x\right)$ (see [1] (Formula (66)):

$$\begin{array}{c}{p}_{xy}^{\prime }\left(t\right)=\hfill \\ p_{xy}\left(t\right)\left(-4p_{yx}\left(t\right)-2\sum _{y^{\prime }\ne y}{p}_{y{y}^{\prime }}\left(t\right)+\frac{4}{D_x}\sum _{y^{\prime }}{p}_{x{y}^{\prime }}\left(t\right)p_{y^{\prime }x}\left(t\right)+\frac{1}{D_x}\sum _{y^{\prime },y^{″}}{p}_{x{y}^{\prime }}\left(t\right)p_{y^{\prime }{y}^{″}}\left(t\right)-p_{yy}\left(t\right)\right)\\ \hfill +\underset{(\ast )}{\underbrace{\sum _{y^{\prime }\ne y}{p}_{x{y}^{\prime }}\left(t\right)p_{y^{\prime }y}\left(t\right)}}.\end{array}$$

(9)

Here, we use $D_x={\sum }_{y^{\prime }}{p}_{xy}\left(t\right)=1-p_{xx}\left(t\right)$. Note that (6) guarantees that $D_x\le 1$ is independent of the time parameter t.

The following theorem collects some fundamental properties of the normalized curvature flow.

Theorem 2

(see Theorem 1.5 and Prop. 1.6 in [1]). Let $(G,P_0)$ be a Markovian weighted graph. Then, the curvature flow ${(G,P\left(t\right))}_{t\ge 0}$ associated with $(G,P_0)$ with normalization (8) is well defined for all $t\ge 0$ and preserves the Markovian property. If $(G,P_0)$ is non-degenerate, then $(G,P(t\left)\right)$ is also non-degenerate for all $t\ge 0$. Moreover, if the flow converges for $t\to \infty$ to $P^{\infty }={lim}_{t\to \infty }P\left(t\right)$, then the weighted graph $(G,P^{\infty })$ is curvature sharp.

We would like to emphasize that, even in the Markovian case, a flow limit $P^{\infty }={lim}_{t\to \infty }P\left(t\right)$ of a non-degenerate weighted graph $(G,P_0)$ is in most cases no longer non-degenerate, despite the fact that all weighting schemes $P\left(t\right)$ for finite $t\ge 0$ are non-degenerate. In other words, some transition probabilities converge to zero under our normalized curvature flow, as time tends to infinity.

1.4. Other Discrete Curvature Flows

We presented already other curvature flows in [1] and briefly recall them here for the readers’ convenience.

To our knowledge, the curvature flow on discrete Markov chains in this paper is the first one that is based on Bakry–Émery curvature. Other curvature flow investigations for discrete Ricci-type curvature are the following:

• Combinatorial Ricci flows on surfaces [25];
• Ollivier Ricci flow ([26,27], Problem N) and [28,29];
• Forman-Ricci curvature flow [30,31];
• Stochastic Ricci flow for balanced Forman curvature [32,33];
• Entropic super Ricci flow [34,35];
• Resistance curvature flow [36].

2. Curvature Flow Examples

In this section, we investigate normalized curvature flows on some unmixed combinatorial graphs $G=(V,E)$. By unmixed we mean that G does not have one-sided edges. We assume in all examples and statements in this section that our graphs $G=(V,E)$ are finite, simple, unmixed, and connected and that our initial weighting schemes $P_0=P\left(0\right)={\left(p_{xy}\left(0\right)\right)}_{x,y\in V}$ are non-degenerate Markovian without laziness (even if we do not mention this). Curvature flow limits $(G,P^{\infty })$ with $P^{\infty }={lim}_{t\to \infty }P\left(t\right)$ are necessarily curvature sharp by Theorem 2 above. We do not know of any initial Markovian weighted graph that does not converge as $t\to \infty$. Moreover, we know that every finite connected graph with at least two vertices admits many curvature sharp Markovian weighting schemes without laziness (see [1] (Theorem 1.10)). These facts give rise to the following conjecture.

Conjecture 1

(see Conjecture 1.7 in [1]). The curvature flow ${(G,P\left(t\right))}_{t\ge 0}$ with normalization (8) converges for any initial condition $(G,P_0)$ as $t\to \infty$.

Let us now address some practical aspects of the curvature flow implementation. The Python code is freely available as an ancillary file of [2] or via github at https://github.com/georgestagg/graph-curvature-server (accessed on 6 June 2023). The solution of the curvature flow is computed numerically by the Runge–Kutta (RK4) method, which is based on a time discretization with time increments $dt>0$. In the following examples, we choose the step sizes $dt=0.1$ and $dt=0.3$. In order to distinguish between the theoretical curvature flow and its implementation, we refer to the latter as the numerical curvature flow. Since a numerical curvature flow cannot run forever, a suitable numerical convergence criterion needs to be introduced. Our convergence criterion is based on the parameter ${\lim }_{\mathrm{tolerance}}>0$. We say that a numerical curvature flow solution ${\left(P\left(t\right)\right)}_{t\ge 0}$ has converged numerically at time t (with respect to the parameter ${\lim }_{\mathrm{tolerance}}$) if all the entries of $P\left(t\right)$ differ from the corresponding entries of $P(t+10)$ and $P(t+20)$ by less than ${\lim }_{\mathrm{tolerance}}$. The numerical flow limit is then defined to be $P\left(t\right)$. In all the following examples, we set ${\lim }_{\mathrm{tolerance}}=0.001$.

In this section, we are particularly interested in numerical flow limits that are not numerically totally degenerate. An unmixed weighted graph $(G,P)$ is called numerically totally degenerate if there are no edges with numerical non-zero transition rates in both directions, where we consider a transition rate $p_{xy}$ as numerically non-zero (with respect to a parameter $\mathrm{threshold}>0$), if and only if $p_{xy}\ge \mathrm{threshold}$. In all the following examples we set $\mathrm{threshold}=0.001$.

2.1. Random Graphs

In this subsection, we investigate the numerical curvature flow for random weighted graphs $(G,P_0)$ with vertex set V. The relevance of this example is that it provides insights into flow limits for general graphs without any special symmetries. At the same time, readers are familiarized with our software tool, enabling them to run examples of their own interest. The edge set E of the random graph G is generated by an Erdös-Rényi process [37], that is, any pair of vertices is independently and randomly connected by a two-sided edge with a probability $p\in [0,1]$. Similarly, we choose a random initial weighting scheme $P_0$ with the property that all non-zero transition rates $p_{v{v}^{\prime }}\left(0\right)$ lie in the interval $[\mathrm{threshold},1]$, for some positive parameter $\mathrm{threshold}>0$. We are interested in the properties of the numerical flow limits of these graphs. If the parameter $p>0$ in the Erdös-Rényi process is too small, these flow limits are always numerically totally degenerate. A reasonable choice to obtain not numerically totally degenerate flow limits for such random graphs with, say, 10 vertices in roughly half of the cases, is $p=0.7$.

Example 1

(A random graph with 10 vertices). Let $(G,P_0)$ be the unmixed weighted Markovian graph with vertex set $V=\{v_0,\dots ,v_9\}$ and

$$P_0=\left(\begin{array}{cccccccccc}0& 0.1& 0.08& 0.17& 0& 0.28& 0.21& 0.08& 0& 0.08\\ 0.08& 0& 0& 0.16& 0& 0.2& 0.07& 0.3& 0.04& 0.15\\ 0.27& 0& 0& 0& 0& 0& 0.3& 0& 0& 0.43\\ 0.02& 0.19& 0& 0& 0.17& 0.17& 0.11& 0.34& 0& 0\\ 0& 0& 0& 1& 0& 0& 0& 0& 0& 0\\ 0.04& 0.21& 0& 0.41& 0& 0& 0.34& 0& 0& 0\\ 0.06& 0.29& 0.14& 0.12& 0& 0.3& 0& 0.09& 0& 0\\ 0.08& 0.31& 0& 0.19& 0& 0& 0.23& 0& 0.19& 0\\ 0& 0.13& 0& 0& 0& 0& 0& 0.25& 0& 0.62\\ 0.1& 0.33& 0.38& 0& 0& 0& 0& 0& 0.19& 0\end{array}\right)$$

This randomly generated initial graph is illustrated on the left-hand side of Figure 1. The numerical curvature flow of $(G,P_0)$ has numerical convergence time $t_{\max }=20.7$ (with respect to ${\lim }_{\mathrm{tolerance}}=0.001$). The numerical flow limit $(G,P\left(t_{\max }\right))$ is presented on the right hand side of Figure 1. Let us briefly explain the illustration of the edges of this flow limit: Edges with numerical non-zero transition rates in both directions are displayed in green. Edges with only one-sided numerical non-zero transition rates are displayed as dashed red lines with arrows. Edges whose transition rates shrink in both directions numerically to zero under the curvature flow are displayed as dotted black lines. The corresponding non-zero transition rates are written along these edges. The vertices of the green edges of the flow limit in Figure 1 are given by

$$W=\{v_0,v_1,v_3,v_5,v_6,v_7\}.$$

Since there are no numerical non-zero transition rates from these vertices to the other vertices $v_2,v_4,v_8,v_9\in V\backslash W$, the vertex set W together with the green edges and their transition rates represent a highly connected non-degenerate Markovian weighted subgraph $(G_W,P_W)$. The combinatorial graph $G_W$ with vertex set W can be viewed as a double cone over the complete graph $K_4$ of the four vertices $v_0,v_1,v_3,v_6$ with the vertices $v_5,v_7$ as its cone tips. The transition rates of the weighting scheme $P_W$ towards all vertices in $K_4$ are $0.25=1/4$, all transition rates towards $v_5$ are $x=0.05$ and all transition rates towards $v_7$ are $y=0.2$. Note that such a weighted double cone over $K_4$ with these transition rates for any choice of $x,y>0$ satisfying $x+y=1/4$ is curvature sharp (this follows readily from the geometric criterion in [1] (Theorem 3.15), since the weighting scheme is volume homogeneous in all vertices and reversible with $\pi \left(v_5\right)=4x/5$, $\pi \left(v_7\right)=4y/5$ and ${\pi |}_{K_4}\equiv 1/5$). Therefore, not only is the flow limit itself curvature sharp in Theorem 2 but also the highly connected non-degenerate weighted subgraph $(G_W,P_W)$.

Figure 2 presents the transition rates of the vertices $v_3,v_7$ under the curvature flow as functions over the interval $[0,t_{\max }]$. While most transition rates converge to strictly positive limits, the transition rates $p_{v_3{v}_4}\left(t\right)$ and $p_{v_8{v}_9}\left(t\right)$ shrink to zero. Consequently, the corresponding edges on the right-hand side of Figure 1 are represented by a dashed red line and a black dotted line, respectively.

Let us finally consider the curvatures $t\mapsto K_{P\left(t\right),N}\left(v_j\right)$ of the vertices $v_j$ under the curvature flow. We focus primarily on the vertices $v_1$ and $v_9$ and the dimension parameter $N=\infty$. Figure 3 presents the ∞-curvature (in blue) and upper curvature bound $K_{P\left(t\right),\infty }^{d(v,·)}\left(v\right)$ (in orange) of $v\in \{v_1,v_9\}$ as functions over the interval $[0,t_{\max }]$. Note that the absence of laziness implies that the upper curvature bounds $K_{P\left(t\right),\infty }^{d(v,·)}\left(v\right)$ are all $\ge 0$ (see [1] (19)). For both vertices, the curvature and upper curvature bound functions are numerically asymptotic as $t\to t_{\max }$, indicating that $v_1$ and $v_9$ of the flow limit are ∞-curvature sharp. (In fact, all vertices in this example are ∞-curvature sharp with respect to the flow limit.) This is not always the case, but Theorem 1 confirms that all vertices of a flow limit $(G,P^{\infty })$ are at least N-curvature sharp for dimension $N=2$. Moreover, the initial and final ∞-curvatures of all vertices under the curvature flow are presented in

Table 1. Vertex curvatures for the dimension ∞ of a random graph with 10 vertices at the beginning and the end of the curvature flow.

j	${\mathit{K}}_{\mathit{P}\left(0\right),\mathit{\infty }}\left({\mathit{v}}_{\mathit{j}}\right)$	${\mathit{K}}_{\mathit{P}\left({\mathit{t}}_{\mathbf{max}}\right),\mathit{\infty }}\left({\mathit{v}}_{\mathit{j}}\right)$
0	$0.406$	$0.875$
1	$0.293$	$0.875$
2	$0.168$	$0.125$
3	$0.346$	$0.875$
4	$0.34$	0
5	$0.527$	$0.476$
6	$0.404$	$0.875$
7	$0.202$	$0.773$
8	$0.246$	$0.11$
9	$0.236$	$0.125$

. The final curvatures of all vertices in $K_4$ assume the highest values $0.875$, followed by the final curvature values $0.773$ and $0.476$ of the vertices $v_7$ and $v_5$, respectively. All other vertices in $V\backslash V_0$ have much lower final curvatures with values $\le 0.125$.

Example 1 and its illustrations were generated by the following code.

Let us explain the commands of this program in detail. After initializing various parameters in lines 1–11, a random graph $G=(V,E)$ is generated in line 13 via an Erdös-Renyi process with probability $p=0.7$. A corresponding random non-degenerate Markovian weighting scheme without laziness is generated in line 14, with each non-zero transition rate satisfying $p_{v{v}^{\prime }}\left(0\right)\in [\mathrm{threshold},1]$. The numerical convergence time $t_{\max }\ge 0$ is determined in lines 16–18. In most cases, the convergence time does not exceed 100. (If convergence is not achieved by $t_{\lim }=10000$, the curvature flow computation stops at that time and notifies the user.) The numerical curvature flow on the interval $[0,t_{\max }]$ is solved again in line 21. Initial and final weighting schemes are displayed by the commands in lines 23 and 24, respectively, providing the illustrations given in Figure 1. The transition rates of the vertices $v_3$ and $v_8$ are displayed by the command in line 25, providing the illustrations given in Figure 2. The curvatures and upper curvature bounds at time steps $j·\frac{k·t_{\max }}{dt}$, $j=0,1,\dots$, are computed in lines 27 and 28, respectively, with the choice $k=1$. They are displayed for the vertices $v_1$ and $v_9$ via the command in line 29, providing the illustrations given in Figure 3.

When running this program, users may be faced with the following message:

‘norm_tolerance’ has been exceeded at one or more vertices, at time t = … Would you like to:

A = Stop calculation and return list of P-matrices so far

B = Apply manual normalization now, and apply it again when necessary without asking (you will still be notified when it is applied)

C = Apply manual normalization now, and ask again before reapplying it

Please enter A, B or C here:

The reason behind this message is the following. The computation of the numerical curvature flow is based on a time discretization. Therefore, the solution will increasingly depart from the Markovian property after each time increment $dt=0.1$. If the sum of entries of one of the rows of $P\left(t\right)$ at time t differs from one by more than ${\mathrm{norm}}_{\mathrm{tolerance}}=0.001$, the program informs the user that a normalization of the weighting scheme is needed for the continuation of the flow calculations. After choosing the option ‘B’, the program will continue with its flow calculations without further interruptions, and the user is simply notified about the times at which the program applies further artificial normalizations of the transition rates. The user can suppress this message entirely by changing line 2 of the program into “stoch_corr = True”, in which case the program applies stochastic corrections automatically, each time with the message

Transition rates have been artificially normalized at time t = …

The numerical observations of Example 1 suggest that similar properties may also hold in general in the theoretical setting. Firstly, we call an edge $\{x,y\}\in E$ in an unmixed weighted graph $(G,P)$non-degenerate if $p_{xy},p_{yx}>0$. An unmixed weighted graph $(G,P)$ is called totally degenerate if it does not have any non-degenerate edges. Secondly, we denote the vertices of all non-degenerate edges of a flow limit $(G,P^{\infty })$ by $W\subset V$, and $G_W$ denotes the subgraph of G consisting of the vertices W and all non-degenerate edges of $(G,P^{\infty })$. Moreover, $P_W$ denotes the restriction of the weighting scheme $P^{\infty }$ to the vertex set W. We call $G_W$ the non-degenerate subgraph of $(G,P^{\infty })$ and we conjecture the following:

Conjecture 2.

If the normalized curvature flow of a nondegenerate unmixed Markovian weighted graph $(G,P_0)$ converges to a not totally degenerate limit $(G,P^{\infty })$, then the non-degenerate subgraph $G_W$ coincides with the induced subgraph (of G) of the subset $W\subset V$ and all transition rates from W to $V\backslash W$ are zero. $(G_W,P_W)$ is a non-degenerate Markovian weighted graph, which is itself curvature sharp.

Figuratively speaking, the curvature flow converges towards the non-degenerate Markovian subgraph $(G_W,P_W)$ consisting of highly connected components. Moreover, each vertex of $V\backslash W$ is usually connected to the set W by a sequence of edges with one-sided non-zero transition rates pointing towards the set W, and we generally expect that the ∞-curvature values of the vertex set W in $(G,P^{\infty })$ are significantly larger than the ∞-curvature values of the set $V\backslash W$ in $(G,P^{\infty })$.

In Example 1, the weighted Markovian subgraph $(G_W,P_W)$ has only one connected component, but we will see in the next subsection in the case of paths and cycles that $(G_W,P_W)$ may be composed of more than just one connected component.

2.2. Paths and Cycles

Paths and cycles are the easiest examples of graph families and are therefore natural examples to study. We will see that both families have flow limits of very similar types. Let $G=(V,E)$ be a path of length $N\ge 2$, that is, $V=\{v_0,\dots ,v_{N-1}\}$ with a two-sided edge between $v_i$ and $v_j$ if and only if $|i-j|=1$. If $N=2$, G is a trivial case of a star graph and any Markovian weighting scheme satisfying $p_{01}=p_{21}=1$ and $p_{10}+p_{12}=1$ is curvature sharp (see [1], Example 4.3). For that reason, we consider only paths of lengths $N\ge 3$. A cycle of length 3 is the complete graph $K_3$, and a full list of all curvature sharp weighting schemes was given in [1], Prop. 1.9, so we consider only cycles of length $N\ge 4$. The following example provides some insights into some features of curvature flow limits of weighted paths and cycles.

Example 2

(A path and a cycle of length 12). Figure 4 presents numerical curvature flow limits of a Markovian weighted path with vertices 12 vertices (left-hand side) and of a weighted cycle with 12 vertices (right-hand side). The non-zero transition rates of the initial weighting scheme $P_0={\left(p_{ji}\right)}_{0\le i,j\le 11}$ for the path limit in Figure 4 were chosen as follows:

$p_{0,1}$	$p_{1,2}$	$p_{2,3}$	$p_{3,4}$	$p_{4,5}$	$p_{5,6}$	$p_{6,7}$	$p_{7,8}$	$p_{8,9}$	$p_{9,10}$	$p_{10,11}$
1	$0.25$	$0.72$	$0.46$	$0.23$	019	$0.84$	$0.71$	$0.62$	$0.9$	$0.55$
$p_{1,0}$	$p_{2,1}$	$p_{3,2}$	$p_{4,3}$	$p_{5,4}$	$p_{6,5}$	$p_{7,6}$	$p_{8,7}$	$p_{9,8}$	$p_{10,9}$	$p_{11,10}$
$0.75$	$0.28$	$0.54$	$0.77$	$0.81$	$0.16$	$0.29$	$0.38$	$0.1$	$0.45$	1

The non-degenerate subgraph $G_W$ of the path limit consists of $W=\{v_1,v_2,v_3,v_9,v_{10},v_{11}\}$ as its vertex set together with the four green edges. Moreover, $G_W$ has two paths of length 2 as its connected components. For each vertex $v\in V\backslash W$, there exists a directed path to one of the components of $G_W$ via a sequence of one-sided transition rates. For example, the vertices $v_5,v_6,v_7,v_8$ are connected to the vertex $v_9\in W$ via such one-sided paths, and $v_4,v_5$ are connected to the vertex $v_3$ via such one-sided paths.

The cycle limit on the right-hand side of Figure 4 is totally degenerate with no green edges, and all one-sided non-zero transition rates are oriented in a clockwise direction. The experiments show that Figure 4 exhibits generic limit properties: flow limits of weighted paths are never totally degenerate and their non-degenerate subgraphs $G_W$ consist of disjoint paths of length $\le 2$. Such types of limits also appear in the case of weighted cycles. However, in contrast to the path case, sometimes a cycle limit is totally degenerate with all its one-sided transition rates oriented either clockwise or anti-clockwise. The code for running the curvature flow for paths and cycles of length n reads as follows:

Before we provide the following result about path limits, let us first introduce the notion of a two-sided degenerate edge of weighted graphs $(G,P)$: An edge $\{x,y\}\in E$ of G is called two-sided degenerate if its transition rates vanish in both directions, that is, $p_{xy}=p_{yx}=0$. Similarly, an edge is called numerically two-sided degenerate, if its transition rates in both directions are $<\mathrm{threshold}$. Such edges are displayed in the routine display_weighted_graph as dotted black lines (see, e.g., the edges $\{v_2,v_9\}$ and $\{v_8,v_9\}$ on the right hand side of Figure 1).

Proposition 1

(Flow limits of paths of length $\ge 3$). Let $(G,P_0)$ be a weighted path of length $N\ge 3$ with consecutive vertices $v_0,\dots ,v_{N-1}$. Let $P\left(t\right)$ be its corresponding curvature flow converging to a limit $P^{\infty }={lim}_{t\to \infty }P\left(t\right)$, such that $(G,P^{\infty })$ does not have two-sided degenerate edges. Then, this limit is neither totally degenerate nor is it not non-degenerate (that is, it contains both green and dashed red edges). Moreover, the components of the non-degenerate subgraph $G_W$ are paths of lengths $\le 2$, and they are separated from each other by at least two degenerate edges. If a component of $G_W$ is a path of length 1, that is, just one edge $\{v_j,v_{j+1}\}$, then we have $p_{j,j+1}^{\infty }=p_{j+1,j}^{\infty }=1$.

Proof. 

By the last statement of Theorem 2, we only need to prove these statements for curvature sharp weighting schemes $P={\left(p_{ij}\right)}_{0\le i,j\le N}$ of paths of length $N\ge 3$. Such weighting schemes can never be non-degenerate by [1] (Prop. 1.11).

For the proof that curvature sharp weighting schemes can never be totally degenerate, we note that we cannot have two consecutive degenerate edges $\{v_j,v_{j+1}\},\{v_{j+1},v_{j+2}\}\in E$ with $p_{j,j+1}=1=p_{j+2,j+1}$ (since this would imply $p_{j+1,j}=p_{j+1,j+2}=0$, contradicting to $p_{j+1,j}+p_{j+1,j+2}=1$). Therefore, since $p_{01}=1$, any totally degenerate weighting scheme would require $p_{j,j+1}=1$ for any $j\ge 0$, in particular, $p_{N-2,N-1}=1$, but this contradicts the fact that we also have $p_{N-1,N-2}=1$ and that the edge $\{v_{N-1},v_{N-2}\}$ needs to be degenerate.

A curvature sharp weighting scheme cannot have more than two consecutive non-degenerate edges. This can be seen as follows: At any vertex $v_j$ with $p_{j,j-1},p_{j,j+1}>0$ we must have $p_{j-1,j}=p_{j+1,j}$. This follows from the arguments in the proof of [1] (Lemma 4.1) (namely, since the vertex $v_j$ is not contained in any triangle, we have $p_{j-1,j}=p_{y+1,j}=p_{j,j-1}{p}_{j-1,j}+p_{j,j+1}{p}_{j+1,j}$.) If $0<p_{j-1,j}=p_{j+1,j}<1$, we could iterate this argument backward and forward and would end up with the fact that all entries of the matrix P above the diagonal and below the diagonal would lie in $(0,1)$, which is a contradiction to $p_{01}=1$. So, we must have $p_{j-1,j}=p_{j+1,j}\in \{0,1\}$. Therefore, we cannot have two consecutive indices $j\in \{0,\dots ,N-1\}$ with $0<p_{j-1,j}=p_{j+1,j}<1$, which would exist in the case of three consecutive non-degenerate edges. So, the components of $G_W$ are paths of length $\le 2$.

Moreover, any gap between two consecutive non-degenerate edges must be at least two edges: this follows from the fact that if a non-degenerate edge $\{v_j,v_{j+1}\}$ is followed by a degenerate edge $\{v_{j+1},v_{j+2}\}$, then we must have $p_{j+1,j+2}=0$: if we had $p_{j+1,j+2}>0$, then we had $p_{j+1,j},p_{j+1,j+2}>0$ and, therefore, $0<p_{j+1,j}=p_{j+2,j+1}$ and $\{v_{j+1},v_{j+2}\}$ would be non-degenerate, which is a contradiction. Similarly, if a degenerate edge $\{v_k,v_{k+1}\}$ is followed by a non-degenerate edge $\{v_{k+1},v_{k+2}\}$, we must have $p_{k+1,k}=0$. Combining both facts implies that there cannot be a single degenerate edge separating two components of $G_W$. These arguments also show that components of $G_W$ that are single edges $\{v_j,v_{j+1}\}$ must satisfy $p_{j,j+1}=p_{j+1,j}=1$ since the adjacent degenerate edges have one-sided transition rates pointing towards this component.    □

Similar arguments to the above can be used to prove that for cycles, the components of any flow limit of a weighted cycle are again paths of length $\le 2$, unless the limit is non-degenerate. Such non-degenerate limits exist for cycles, namely the simple random walks, but experiments show that simple random walks are very unstable stationary solutions of the curvature flow. Small perturbations of simple random walks do not converge back to the simple random walk (unless our cycle is $K_3$) but converge usually to a degenerate limit. Finally, if a totally degenerate cycle limit does not have two-sided degenerate edges, then its transition rates must all be either oriented clockwise or anti-clockwise; otherwise, we would necessarily have a vertex with transition rates of both incident degenerate edges pointing towards this vertex. This would fail to satisfy the Markovian condition at this vertex.

Paths and cycles of length $N\ge 3$ have the property that no edge is contained in a triangle. We would like to finish the section with a general statement about the curvature flow for edges not contained in triangles.

Proposition 2.

Let $(G,P_0)$ be a weighted Markovian graph without laziness and ${\left(P\left(t\right)\right)}_{t\ge 0}$ be its associated normalized curvature flow. If we have, for some $t_0\ge 0$ and an edge $e=\{x,y\}\in E$, $p_{xy}\left(t_0\right)=0$ and e is not contained in a triangle of G, then we have

$$p_{xy}\left(t\right)=0\mathit{for}\mathit{all}t\ge t_0.$$

Proof. 

This proposition is an easy consequence of the flow Equation (9). Since we assume no laziness, we have $p_{yy}\left(t\right)=0$, and since e is not contained in a triangle, the last term of (9), denoted by $(\ast )$, is zero and the statement follows now from the uniqueness of the solution satisfying $p_{xy}\left(t_0\right)=0$.    □

2.3. Complete Graphs

Complete graphs are the natural choice for the study of the behavior of curvature flows of highly connected graphs. For these graphs, the simple random walks turn out to be the flow limits of non-degenerate initial weighting schemes. The simple random walk on any complete graph is a non-degenerate curvature sharp Markovian weighting scheme. Our experiments show that any non-degenerate initial weighting scheme $P_0$ on $K_n$ converges to the simple random walk, that is $p_{jk}^{\infty }=\frac{1}{n-1}$. The following example shows that convergence to the simple random walk appears even if the initial weighting scheme is degenerate. This is not in contradiction to Proposition 2 since, in a complete graph $K_n$, $n\ge 3$, every edge is contained in a triangle. Example 3 is the only exception in this section where we allow an initial weighting scheme to have degenerate edges.

Example 3

(A degenerate weighted complete graph with six vertices). Let $(G,P_0)$ be the complete weighted Markovian graph with vertex set $V=\{v_0,\dots ,v_5\}$ and

$$P_0=\left(\begin{array}{cccccc}0& 0.2& 0.1& 0.2& 0& 0.5\\ 0.1& 0& 0.3& 0.25& 0.25& 0.1\\ 0.2& 0& 0& 0.3& 0.15& 0.35\\ 0.3& 0.5& 0.1& 0& 0.1& 0\\ 0.2& 0.3& 0.3& 0.2& 0& 0\\ 0.6& 0.1& 0.1& 0.2& 0& 0\end{array}\right),$$

as illustrated in Figure 5. Note that the edges $\{v_0,v_4\},\{v_1,v_2\},\{v_3,v_5\}$ and $\{v_4,v_5\}$ of this initial weighting scheme are degenerate, with the latter being two-sided degenerate. The numerical curvature flow has numerical convergence time $t_{\max }=18.5$ (with respect to ${\lim }_{\mathrm{tolerance}}=0.001$) with the simple random walk as its numerical flow limit. Figure 6 presents the transition rates of the vertices $v_2$ and $v_5$. Of particular interesting are the functions $p_{2,1}\left(t\right)$ and $p_{5,4}\left(t\right)$ for $t\in [0,t_{\max }]$, since their initial values are zero.

Based on our experiments, we conjecture the following, which is a strengthening of [1] (Conjecture 1.8).

Conjecture 3

(Curvature flow of complete graphs). Let $P_0$ be a Markovian weighting scheme without laziness on a complete graph $K_n=(V,E)$ with $n\ge 2$ such that, for every proper subset $W\subset V$, there exist $x,x^{\prime }\in W$ and $y,y^{\prime }\in \in V\backslash W$ with $p_{xy}>0$ and $p_{y^{\prime }{x}^{\prime }}>0$. Then, the curvature flow has a limit $P^{\infty }$, which is the simple random walk.

2.4. Wedge Sums of Complete Graphs

Let $G_1=(V_1,E_1)$ and $G_2=(V_2,E_2)$ be two combinatorial graphs and $x_1\in V_1$ and $x_2\in V_2$. By merging the vertices $x_1$ and $x_2$ into one new vertex x, which inherits the incident edges of both vertices $x_1$ and $x_2$, we obtain a new combinatorial graph in which we denote the wedge sum of $G_1$ and $G_2$. In this subsection, we consider the flow limits of wedge sums of complete graphs. The study of these graphs leads to an interesting dynamical aspect of our curvature flow, namely that its limits concentrate on one of the components of these wedge sums.

Example 4

(A wedge sum of a $K_4$, $K_5$, $K_2$, and $K_3$). We consider the curvature flow on the wedge sum $G=(V,E)$ of complete graphs presented in Figure 7 with random initial weighting schemes $P=P_0$. The adjacency matrix A of this graph is generated with our code in the following way:

Our experiments show that, depending on the initial weighting scheme $P_0$, the curvature flow converges to a limit that is concentrated in one of the complete graphs. More precisely, the limit weighting scheme $P^{\infty }$ represents a simple random walk on one of the complete graphs while, from all other vertices, there is a directed path of $\{0,1\}$ transition rates towards this particular complete graph. Figure 8, Figure 9 and Figure 10 show the numerical curvature flow limits of various random initial weighting schemes. We carried out several runs of 100,000 numerical curvature flows with random initial weighting schemes to describe the limit behavior of these flows quantitatively. The results are presented in

Table 2. Statistics of flow limit concentrations on the complete subgraphs $K_4,K_5,K_2$ and $K_3$ of their wedge sum, together with mean convergence times.

	${\mathit{K}}_4$	${\mathit{K}}_5$	${\mathit{K}}_2$	${\mathit{K}}_3$
concentration of flow limits	$8.7\%$	$80.5\%$	$0\%$	$10.8\%$
mean convergence time	$18.9$	$24.7$	–	$17.7$

. While more than $80\%$ of limits concentrate on the largest clique $K_5$, it is somewhat surprising that more limits concentrate on $K_3$ than on the larger subgraph $K_4$. Not a single flow limit ended up concentrating on $K_2$. The mean numerical convergence time is shortest for $K_3$, followed by $K_4$ and $K_5$ (with respect to ${lim}_{\mathrm{tolerance}}=0.001$). While most convergence times are below 100, there were maximal numerical convergence times well above 500.

The limit behavior described in the above example seems to be common for many wedge sums of complete graphs. It is, however, not always true that flow limits concentrate on just one of the constituents of a wedge sum. A path of length $\ge 2$ can be viewed as a wedge sum of consecutive $K_2$’s, and we have seen in Section 2.2 that flow limits will concentrate on more than only one of these $K_2$’s (see left-hand side of Figure 4). Another special case of a wedge sum is a dumbbell which is our next example.

Example 5

(A symmetric weighted dumbbell). The weighted graph $(G,P_0)$ in this example is a wedge sum of a $K_5$, $K_2$ and another $K_5$, together with a simple random walk as initial weighting scheme (see line 4 in the code). This situation can be set up by the following code:

The numerical convergence time is 79 and the limit of the numerical curvature flow concentrates on the “bridge” $K_2$ between the two $K_5$’s, as illustrated in Figure 11. (To obtain the flow limit illustrated at the right-hand side of this figure, users should choose ${lim}_{\mathrm{tolerance}}=0.0001$.)

This limit could have been predicted assuming that the initial weighted graph symmetry across the “bridge” is preserved under the curvature flow and that the limit concentrates on only one of the complete graphs $K_5,K_2,K_5$. If instead of the simple random walk, a random initial weighting scheme would have been chosen on the dumbbell G, the limit would have usually concentrated on one of the two $K_5$’s.

2.5. Cartesian Products of Complete Graphs

It is tempting to assume that simple random walks are always the preferred curvature flow limits in relation to complete graphs. However, this is not always the case as the following example shows. Namely, we have the following result for Cartesian products of complete graphs:

Theorem 3.

Let G be the Cartesian product of two complete graphs $K_{n+1}$ and $K_{m+1}$ with the non-lazy simple random walk P as the initial weighting scheme. Then, the curvature flow $P\left(t\right)$ converges to a limit as $t\to \infty$ with limit transition rates

$$\begin{array}{ccc}\hfill a^{\infty }& =& \frac{m+3}{2nm+3n+3m},\hfill \\ \hfill b^{\infty }& =& \frac{n+3}{2nm+3n+3m},\hfill \end{array}$$

where $a^{\infty }$ are the transition rates along edges between $(x,x^{\prime })$ and $(y,x^{\prime })$ with $x,y\in K_{n+1}$, $x\sim y$ and $x^{\prime }\in K_{m+1}$ (“horizontal edges”), and $b^{\infty }$ are the transition rates along edges between $(x,x^{\prime })$ and $(x,y^{\prime })$ with $x\in K_{n+1}$ and $x^{\prime },y^{\prime }\in K_{m+1}$, $x^{\prime }\sim y^{\prime }$ (“vertical edges”).

Proof. 

By symmetry of the configuration, we have only two types $a\left(t\right)$ and $b\left(t\right)$ of transition rates for the curvature flow at time $t\ge 0$ (for horizontal and vertical edges, respectively) with

$$a\left(0\right)=b\left(0\right)=\frac{1}{n+m}$$

and

$$na\left(t\right)+mb\left(t\right)=1,$$

due to the Markovian property. This implies that $b\left(t\right)=\frac{1-na\left(t\right)}{m}$ and $b^{\prime }\left(t\right)=-\frac{-n{a}^{\prime }\left(t\right)}{m}$, and we only need to consider an ordinary differential equation for a with initial condition $a\left(0\right)=\frac{1}{n+m}$. We can also assume without loss of generality that $n\le m$. We derive from the explicit description (9) of the curvature flow that

$$\begin{array}{c}{a}^{\prime }\left(t\right)=\hfill \\ a\left(t\right)\left(-4a\left(t\right)-2(n-1)a\left(t\right)+4(n{a}^2\left(t\right)+m{b}^2\left(t\right))+n(n-1)a^2\left(t\right)+m(m-1)b^2\left(t\right)\right)\\ +(n-1)a^2\left(t\right)=n(n+3)a^3\left(t\right)+m(m+3)a\left(t\right){\left(\frac{1-na\left(t\right)}{m}\right)}^2-(n+3)a^2\left(t\right)=\\ \hfill a\left(t\right)\left(n\left(2n+3+\frac{3n}{m}\right)a^2\left(t\right)-\left(3n+3+\frac{6n}{m}\right)a\left(t\right)+1+\frac{3}{m}\right).\end{array}$$

(10)

If this differential equation converges as $t\to \infty$, its limit $a^{\infty }$ must satisfy

$$a^{\infty }\in \left\{0,\frac{1}{n},\frac{m+3}{2nm+3n+3m}\right\},$$

that is

$$(a^{\infty },b^{\infty })\in \left\{\left(0,\frac{1}{m}\right),\left(\frac{1}{n},0\right),\left(\frac{m+3}{2nm+3n+3m},\frac{n+3}{2nm+3n+3m}\right)\right\},$$

with

$$b^{\infty }=\underset{t\to \infty }{lim}b\left(t\right)=\underset{t\to \infty }{lim}\frac{1-na\left(t\right)}{m}=\frac{1-n{a}^{\infty }}{m}.$$

Our assumption $n\le m$ implies that we have

$$\frac{1}{n+m}\le \frac{m+3}{2nm+3n+3m}<\frac{1}{n}$$

and that the right hand side of (10) is strictly positive for $a\left(t\right)$ in the interval

$$\left[\frac{1}{n+m},\frac{m+3}{2nm+3n+3m}\right),$$

zero at $a\left(t\right)=\frac{m+3}{2nm+3n+3m}$, and strictly negative for $a\left(t\right)$ on the interval

$$\left(\frac{m+3}{2nm+3n+3m},\frac{1}{n}\right).$$

These monotonicity properties force the function $a\left(t\right)$ to converge to the limit $a^{\infty }=\frac{m+3}{2nm+3n+3m}$.    □

Example 6

(Flow limit of $K_3\times K_4$ with simple random walk). The following code computes the numerical flow limit for $K_3\times K_4$ with the simple random walk as initial weighting scheme.

The initial transition rates are all equal to $1/5=0.2$, and the transition rates of the numerical limit are $0.22\cdots \approx 2/9$ and $0.19\cdots \approx 5/27$, as predicted by the above theorem.

2.6. Totally Degenerate Flow Limits

As mentioned before, many curvature limits are totally degenerate. Therefore, it is worth investigating the properties of those particular limits.

Example 7

(The octahedron with a totally degenerate flow limit). Let $G=(V,E)$ with $V=\{v_0,\dots ,v_5\}$ be the unmixed graph representing the octahedron, as illustrated in Figure 12 (left hand side). While the simple random walk without laziness is a stationary solution of the normalized curvature flow, any small perturbation of this initial weighting scheme leads to another curvature sharp limit, which is totally degenerate. For example, the initial weighting scheme

$$P_0=\left(\begin{array}{cccccc}0& 0.26& 0& 0.24& 0.25& 0.25\\ 0.25& 0& 0.25& 0& 0.25& 0.25\\ 0& 0.25& 0& 0.25& 0.25& 0.25\\ 0.25& 0& 0.25& 0& 0.25& 0.25\\ 0.25& 0.25& 0.25& 0.25& 0& 0\\ 0.25& 0.25& 0.25& 0.25& 0& 0\end{array}\right)$$

(11)

converges to the totally degenerate limit illustrated in Figure 12 (right-hand side) under the numerical curvature flow.

The following considerations show that the limit in Example 7 is essentially the only totally degenerate curvature sharp weighting scheme without two-sided degenerate edges for the octahedron (see Proposition 3 below). We start with a general unmixed combinatorial graph $G=(V,E)$ without isolated vertices. A totally degenerate weighting scheme P assigns to each edge $\{x,y\}\in E$ either a direction ($x\to y$ if $p_{xy}>0$ and $y\to x$ if $p_{yx}>0$, illustrated by a red dashed line with an arrow) or the edge is two-sided degenerate (that is $p_{xy}=p_{yx}=0$, illustrated by a black dotted line). A first observation is that none of the vertices $x\in V$ can be a sink: the Markovian property requires that at least one edge incident to x must be outward directed. The following lemma presents a useful property of triangles.

Lemma 1.

Let $(G,P)$ be a totally degenerate curvature sharp Markovian weighted graph without laziness. Then, the directions of a triangle $T=\{x,y,z\}\subset V$ without two-sided degenerate edges cannot be oriented, that is, its edges cannot have the orientations $x\to y\to z\to x$ or $x\to z\to y\to x$.

Proof. 

Assume that a totally degenerate curvature sharp Markovian weighting scheme without laziness contains a triangle $T=\{x,y,z\}$ with $p_{xz},p_{zy},p_{yx}>0$, that is, we have an orientation $x\to z\to y\to x$. This means, in particular, that $p_{xy}=0$. We can read off the flow equation (9) that curvature sharpness of a totally degenerate Markovian weighting scheme means

$$p_{xy}(-2\sum _{y^{\prime }\ne y}{p}_{y{y}^{\prime }}+\sum _{y^{\prime },y^{″}}{p}_{x{y}^{\prime }}{p}_{y^{\prime }{y}^{″}})+\sum _{y^{\prime }\ne y}{p}_{x{y}^{\prime }}{p}_{y^{\prime }y}=0.$$

(12)

Since we have $p_{xy}=0$, this means

$$0\le p_{xz}{p}_{zy}\le \sum _{y^{\prime }\ne y}{p}_{x{y}^{\prime }}{p}_{y^{\prime }y}=0,$$

in contradiction to $p_{xz},p_{zy}>0$. The orientation $x\to y\to z\to x$ can be ruled out similarly.    □

The main tool in the proof of the following proposition can be found in [38] (Exercise 25.14). Assume a tessellation of the two-dimensional sphere carries an orientation along all its edges. For every vertex v of the tessellation, let $\mathrm{ind}\left(v\right)=1-c\left(v\right)/2$, where $c\left(v\right)$ is the number of changes in the orientation of edges adjacent to v (in cyclic order). For a face f, let $\mathrm{ind}\left(f\right)=1-c\left(f\right)/2$, where $c\left(f\right)$ is the number of changes in the orientation (clockwise vs. anti-clockwise) of edges of f. Then, we have

$$\sum _v\mathrm{ind}\left(v\right)+\sum _f\mathrm{ind}\left(f\right)=2.$$

Proposition 3.

Let $G=(V,E)$ be the octahedron, as illustrated in Figure 12 (left-hand side). Then, G has essentially only one totally degenerate non-lazy curvature sharp Markovian weighting scheme without two-sided degenerate edges, namely, we have (up to a permutation of the vertices corresponding to a graph automorphism)

$$p_{v_0,v_1}=p_{v_1,v_2}=p_{v_2,v_3}=p_{v_3,v_0}=1$$

and

$$p_{v_4,v_0}=p_{v_4,v_1}=p_{v_4,v_2}=p_{v_4,v_3}=1/4,p_{v_5,v_0}=p_{v_5,v_1}=p_{v_5,v_2}=p_{v_5,v_3}=1/4.$$

Proof. 

We can think of the octahedron as a tessellation of the sphere by eight triangles, all of them not oriented, by Lemma 1. This means that $\mathrm{ind}\left(f\right)=0$ for all faces of the octahedron. Therefore, we must have

$$\sum _v\mathrm{ind}\left(v\right)=2,$$

that is, at least two vertices of the octahedron must have $\mathrm{ind}\left(v\right)=1$, which means that each of them must be a source or a sink. The Markovian property rules out sinks, and two sources cannot be adjacent (otherwise they would be connected by a non-degenerate edge). Therefore, the octahedron must have two sources at distance 2, which we denote by $v_4,v_5$. Their edges are all directed towards a cycle of length 4. The directions of the edges in this cycle must be directed, for otherwise there would be a sink in this cycle, which is not possible. We denote this oriented cycle by $v_0\to v_1\to v_2\to v_3\to v_0$, and we must have

$$p_{01}=p_{12}=p_{2,3}=p_{3,0}=1,$$

where we used the notation $p_{ij}=p_{v_i,v_j}$ for simplicity. Applying formula (12) to $(x,y)=(v_4,v_0)$ yields

$$p_{40}(-2(p_{01}+p_{03})+\sum _{i=0}^3{p}_{4i})+p_{41}{p}_{10}+p_{43}{p}_{30}=p_{40}(-2+1)+p_{43}=0,$$

that is $p_{40}=p_{43}$. Similarly, we can show that all transition rates $p_{4i}$ must coincide and, therefore, $p_{40}=p_{41}=p_{42}=p_{43}=1/4$. The same arguments apply to the vertex $v_5$, finishing the proof of the proposition.    □

Conjecture 4

(Flow limits of the octahedron). The normalized curvature flow on the octahedron for any non-degenerate initial weighting scheme without laziness different from the simple random walk always converges to a limit described in Proposition 3.

3. Asymptotically Stable and Unstable Curvature Sharp Markovian Weighting Schemes

Our curvature flow can be viewed as a continuous dynamical system with its limits as equilibrium states. This leads naturally to stability questions at these equilibria, which is the concern of this section. Further information about the relevance and the background theory of stability theory for dynamical systems can be found, for example, in [39] (Chapter 9).

Since every curvature sharp Markovian weighting scheme on a given combinatorial graph is a stationary solution of the normalized curvature flow, the question arises whether such a stationary solution $P^s$ is asymptotically stable, that is, whether any closeby Markovian weigthing scheme P converges back to this equilibrium $P^s$ as $t\to \infty$. This can be decided via the linearization of the curvature flow equations around such an equilibrium. In the next subsection, we will describe this linearization in full detail before we consider various examples in the following subsection.

3.1. Linearization of the Curvature Flow Equations at Equilibria

Let $P^s$ be a curvature sharp Markovian weighting scheme of a finite simple mixed combinatorial graph $G=(V,E)$. We consider Markovian weighting schemes P near $P^s$ with the same laziness, that is, $p_{xx}=p_{xx}^s$ for all $x\in V$. Let $E^{\mathrm{dir}}=\{(x,y)\in V\times V:d(x,y)=1\}$. The linearization of F at the equilibrium $P^s$ of the normalized curvature flow is given for each component function $F_{xy}$, $(x,y)\in E^{\mathrm{dir}}$, by

$$\begin{array}{c}D{F}_{xy}\left(P^s\right){\left(q_{uv}\right)}_{(u,v)\in E^{\mathrm{dir}}})=\hfill \\ \left(-4p_{yx}^s-2\sum _{y^{\prime }\ne y}{p}_{y{y}^{\prime }}^s+\frac{4}{D_x}\sum _{y^{\prime }}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }x}^s+\frac{1}{D_x}\sum _{y^{\prime },y^{″}}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }{y}^{″}}^s-p_{yy}^s\right)q_{xy}\\ +p_{xy}^s\left(-4q_{yx}-2\sum _{y^{\prime }\ne y}{q}_{y{y}^{\prime }}+\frac{4}{D_x}\sum _{y^{\prime }}\left(p_{x{y}^{\prime }}^s{q}_{y^{\prime }x}+p_{y^{\prime }x}^s{q}_{x{y}^{\prime }}\right)+\frac{1}{D_x}\sum _{y^{\prime },y^{″}}\left(p_{x{y}^{\prime }}^s{q}_{y^{\prime }{y}^{″}}+p_{y^{\prime }{y}^{″}}^s{q}_{x{y}^{\prime }}\right)\right)\\ \hfill =\sum _{y^{\prime }\in S_1\left(x\right)}{B}_{xy}\left(x{y}^{\prime }\right)q_{x{y}^{\prime }}+\sum _{y^{\prime }\in S_1\left(x\right)}\sum _{z\in B_1\left(x\right)}{B}_{xy}\left(y^{\prime }z\right)q_{y^{\prime }z}.\end{array}$$

Here, $y^{\prime },y^{\prime }$ are vertices in $S_1\left(x\right)$ and the potentially non-zero B-coefficients are given by

$$\begin{array}{cc}\hfill B_{xy}\left(xy\right)=& -4p_{yx}^s-p_{yy}^s-2\sum _{y^{\prime }\ne y}{p}_{y{y}^{\prime }}^s\hfill \\ & +\frac{1}{D_x}\left(4p_{xy}^s{p}_{yx}^s+4\sum _{y^{\prime }}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }x}^s+p_{xy}^s\sum _{y^{\prime }}{p}_{y{y}^{\prime }}^s+\sum _{y^{\prime },y^{″}}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }{y}^{″}}^s\right),\hfill \end{array}$$

(13)

$$\begin{array}{ccc}\hfill B_{xy}\left(x{y}^{\prime }\right)& =& \frac{4}{D_x}{p}_{xy}^s{p}_{y^{\prime }x}^s+\frac{1}{D_x}{p}_{xy}^s\sum _{y^{″}}{p}_{y^{\prime }{y}^{″}}^s+p_{y^{\prime }y}^s\mathrm{if}{y}^{\prime }\ne y,\hfill \end{array}$$

(14)

$$\begin{array}{ccc}\hfill B_{xy}\left(yx\right)& =& 4p_{xy}^s\left(\frac{1}{D_x}{p}_{xy}^s-1\right),\hfill \end{array}$$

(15)

$$\begin{array}{ccc}\hfill B_{xy}\left(y^{\prime }x\right)& =& \frac{4}{D_x}{p}_{xy}^s{p}_{x{y}^{\prime }}^s\mathrm{if}{y}^{\prime }\ne y,\hfill \end{array}$$

(16)

$$\begin{array}{ccc}\hfill B_{xy}\left(y{y}^{\prime }\right)& =& p_{xy}^s\left(\frac{1}{D_x}{p}_{xy}^s-2\right)\mathrm{if}{y}^{\prime }\ne y,\hfill \end{array}$$

(17)

$$\begin{array}{ccc}\hfill B_{xy}\left(y^{\prime }y\right)& =& p_{x{y}^{\prime }}^s\left(\frac{1}{D_x}{p}_{xy}^s+1\right)\mathrm{if}{y}^{\prime }\ne y,\hfill \end{array}$$

(18)

$$\begin{array}{ccc}\hfill B_{xy}\left(y^{\prime }{y}^{″}\right)& =& \frac{1}{D_x}{p}_{xy}^s{p}_{x{y}^{\prime }}^s\mathrm{if}{y}^{\prime }\ne y\mathrm{and}{y}^{″}\ne y,y^{\prime }.\hfill \end{array}$$

(19)

All other B-coefficients are chosen to be zero. Note, however, that the transition probabilities ${\left(p_{xy}\right)}_{y\in S_1\left(x\right)}$ are not independent, and therefore, the choice of the B-coefficients is not unique, as explained in the following remark.

Remark 1.

Since we have ${\sum }_{v\in S_1\left(u\right)}{q}_{uv}=0$ for all $u\in V$, there is a degree of freedom in the choice of the B-coefficients. For example, in the case that G is an unmixed complete graph, we can replace $B_{xy}(u,v)$ by $B_{xy}^{\prime }(u,v)=C_u+B_{xy}(u,v)$ with arbitrary constants $C_u$. This allows us to modify the B-coefficients $B_{xy}\left(y{y}^{\prime }\right)$ and $B_{xy}\left(y^{\prime }{y}^{″}\right)$ in (17) and (19) to vanish, and we can use instead

$$D{F}_{xy}\left(P^s\right)\left({\left(q_{uv}\right)}_{(u,v)\in E^{\mathrm{dir}}}\right)=\sum _{y^{\prime }\ne x}{B}_{xy}^{\prime }\left(x{y}^{\prime }\right)q_{x{y}^{\prime }}+\sum _{y^{\prime }\ne x}{B}_{xy}^{\prime }\left(y^{\prime }x\right)q_{y^{\prime }x}+\sum _{y^{\prime }\ne x,y}{B}_{xy}^{\prime }\left(y^{\prime }y\right)q_{y^{\prime }y}$$

with

$$\begin{array}{ccc}\hfill B_{xy}^{\prime }\left(xy\right)& =& -4p_{yx}^s-p_{yy}^s-2\sum _{y^{\prime }\ne y}{p}_{y{y}^{\prime }}^s\hfill \\ & & +\frac{1}{D_x}\left(4p_{xy}^s{p}_{yx}^s+4\sum _{y^{\prime }}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }x}^s+p_{xy}^s\sum _{y^{\prime }}{p}_{y{y}^{\prime }}^s+\sum _{y^{\prime },y^{″}}{p}_{x{y}^{\prime }}^s{p}_{y^{\prime }{y}^{″}}^s\right),\hfill \\ \hfill B_{xy}^{\prime }\left(x{y}^{\prime }\right)& =& \frac{4}{D_x}{p}_{xy}^s{p}_{y^{\prime }x}^s+\frac{1}{D_x}{p}_{xy}^s\sum _{y^{″}}{p}_{y^{\prime }{y}^{″}}^s+p_{y^{\prime }y}^s\mathrm{if}{y}^{\prime }\ne x,y,\hfill \\ \hfill B_{xy}^{\prime }\left(yx\right)& =& p_{xy}^s\left(\frac{3}{D_x}{p}_{xy}^s-2\right),\hfill \\ \hfill B_{xy}^{\prime }\left(y^{\prime }x\right)& =& \frac{3}{D_x}{p}_{xy}^s{p}_{x{y}^{\prime }}^s\mathrm{if}{y}^{\prime }\ne x,y,\hfill \\ \hfill B_{xy}^{\prime }\left(y^{\prime }y\right)& =& p_{x{y}^{\prime }}^s\mathrm{if}{y}^{\prime }\ne x,y.\hfill \end{array}$$

Note that in the case of the unmixed complete graph, we have $S_1\left(x\right)=V\backslash \left\{x\right\}$ and, in the formulas for the $B^{\prime }$-coefficients, $y^{\prime }{y}^{″}$ represent all vertices different from x, as before.

To end up with a uniquely defined Jacobi matrix of F, we need to restrict transition probabilities that are independent. For that, we introduce the subset $E^{\mathrm{ess}}\subset E^{\mathrm{dir}}$ of “essential” transition probabilities by removing, for each $x\in V$ with outgoing directed edges, that is $S_1\left(x\right)\ne \varnothing$, one pair $(x,y)$ from $E^{\mathrm{dir}}$. The cardinality of $E^{\mathrm{ess}}$ is $M:=|E^1|+2|E^2|-|V_0|$, where $V_0\subset V$ is the subset of vertices $x\in V$ for which we have $S_1\left(x\right)\ne \varnothing$. Any choice ${\left(p_{uv}\right)}_{(u,v)\in E^{\mathrm{ess}}}$ determines then a weighting scheme P by setting $p_{uv}=D_u-{\sum }_{v^{\prime }\in S_1\left(u\right)\backslash \left\{v\right\}}{p}_{u{v}^{\prime }}$ for the directed edge $(u,v)\notin E^{\mathrm{ess}}$. Similarly, any choice ${\left(q_{uv}\right)}_{(u,v)\in E^{\mathrm{ess}}}$ also determines the parameters $q_{uv}$ with $(u,v)\notin E^{\mathrm{ess}}$ by setting $q_{uv}=-{\sum }_{v^{\prime }\in S_1\left(u\right)\backslash \left\{v\right\}}{q}_{u{v}^{\prime }}$. Then, $DF\left({\left(p_{uv}^s\right)}_{(u,v)\in E^{\mathrm{ess}}}\right)$ is a square matrix of size M, and the weighting scheme $P^s$ corresponding to ${\left(p_{uv}^s\right)}_{(u,v)\in E^{\mathrm{ess}}}$ is asymptotically stable if and only if the real parts of all eigenvalues of this square matrix are negative, and the weighting scheme $P^s$ is unstable if and only if at least one of the real parts of these eigenvalues is positive.

Let us reformulate this restriction in terms of matrix multiplications. We start by enumerating the vertices of the graph $G=(V,E)$: $V=\{v_0,\dots ,v_N\}$. We also introduce the following enumeration on the directed edges in $E^{\mathrm{dir}}$: Let $1=j_0$ and $a_{j_0},\dots ,a_{k_0}$ be the edges of the type $(v_0,\ast )$ (where second vertices are chosen with increasing indices) in $E^{\mathrm{dir}}$, $j_1=k_0+1$ and $a_{j_1},\dots ,a_{k_1}$ be the edges of the type $(v_1,\ast )$ in $E^{\mathrm{dir}}$, and so on. For all vertices $v_l\in V$ with $S_1\left(v_l\right)=\varnothing$, we set $j_l=k_{l-1}+1$ and $k_l=k_{l-1}$. We remove the edges $a_{k_0},a_{k_1},\dots ,a_{k_N}$ from $E^{\mathrm{dir}}$ to obtain $E^{\mathrm{ess}}$. For simplicity, we use the notation $p_j^s$ and $q_j$ for $p_{a_j}^s$ and $q_{a_j}$, and we can write for all $a_j\in E^{\mathrm{ess}}$,

$$D{F}_{a_j}\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)\left({\left(q_k\right)}_{a_k\in E^{\mathrm{ess}}}\right)=\sum _{l=0}^N\sum _{k=j_l}^{k_l}{B}_{a_j}\left(a_k\right)q_k=\sum _{l=0}^N\sum _{k=j_l}^{k_l-1}(B_{a_j}\left(a_k\right)-B_{a_j}\left(a_{k_l}\right))q_k,$$

where the last expression involves only parameters $q_k$ corresponding to essential directed edges $a_k\in E^{\mathrm{ess}}$. Consequently, the Jacobi matrix $DF\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)\left({\left(q_k\right)}_{a_k\in E^{\mathrm{ess}}}\right)$ can be written as

$$DF\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)=P_{1}B{P}_2$$

(20)

where B is the square matrix of size $k_N$ with $B_{jk}=B_{a_j}\left(a_k\right)$ for $a_j,a_k\in E^{\mathrm{dir}}$, $P_1$ is obtained from the identity matrix $I_{k_N}$ by removing the rows $k_0,k_1,\dots ,k_N$, and $P_2=P_1^{\top }-P_3$ with $P_3$ a $k_N\times M$ matrix whose first $k_0-1$ columns are all the standard basis vector $e_{k_0}$, the next $k_1-j_1$ columns are all the standard basis vector $e_{k_1}$, and so on.

3.2. Examples of Asymptotically Stable and Unstable Equilibria

In this subsection, we investigate curvature sharp weighting schemes of various examples of unmixed combinatorial graphs.

Example 8

(Curvature sharp weighting schemes on a cycle). Let $C_N=(V,E)$ be a cycle of length $N\ge 4$, that is, $V=\{v_0,v_1,\dots ,v_{N-1}\}$ and $v_i\sim v_{i+1}$ with indices i modulo N. For simplicity, we refer to vertex $v_i$ henceforth as i. We assume $P^s$ to be a non-lazy curvature sharp weighting scheme on $C_N$, and we remove the directed edges $(i,i+1)\in E^{\mathrm{dir}}$ to obtain $E^{\mathrm{ess}}$, and the only coefficients $B_{a_j}\left(a_k\right)$ with $a_j\in E^{\mathrm{ess}}$ and $a_k\in E^{\mathrm{dir}}$, which may be potentially non-zero are (see (13)–(16)),

$$\begin{array}{ccc}\hfill B_{i,i-1}(i,i-1)& =& -4p_{i,i-1}^s+8p_{i,i-1}^s{p}_{i-1,i}^s+4p_{i,i+1}^s{p}_{i+1,i}^s,\hfill \\ \hfill B_{i,i-1}(i,i+1)& =& 4p_{i,i-1}^s{p}_{i+1,i}^s,\hfill \\ \hfill B_{i,i-1}(i-1,i)& =& 4{\left(p_{i,i-1}^s\right)}^2-4p_{i,i-1}^s,\hfill \\ \hfill B_{i,i-1}(i+1,i)& =& 4p_{i,i-1}^s{p}_{i,i+1}^s.\hfill \end{array}$$

This implies

$$\begin{array}{c}D{F}_{i,i-1}\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)\left({\left(q_k\right)}_{a_k\in E^{\mathrm{ess}}}\right)\hfill \\ =(B_{i,i-1}(i,i-1)-B_{i,i-1}(i,i+1))q_{i,i-1}+B_{i,i-1}(i+1,i)q_{i+1,i}-B_{i,i-1}(i-1,i)q_{i-1,i-2}\\ =4({\left(p_{i,i-1}^s\right)}^2-p_{i,i-1}^s)q_{i-1,i-2}+4(-p_{i,i-1}^s+2p_{i,i-1}^s{p}_{i-1,i}^s+p_{i,i+1}^s{p}_{i+1,i}^s-p_{i,i-1}^s{p}_{i+1,i}^s)q_{i,i-1}\\ \hfill +4p_{i,i-1}^s{p}_{i,i+1}^s{q}_{i+1,i}.\end{array}$$

In the case of the simple random walk $p_{i,i-1}^s=p_{i,i+1}^s=1/2$, $i\in \{0,\dots ,N-1\}$, this simplifies to

$$D{F}_{i,i-1}\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)\left({\left(q_k\right)}_{a_k\in E^{\mathrm{ess}}}\right)=q_{i-1,i-2}+q_{i+1,i},$$

and the corresponding matrix $DF\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)$ coincides with the adjacency matrix of $C_N$ whose largest eigenvalue is 2. Therefore, the simple random walk on $C_N$, $N\ge 4$, is an unstable equilibrium. This is in contrast to the simple random walk on $C_3=K_3$, which is asymptotically stable, as we will see in Example 10.

In the case of the totally degenerate clockwise weighting scheme $p_{i,i-1}=1$ and $p_{i,i+1}=0$ (see right-hand side of Figure 4), we have

$$D{F}_{i,i-1}\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)\left({\left(q_k\right)}_{a_k\in E^{\mathrm{ess}}}\right)=-4q_{i,i-1},$$

that is, $DF\left({\left(p_k^s\right)}_{a_k\in E^{\mathrm{ess}}}\right)=-4{\mathrm{Id}}_N$, and this curvature sharp Markovian weighting scheme is asymptotically stable. This agrees with the fact that many initial weighting schemes end up in this limit under the curvature flow. The same holds true for the corresponding totally degenerate anti-clockwise weighting scheme with $p_{i,i+1}=1$ and $p_{i,i-1}=0$.

Example 9

(Flow limits of the octahedron). We know from Example 7 that stationary solutions of the normalized curvature flow on the octahedron are the simple random walk without laziness as well as the totally degenerate weighting scheme given as matrix $P^s$ in the following code (see also the right-hand side of Figure 12). The linearization $DF\left(P^s\right)$ is analyzed in line 10 of the program, and the function equilibrium_type with the parameters chosen in lines 6 and 7 returns one of the values $-1,0,1$ (corresponding to “asymptotically stable”, ”undecided”, “unstable”, respectively), followed by a list of its eigenvalues. That is, after execution of line 10, result[0] is one of the values $-1,0,1$ and result[1] is a list of the 18 eigenvalues ${\lambda }_j$ of $DF\left(P^s\right)$.

The program provides us with the following list of complex eigenvalues:

${\lambda }_j$	$-1$	$-1+i$	$-1-i$	$-2$	$-2+i$	$-2-i$	$-3$	$-4$
multiplicity	2	2	2	2	2	2	2	4

This shows that the totally degenerate curvature sharp weighting scheme $P^s$ on the octahedron is asymptotically stable. This is expected since the initial weighting scheme $P_0$ in (11) of Example 7 converges to this limit under the normalized numerical curvature flow.

Running the same code for the simple random walk without laziness instead (by uncommenting line 3 and commenting out lines 4 and 5 in the above code) shows that this second curvature sharp weighting scheme is unstable. The eigenvalues, in this case, are all real-valued, one of them $0.5$, and given as follows:

${\lambda }_j$	$0.5$	0	$-0.75$	$-1$	$-1.5$
multiplicity	3	3	2	6	4

Example 10

(Simple random walk on a complete graph). Let $K_{n+1}=(V,E)$ be the complete unmixed graph with $n+1\ge 3$ vertices. Instead of the B-coefficients, we make use of the $B^{\prime }$-coefficients introduced in Remark 1. The simple random walk $p_{xy}=\frac{1}{n}$ for $x\ne y$ is a curvature sharp Markovian weighting scheme, and the non-zero $B^{\prime }$-coefficients are then given by

$$B_{xy}^{\prime }\left(xy\right)=\frac{(n+1)(3-n)}{n^2},B_{xy}^{\prime }\left(x{y}^{\prime }\right)=\frac{2n+3}{n^2},B_{xy}^{\prime }\left(yx\right)=\frac{3-2n}{n^2},B_{xy}^{\prime }\left(y^{\prime }x\right)=\frac{3}{n^2},B_{xy}^{\prime }\left(y^{\prime }y\right)=\frac{1}{n},$$

where $y^{\prime }\in V$ is an arbitrary vertex different from $x,y$. Let $\{v_0,v_1,\dots ,v_n\}$ be the vertex set of $K_{n+1}$ and, as in the previous example, we refer to vertex $v_i$ as i, for simplicity.

Let us now consider the case $n=2$. The process of removing edges from $E^{\mathrm{dir}}$ described earlier leads to the following remaining edges in $E^{\mathrm{ess}}$:

$$a_1=(0,1),a_3=(1,0),a_5=(2,0).$$

Choosing the simple random walk $p_{01}^s=p_{10}^s=p_{20}^s=1/2$, we have

$$\begin{array}{ll}DF& (p_{01}^s,p_{10}^s,p_{20}^s)=\left(\begin{array}{c}{e}_1^{\top }\\ e_3^{\top }\\ e_5^{\top }\end{array}\right)B^{\prime }\left(\begin{array}{ccc}{e}_1-e_2& e_3-e_4& e_5-e_6\end{array}\right)\\ & =\left(\begin{array}{ccc}{B}_{01}^{\prime }\left(01\right)-B_{01}^{\prime }\left(02\right)& B_{01}^{\prime }\left(10\right)-B_{01}^{\prime }\left(12\right)& B_{01}^{\prime }\left(20\right)-B_{01}^{\prime }\left(21\right)\\ B_{10}^{\prime }\left(01\right)-B_{10}^{\prime }\left(02\right)& B_{10}^{\prime }\left(10\right)-B_{10}^{\prime }\left(12\right)& B_{10}^{\prime }\left(20\right)-B_{10}^{\prime }\left(21\right)\\ B_{20}^{\prime }\left(01\right)-B_{20}^{\prime }\left(02\right)& B_{20}^{\prime }\left(10\right)-B_{20}^{\prime }\left(12\right)& B_{20}^{\prime }\left(20\right)-B_{20}^{\prime }\left(21\right)\end{array}\right)=\left(\begin{array}{ccc}-1& \frac{1}{4}& \frac{1}{4}\\ \frac{1}{4}& -1& \frac{1}{4}\\ \frac{1}{4}& \frac{1}{4}& -1\end{array}\right),\end{array}$$

(21)

which is a negative definite matrix with eigenvalues $-\frac{1}{2}$, $-\frac{5}{4}$, $-\frac{5}{4}$, showing that the simple random walk on $K_3$ is an asymptotically stable equilibrium. This result can be numerically verified via the following code:

The program returns the eigenvalues of the linearized flow matrix of $K_3$ at the simple random walk. Note, however, that the computed matrix is based here on the B-coefficients instead of the $B^{\prime }$-coefficients, so this matrix is slightly different from the one given above, while the eigenvalues are the same.

In the case $n=3$, we have $p_{01}^s=p_{02}^s=p_{10}^s=p_{12}^s=p_{20}^s=p_{21}^s=p_{30}^s=p_{31}^s=1/3$ and the returned matrix by the program (after changing n=2 into n=3 in line 1 of the code) is as follows:

$$DF(p_{01}^s,p_{02}^s,p_{10}^s,p_{12}^s,p_{20}^s,p_{21}^s,p_{30}^s,p_{31}^s)=\left(\begin{array}{cccccccc}-1& 0& -1/3& 0& 1/3& 1/3& 1/3& 1/3\\ 0& -1& 1/3& 1/3& -1/3& 0& 0& -1/3\\ -1/3& 0& -1& 0& 1/3& 1/3& 1/3& 1/3\\ 1/3& 1/3& 0& -1& 0& -1/3& -1/3& 0\\ 0& -1/3& 1/3& 1/3& -1& 0& 0& -1/3\\ 1/3& 1/3& 0& -1/3& 0& -1& -1/3& 0\\ 1/3& 1/3& 0& -1/3& 0& -1/3& -1& 0\\ 0& -1/3& 1/3& 1/3& -1/3& 0& 0& -1\end{array}\right).$$

Its eigenvalues are $-\frac{2}{3}$ with multiplicity 6 and $-2$ with multiplicity 2. This shows that the simple random walk on $K_4$ is again an asymptotically stable equilibrium.

Since numerical experiments show that all non-degenerate Markovian weighting schemes without laziness on $K_{n+1}$ converge to the simple random walk, we expect that the simple random walk on $K_{n+1}$ is an asymptotically stable equilibrium for all $n\ge 2$ (see also our Conjecture 3, which is an even stronger statement). To this end, our code provides the following numerical results: The eigenvalues of $DF\left(P^s\right)$ for the simple random walk without laziness on the complete graph $K_{n+1}$ are real-valued and given by:

• $-\frac{7}{16}$ with multiplicity 4, $-\frac{3}{4}$ with multiplicity 6, and $-\frac{7}{4}$ with multiplicity 5 for $n=4$;
• $-\frac{8}{25}$ with multiplicity 5, $-\frac{4}{5}$ with multiplicity 10, and  $-\frac{8}{5}$ with multiplicity 9 for $n=5$;
• $-\frac{1}{4}$ with multiplicity 6, $-\frac{5}{6}$ with multiplicity 15, and $-\frac{3}{2}$ with multiplicity 14 for $n=6$;
• $-\frac{10}{49}$ with multiplicity 7, $-\frac{6}{7}$ with multiplicity 21, and $-\frac{10}{7}$ with multiplicity 20 for $n=7$.

These results give rise to the following conjecture:

Conjecture 5.

The non-lazy simple random walk $P^s$ on the unmixed complete graph $K_{n+1}$, $n\ge 2$, is asymptotically stable and the eigenvalues of $DF\left(P^s\right)$ are given by

$$\begin{array}{ccc}\hfill -\frac{n-1}{n}& & \mathrm{with}\mathrm{multiplicity}\left(\genfrac{}{}{0pt}{}{n}{2}\right),\hfill \\ \hfill -\frac{n+3}{n^2}& & \mathrm{with}\mathrm{multiplicity}n,\hfill \\ \hfill -\frac{n+3}{n}& & \mathrm{with}\mathrm{multiplicity}\left(\genfrac{}{}{0pt}{}{n}{2}\right)-1.\hfill \end{array}$$

As $n\to \infty$, we have $-\frac{n-1}{n}\to -1$, $-\frac{n+3}{n^2}\approx -\frac{1}{n}\to 0$ and $-\frac{n+3}{n}\to -1$.

Example 11

(Simple random walks on hypercubes). We know from [1] (Corollary 1.14) that the simple random walk without laziness on a hypercube $Q^d={\left(K_2\right)}^d$ is curvature sharp, and it is the only non-degenerate curvature sharp weighting schemes without laziness if and only if d is odd. For $d=2$, a non-degenerate curvature sharp weighting scheme on $Q^2$ different from the simple random walk is given in line 4 of the following code:

All eigenvalues ${\lambda }_j$ of the corresponding linearized flow matrix provided by this code via line 6 seem to be real. Two of them are numerically zero and the other two are given by $\pm \lambda$ with a non-zero real λ. So this equilibrium is unstable.

Our experiments for arbitrary weighting schemes on $Q^d$, $d\ge 2$, close to the simple random walk show that the numerical flow converges usually to a degenerate limit. This agrees with our observation that the linearized flow matrix at the simple random walk always seems to have real eigenvalues with some of them being positive. These eigenvalues can be obtained via the following code:

Let us finish this section with a conjecture that we verified numerically for $d=2,3,\dots ,9$:

Conjecture 6.

The non-lazy simple random walk $P^s$ on the hypercube $Q^d$, $d\ge 2$, is unstable and the eigenvalues of the corresponding linearized curvature flow matrix are all real with the largest eigenvalue ${\lambda }_{max}$ equals

$${\lambda }_{max}=\frac{4}{d}.$$

4. Conclusions and Open Questions

In Section 2, we investigated curvature flows of various concrete examples. Our findings and possible directions for further research are as follows:

(a)

Random graphs: We observed an edge density threshold for Erdös-Rényi graphs. For values below this threshold, the flow usually converges to a totally degenerate limit, while for values above this threshold, the flow limit concentrates on a non-degenerate cluster (see Conjecture 2). Further research could encompass curvature flows of random regular graphs with simple random walks as their initial weighting schemes.

(b)

Paths and cycles: Our observed limits for these graphs are usually directed graphs with isolated clusters of three vertices (see Figure 4 and Proposition 1).

(c)

Complete graphs and their wedge sums and Cartesian products: In the case of complete graphs, non-degenerate initial weighting schemes seem to converge always to simple random walks (see Conjecture 3). In wedge sums of complete graphs, curvature flow limits usually concentrate totally on one of its components, which is most often the largest component (see Table 2). This observation might motivate further research on the largest clusters of general highly connected graphs. We also describe the flow limit of Cartesian products of complete graphs with the simple random walk as its initial weighting scheme (see Theorem 3). A natural question here is the limit behavior of Cartesian products of graphs with curvature sharp initial weighting schemes.

(d)

Octahedron: In the case of the octahedron, there is a unique totally degenerate flow limit (up to symmetries), which has two sources and a directed 4-cycle (see Figure 12 and Proposition 3). The shape of totally degenerate flow limits for general graphs is obscure and would benefit from further investigations.

As mentioned already in [1], curvature flow limits are always curvature sharp. A rough interpretation of curvature sharpness at a vertex is that its Ricci curvature is constant in all directions (see [1] Equation (8)). Section 3 was concerned with a dynamical aspect of these curvature flow limits, namely their stability properties. Our findings are the following:

(e)

Cycles: Simple random walks on cycles $C_n$ for $n\ge 4$ are unstable, and totally degenerate directed weighting schemes are asymptotically stable (see Example 8).

(f)

Complete graphs: We have strong numerical evidence that simple random walks on complete graphs $K_n$ are asymptotically stable (numerically verified for $n=3,\dots ,8$); see Example 10. Moreover, our numerics provide explicit information about all eigenvalues of the linearized curvature flow as stated in Conjecture 5.

(g)

Hypercubes: Simple random walks of hypercubes $Q^d$ for $d\ge 2$ seem to be unstable and our experiments give rise to a conjecture about the largest eigenvalue of the linearized curvature flow (see Conjecture 6).

(h)

Octahedron: The totally degenerate flow limit of the octahedron is asymptotically stable (see Example 9).

(i)

Edge-transitive graphs: In view of (e), (f), and (g), an interesting research direction would be to investigate the stability of the simple random walk for general edge-transitive graphs.

Let us end this section with some final reflections: After having introduced a discrete curvature flow based on the Bakry–Émery calculus in [1], the main concern of this paper was to gain further insights into this flow from a dynamical viewpoint using a specially developed Python software. The most pressing and still open question is whether our curvature flow is always convergent as time tends to infinity (see Conjecture 1), or whether other limiting behavior may occur such as periodic orbits. Moreover, this could be just the starting point of various further questions, for example, whether every curvature flow allows for asymptotically stable limits—a natural technique here might be to use particular Lyapunov functions (see, e.g., [39] (Section 9.3)). The dynamical viewpoint also gives rise to other questions such as the existence of heteroclinic orbits between equilibria. All these questions may lead to further research in the future.

5. Implementation of the Curvature Flow and Useful Tools

This section provides a short description of a program, written in Python, which allows users to carry out their own curvature flow experiments. This program is provided as an ancillary file of [2] or via github at https://github.com/georgestagg/graph-curvature-server (accessed on 6 June 2023) and can be used in combination with the “Graph Curvature Calculator”, which can be freely accessed at http://www.mas.ncl.ac.uk/graph-curvature (accessed on 6 June 2023).

The Graph Curvature Calculator is a powerful but easy-to-use interactive Web tool to draw graphs and to compute various types of curvatures such as Bakry–Émery curvature on its vertices or Ollivier Ricci curvature on its edges (for details see [40]). This Web tool allows users to obtain the adjacency matrix of the graph under consideration, which can then be used as input for the curvature flow code.

The curvature flow code provides functions and routines, which can be divided into six categories:

• Functions related to graphs and their adjacency matrices;
• Functions related to weighting schemes;
• Functions testing combinatorial and weighted graphs;
• Curvature flow computation routines;
• Curvature computation routines;
• Display routines.

Each category is discussed in one of the following subsections.

5.1. Functions Related to Graphs and Their Adjacency Matrices

Recall that the topology of a given Markov chain $(G,P)$ is contained in the combinatorial graph $G=(V,E)$. This information is provided by the adjacency matrix $A=A_G$ of the graph. Firstly, we go through some useful functions performing operations with these adjacency matrices. Every graph generated by one of the functions below is returned as such an adjacency matrix; in particular, as a NumPy array. The names of these functions are as follows:

• rand_adj_mat(n, p, connected=False);
• complete(n);
• path(n);
• cycle(n);
• wedge_sum(A, B, i, j);
• bridge_at(A, B, i, j);
• hypercube(n);
• cart_prod(A, B);
• onespheres(A).

rand_adjmat(n, p) returns a random adjacency matrix with n vertices, where p is the probability that an edge exists between any two vertices. Therefore, higher values of p usually lead to better-connected graphs. Choosing connected=True instead guarantees that the returned adjacency matrix provides a connected combinatorial graph.

The functions complete(n), path(n) and cycle(n) return the adjacency matrix of a complete graph with n vertices, a path of length n and a cycle of length n, respectively.

wedge_sum(A, B, i, j) produces an adjacency matrix of a graph that is the wedge sum of two pointed graphs represented by the adjacency matrices A and B and the ith vertex of A and the jth vertex of B, that is, the new graph is obtained from the disjoint union of these two graphs by only identifying these two vertices and keeping all other edges and vertices disjoint.

There is also bridge_at(A, B, i, j), which returns an adjacency matrix of the graph formed by connecting the graphs represented by A and B by a single edge at the ith vertex of A and the jth vertex of B.

The function hypercube(n) returns the adjacency matrix of the n-dimensional hypercube.

The Cartesian product of two graphs represented by adjacency matrices A and B is returned by the function cart_prod(A, B). (The n-dimensional hypercube can also be alternatively generated as the Cartesian product of n complete graphs $K_2$.)

Finally, for a graph $G=(V,E)$ given by the adjacency matrix A and the vertex set $V=\{0,1,\dots ,n-1\}$, the function onespheres(A) returns a list of lists whose i-th entry is a list of all neighbors of vertex $i\in V$, and whose n-th entry is a list of the combinatorial degrees of the vertices in V. This function is mainly used in the execution of other functions.

5.2. Functions Related to Weighting Schemes

The weighting scheme of a Markov chain $(G,P)$ is provided via the weighted matrix $P=P_G$. The functions in this category are the following:

• randomizer(A, threshold=0.001, laziness=False);
• srw(A, laziness=False);
• cart_prod_prob(P, Q, p, q).

randomizer(A) and srw(A) are two useful functions that can be used to give weighted matrices from an adjacency matrix.

The function randomizer(A) returns a weighting scheme for the graph G with random, numerically non-degenerate transition rates, that is, no transition rate is chosen to be below the parameter threshold. Usually, the returned weighting schemes are without laziness, but choosing laziness=True returns weighting schemes with all vertices having laziness ≥ threshold.

srw(A) returns the weighting scheme corresponding to the non-lazy simple random walk on the graph $G=(V,E)$ represented by the adjacency matrix A. Choosing laziness=True, the transition rates of a vertex $v\in V$ with degree n to its neighbors are chosen to be $\frac{1}{n+1}$ and its laziness is also chosen to be $\frac{1}{n+1}$.

The function cart_prod_prob(P, Q, p, q) is a “weighted” analog of the function cart_prod from the previous subsection for the Cartesian product of two weighting schemes P, Q with weights p, q, with $\mathtt{p}+\mathtt{q}=1$. If P and Q are of size n and m, respectively, this function returns the matrix $\mathtt{p}\mathtt{P}\otimes I_n+\mathtt{q}{I}_m\otimes \mathtt{Q}$ of size $nm$, where $I_n$ is the identity matrix of size n and $A\otimes B$ is the Kronecker product of A and B.

5.3. Functions Testing Combinational and Weighted Graphs

For combinatorial graphs given by their adjacency matrices A and weighted graphs given additionally by their weighting scheme P, we have the following test functions:

• is_connected(A);
• is_weakly_connected(A, threshold=0.001);
• is_totally_degenerate(A, P, threshold=0.001);
• is_markovian(P, norm_tolerance=0.001);
• is_curvature_sharp(A, P, norm_tolerance=0.001, threshold=0.001);
• equilibrium_type(A, P, eigenvalues=False, jacobi_matrix=False,
  norm_tolerance=0.001, threshold=0.001).

The function is_connected(A) returns True if and only if the adjacency matrix A represents a connected graph.

The function is_weakly_connected(P) is a “weighted” analog, which returns True if and only if the weighted matrix P represents a weakly connected graph. It does this by forming an adjacency matrix with a 1 in the $(\mathtt{i},\mathtt{j})$th entry if and only if P[i, j] > threshold or P[j, i] > threshold and then testing for connectedness.

Recall that a weighted graph $(G,P)$ is called numerically totally degenerate if there are no two-sided edges with numerical non-zero transition rates in both directions and no one-sided edges with numerical non-zero transition rate, where we consider a transition rate $p_{xy}$ as numerically non-zero if and only if $p_{xy}\ge$ threshold. The function is_totally_degenerate(A,P) tests this property.

The function is_markovian(P) tests whether the entries of each of the columns of P add up numerically to 1 up to an error $\le \mathtt{norm}\_\mathtt{tolerance}$.

Numerical curvature sharpness (up to an error $\le \mathtt{threshold}$) of Markovian weighted graphs given by (A,P) is tested by is_curvature_sharp(A,P). If (A,P) fails to be Markovian (with respect to norm_tolerance), this function returns NONE and gives notice to the user by an error message.

For dynamical investigations of curvature flow equilibria, we have the function equilibrium_type(A,P), which always returns a list of length three. The function checks first whether (A,P) satisfies the Markovian property and is numerically curvature sharp. If this is not the case, it returns a list of three NONE values. Otherwise, the function investigates the real parts of the eigenvalues ${\lambda }_j$ of the linearized curvature flow matrix at the equilibrium P. The first entry of the return list is $-1,0$ or 1 depending on the maximum ${max}_j\mathrm{Re}\left({\lambda }_j\right)$. If this maximum is $\ge \mathtt{threshold}$, the return value is 1 (for “unstable”) and if this maximum is $\le -\mathtt{threshold}$, the return value is $-1$ (for “asymptotically stable”). Otherwise, the dynamical nature of the equilibrium cannot be numerically decided and the function returns the value 0. The following two entries of the return list are usually NONE unless the user made the choices eigenvalues=True or jacobi_matrix=True. In the first case, the second entry of the return list is a list of all eigenvalues of the linearized curvature flow matrix, and in the second case the third entry of the return is the linearized curvature flow matrix itself.

5.4. Curvature Flow Computation Routines

At the heart of the program are the curvature flow routines solving the initial value ordinary differential Equations (6) and (7). The relevant routines are the following:

• curv_flow(A, P, t_max, dt=0.3, C=zeroes),
• norm_curv_flow(A, P, t_max, dt=0.3, stoch_corr=True, norm_tolerance=0.001),
• norm_curv_flow_lim(A, P, dt=0.3, stoch_corr=True, norm_tolerance=0.001, lim_tolerance=0.001, t_lim=10000).

The initial Markov chain $(G,P_0)$ with $G=(V,E)$ is entered by the adjacency matrix A describing the topology of the graph G and the weighting scheme P containing the initial probability transitions $p_{xy}\left(0\right)$.

The first routine curv_flow(A,P) computes the non-normalized numerical curvature flow with coefficients $C_x\left(t\right)=0$ in (7). If users decide to choose other coefficient functions, they need to modify the input parameter C=zeroes. Note that zeroes(A,P) is a function returning simply a list of zeroes of length $\left|V\right|$. Users can investigate modifications of the curvature flow by choosing their own coefficient functions with input values (A,P) and returning a list of length $\left|V\right|$ of real values. Using the discretization parameter dt=0.3 for the discrete time steps starting at t=0, the curvature flow routine creates a list P_list of weighting schemes (represented by NumPy arrays of size $\left|V\right|\times \left|V\right|$) at each time increment using the Runge–Kutta algorithms RK4. There are two internal subroutines involved, which we would like to mention briefly. Pvecs_to_P translates a weighting scheme given by a list of lists (where each inner list contains the transition rates of the corresponding vertex) into the corresponding NumPy array. Pvecs_prime computes, for a given weighting scheme P, the right-hand side of the ordinary differential equations describing the curvature flow. The representation of this right-hand side is again a list of lists, as described before. The computation stops just before the discrete time steps exceed the limit time t_max and the routine returns P_list. Note that in this general setting, transition rates can assume arbitrary values and even negative ones or diverge to infinity in finite time which may lead to system error messages. Users need to be aware of this possibility.

The normalized curvature flow, using the coefficients $C_x\left(t\right)=K_{P\left(t\right),\infty }^{d(x,·)}\left(x\right)$ (see (8)), is numerically computed by the routine norm_curv_flow(A,P). This special flow is the main focus of this paper and could also be mimicked by choosing C=K_inf in curv_flow. The function K_inf(A,P) returns a list of upper curvature bounds for all vertices of the Markovian weighted graph represented by (A,P) (see [1] (Formula (66)) for an explicit expression of this function in terms of transition rates). During the numerical computations of subsequent time steps, the Markovian property of the corresponding weighting schemes may be slightly violated. If this violation exceeds the threshold norm_tolerance, the routine prepares for a potential correction according to the Boolean variable stoch_corr. If stoch_corr=False, the program stops with a message to the user as discussed in Example 1. Otherwise, the program carries out the following automatic Markovian renormalization of the currently considered weighting scheme: while the diagonal entries (the laziness values) are unchanged, the off-diagonal entries in each row are rescaled by a suitable factor to guarantee that the resulting matrix becomes stochastic again. As before, this routine returns a list P_list of consecutive weighting schemes up to the time limit t_max.

While the user needs to specify the time limit t_max in the above two routines, the third routine norm_curv_flow_lim(A,P) continues computing the normalized numerical curvature flow until a numerical flow limit is reached. This limit is determined by the parameter lim_tolerance. The details for this numerical limit are explained in the introductory part of Section 2. This routine returns a list of length two: the limiting weighting scheme as a NumPy array followed by the numerical convergence time. Since it may happen that a normalized numerical flow does not converge at all (even though we are not aware of any such example), the parameter t_lim provides an upper time limit beyond which the routine will not continue. The parameters stoch_corr and norm_tolerance play the same role as in the routine norm_curv_flow.

5.5. Curvature Computation Routines

The functions in this section calculate Bakry–Émery curvatures and curvature upper bounds of graphs with given weighting schemes at all vertices.

• curvatures(A, P, N=inf, onesps=[], q=None);
• calc_curvatures(A, P_list, N=inf, k=1);
• K_inf(A, P);
• calc_curv_upper_bound(A, P_list, N=inf, k=1).

The routine curvatures(A,P) computes, for a weighted graph $(G,P)$ with $G=(V,E)$ and represented by (A,P), the curvatures of all vertices for dimension N$=\infty$ and returns them as a list of length $\left|V\right|$. If users are interested in curvatures for other dimensions, they need to change the parameter N=inf. There are two other inputs that can speed up the curvature calculations: if the number q of vertices in V is given, it can be specified to avoid its repeated recalculation, for example, during a curvature flow process. Similarly, if onespheres(A) has already been calculated earlier, this information can be communicated to the routine via the input variable onesps.

After a curvature flow computation with corresponding list P_list of consecutive weighting schemes, the routine calc_curvatures(A,P_list) computes the corresponding evolution of vertex curvatures by calling curvatures(A,P_list[j]) and returns it as a list of lists. Here, the j-th inner list contains the curvature evolution of the j-th vertex of the graph $G=(V,E)$ represented by A. The dimension parameter N=inf plays the same role as before. calc_curvatures computes curvatures only of each k-th weighting scheme provided by P_list. Where appropriate, this can help to reduce computation time.

The routine K_inf(A,P) was already discussed in the previous subsection and provides a list of upper curvature bounds for all vertices of the weighted graph $(G,P)$ represented by (A,P).

calc_curv_upper_bound is completely analogous to calc_curvatures, but it calls K_inf instead of curvatures.

5.6. Display Routines

The main display routines for users are the evolution of curvatures at various vertices during the curvature flow, the evolution of transition rates of edges emanating from vertices, and the display of individual weighted graphs with vertices arranged in a circle. The relevant routines are the following:

• display_curvatures(curv, dt=0.3, is_Markovian=True, N=inf, k=1,
  curv_bound=[], vertex_list=[]);
• display_trans_rates(A, P_list, dt=0.3, vertex_list=[]);
• display_weighted_graph(A, P, title=None, threshold=10**(-3),
  display_options=[10, True, 2, []], laziness=False).

Given the evolution of vertex curvatures during a curvature flow process via a list of lists, where the j-th inner list is the curvature evolution of the j-th vertex, display_curvatures(curv) displays the curvature evolution for each consecutive vertex separately, as illustated, for example, in Figure 3. If this information should be only given for specific vertices, this can be specified by the input parameter vertex_list. The time step dt=0.3 and the value of k together determine the labeling of the horizontal time axis. For the role of k, we refer readers to our explanation about the routine calc_curvatures. Upper curvature bounds can be inserted into the displays by the input parameter curv_bound, which needs to be given in the same format as the vertex curvatures. Constant lower and upper curvature bounds $-1$ and 2 are plotted alongside if the Boolean is_Markovian is chosen to be True and if the dimension parameter N is $\ge 2$. These bounds appear, for example, in the illustrations given in Figure 3.

Given the evolution of weighting schemes during a curvature flow process on a graph with adjacency matrix A by a list P_list of NumPy arrays, display_trans_rates(A,P_list) displays the evolution of transition rates of emanating edges for each consecutive vertex separately, as illustrated, for example, in Figure 2. The input parameters dt and vertex_list play the same role as in the previous routine.

Finally, there is the routine display_graph(A,P) with input parameters A and P representing a weighted graph $(G,P)$. This routine produces a MatPlotLib plot of this weighted graph, with the vertices arranged counter-clockwise in a circle. The plot uses the following convention to illustrate different types of edges:

• Green, solid lines represent numerically non-degenerate edges, that is $\{x,y\}\in E$ with both $p_{xy},p_{yx}\ge$ threshold;
• Red, dashed lines with an arrow represent numerically degenerate edges, that is, $\{x,y\}\in E$ with exactly one of $p_{xy}$ and $p_{yx}$ strictly less than threshold;
• Black, dotted lines represent edges with numerically vanishing transition rates in both directions, that is $\{x,y\}\in E$ with both $p_{xy},p_{yx}$ strictly less than threshold.

Users can add a title to the display by specifying the input parameter title. For Markovian weighted graphs with non-vanishing laziness, the option laziness=True labels each vertex with its corresponding laziness.

The input parameter display_options of the display_graph routine remains to be discussed. This parameter is a list of four entries. The first entry determines the size of the plot. Usually, the transition rates are printed above the edges, but if the second entry is chosen to be False, this information about the transition rates is omitted. Otherwise, the transition rates are given to a number of decimal places determined by the third entry. The default positions of these transition rates are $\frac{1}{6}$ of the way along the edges, with the number closest to the vertex x in the edge $(x,y)$ being $p_{xy}$, but this can be altered manually to avoid overlapping by specification in the fourth entry of display_options. For example, if one wishes the $p_{45}$ label to be moved to a position $\frac{1}{4}$ of the way from vertex 4 to vertex 5 and the $p_{62}$ label to be moved to a position $\frac{1}{5}$ of the way from vertex 6 to vertex 2, this fourth entry should be chosen to be [[4, 5, 1/4, 1/6], [6, 2, 1/5, 1/6]].

This completes the description of the functions and routines in the accompanying Python program to this article.