Analysis of the stability and the bifurcations of two heterogeneous triopoly games with an isoelastic demand

1.   Introduction

Unlike a competitive market with a large number of relatively small companies producing substitute products and competing with each other, an oligopoly is a market supplied by only a few firms. It is well known that Cournot developed the first formal theory of oligopoly in ^[10], where a static duopoly model is investigated and each firm is supposed to have perfect information about its rival's strategic behavior. However, a dynamic oligopolistic model is much more complex than a competitive model although only a few firms are involved, since a firm's decision, as well as the reactions of its competitors, can change the price and further influence their profits.

Puu ^[37] followed Cournot's seminal work but introduced an isoelastic demand function by setting the price simply the reciprocal of the total supply of the market. In his game, where the players have naive expectations, it was shown that the evolution of the dynamics can develop a period-doubling bifurcation and finally reach chaos. It was pointed out in ^[17] that expectations and behavioral rationalities play ambitious roles in modeling economic phenomena. To our knowledge, there are mainly four types of expectations in the current literature, i.e., the naive ^[37], the boundedly rational ^[12], the adaptive ^[6], and the local monopolistic approximation (or LMA for short) ^[43].

The monopoly, a market ruled by one single firm, is the simplest oligopoly but may exhibit quite complex dynamic behavior such as quasi-period orbits and chaos ^[3,31,32,38]. In the real world, however, it is more possible that complex dynamic phenomena appear in a market supplied by two, three, or even four companies, which is called a duopoly, a triopoly, or a quadropoly, respectively. Among them, the duopoly is most intensively studied, see, e.g., ^{[1,2,5,8,14,20,24,28,36,39,48]}. Furthermore, the reader can refer to ^{[4,12,27,29,41,42]} for the explorations of triopoly games, which are also of concern in this paper. However, there are relatively few studies on quadropoly models, and the reader may refer to ^[13,30].

Oligopoly games with heterogeneous players are more realistic than those with homogeneous players. The reason is that we rarely find different firms that behave according to identical rules among a large number of possible strategies. Furthermore, in ecology, one knows that it is impossible for species with identical niches can coexist indefinitely according to the competitive exclusion principle. This principle also applies to economics, since in the real economy we often observe that companies with different business strategies coexist in the same industry. However, homogeneous firms in many industries, e.g., the Internet industry, can just coexist temporarily before the state of the economic system reaches an equilibrium.

In this work, we investigate two distinct heterogeneous triopoly games with the same isoelastic demand function proposed by Puu ^[37]. One game was first introduced in ^[41], where an adaptive firm, a naive firm, and a boundedly rational firm compete with each other. By numerical simulations, it was observed that this game could exhibit both the period-doubling bifurcation and the Neimark-Sacker bifurcation as the parameters vary. The other game ^[34] is similar but with the first adaptive player replaced by an LMA player. For this game, it was also discovered that there exist two possible routes to chaos, through the period-doubling bifurcation and the Neimark-Sacker bifurcation.

We try to apply the symbolic approach, proposed by the author of this paper and his coworker ^[26], to algebraically analyze the stability and the bifurcations of the above two models. Consider the parametric system of the following form.

where ${\boldsymbol{u}} = (u_1, u_2, \ldots, u_d)$ are the parameters, ${\boldsymbol{x}} = (x_1, x_2, \ldots, x_n)$ are the variables, and $F_1, \ldots, F_n$, $N_1, \ldots, N_m$, $P_1, \ldots, P_{a+b}$ are all polynomials. The symbolic approach ^[26] permits us to obtain the necessary and sufficient condition that system (1.1) has a given number of distinct real solutions. The basic idea of this approach is to first transform system (1.1) into a new system in one unique variable and then analyze the real solutions of the resulting system. Different from those based on numerical computations, this approach can establish exact and rigorous results and thus is more adequate for the theoretical study of economic models.

Take the following system as an illustrative example.

where $u, s$ are parameters. In the first step, using the triangular decomposition method, one can transform system (1.2) into a univariate system, i.e., (1.3), without changing the roots. The triangular decomposition method, used here, can be viewed as an extension of the Gaussian elimination method. The main functions of both methods are to transform a system into a triangular form. However, the triangular decomposition method is for polynomial systems, while the Gaussian elimination method is for linear systems. The reader may refer to ^{[19,25,44,47]} and references therein for the triangular decomposition.

where

By virtue of the concept of border polynomial (see ^[26] for its definition), we can then decompose the parameter space $\{(u, s)\, |\, u\in\mathbb{R}, s\in\mathbb{R}\}$ into regions such that in each region the number of distinct real solutions of (1.3) or (1.2) is invariant. Therefore, the number of real solutions in each region can be determined by counting it at a selected sample point. Here, we neglect the rest of the tedious minor steps. Finally, our computational results show that provided

where $R = s^6-3\, s^4-8\, us^2+3\, s^2-1$, the number of distinct real solutions of system (1.2) is

