# Complexity Analysis of a 2D-Piecewise Smooth Duopoly Model: New Products versus Remanufactured Products

Authors: S. S. Askar; A. Al-khedhairi
Journal: Complexity
 (2022)

---

## Abstract

Recent studies on remanufacturing duopoly games have handled them as smooth maps and have observed that the bifurcation types that occurred in such maps belong to generic classes like period-doubling or Neimark-Sacker bifurcations. Since those games yield piecewise smooth maps, their bifurcations belong to the so-called border-collision bifurcations, which occur when the map’s fixed points cross the borderline between the smooth regions in the phase space. In the current paper, we present a proper systematic analysis of the local stability of the map’s fixed points both analytically and numerically. This includes studying the border-collision bifurcation depending on the map’s parameters. We present different multistability scenarios of the dynamics of the game’s map and show different types of periodic cycles and chaotic attractors that jump from one region to another or just cross the borderline in the phase space.


---

## Body

## 1. Introduction

Recently, many countries have adopted the remanufacturing process in order to promote their economy. This process requires remanufacturing of used or default products as remanufacturing is considered a friendly environmental process. There are many popular firms that have commenced this process, such as Ford, Xerox, and Caterpillar [1, 2]. Indeed, the remanufacturing process may be beneficial for some companies to increase their profits, but it may be worse for other companies. It may not be useful for companies which produce new products because it may deprive sales of new products and therefore hurt their expected profits. For these reasons, companies that produce new or original products prefer not to follow remanufacturing and designate this strategy to third-party companies. Based on environmental conditions, those third-party companies accumulate used or default products so that they remanufacture them again and then send them to the market for selling [3–6]. According to USA economy magazines, the industry of remanufacturing in the USA is worth $53 billion [7, 8]. It reveals the fact that original products manufacturers face a competitive threat from remanufacturer companies. Such competition between those companies and its complex dynamic characteristics can be described and investigated by duopoly games.Competition in duopoly games includes only two competed firms whose strategies may be quantities (as in Cournot) or prices (as in Bertrand). Duopoly games and their complex dynamic characteristics have been analyzed by several authors in the literature. For instance, the duality of quantities and prices in a differentiated duopoly game proposed by Dixit has been analyzed in [9]. In [10], the dynamics of a Cournot duopoly game whose players are boundedly rational and adopt a gradient adjustment mechanism have been investigated. A duopoly game with firms that adopted two different adjustment mechanisms, local monopolistic approximation and gradient-based, has been studied in [11]. Price competition of a nonlinear duopoly game has been illustrated when price demand elasticity varies in [12]. For more information on the adjustment mechanisms, readers are advised to see [13]. In [14], a Bertrand duopoly game whose players want to maximize their relative profits has been introduced and discussed. Different duopoly models based on Cobb-Douglas utility have been analyzed in [15]. A duopoly model of technological innovation based on constant conjectural variation and boundedly rational players has been studied in [16]. In [17], the dynamic properties of a Bertrand duopoly game with two-stage delay have been investigated. For more detailed studies and information on the complex dynamic characteristics of duopoly games, we refer to the literature [18–22].There are few studies on the remanufacturing of duopoly games in literature. Those studies have introduced such games based on smooth models and hence they have investigated the stability conditions of the fixed points of that models. For example, the dynamic characteristics of the fixed points of this game have been introduced and discussed in [3, 23]. In spite of the stability investigations given in those previous references, they have improperly handled their games. They have studied their games as a game described by a smooth discrete dynamical system. The current paper is motivated by the game given in [3]. Our motivation gives rise to a proper investigation of the game and in addition, introduces a rich analysis of the stability of the game’s fixed point by analyzing the border-collision bifurcation by which the fixed point becomes unstable. Indeed, such games should be described by a piecewise smooth map as given in the current manuscript. Piecewise-smooth map means that the map is defined by many parts and the phase plane is divided into regions separated by a border curve. Consequently, the dynamic of map in these regions should not be studied separately as given in [3, 23], simply because the dynamic of destabilization of the map’s fixed point is usually never related to only one region. The dynamics of the map should be studied based on parameter values belonging to these regions where the map is defined.The current paper gives a proper description of the map used to define this game in [3]. Before that, we recall some important aspects of piecewise smooth maps. There are several applications of different disciplines in the literature that were modeled using piecewise-smooth maps. Such applications include physiological and economic systems, mechanical systems, and switching circuits. They have reported many interesting observations regarding the kinds of bifurcations that destabilize such maps. They have observed that the bifurcation that occurred in these piecewise smooth maps does not belong to a generic class of bifurcation such as period doubling or Neimark-Sacker. Instead, new types of bifurcations have emerged from these maps when the fixed point crosses the border curve separating the smooth regions on where the maps are defined. Now, we report some important studies on the complex dynamic characteristics of piecewise smooth maps. For example, in [24], a general three-dimensional piecewise smooth map and its multiple attractors have been investigated. A global analysis of a piecewise-smooth Cournot duopoly model derived from an isoelastic demand function has been performed in [25]. In [26], multistability investigations and border bifurcations have been discussed for a 2D piecewise-smooth map. In [27], a 2D map that is continuous, noninvertible, and piecewise smooth has been used to characterize the innovation activities in the model of trade and product innovation of two countries. The global dynamics of the piecewise-smooth map that describes the spatial distribution of long-run industrial activity for the economic geography model have been analyzed in [28]. Interesting readers of piecewise smooth maps and their dynamic characteristics may be directed to see literature [29–34].Further, we have to highlight other investigations in literature that uphold this research direction and raise new directions in future work. In [35], a review report on new research directions on evolutionary games and their rules which are adopted to exploit the benefits of cooperation that might be occurred among players. The characteristics of cooperation under the synergies between evolutionary game and networks have been reviewed in [36]. Discussing the factors that might affect the overall performance of a firm has been determined and analyzed in [37]. Recently, in [38], the authors have identified the most influential invaders in cooperative communities based decomposition of the weighted-degree mechanism.After the above introduction, the current paper is briefly organized as follows. In Section2, we give a proper version of the map used to describe the remanufacturing game presented in [3]. In Section 3, local analysis of the map’s fixed point on the smooth regions where the map is defined is analytically investigated. Section 3 presents rich numerical simulation experiments which uphold obtained results in Section 2. In Section 5, we give a brief conclusion.


