The Center Problem for the Lotka Reactions with Generalized Mass-Action Kinetics

Chemical reaction networks with generalized mass-action kinetics lead to power-law dynamical systems. As a simple example, we consider the Lotka reactions and the resulting planar ODE. We characterize the parameters (positive coefficients and real exponents) for which the unique positive equilibrium is a center.


Introduction
Lotka [7] considered a series of three chemical reactions, transforming a substrate into a product via two intermediates, X and Y. If the reactions producing X and Y, respectively, are assumed to be autocatalytic, then the resulting ODE is the classical Lotka-Volterra predator-prey system [8,9].
Farkas and Noszticzius [4] and Dancsó et al. [2] considered generalized Lotka-Volterra schemes, arising from the Lotka reactions with power-law kinetics. They studied the ODEẋ = k 1 xp − k 2 x p y q , y = k 3 x p y q − k 4 yq (1) with positive coefficients k 1 , k 2 , k 3 , k 4 > 0 and real exponents p, q,p,q ≥ 1. (The special case p = q =p =q = 1 is the classical Lotka-Volterra system.) Dancsó et al. [2] provided a local stability and bifurcation analysis. In particular, by finding first integrals, they determined four cases where the ODE admits a center.
In this work, we allow arbitrary real exponents p, q,p,q ∈ R in the ODE (1) and study the dynamics on the positive quadrant. In addition to the four known cases, we identify two new cases of centers, by showing that they correspond to reversible systems. Moreover, we prove that centers are characterized by these six cases. In a complementary work [1], we provide a global stability analysis for the Lotka reactions with generalized mass-action kinetics.
The paper is organized as follows. In Sect. 2, we elaborate on the chemical motivation of the ODE under study, and in Sect. 3, we present our main result.

The Lotka Reactions with Generalized Mass-Action Kinetics
As in the original work by Lotka [7], we start by considering a series of net reactions, S → X, X → Y, and Y → P, which transform a substrate into a product. We are interested in the dynamics of X and Y only, in particular, we assume that the substrate is present in constant amount and that the product does not affect the dynamics. As a consequence, we omit substrate and product from consideration and arrive at the simplified reactions To obtain a classical Lotka-Volterra system as in [8,9], one assumes the first and the second reaction to be autocatalytic, in particular, one defines the kinetics of the with rate constants k 0→X , k X→Y , k Y→0 > 0 and concentrations [X], [Y] ≥ 0. In this work, we consider the Lotka reactions with arbitrary power-law kinetics. In terms of chemical reaction network theory, we assume generalized mass-action kinetics [10,11], that is, with arbitrary real exponents α 1 , β 1 , α 2 , β 2 , α 3 , β 3 ∈ R. The resulting ODE for the concentrations x = [X] and y = [Y] amounts to where k 1 = k 0→X , k 2 = k 3 = k X→Y , and k 4 = k Y→0 . Since we allow real exponents, we consider the dynamics on the positive quadrant. In fact, we study an ODE which is orbitally equivalent to (2) on the positive quadrant and has two exponents less, where Further, we assume that the ODE admits a positive equilibrium (x * , y * ) and use the equilibrium to scale the ODE (3). We introduce K = k 3 k 2 x * y * > 0 and obtaiṅ Clearly, the ODE (4) admits the equilibrium (1, 1) which is not necessarily unique, and the Jacobian matrix at (1, 1) is given by Dancsó et al. [2] studied the ODE (4) in the orbitally equivalent forṁ They stated four cases where the equilibrium (1, 1) is a center and provided first integrals. In this work, we identify two new cases and show that they correspond to reversible systems. Moreover, we prove that every center belongs to one of the six cases.

Main Result
An equilibrium is a center if all nearby orbits are closed.

Theorem 1
The following statements are equivalent.
The eigenvalues of the Jacobian matrix at (1, 1) are purely imaginary, that is, tr J = 0 and det J > 0, and the first two focal values vanish. 3. The parameter values a 1 , b 1 , a 3 , b 3 ∈ R, and K > 0 belong to one of the six cases in Table 1.
The first four cases were already stated in Dancsó et al [2], where first integrals have been given. The last two cases correspond to reversible systems  Table 1, first integrals have been given by Dancsó et al. [2]. In fact, they determined all the cases for which a first integral can be found by using an integrating factor of the form x A y B . See Table 2 and [2, p. 122, Table I]. Case (i) includes the classical Lotka-Volterra systems; the corresponding first integral is of the type of separated variables and was already stated in Farkas and Noszticzius [4]. In case (iv), there is a typo in [2]; the correct formula is C = 1 q−1 . The remaining cases, (r1) and (r2), are reversible systems. See Sect. 3.2.