● $0$ if and only if $R < 0$ and $s < 0$;

● $1$ if and only if $R > 0$;

● $2$ if and only if $R < 0$ and $s > 0$.

The reader can refer to ^[26] for more technical details of the above approach.

The novelty of our work is twofold: the methodology and the theoretical results. On one hand, we introduce the symbolic approach proposed by the author into the study of the dynamic behavior of oligopoly games. We show that this approach can establish exact results and might be a powerful tool in theoretical economics, such as theorem proving (see, e.g., the proof of Theorems 1 and 2). In particular, for the triopoly models of ^[34,41], we obtain the analytical conditions of the local stability and prove the existence of double routes to chaos, i.e., through the period-doubling bifurcation and the Neimark-Sacker bifurcation. It is worth noting that, for the two games considered in this paper, the existing investigations are just built on observations through numerical simulations. Our analytical results not only confirm the numerical results in the current literature but also provide a solid foundation for the theoretical development of these games. On the other hand, for the special case of the involved firms having identical marginal costs, we acquire the conditions of the local stability for both models. By further analyzing these conditions, we confirm the statement in ^[7,9,33] that the presence of the LMA mechanism might have a stabilizing effect. Furthermore, we find that there may exist a connection between the two models by showing the stable regions of them would be the same if the proportion parameter of the first model takes the value of $2/3$.

The rest of this paper is structured as follows. In Section 2, we revisit the construction of the two models. In Section 3, the local stability is thoroughly investigated, and the stable regions of the two games are compared analytically in the special case of identical marginal costs. In Section 4, employing the algebraic criteria proposed by Wen and others ^[45,46], we obtain the complete conditions of the occurrence of the period-doubling bifurcation and the Neimark-Sacker bifurcation. Moreover, by the bifurcation continuation, we provide some preliminary numerical results of other observed bifurcations. The paper is concluded with some remarks in Section 5.

2.   Models

Let us consider a market served by three firms producing substitute products. We use $x(t)$, $y(t)$, and $z(t)$ to denote the outputs of the three firms at period $t$, respectively. Assume that the demand function of the market is isoelastic, which is founded on the hypothesis that the consumers have the Cobb-Douglas utility function. Specifically, the price of the product (the inverse demand function) is set to be

where $Q = x+y+z$ is the total supply. Moreover, the marginal cost of each firm is supposed to be constant. To be exact, the cost functions of the three firms are set to be $c_1x$, $c_2y$, and $c_3z$ ($c_1, c_2, c_3 > 0$), respectively.

Under the above assumptions, the first firm should have the profit function

Evidently, its marginal profit would be

The profits $\Pi_2, \Pi_3$, and the marginal profits $\phi_2, \phi_3$ of the second and the third firm can be similarly obtained.

Suppose that the first firm is an adaptive player. At each period, the first firm decides the output $x(t+1)$ according to its previous output $x(t)$ and its expectations of other competitors' outputs at period $t+1$. Specifically, firm $1$ naively expects that at period $t+1$ firms $2$ and $3$ would produce the same quantity as at period $t$. Therefore, the first firm could calculate the optimal output $x_{opt}$ to maximize its expected profit of period $t+1$ such that the first order condition $\phi_1(x_{opt}, y(t), z(t)) = 0$ is fulfilled. That is

The adaptive decision mechanism for firm $1$ is that it chooses the output $x(t+1)$ proportionally to be

where $l\in[0, 1]$ is a parameter controlling the proportion.

The second firm is simply a naive player, who expects that the outputs of its rivals would be equal to those at the last period. Thus, to gain the maximal expected profit, it would choose the best response to produce

Our third firm is assumed to be a boundedly rational player, that adopts the so-called gradient adjustment mechanism, i.e.,

where $k > 0$ is a parameter controlling the adjustment speed. In other words, the third firm increases/decreases its output according to the information given by its marginal profit of the last period. It is worth noting that the adjustment speed depends upon not only the parameter $k$ but also the size of the firm $z(t)$.

In short, the above triopoly game could be described by the following $3$-dimensional discrete dynamic system denoted by $T_{ANB}(x, y, z)$.

The above game was first explored by Tramontana and Elsadany ^[41]. Motivated by their work, a distinct heterogeneous triopoly game was investigated by Naimzada and Tramontana ^[34], where the first firm is replaced with a firm that has further reduced rationality and adopts the so-called local monopolistic approximation (LMA) mechanism ^[43]. The LMA mechanism is also a boundedly rational process, where the player has limited knowledge of the demand function. Specifically, the firm can just observe the current market price $p(t)$ and the corresponding total supply $Q(t)$ and is able to correctly estimate the slope $p'(Q(t))$ of the price function around the point $(p(t), Q(t))$. Then the firm uses such information to conjecture the demand function and expects the price at period $t+1$ to be