## 2. The Model

Let us consider an economic market populated by two competed firms. The first firm is called a manufacturer and supports the market with new products, while the second firm is called a third-party remanufacturer and supports the market with differentiated remanufactured products. The customers are willing to pay for productions by those two firms. The demand productions sent to the market by firms are denoted byx1 and x2. The competition between these firms is carried out in discrete time periods, t=0,1,2,…. The first firm sends the new quantity x1 for selling in the market at time t, while the second firm can receive returned quantity x2 for remanufacturing and sells it again in the market at time t+1. According to the literature [3, 23], this discussion can be described by the following inverse demand functions (The inverse demand functions areas in [1]). But the authors have defined their game’s map as a smooth one which is incorrect. They have studied each part of the map separately; however, it had to be studied as a piecewise-smooth map. In a addition, they had to investigate the fixed points as a function of the map’s parameters when they belong to two different regions separated by border line on where the piecewise map was defined. This is for simple reason that is the dynamics of destabilization of the fixed points are usually never related to only one region. Thus the map’s dynamics had to be really analyzed as a piecewise one with points jumping from one region to the other:(1)p1=1−x1−δx2,p2=δ1−x1−x2.The parameterδ has an important meaning in this game. If δ=1, it means that the customers are prepared to pay the same price for new and remanufactured products. From an economic perspective, this may not be approved. If δ=0, the customers will not pay any price for remanufactured products. Because of the variety of customers, we restrict this parameter to δ∈0,1. Assuming that Cixi refers to the cost of the quantity xi and is given by the following linear form:(2)Cixi=cixi,i=1,2,where c1 and c2 refer to the marginal costs for both firms, respectively. Therefore, the profits of both firms are given as follows:(3)π1=1−x1−δx2−c1x1,π2=δ1−x1−x2−c2x2.And their marginal functions become as follows:(4)∂π1∂x1=1−2x1−δx2−c1,∂π2∂x2=δ1−x1−2x2−c2.Due to incomplete information both firms can get from the market, they update their productions based on partial information. So we assume that the two firms are heterogeneous and adopt different adjustment mechanisms in order to update their productions. We presume that the first firm will behave as a bounded rational firm and consequently will maximize its profit. It is easy to see that forδ∈0,1 and c1>c2, the marginal profit ∂π1/∂x1 is always positive and lies within the first quadrant. Therefore, the first firm will increase its output to get the maximum of its profit. This mechanism is called bounded rationality and has been intensively adopted in the literature ([10–15]) to describe such firm’s behavior. Based on this reasoning the first firm will update its output at period t+1 according to the following form:(5)x1,t+1=x1,t+kx1,t∂π1,t∂x1,t,where k is a positive parameter. On the other hand, we assume that the second firm seeks to share the market with a certain profit. The weighted sum method [39] is used to handle this behavior. It starts with assuming that it seeks a complete market share maximization. This makes π2=0 and then its optimum output becomes as follows:(6)x¯2=1−x1−c2δ.But when it completely seeks profit maximization, its marginal profit will vanish and then we have the following:(7)x^2=121−x1−c2δ.According to some weights, the second firm will be traded off between market share and profit as follows:(8)x˜2=ωx¯2+1−ωx^2=1+ω21−x1−c2δ,where ω∈0,1. When ω=0, it means the second firm seeks profit maximization only while ω=1 means it seeks market share maximization. But as it trades off between both, we have restricted the parameter ω to the interval 0,1. Now we assume this firm updates its output according to the following adaptive mechanism:(9)x2,t+1=1−βx2,t+βx˜2,where β is a positive parameter and is restricted to the interval β∈0,1. Since the output of the second firm at time t+1 should be less than or equal to those of the first firm at time t, x1,t+1≤x1,t, therefore, (9) will be modified as follows:(10)x2,t+1=min1−βx2,t+βx˜2,x1,t.Using (5) and (10), one gets the map that describes this game as follows:(11)x1,t+1=x1,t+kx1,t1−2x1,t−δx2,t−c1,x2,t+1=min1−βx2,t+βx˜2,x1,t,The map (11) is a two-dimensional piecewise smooth map and is constructed to describe the proposed duopoly game in this paper. In order to study the stability of its fixed points we should study it as a piecewise-smooth map not as in [3]. It is worthwhile to mention that the authors in [3] have improperly studied their map as a smooth map. Here, in this paper we introduce a proper investigation for the map (11). The stability of this map depends on the second equation where the function Ϝ=1−βx2,t+βx˜2 exists. It means that any order pair x1,x2 belongs to Ϝ<x1 gives only one part of the map to be dynamically studied, while if Ϝ>x1 the other part of the map should be studied. Hence, we get the fact that there is a borderline Ϝ=x1 where the map is continuous. Furthermore, the map’s phase plane will be divided by the border line into two regions (Left region, Rℓ and right region Rr). This borderline takes the following form:(12)x1=12+β1+ω21−βx2+βδ1+ωδ−c2≕gx2.So the map (11) is modified to the following:(13)x1,t+1=x1,t+kx1,t1−2x1,t−δx2,t−c1,x2,t+1=x2,t+1=1−βx2+12β1+ω1−x1−c2δ,ifx1,t≥gx2,x1,t,ifx1,t≤gx2.or, equivalently, it becomes as follows:(14)x1,t+1,x2,t+1=x1,t+kx1,t1−2x1,t−δx2,t−c1,x2,t+1=1−βx2+12β1+ω1−x1−c2δ,ifx1,t≥gq2regionRr,x1,t+1,x2,t+1=x1,t+kx1,t1−2x1,t−δx2,t−c1,x1,tifx1,t≤gq2regionRℓ.Now, we come to the fact that the stability of the map (14) depends on points located in region Rr from the borderline or on points located in the region Rℓ. Furthermore, other attracting sets may appear at the same parameter values. But when the map’s fixed points are destabilized, this means that their dynamics are usually never related to only one region in the phase space. As a result, the dynamic of the map (14) must be studied as a piecewise smooth map with points jumping from region Rr to region Rℓ and vice versa.