Case Distinction
Using tr J = 0, that is, a 1 = K b 3 by Eq. (5), we compute det J and the first two focal values, L 1 and L 2 . We find and note that det J > 0 implies a 3 , b 1 = 0. Further, using the Maple program in [6], we find If α is zero in x α α (in a first integral), replace x α α by ln x Expressions for L 2 (in case L 1 = 0) will be given below. We show that all parameters a 1 , b 1 , a 3 , b 3 ∈ R and K > 0 in the ODE (4) for which tr J = L 1 = L 2 = 0 and det J > 0 belong to one of the six cases in Table 1.
To begin with, L 1 = 0 implies either In case (a), where b 3 = 0 (and hence a 1 = 0), we find det J = K a 3 b 1 . Hence, the situation is covered by case (i) in Table 1.
In case (b), where D = 0, we find b 3 = 1 (otherwise b 1 = 0) and, using the Maple program in [6], Now, L 2 = 0 implies that at least one of six factors is zero. The first subcase a 3 +b 3 = 0 implies D = 1 − b 3 and hence b 1 = a 3 K and det J = 0. As shown above, the subcase b 3 = 0 is covered by case (i) in Table 1. a 1 − b 3 ) which corresponds to case (ii). The subcase b 3 K = 1 (and hence a 1 = 1) implies D = a 3 (1 − K ) and hence b 1 = −1. Now, det J = −K (a 3 + b 3 ), and the situation is covered by case (iv). The subcase K = 1 (and hence a 1 = b 3 ) implies D = 1 − b 3 and hence b 1 = a 3 . Now, det J = b 2 1 − a 2 1 , and the situation is covered by case (r1). Finally, the subcase 1 + a 3 In case (c1), where b 3 = 1 and a 3 = −1 (and hence a 1 = K ), we find det J = −K (a 1 + b 1 ). Hence, the situation is covered by case (iii) in Table 1. In case (c2), where b 3 = 1 and K = 1, we find Now, L 2 = 0 implies that at least one of four factors is zero. The first subcase a 3 = 0 implies det J < 0. As shown above, the subcase a 3 = −1 is covered by case (iii). Finally, the subcase b 1 = −1 is covered by case (iv), and the subcase a 3 = b 1 is covered by case (r1).

Reversible Systems
Let R : R 2 → R 2 be a reflection along a line. A vector field F : R 2 → R 2 (and the resulting dynamical system) is called reversible w.r.t. R if It is easy to see that, for any function f : R 2 → R, the systeṁ is reversible w.r.t. the reflection R : (x, y) → (y, x). The following is a well-known fact, see e.g. [12, 4.6571] or, more generally, [3,Theorem 8.1]. An equilibrium of a reversible system which has purely imaginary eigenvalues and lies on the symmetry line of R is a center. Now we are in a position to deal with the last two cases in Table 1.
This vector field is of the form (8), and hence it is reversible.
We apply the coordinate transformation u = x K , v = y −1 and obtainu Finally, we multiply the vector field with the positive function K −1 v b 1 and obtaiṅ Since 1− 1 K = 2+b 1 −b 3 , this vector field is of the form (8), and hence it is reversible. Since (r1) and (r2) lead to centers, analytic first integrals must exist. However, it seems difficult to find them. So far we succeeded only in the intersection of (r1) and (r2), that is, the case where a 1 = b 3 = b 1 + 2, a 3 = b 1 , K = 1, and b 1 < −1 (a one parameter family). See Table 2.

Limit Cycles
Dancsó et al. [2] and Boros et al. [1] identified ODEs (1) and (2) that admit one limit cycle (via a super-or sub-critical Hopf bifurcation). As a simple consequence of our characterization of the center variety, we can construct ODEs (4) with two limit cycles via a degenerate Hopf or Bautin bifurcation, see [5,Sect. 8.3]. We pick a system with tr J = L 1 = 0 and L 2 = 0, in particular, we consider case (c2) of our case distinction: we take b 3 = a 1 = 1, K = 1 and hence tr J = L 1 = 0 and choose b 1 and a 3 such that L 2 < 0 with L 2 given by Eq. (7); for example, b 1 = −2, a 3 = −3 or b 1 = 2, a 3 = 1. If we now slightly perturb K (keeping a 1 = K ) such that L 1 > 0 (and tr J = 0), the resulting system has a stable limit cycle. Finally, if we slightly change a 1 such that tr J < 0, we create a small unstable limit cycle via a subcritical Hopf bifurcation.
It remains open, if the ODE (4) admits more than two limit cycles. For a computational algebra approach to this question, see [13].

Supplementary Material
We provide a Maple worksheet containing (i) the program from [6] for the computation of the first two focal values and (ii) the case distinction described in Sect. 3.1. Further we provide a CDF file (created with Mathematica) containing 3-dimensional visualizations of the center variety. (Thereby, we start from the 5 parameters a 1 , b 1 , a 3 , b 3 , and K , use tr J = 0, that is, a 1 = K b 3 , and fix K . As a result, we obtain plots in the 3 parameters a 1 , b 1 , and a 3 .) The material is available at http://gregensburger.com/softw/center/.