where $Q^e(t+1)$ represents the expected aggregate production at period $t+1$. Moreover, the new firm is also assumed to use the naive expectations of its rivals, i.e., $Q^e(t+1) = x(t+1)+y(t)+z(t)$. Thus, we have that

The expected profit of firm 1 in the new model is

To maximize the expected profit, according to the first order condition, the firm chooses its output at period $t+1$ to be

Therefore, the new model can be described by the following $3$-dimensional discrete dynamic system denoted by $T_{LNB}(x, y, z)$.

3.   Local stability

In order to compute the equilibria of the two models introduced in Section 2, we set $x(t+1) = x(t) = x$, $y(t+1) = y(t) = y$ and $z(t+1) = z(t) = z$. Then for model $T_{ANB}$,

and for model $T_{LNB}$,

Both of them can be solved by two equilibria:

It should be noted that the origin $(0, 0, 0)$ is not an equilibrium since it is not defined for maps (2.1) and (2.2). From an economic point of view, only the equilibrium $E_2$ is of our concern because the equilibrium $E_1$, where the third firm is out of the market forever, is unstable as mentioned in ^[41].

The Jacobian matrix takes the form

We use $J_{ABG}$ and $J_{LBG}$ to denote the Jacobian matrices of our two models, respectively. Then we have

and

One can linearize a dynamic discrete system around the equilibrium and then determine the stability conditions via the eigenvalues of the Jacobian matrix of the system. The following proposition follows.

Proposition 1. Let $E$ be an equilibrium of a discrete dynamic system. If all the eigenvalues $\lambda_i$ of the Jacobian matrix at $E$ lie in the open unit disk, i.e., $|\lambda_i| < 1$ for all $i$, then $E$ is locally stable. Moreover, if the Jacobian matrix at $E$ has at least one eigenvalue $\lambda_0$ outside the open unit disk, i.e., $|\lambda_0| > 1$, then $E$ is unstable.

Let $A$ be the characteristic polynomial of the Jacobian matrix $J$. The eigenvalues of $J$ are simply the roots of $A$. So the problem of stability analysis can be reduced to that of determining whether all the roots of $A$ lie in the open unit disk $|\lambda| < 1$. To our knowledge, in addition to the Routh-Hurwitz criterion ^[35] generalized from the corresponding criterion for continuous systems, there are two other criteria, i.e., the Schur-Cohn criterion [11,pp. 246–248] and the Jury criterion ^[18], feasible for discrete dynamic systems. In what follows, we provide a short review of the Schur-Cohn criterion.

Proposition 2 (Schur-Cohn Criterion). For an $n$-dimensional discrete dynamic system, assume that the characteristic polynomial of its Jacobian matrix is

Consider the sequence of determinants $D^\pm_1$, $D^\pm_2$, $\ldots$, $D^\pm_n$, where

The characteristic polynomial $A$ has all its roots inside the unit open disk if and only if

1) $A(1) > 0$ and $(-1)^nA(-1) > 0$,

2) $D^\pm_1 > 0, D^\pm_3 > 0, \ldots, D^\pm_{n-3} > 0, D^\pm_{n-1} > 0$ (when $n$ is even), or $D^\pm_2 > 0, D^\pm_4 > 0, \ldots, D^\pm_{n-3} > 0, D^\pm_{n-1} > 0$ (when $n$ is odd).

Consider a $3$-dimensional discrete dynamic system with its characteristic polynomial being of the form

According to the Schur-Cohn criterion, an equilibrium $E$ is locally stable if the following inequalities are satisfied at $E$.

Consequently, the problem of determining the stability of model $T_{ANB}$ can be transformed into determining the existence of real solutions of

where $x > 0$, $y > 0$, $z > 0$ are used to let us focus on $E_2$ rather than $E_1$. It should be noted that $a_0$, $a_1$ and $a_2$ are rational functions, which are extremely complex and involve the radical expressions $\sqrt{\frac{y+z}{c_1}}$ and $\sqrt{\frac{x+z}{c_2}}$.

The approach illustrated in Section 1 can be used to count real solutions but are only feasible for polynomial systems. Here, the obvious obstacle is that system (3.4) involves radical expressions, and thus can not be directly tackled with this approach. However, the following substitutions are helpful.

We use $\bar a_0$, $\bar a_1$ and $\bar a_2$ to denote the results of $a_0$, $a_1$ and $a_2$ after the above substitutions. Evidently, system (3.4) is equivalent to the following system.

where ${\rm{Numer}}(\cdot)$ and ${\rm{Denom}}(\cdot)$ represent the numerator and the denominator, respectively. System (3.6) can be viewed as a pure polynomial system with the variables $x, y, z, u, v$ and the parameters $k, l, c_1, c_2, c_3$.