## 3. Local Analysis

It is clear that the map (14) is defined in two different regions separated by a borderline. This gives rise to two fixed points in the two regions. We denote these two points Er (for the one in region Rr) and Eℓ (for the point in region Rℓ), where(15)Er=21−c1−1+ωδ−c24−1+ωδ,1+ω1+c1δ−2c24−1+ωδδ,Eℓ=1−c12+δ,1−c12+δ.These fixed points are positive under the conditions,δ/21−c11+ω<1+ωδ−c2<21−c1 and c1<1. The Jacobian matrices at (15) become as follows:(16)Jr=1−41−c1−21+ωδ−c24−1+ωδk−21−c1−1+ωδ−c24−1+ωδδk−1+ωβ21−β,Jℓ=1−21−c1δ+2k−δ1−c1δ+2k10.The eigenvalues for the above Jacobian matrices take the form,λ1,2=1/2τ±τ2−4Δ, where τ and Δ refer to the trace and determinant respectively. It should be noted that τ and Δ may be taken as τr and Δr (or τℓ and Δℓ) depending on whether the fixed point lies within region Rr (or region Rℓ). Those traces and determinants in both regions take the following form:(17)τr=2−β−221−c1−1+ωδ−c24−1+ωδk,τℓ=1−21−c1δ+2k,Δr=1−β+β1+ω2δ2+δ1+ω4−c2β1+ω−2β3−c1−41−β21−c1+c21+ω24−1+ωδk,Δℓ=δ1−c1δ+2k.and the eigenvalues are as follows:(18)λ1r,2r=1−β2−A1k±β24−β2A2k−A12k2,λ1ℓ,2ℓ=12−1−c12+δk±14−1−c11+δ2+δk+1−c12+δ2k2,where A1=21−c1−1+ωδ−c2/4−1+ωδ and B=2−1+ωδ21−c1−1+ωδ−c2/4−1+ωδ. Or equivalently,(19)λ1r,2r=12τr±τr2−4Δr,λ1ℓ,2ℓ=12τℓ±τℓ2−4Δℓ.We should highlight that the map (14) is continuous, and it also has continuous derivatives in each region. Its derivative is discontinuous at the borderline Ϝ=x1. Moreover, the map may have no fixed point in half of the phase space. It means that the location Eℓ (or Eℓ) may turn out to be in region Rr (or region Rℓ) and consequently, the character of a virtual fixed point uprises. Furthermore, the triangle of stability for each fixed point is given by the following:(20)Sr=τr,Δr:1+τr+Δr>0,1−τr+Δr>0,1−Δr>0,Sℓ=τℓ,Δℓ:1+τℓ+Δℓ>0,1−τℓ+Δℓ>0,1−Δℓ>0,where(21)1+τr+Δr=8−4β−21−c1−1+ωδ−c21+ωδβ+42−β4−1+ωδk,1−τr+Δr=β221−c1−1+ωδ−c2k,1−Δr=β−21−c1−1+ωδ−c21+ωδβ+41−β24−1+ωδk,1+τℓ+Δℓ=2−1−c12−δ2+δk,1−τℓ+Δℓ=1−c1k,1−Δℓ=1−δ1−c12+δk.