In ^[41], numerical simulations rather than symbolic deductions are used to explore the magnitude of the stability region of model $T_{ANB}$. The reason may be that the expressions of the eigenvalues of the Jacobian matrix are so complicated that any rigorous analysis seems to be impossible. However, by the Schur-Cohn criterion, we obtain several analytical results. It is impossible to deal with the complex computation results if we keep all $5$ parameters. For comparison, as in ^[41], we just keep $c_3$ and $k$ as free parameters and set $c_1 = 1.63$, $c_2 = 2.1$ and $l = 0.6$, then we acquire the following theorem using our symbolic approach illustrated in Section 1.

Theorem 1. Assume that $c_1 = 163/100$, $c_2 = 21/10$ and $l = 3/5$^*. The equilibrium $E_2$ of model $T_{ANB}$ is locally stable if

^*Here, we use fractions instead of floats to underline that the inputs of our algorithms are rigorous expressions.

where

and $LS_{ANB}^3$, $LS_{ANB}^4$ are complicated polynomials listed in Appendix.

Remark 1. Let $LS_{ANB}^3 = 0$. One can obtain

When $c_3 = 373/100$, we have

and

Therefore, $k\rightarrow +\infty$ as $c_3$ tends to $373/100$ from the left, i.e., $c_3\rightarrow (373/100)^-$. Consequently, the curves of $LS_{ANB}^2 = 0$ and $LS_{ANB}^3 = 0$ approach arbitrarily close to each other around $c_3 = 373/100$. From an economic point of view, it is known that for any adjustment speed $k$ of firm 3 (even quite large), there always exist some values of the marginal cost $c_3$ of firm 3 such that model $T_{ANB}$ is locally stable. This finding is interesting because for oligopolies with boundedly rational players it is common that the equilibrium would lose its stability if the adjustment speed is fast enough.

The $(c_3, k)$ parameter plane of model $T_{ANB}$ is depicted in Figure 1. One can verify that the equilibrium $E_2$ is locally stable in the grey region surrounded by the green curve, the blue curve, the two black vertical lines, and the horizontal axis. For example, at $(35/10, 5)$ (the red point in Figure 1), we have that

which satisfies the stability conditions given in Theorem 1.

For model $T_{LNB}$, by setting $c_1 = 0.5$ and $c_2 = 0.55$ as in ^[34], we obtain the following theorem.

Theorem 2. Assume that $c_1 = 1/2$ and $c_2 = 11/20$^†. The equilibrium $E_2$ of model $T_{LNB}$ is locally stable if

^†Here, we use fractions instead of floats to underline that inputs of our algorithms are rigorous expressions.

where

and the expressions of $LS_{LNB}^3$, $LS_{LNB}^4$ can be found in Appendix.

Remark 2. If $LS_{LNB}^1 = 0$, then

Furthermore, when $c_3 = 21/20$, we have

Accordingly, $k\rightarrow +\infty$ as $c_3\rightarrow (21/20)^-$, which means that the curve of $LS_{LNB}^1 = 0$ approaches arbitrarily close to $LS_{LNB}^3$ as $c_3\rightarrow (21/20)^-$. Consequently, it is known that for model $T_{LNB}$, even if the adjustment speed $k$ is extremely large, there exist values of the marginal cost $c_3$ such that the equilibrium $E_2$ is locally stable. Therefore, $k$ is large enough and is not a sufficient condition of instability.

The $(c_3, k)$ parameter plane of model $T_{LNB}$ is displayed in Figure 2. One can see that the equilibrium $E_2$ is locally stable in the grey region surrounded by the green curve, the blue curve, the two black vertical lines, and the horizontal axis. For instance, at $(1, 10)$, the red point in Figure 2,

which satisfies the stability conditions of Theorem 2.

By the analytical results of Theorems 1 and 2, we confirm the observations by numerical simulations of ^[34,41] (see Figures 1 and 2). There are several benefits of acquiring analytical, instead of simulated, results. Simulations can only be performed in a bounded area of the parameter space and provide limited information on stability. Analytical results, however, may unveil more features of models' stability. For example, according to Theorem 1 or 2, we know that there exists one single piece of the stable region (the grey area in Figure 1 or 2), and other ones are impossible. Moreover, analytical, rather than simulated, results could be reused to produce further theoretical conclusions such as Remarks 1 and 2.

Bischi and others ^[7] proved that for an oligopoly with an isoelastic demand and identical firms adopting the LMA mechanism, instability arises if there exist at least $5$ firms. While for an oligopoly with a linear demand and identical firms adopting the best response mechanism, the equilibrium goes unstable if at least $4$ firms compete with each other ^[40]. Accordingly, the LMA mechanism may improve the model stability compared to the best response mechanism for games with homogeneous players. In ^[33], it was proved that the equilibrium is globally stable in a homogeneous duopoly game if both firms adopt the LMA mechanism. Moreover, Cavalli and Naimzada ^[9] investigated several heterogeneous duopoly games and found the similar phenomenon that the LMA mechanism may enlarge the size of the stable region, which suggests that the stability is influenced not only by the rationality level of the firms but also by the reaction way the least rational firm uses.

Therefore, it would be interesting to investigate how the LMA mechanism is connected to the local stability of heterogenous triopoly games. The comparison of the two games considered in this paper may expose this connection. To obtain results possible for comparison, we simplify the situation by setting the marginal costs of the three firms identical. The following theorem is finally acquired.

Theorem 3. Suppose that $c_1 = c_2 = c_3 = c$. For model $T_{ANB}$, the equilibrium $E_2$ is locally stable if

Moreover, for model $T_{LNB}$, the equilibrium $E_2$ is locally stable if

The $(c, k)$ parameter plane of the two games is depicted in Figure 3. The stable regions are surrounded by the horizontal axis, the vertical axis, and the corresponding curves. The red curve represents the border of the region described by (3.10). The green and blue curves are the borders of (3.9) with $l = 1$ and $l = 0$, respectively.

According to Theorem 3, if $l = 2/3$, the stable regions of these two models would be the same. This fact is informative and unveils that there might be a connection between an adaptive player and an LMA player^‡. Furthermore, if $0\leq l < 2/3$, the stable region of model $T_{ANB}$ is larger than that of model $T_{LNB}$; if $2/3 < l\leq 1$, the conclusion is the opposite. In the special case when $l = 1$, the first firm in model $T_{ANB}$ degenerates to a naive player adopting the best response mechanism. This confirms the observation in ^[7,9,33] that the LMA mechanism has an effect of stability enhancement and suggests that we may extend a similar conclusion to the case of heterogeneous oligopolies with $n$ players. However, if $0 < l < 1$, a general adaptive decision mechanism is used by the first firm in model $T_{ANB}$. As the proportion parameter $l$ decreases from $1$ to $0$, the first firm adjusts its output with more caution and the stability region becomes larger. The LMA mechanism might lead to a global underestimation of the prices and result in more careful adjustments. Therefore, the stability of an oligopoly may be essentially connected to the adjustment speed of the output.

^‡The two models considered in this paper are differentiated only by the first firm, which is an adaptive player in $T_{ANB}$ but an LMA player in $T_{LNB}$.

4.   Bifurcations

The technique of linearization used for stability analysis may fail at bifurcation points, for near such points the dynamic behavior of the system may differ qualitatively from that of its linearization. This means that a small smooth change of the bifurcation parameters may cause a sudden qualitative change to the dynamic behavior of the system. In this paper, we analytically and numerically investigate serval types of bifurcations, such as period-doubling, Neimark-Sacker, and resonance.

4.1.   Period-doubling bifurcation

The numerical analysis of the period-doubling and the Neimark-Sacker bifurcations has been carried out in ^[34,41]. Instead, we focus on the rigorous proof of the existence of these two bifurcations. First, consider the period-doubling bifurcation. We need the following formal definition, which could also be found in ^[21].

Definition 1 (Period-doubling bifurcation). Consider an $n$-dimensional iteration map $F_{ {\boldsymbol{u}}}:\, {\boldsymbol{x}}(t)\mapsto {\boldsymbol{x}}(t+1)$ with an equilibrium $E$, where ${\boldsymbol{x}}(t) = (x_1(t), x_2(t), \ldots, x_n(t))$ is the state vector at period $t$ and ${\boldsymbol{u}} = (u_1, u_2, \ldots, u_d)$ stands for all the parameters. A period-doubling bifurcation, also called flip bifurcation, takes place at a bifurcation parameter point ${\boldsymbol{u}} = \overline {\boldsymbol{u}}$ if and only if the following two conditions are satisfied.

$1)$ Eigenvalue assignment: the Jacobian matrix at $E$ has one real eigenvalue $\lambda_1(\overline {\boldsymbol{u}}) = 1$, and $|\lambda_j(\overline {\boldsymbol{u}})| < 1$ for all the rest eigenvalues $\lambda_j$, $j\neq 1$.

$2)$ Transversality condition: for any parameter $u_j\; (1\leq j \leq d)$,

It is known that ranges of parameter values exist where the dynamics are periodic if the period-doubling bifurcation occurs. From an economic point of view, it is quite realistic to assume that a boundedly rational firm can not learn the pattern behind quantities and profits if a long period of periodic dynamics takes place. Wen and others ^[46] generalized the Schur-Cohn criterion and proposed a new criterion for the detection of the period-doubling bifurcation, which is as follows.

Proposition 3. Let $A = \lambda^n + a_{n-1}\lambda^{n - 1} + \cdots + a_0$ be the characteristic polynomial of the Jacobian matrix of an $n$-dimensional iteration map $F_{ {\boldsymbol{u}}}$ and $D^\pm_i$ be the same as in the Schur-Cohn criterion. A period-doubling bifurcation appears if and only if the following conditions are satisfied.

$1)$ Eigenvalue assignment: at the bifurcation point ${\boldsymbol{u}} = \overline {\boldsymbol{u}}$,

(a) $A(1) > 0$ and $A(-1) = 0$,