## 4. Simulation

We carry out in this section several simulation experiments to get more insights into the local and global analysis of the map (14). This simulation will deeply investigate different qualitative behaviors of that map. We analyze the influences of map’s parameters on the stability of the fixed points using some tools of numerical simulations such as 1D and 2D bifurcation diagrams, largest Lyapunov exponent, time series plot, phase diagram, and basin of attraction. We start our numerical experiments with the initial datum x1,0,x2,0=0.19,0.15 with the parameter values, ω=0.5,c1=0.5,c2=0.3, and δ=0.5525470489. At these parameter values we get Er=Eℓ=0.19588,0.0.19588. This means that we have two coincided fixed points born in the border line Ϝ=x1. At the same time, the other two parameters k and β are considered the bifurcation parameters. Both the Jacobians defined in (16) depend on those parameters and to calculate them we assume β=0.5 and k=5.72. This gives the following:(22)Jℓ=−1.2409−0.6191010,Jr=−1.2409−0.61910−0.3750.5,where Δr=−0.85261,τr=−0.7409,Δℓ=0.61910, and τℓ=−1.2409 with λ1ℓ,2ℓ=−0.62045±0.48388i and λ1r,λ2r=0.62446,−1.3654. It is clear that Δr<1 and Δℓ<1 that means the map’s dynamics in both regions are dissipative. Furthermore, we have −1<Δr<0,0<λ1r<1 and λ2r<−1 and hence a transverse homoclinic intersection may not exist and the obtained attractor may not be chaotic. But since 1+τrτℓ−Δr−Δℓ+ΔrΔℓ>0, then the attractor will be a cycle of 2-period only. Figure 1(a) shows the bifurcation diagram at those parameter values when varying k. The figure gives rise to a bifurcation diagram that stops at the cycle of 2-period only for any increase of that parameter. That is confirmed in the corresponding Largest Lyapunov exponent (or LLE) given in Figure 1(b). In Figure 1(c), one can see that there are three different sections and the border line crosses them. Those sections are colored by cyan, yellow, and grey. The cyan and yellow colors show the basin of attraction of the locally stable period 2-cycle while the grey one depicts the basin of attraction of diverging trajectories. Moreover, the border line contains the two coincided fixed points born on it. Fixing the previous parameter values and increasing β above gives rise to higher periodic cycles with respect to k. For example, at β=0.8 a stable period 4-cycle exists when varying the parameter k. Figure 2(c) presents a stable period 4-cycle that is born due to period-doubling bifurcation. The basin of attraction of this cycle is given in Figure 1(e). Increasing β further makes periodic cycles such as period-8, period-16, and higher period cycles appear till β approaches close to 1 on, where the dynamics of the map evolve in a chaotic attractor. Figure 1(e) shows four unconnected chaotic areas in both regions. Now, we investigate the influence of the parameter k when β takes values close to 0. At the parametr values, ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.11 and k=5 we get the following:(23)Jℓ=−0.95883−0.5411710,Jr=−0.95883−0.54117−0.08250.89,where Δr=−0.8980,τr=−0.06883,Δℓ=0.54117, and τℓ=−0.95882 with λ1ℓ,2ℓ=−0.47941±0.55797i and λ1r,λ2r=0.91384,−0.98267. It is clear that Δr<1 and Δℓ<1 that means the map’s dynamics in both regions are also dissipative. Simulation shows that at small values of the parameter β, the influence of the parameter k becomes very bad. Figure 2(a) shows the bifurcation diagram when varying the parameter k at the same parameter values. At k=5.1 the dynamic of the map gives a stable period 2-cycle. Its basin of attraction is given in Figure 2(b) where there are two periodic points separated by the border line. Increasing the parameter k more gives unstable periodic cycle and routes to chaotic attractors are obtained as shown in Figures 2(c) and 2(d). Both figures present chaotic attractors that are separated by the border line or cross it.Figure 1
(a) The stable period 2-cycle when varying the parameterk at the parameter values: ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489, and β=0.5. (b) Largest Lyapunov exponent corresponds to the bifurcation diagram in Figure 1(a). (c) Basin of attraction of the stable period 2-cycle at ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.5 and k=5.72. (d) The stable period 2-cycle when varying the parameter k at the parameter values: ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489 and β=0.8. (e) Basin of attraction of the stable period 4-cycle at ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.8 and k=6.2. (f) Phase plane for the chaotic attractor at: ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.999 and k=6.2.
(a)(b)(c)(d)(e)(f)Figure 2
(a) Bifurcation diagram at varying the parameterk. (b) Basin of attraction of the period 5-cycle at k2=5.1. (c) Phase plane for the chaotic attractor at ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.11 and k=5.78. (d) Phase plane for the chaotic attractor at ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.11 and k=5.95.
(a)(b)(c)(d)Now, we analyze the impact of the parameterβ on the map’s dynamic while keeping the other parameter values fixed. Let us assume ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489, and k=5. Figure 3(a) presents the bifurcation diagram when varying β. It is clear that the two fixed points born on the border line are locally stable and above β=0.1473 a period 2-cycle arises. Numerical experiments show that the influence of the parameter β changes as k increases further. For example, at k=5.09 and keeping the other parameter values fixed, a bifurcation diagram is given in Figure 3(b). Another bad impact of the parameter β on the map’s dynamics is given in the bifurcation diagram in Figure 3(c) at k=6. This makes us investigate more the dynamic behavior of the map when assuming different values of those two parameters. This can be done by plotting the 2D bifurcation diagram for them. It is depicted in Figure 3(d) on where different periodic cycles can be obtained. The numbers from 1 to 10 show different types of periodic cycles. For the parameter values ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489,β=0.9244, and k=6.1133, the basin of attraction of stable period 6-cycle is depicted in Figure 3(e). At the same set of parameter values used in Figure 3(e) but for β=0.182 and k=5.989, we get a chaotic behavior of the map’s dynamic which after that is turned into a stable period 7-cycle. We simulate this behavior of the map by the time series given in Figure 3(f).Figure 3
(a) Bifurcation diagram at varying the parameterβ. (b) Bifurcation diagram at varying the parameter β and k=5.09. (c) Bifurcation diagram at varying the parameter β and k=6. (d) 2D Bifurcation diagram at ω=0.5,c1=0.5,c2=0.3,δ=0.5525470489. (e) The basin of attraction of stable period 6-cycle. (f) The time series for the periodic 7-cycle.
(a)(b)(c)(d)(e)(f)At the same set of parameters and forβ=0.158 and k=5.491 a stable period 8-cycle arises and is distributed in both regions from the border line. Another stable period 10-cycle is given in Figure 4(b) at the same set of parameter values but for β=0.158 and k=5.413. Carrying out more numerical experiments about the dynamic of the map (14) gives rise to more complicated behavior of it. For example, at the same set of parameter values used in Figures 4(a) and 4(b) and for β=0.3088 and k=5.9733, a one-piece chaotic attractor that crosses the border line is given in Figure 4(c). From the above discussions that are obtained when the second firm uses symmetric weights ω=0.5, we conclude that the smooth-piecewise map (14) describing the heterogeneous duopoly is characterized by higher degree of unpredictability. In order to end our discussion in this paper, we give two examples for the asymmetric case (when ω≠0.5). Assuming the following parameter values, ω=0.33,c1=0.5,c2=0.3,δ=0.5822914793,β=0.3088 and k=5.7 a stable period of 18-cycle arises. This cycle is plotted with the coincided fixed points born on the border line in Figure 4(d). Assuming the following parameter values, ω=0.65,c1=0.5,c2=0.3,δ=0.5326030850,β=0.3088 and k=6, Figure 4(e) presents a chaotic behavior of the map that is changed into a stable period 18-cycle. The time series given in Figure 4(f) simulates the behavior given in Figure 4(e).Figure 4
(a) Basin of attraction of the period 8-cycle. (b) Basin of attraction of the period 10-cycle. (c) Phase space of chaotic attractor. (d) Basin of attraction of the period 18-cycle. (e) Chaotic attractor that is changed into stable period 18-cycle. (f) Times series for the chaotic attractor given in Figure4(e).
(a)(b)(c)(d)(e)(f)


## 5. Conclusion

We have analyzed in this paper the dynamic behavior of a remanufactured duopoly game. Previous works in literature have studied and discussed such a game as a smooth map. In the present paper, we have given a proper investigation of the competition carried out between firms in this game. Our investigation and analysis are entirely different than those given in [3, 23] for the same game. Analytically, we have analyzed the piecewise-smooth map describing such a game and illustrated the local stability conditions of its fixed points in both regions where the map was defined. Numerically, we have enriched this paper with intensive numerical simulation experiments about the global analysis of the map’s fixed points that have shown how the map’s dynamics may become quite complex. To the extent of our knowledge, the obtained results in this paper have provided new analytical results and proper investigations of the model investigated by the authors in [3]. This includes studying the game’s model as a piecewise smooth map. Our interesting results have detected different scenarios of multistability situations and period cycles points jumping from one region to another or passing through the borderline. Our future studies will be directed to the importance of networks in dealing with such remanufacturing games.



---

*Source: 8856009-2022-04-30.xml*