(b) $D^\pm_1 > 0, D^\pm_3 > 0, \ldots, D^\pm_{n-3} > 0, D^\pm_{n-1} > 0$ (when $n$ is even), or

$D^\pm_2 > 0, D^\pm_4 > 0, \ldots, D^\pm_{n-3} > 0, D^\pm_{n-1} > 0$ (when $n$ is odd).

2) Transversality condition: for any parameter $u_j\; (1\leq j\leq d)$,

where $\partial a_i/ \partial u_j$ is the partial derivative of $a_i({\boldsymbol{u}})$ with respective to $u_j$.

For a $3$-dimensional iteration map as (2.1) and (2.2), the criterion becomes simple and is restated in the following corollary.

Corollary 1. Consider a $3$-dimensional discrete dynamic system. Let $A = \lambda^3 + a_2\lambda^2 + a_1\lambda + a_0$ be the characteristic polynomial of its Jacobian matrix. There is a period-doubling bifurcation if and only if

For model $T_{ANB}$, in order to calculate the bifurcation parameter surface determined by $-1+a_2-a_1+a_0 = 0$, we consider

and then transform it into a polynomial system by applying the substitutions (3.5). Using the triangular decomposition method, we obtain three triangular sets, i.e., (4.3). The union of the zeros of the resulting triangular polynomial systems is the same as the zeros of the system generated by substituting (3.5) into (4.2).

where $PD_{ANB}$ is complicated and can be found in Appendix. Obviously, $\mathcal{T}_1$ is corresponding to $(x, y, z) = (0, 0, 0)$, where map (2.1) is not defined. The triangular set $\mathcal{T}_3$ is corresponding to $E_1$. However, the equilibrium $E_1$ is not of our concern, so we focus only on $\mathcal{T}_2$. It is clear that in $\mathcal{T}_2$ the last polynomial $PD_{ANB}$ involves only $c_1$, $c_2$, $c_3$, $l$ and $k$ (all the parameters) but none of $u$, $v$, $x$, $y$ and $z$, thus the bifurcation parameter surface of the period-doubling bifurcation for model $T_{ANB}$ is $PD_{ANB} = 0$. In addition, if the bifurcation could take place, it is also required the transversality condition is satisfied. Consequently, we obtain the following result.

Theorem 4. For model $T_{ANB}$, the bifurcation parameter surface of the period-doubling bifurcation is $PD_{ANB} = 0$. Furthermore, if a period-doubling bifurcation takes place, the following transversality condition should be satisfied.

See Appendix for the expressions.

Similarly, we have the following theorem for model $T_{LNB}$.

Theorem 5. For model $T_{LNB}$, the bifurcation parameter surface of the period-doubling bifurcation is $PD_{LNB} = 0$. Furthermore, if a period-doubling bifurcation takes place, the following transversality condition should be satisfied.

See Appendix for the expressions.

It is easy to verify that $PD_{ANB}$ becomes $LS_{ANB}^3$ if setting $c_1 = 1.63$, $c_2 = 2.1$, $l = 0.6$ as in Theorem 1. Moreover, $PD_{LNB}$ becomes $LS_{LNB}^3$ if setting $c_1 = 0.5$, $c_2 = 0.55$ as in Theorem 2. Thus, one knows that the equilibrium $E_2$ would lose its stability through the period-doubling bifurcation when a parameter point in the stable region goes across the blue curve of Figure 1 or 2. Furthermore, for we keep only $c_3$ and $k$ as the parameters in Section 3, the transversality conditions with respect to $c_3$ and $k$ are required, which are $PT^k_{ANB}\neq 0, PT^{c_3}_{ANB}\neq 0$ for model $T_{ANB}$ and $PT^k_{LNB}\neq 0, PT^{c_3}_{LNB}\neq 0$ for model $T_{LNB}$, respectively. These conditions are shown in Figures 1 and 2 with the red curves. One can see that the red curves do not pass through the segments of the blue curves where the period-doubling bifurcation takes place for both models. Therefore, the problem of violating the transversality condition does not occur in our settings.

Furthermore, in the case of the identical marginal costs $c_1 = c_2 = c_3 = c$, $PD_{ANB}$ would be $(33\, ckl-68\, ck-102\, l+192)c^3$ corresponding to (3.9), and $PD_{ANB}$ would be $(23\, ck-62)c^2$ corresponding to (3.10). This means that for both models the equilibrium $E_2$ can lose its stability via the period-doubling bifurcation when the marginal costs of the firms are identical. The bifurcation diagrams with respect to $k$ are given in Figure 4, where the bifurcations against $x$, $y$, and $z$ are marked in red, blue, and green, respectively.

4.2.   Neimark-Sacker bifurcation

Quasi-periodic solutions are more difficult to deal with than periodic solutions. From an economic point of view, a quasi-periodic trajectory looks like chaos to boundedly rational firms with incomplete information and limited capability of calculation. Furthermore, if quasi-periodic dynamics appear, the pattern behind productions and profits is also hard to learn even for completely rational players. Thus, it is extremely challenging for a firm to handle a quasi-periodic economy, where no market rules could be discovered and followed. When a Neimark-Sacker bifurcation takes place, dynamics suddenly become quasi-periodic. The reader may refer to ^[16] for the following formal definition of the Neimark-Sacker bifurcation.

Definition 2 (Neimark-Sacker bifurcation). Consider an $n$-dimensional iteration map $F_{ {\boldsymbol{u}}}:\, {\boldsymbol{x}}(t)\mapsto {\boldsymbol{x}}(t+1)$ with an equilibrium $E$, where ${\boldsymbol{x}}(t) = (x_1(t), x_2(t), \ldots, x_n(t))$ is the state vector at period $t$ and ${\boldsymbol{u}} = (u_1, u_2, \ldots, u_d)$ stands for all the parameters. A Neimark-Sacker bifurcation takes place at a bifurcation parameter point ${\boldsymbol{u}} = \overline {\boldsymbol{u}}$ if and only if the following two conditions are satisfied.

$1)$ Eigenvalue assignment: the Jacobian matrix at $E$ has a pair of complex conjugate eigenvalues $\lambda_1(\overline {\boldsymbol{u}})$, $\lambda_2(\overline {\boldsymbol{u}})$ on the unit circle, i.e., $|\lambda_1(\overline {\boldsymbol{u}})| = |\lambda_2(\overline {\boldsymbol{u}})| = 1$, and $|\lambda_j(\overline {\boldsymbol{u}})| < 1$ for all the rest eigenvalues $\lambda_j$, $j\neq 1, 2$.

$2)$ Transversality condition: for any parameter $u_j\; (1\leq j \leq d)$,

Based on the Schur-Cohn criterion, a criterion for the identification of the Neimark-Sacker bifurcation was proposed in ^[45], which is restated below.

Proposition 4. Let $A = \lambda^n + a_{n-1}\lambda^{n - 1} + \cdots + a_0$ be the characteristic polynomial of the Jacobian matrix of an $n$-dimensional iteration map $F_{ {\boldsymbol{u}}}$ and $D^\pm_i$ be the same as in the Schur-Cohn criterion. A Neimark-Sacker bifurcation occurs if and only if the following conditions are satisfied.

1) Eigenvalue assignment: at the bifurcation point ${\boldsymbol{u}} = \overline {\boldsymbol{u}}$,

(a) $A(1) > 0$ and $(-1)^nA(-1) > 0$,

(b) $D^\pm_1 > 0, D^\pm_3 > 0, \ldots, D^\pm_{n-3} > 0, D^+_{n-1} > 0, D^-_{n-1} = 0$ (when $n$ is even), or $D^\pm_2 > 0, D^\pm_4 > 0, \ldots, D^\pm_{n-3} > 0, D^+_{n-1} > 0, D^-_{n-1} = 0$ (when $n$ is odd).

2) Transversality condition: for any parameter $u_j\; (1\leq j\leq d)$,

Corollary 2. Consider a $3$-dimensional discrete dynamic system. Let $A = \lambda^3 + a_2\lambda^2 + a_1\lambda + a_0$ be the characteristic polynomial of its Jacobian matrix. The Neimark-Sacker bifurcation occurs if and only if

Through a series of computations, the following two theorems are obtained.

Theorem 6. For model $T_{ANB}$, the bifurcation parameter surface of the Neimark-Sacker bifurcation is $NS_{ANB} = 0$. Furthermore, if a Neimark-Sacker bifurcation appears, the following transversality condition should be satisfied.

See Appendix for these expressions.

Theorem 7. For model $T_{LNB}$, the bifurcation parameter surface of the Neimark-Sacker bifurcation is $NS_{LNB} = 0$. Furthermore, if a Neimark-Sacker bifurcation occurs, the following transversality condition should be satisfied.

See Appendix for these expressions.

We could see that $NS_{ANB}$ and $NS_{LNB}$ become $LS_{ANB}^4$ and $LS_{LNB}^4$, respectively, if setting $c_1 = 1.63$, $c_2 = 2.1$, $l = 0.6$ as in Theorem 1 and setting $c_1 = 0.5$, $c_2 = 0.55$ as in Theorem 2. This means that across the green curve in Figure 1 or 2 as the parameters vary, the equilibrium $E_2$ would lose its stability through the Neimark-Sacker bifurcation. Moreover, since we only keep $c_3$ and $k$ as free parameters, the transversality conditions with respect to $c_3$ and $k$ are needed, which are $NT^k_{ANB}\neq 0, NT^{c_3}_{ANB}\neq 0$ for model $T_{ANB}$ and $NT^k_{LNB}\neq 0, NT^{c_3}_{LNB}\neq 0$ for model $T_{LNB}$. These conditions are plotted as magenta curves in Figures 1 and 2. It is obvious that the magenta curves do not pass through the segments of the green curves where the Neimark-Sacker bifurcation is developed for both models. Therefore, in our settings, we do not need to worry about the possibility that the transversality condition might be violated. Furthermore, in the case of identical marginal costs, plugging $c_1 = c_2 = c_3 = c$ into $NS_{ANB}$ and $NS_{LBN}$, one can check that it is impossible that $NS_{ANB}$ or $NS_{LBN}$ vanishes, which means the Neimark-Sacker bifurcation can not occur for our two models.

4.3.   Other observed bifurcations

In what follows, we use continuation methods^§ to explore other types of bifurcations as complements to the current literature on numerical investigations of the two models considered in this paper. For this purpose, we employ $\textsf{MatContM}$, a powerful numerical tool implemented by Govaerts and others ^[15]. Specifically, $\textsf{MatContM}$ is a MATLAB toolbox, which can be applied to detect the dynamic behavior of an iteration map such as equilibria, periodic orbits, homoclinic orbits, heteroclinic orbits, and bifurcation points.

^§The main idea of continuation methods is to trace solution manifolds of various objects while we continuously change the parameters of the considered iteration map.

Let $k$ and $c_1$ be free parameters. Using $\textsf{MatContM}$, we observe other bifurcations such as resonance, generalized period-doubling, and Chenciner bifurcations. For model $T_{ANB}$ or map (2.1), we find three kinds of resonance, i.e., the 1:2 resonance, the 1:3 resonance, and the 1:4 resonance labeled R2, R3, and R4 in Figure 5 (a), respectively. The 1:2 resonant point is the intersection of the Neimark-Sacker bifurcation curve and the period-doubling bifurcation curve. The reader may refer to ^[22,23] for more information on resonance bifurcations. The Chenciner bifurcation is also identified (labeled CH in Figure 5 (a)) to provide insights into the transition between different closed invariant circles. The Chenciner bifurcation may introduce much more complex dynamic behavior including the coexistence of two stable attractors and chaos. The label GPD in Figure 5 (a) is used to locate the generalized period-doubling point generated from the 1:2 resonant point.

Furthermore, Figure 5 (b) is plotted to display the numerical results of model $T_{LNB}$ or map (2.2) discovered by $\textsf{MatContM}$. Similarly, we find three different resonance bifurcations, i.e., the 1:2 resonance (R2), the 1:3 resonance (R3), and the 1:4 resonance (R3). In addition, two generalized period-doubling (GPD) points are detected. One may notice that in Figure 5 (b) the left segment of the Neimark-Sacker bifurcation curve marked in magenta is incorrect. See Figure 2, based on analytical conditions of Theorem 2, for the correct Neimark-Sacker bifurcation curve. This means that numerical results may be unreliable in some cases and need to be confirmed by symbolic results.

In the above discussion, we simply provide the results by the bifurcation continuation and leave the detailed bifurcation analysis to our future research.

5.   Conclusions

In this paper, we have analyzed the dynamics of two different heterogeneous triopoly games under the assumption that the demand function of the market is isoelastic. Nonlinearities are introduced both in the demand function and in the decisional mechanisms adopted by the firms. We have considered several decisional mechanisms including the naive (best response) mechanism, the boundedly rational (gradient adjustment) mechanism, the adaptive mechanism, and the local monopolistic approximation (LMA) mechanism. In particular, in both games, the second firm is a naive player, while the third firm is a boundedly rational player. The only difference is that firm 1 is an adaptive player in model $T_{ANB}$, but an LMA player in model $T_{LNB}$.

The existing results of triopoly games in the literature are usually based on observations through numerical simulations. In comparison, using the approach proposed in ^[26], we have obtained several analytical results for the two models considered in this paper. In particular, we have rigorously proved the existence of double routes, through the period-doubling bifurcation and the Neimark-Sacker bifurcation, from the stable state to the chaotic behavior. In addition, we have considered the special case that the involved firms have identical marginal costs and acquired the analytical conditions of the local stability for each model. It has been found that the presence of the local monopolistic approximation (LMA) mechanism might have a stabilizing effect on heterogeneous triopoly games with the isoelastic demand. Using $\textsf{MatContM}$, a numerical toolbox, we have observed other bifurcations such as resonance, generalized period-doubling, and Chenciner bifurcations.

In our future work, the rigorous bifurcation analysis and the proof of the existence of chaos are theoretically interesting although they have been observed via simulations. Moreover, we are also curious to learn how to control chaos in oligopoly games using symbolic methods.

Acknowledgments

The author wishes to thank Dr. Bo Li for the beneficial discussions and is grateful to the anonymous referees for their helpful comments. This research was funded by the Philosophy and Social Science Foundation of Guangdong under Grant GD21CLJ01 and the Major Research and Cultivation Project of Dongguan City University under Grant 2021YZDYB04Z.

Conflict of interest

The author declares no conflict of interest.

Appendix