Abstract
Transmission of COVID-19 occurs either through living media, such as interaction with a sufferer, or nonliving objects contaminated with the virus. Recovering sufferers and disinfectant spraying prevent interaction between people and virus become the treatment to overcome it. In this research, we formulate a new mathematical model as a three-dimensional ordinary differential equation system representing an interaction between viruses attached in nonliving media, susceptible, and infected subpopulations, including the treatment to investigate its effect. Disease-free, sterile-media endemic, and two nonsterile media endemic equilibriums exist in the model. The nonexistence of sterile-media equilibria interpreting the nonendemic condition is achieved by crossing the branch point bifurcation of the equilibria point as the infected subpopulation recovery rate increases. Continuation of the limit cycle generated at a Hopf bifurcation point as susceptible-coronavirus interaction prevention rate and period increase trigger two saddle-node bifurcations and a branch point bifurcation of cycle. Stable symmetric cycles with decreasing amplitude that make the dynamic of subpopulation easier to control start to be gained at the branch point bifurcation of cycle between the two saddle-node bifurcation points as the prevention rate increases. Some chaotic attractors which describe a complex and unpredictable pattern of the dynamic in the population are also found at inclination flip bifurcation by a continuation of a homoclinic orbit generated near the Bogdanov-Takens bifurcation point as the prevention rate increases while the recovery rate decreases. Increasing the recovery and prevention rate along with avoiding an increase of the prevention rate while the recovery rate decreases becomes the treatment to optimize the effort in overcoming COVID-19 transmission.
1. Introduction
COVID-19 is a deadly disease first detected in Wuhan, Hubei Province, China, in late December 2019 [1, 2]. The causative of this disease was identified as a novel coronavirus (nCoV) and has been named coronavirus or severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) [3]. The COVID-19 pandemic has emerged since March 2020, as stated by WHO. By 23 September 2022, there have been over 611 million COVID-19 positive cases, including more than 6 million deaths worldwide [4].
Coronavirus reproduces in its natural host, i.e., bats [5]. Transmission of COVID-19 occurs through interaction between animals and humans or between humans carrying coronavirus [6, 7]. Interactions between the suspected subpopulation and the infected subpopulation can cause contagion through droplets and close interaction. Things or nonliving media exposed to the virus is also a media of COVID-19 transmission [8] because the virus is stable for several hours to days in aerosols and on surfaces [9].
The experiment found that SARS-CoV-2 remained viable in aerosols throughout the experiment (3 hours) [9]. No viable SARS-CoV-2 was measured on copper after 4 hours [9]. On cardboard, no viable SARS-CoV-2 was measured after 24 hours [9]. Furthermore, SARS-CoV-2 was more stable on plastic and stainless steel than on copper and cardboard, and a viable virus was detected up to 72 hours after application to these surfaces [9]. On different surfaces of various materials, such as paper, glass, PVC, metal, ceramic, and Teflon, the virus can survive up to 5 days [10]. It provides critical information about the stability of SARS-CoV-2, and it is still possible for the virus to infect people after touching contaminated nonliving media.
Recently, many researchers have established numerous models for COVID-19 under the different concepts of fractional calculus, such as the SIR model, to characterize COVID-19 transmission. But the model still focuses on the interaction between humans as the main cause of COVID-19 transmission [11] and has not considered COVID-19 transmission through nonliving media along with the effort to inhibit it. In this research, we construct a three-dimensional ordinary differential equation system as a new mathematical model representing the interaction between virus, susceptible, and infected subpopulations to investigate the impact of inhibiting the growth of the virus attached to nonliving media on COVID-19 transmission.
2. Model Development
Transmission of COVID-19 occurs through interaction between humans, i.e., the susceptible and the infected subpopulation, and between humans touching contaminated nonliving media. We assume that interaction between susceptible and infected subpopulations cannot be avoided, so repressive treatment is the most appropriate treatment to suppress the negative effect of the interaction. The treatment will affect the virus and susceptible subpopulation by preventing interaction between the subpopulations. Another assumption is that the rate of preventing interaction between virus and susceptible subpopulations can be greater than the rate of interaction between susceptible and virus subpopulations.
The increase of the virus on surfaces or nonliving media is caused by the infected subpopulation of respiratory droplets expelled through the air and contaminating objects. Sneezing and coughing are the common cause of droplet expulsion. Surface viruses can survive various times depending on the surface type, and it is still possible for the virus to infect people after touching it. In this research, we assume that the infection rate from a virus attached to nonliving media is constant.
Based on the interaction of virus, susceptible, and infected subpopulation that has been described, several variables and parameters of the model can be defined in order to construct the model. The model variables are defined in Table 1.
Henceforth, will be written as coronavirus to simplify the writing. The model variables are nonnegative since , , and denote subpopulation. Model parameters are summarized in Table 2.
We define as the difference of susceptible-coronavirus interaction and its prevention rate.
The model parameters are nonnegative because they show the level of interaction between virus, susceptible, and infected subpopulations. In addition, the parameter is not zero because the zero value of the parameter means no interaction in the system or the parameter does not affect the model. Also, some parameter values are set (assumed) to get a qualitative result. The increase in the virus subpopulation is assumed to be expressed by the parameter which does not depend on variables other than virus subpopulation to accommodate all possible factors that increase virus subpopulation and simplify the model.
Based on the interactions between virus, susceptible, and infected subpopulations which are described as variables and parameters of the model, we construct a compartment diagram presented in Figure 1.
Based on the compartment diagram in Figure 1, a mathematical model can be formed as a system of nonlinear ordinary differential equations with three-dimensional variables and eight-dimensional parameters. The model is presented as follows:
Equation (1) represents the rate of change of the virus population with respect to time which is denoted by a logistic form with as the coronavirus addition rate and as the ratio between the coronavirus addition rate to the maximum population of the virus. The maximum population of the virus is .
Equation (2) represents the rate of change in the susceptible population with respect to time. The first term is the increase of susceptible population caused by the susceptible population’s natural birth by . The second term is the reduction of the susceptible population due to its interaction with the virus, where denotes the interaction rate and denotes the interaction prevention rate. The third term is the reduction of susceptible population caused by the susceptible population’s natural death with as its rate. The fourth term is the increase of susceptible population caused by the infected population recovery where represents the rate. The fifth term is reducing the susceptible population due to its interaction with the infected population with as its rate.
Equation (3) represents the rate of change in the infected population with respect to time. The first term is the increase of infected population due to the interaction between the susceptible population and the virus where denotes the interaction rate and denotes the interaction prevention rate. The second term is the increase of the infected population caused by its interaction with the susceptible population by . The third term is the reduction of the infected population caused by the infected population recovery, where represents the rate. The fourth term is the reduction of infected population caused by the infected population death due to coronavirus with as its rate.
3. Results and Discussion
3.1. Equilibrium Points
Equilibrium points represent a state which illustrates a steady phenomenon in a very long time condition. We investigate the equilibrium points of the model by solving [14–16] that interprets a static number of each subpopulation over time. Three equilibrium points of the model were found, i.e., the disease-free equilibria, the sterile media endemic equilibria, and the nonsterile media endemic equilibria.
Theorem 1. The disease-free equilibria of the model is . The sterile media endemic equilibria is . Then, the nonsterile media endemic equilibria is and , where and .
Proof. By setting , we have
From (4), we have or . We substitute to (6) so we have or . Then, we substitute and to (5) so we get . Based on the calculation, we get an equilibria
This equilibrium is called disease-free equilibrium because it has a zero number of the virus and infected subpopulation, i.e., there is no positive case of COVID-19 in this condition.
Then, by substituting and to (5), we obtain . Based on the calculation, we get an equilibria
This equilibrium is called the endemic balance of sterile media because it does not contain viral subpopulations attached to nonliving media, but the number of infected subpopulations is not zero.
By eliminating equations (5) and (6), we get . By substituting and to (6), we have where and .
Based on the calculation, we get an equilibria
where and . These equilibrium points are called the first and second nonsterile media endemic equilibria, respectively, because they contain a nonzero number of the virus subpopulation attached to nonliving media and a nonzero number of the infected subpopulation.
3.2. Existence and Biological Feasibility of Equilibrium Points
The equilibrium points exist if every subpopulation is well-defined mathematically, and they are biologically feasible if every subpopulation has a nonnegative or positive value because it interprets a number of subpopulations which cannot possibly be negative. The disease-free equilibria are biologically feasible if the virus, susceptible, and infected subpopulation have a nonnegative value. The endemic equilibrium points are biologically feasible if the virus and susceptible subpopulation have a nonnegative value while the infected subpopulation has a positive value. It should be positive because the endemic equilibrium points represent the existence of the infected subpopulation.
Theorem 2. The disease-free equilibria exists and it is biologically feasible in every condition. The sterile media endemic equilibria, exists in every condition and it is biologically feasible if or . The first nonsterile media endemic equilibria exists if and it is not biologically feasible in any condition. The second nonsterile media endemic equilibria exists if and it is biologically feasible if , where and .
Proof. In disease-free conditions, exists and it is biologically feasible in every condition because the parameter values are positive, including and . Therefore, .
In the endemic condition when the media is sterile, every subpopulation in is well-defined mathematically, because they have nonzero denominators and real numerators. On the other hand, the infected subpopulation expression in should have a positive value to represent an endemic condition, so that it is biologically feasible. It is satisfied if the numerator, i.e., while the denominator, i.e., is always positive based on the condition that all the parameters are positive. Hence, exists in every condition and it is biologically feasible if which is equivalent to or
In the first endemic condition when the media is not sterile, the value of should be nonnegative in order to make the susceptible and infected subpopulation expressions in well-defined mathematically. Hence, which is equivalent to . Thus, exists if . Then, the value of should be positive and should be nonnegative in order to satisfy the condition that requires a positive value of the infected subpopulation and a nonnegative value of the susceptible subpopulation, so that it is biologically feasible. But, if is nonnegative, we obtain which is equivalent to . It implies which equivalent to This condition is contradictive with the requirement to guarantee the value of the susceptible subpopulation that should be nonnegative, i.e., . Thus, is not biologically feasible in any condition.
In the second endemic condition when the media is not sterile, the value of should be nonnegative in order to make the susceptible and infected subpopulation expressions in well-defined mathematically. It implies which is equivalent to . Thus, exists if . Then, the value of should be positive and should be nonnegative in order to satisfy the condition that requires a positive value of the infected subpopulation and a nonnegative value of the susceptible subpopulation, so that it is biologically feasible. Based on the form of and , the value of is always greater than if . Beside that, if is nonnegative, we obtain the value of which is equivalent to . It implies which equivalent to Hence, we found that is biologically feasible if
3.3. Local Stability of Equilibrium Points
Fluctuation of the number of viruses, susceptible, and infected subpopulations around the equilibrium points are represented by their local stability. We predict the dynamic of each subpopulation starting around the equilibrium point by analyzing the local stability using the linearization method [14–16]. We also use MAPLE to compute some results of algebraic operations, so that calculation error is avoided and accuracy of the calculation is obtained.
Theorem 3. Local stability of , , , and are identified as follows. (i) is saddle if (ii)Local stability of is categorized as follows (1)Consider (a)If , then is unstable (node)(b)If , then is saddle(2)Consider (a)If , then is saddle(b)If , then is unstable (node)(3)Consider (a)If , then is saddle-focus (unstable)(b)If , then is focus-node (stable)(iii)Local stability of is categorized as follows (1)Consider (a)If , then is asymptotically stable(b)If then is saddle(2)Consider . If , then there are two cases as follows (a) implies is asymptotically stable(b) implies is saddle(3)Consider . If , then is neutral stable(4)Consider (a)If , then is focus-node (stable)(b)If , then is saddle-focus (unstable)(iv)Local stability of is categorized as follows (1)Consider (a)If , then is asymptotically stable(b)If , then is saddle(2)Consider . If , then is saddle(3)Consider . If , then is saddle
Proof. Linearization of the model yields the Jacobian matrix of the model as follows.
for , , and
Consider as the eigenvalue of the Jacobian Matrix in (10), we obtain the characteristic equation of as follows:
Hence, we obtain , , and Based on the eigenvalues, is saddle if .
Based on the Jacobian matrix, we found that a characteristic equation of is
Hence, we obtain and , where because of biological feasibility condition. Based on the expression of , we categorize local stability conditions as follows:
(1)Consider . We found two cases below
(a)If , then . Hence, is unstable (node)(b)If , then Hence, is saddle(2)Consider . We obtain two cases below
(a)If , then Hence, is saddle(b)If , then . Hence, is unstable (node)(3)Consider We found two cases below
(a)If , then where Hence, is saddle-focus (unstable)(b)If , then where . Hence, is focus-node (stable)Based on the Jacobian matrix, we found a characteristic equation of as follows:
where and , so that and . Based on the expression of , we categorize local stability conditions as follows:
(1)Consider which is equivalent to . We obtain two cases below
(a)If which is equivalent to , then Hence, is asymptotically stable(b)If which is equivalent to , then Hence, is saddle(2)Consider which is equivalent to . If which is equivalent to which is only satisfied if , then there are two implications, i.e.
(a) which is equivalent to yields , so that Hence, is asymptotically stable(b) which is equivalent to yields , so that Hence, is saddleThe inequality will not be satisfied for any value, because there is no value that satisfies .
(3)Consider which is equivalent to . If which is equivalent to , then Hence, is neutral stable. The inequality will not be satisfied for any value, because there is no value that satisfies (4)Consider which is equivalent to . We obtain two cases below
(a)If which is equivalent to , then where Hence, is focus-node (stable)(b)If which is equivalent to , then where Hence, is saddle-focus (unstable)Based on the Jacobian matrix, we found a characteristic equation of as follows:
where and , so that and . Based on the expression of , we categorize local stability conditions are as follows:
(1)Consider which is equivalent to . We obtain two cases below
(a)If which is equivalent to , then Hence, is asymptotically stable(b)If which is equivalent to , then Hence, is saddle(2)Consider which is equivalent to . If which is equivalent to which is only satisfied if , then the value of are always opposite in sign for any value of . Hence, is saddle. The inequality will not be satisfied for any value, because there is no value that satisfies (3)Consider which is equivalent to . If which is equivalent to which is only satisfied if , then the value of are always opposite in sign. Hence, is saddle. The inequality will not be satisfied for any value, because there is no value that satisfies The inequality which is quivalent to will not be satisfied for any value, because the value on the left-hand side is always nonnegative, while the value on the right-hand side is always negative.
3.4. Bifurcation Analysis
We define the difference between susceptible-coronavirus interaction and its prevention rate as , i.e., . We investigate the effect of treatment-related parameters variation on population dynamic, i.e., infected subpopulation recovery rate and the difference between susceptible-coronavirus interaction and its prevention rate . We choose and as the varying parameters because they are related to preventive or repressive actions that can be taken to overcome corona. Based on the bifurcation analysis, we will investigate the effect of varying those parameters on the dynamics of coronavirus, susceptible, and infected subpopulation in order to determine some criteria used to get an optimal treatment for corona.
The characteristic equation in (13) is equivalent to where , , and . We use the characteristic equation in (15) to identify some bifurcations in codimension-one and codimension-two bifurcation. We also use MAPLE again to compute some results of algebraic operations.
3.5. Codimension-One Bifurcation
Codimension one bifurcation is investigated in order to characterize the effect of a varying parameter on the dynamics of the solution that represents the fluctuation of coronavirus, susceptible, and infected subpopulation.
Theorem 4. Branch point bifurcation of equilibria occurs at
Proof. By considering the existing condition of and , we obtain that and exist in every condition. If we set , we obtain that and are two different existing equilibrium points. Then, if we increase , so that , we found that the infected subpopulation in is and the susceptible subpopulation in is , so that . It means that and collide at . Then, if we increase , so that , we found that and are two different existing equilibrium points again. Thus, a branch point bifurcation of equilibria occurs at .
Theorem 5. Fold bifurcation occurs at the following values: and
Proof. Based on the theory in [17] and the characteristic equation in (15), we found that fold bifurcation occurs if which is equivalent to . It is satisfied if . By doing some algebraic calculations in the equation and defining , , , we found that is equivalent to which is satisfied if . The values of that fulfill it are and .
Consider , , , , and .
Theorem 6. Hopf bifurcations occur at the following values , , and .
Proof. Based on the theory in [17] and the characteristic equation in (15), we obtain that Hopf bifurcation occurs when which is equivalent to . It fulfilled if or . The values of that satisfies are and . Besides that, the values of that satisfies are and .
3.6. Codimension-Two Bifurcation
Codimension-two bifurcation is investigated in order to characterize the effect of two varying parameters on the dynamics of the solution that represent the fluctuation of coronavirus, susceptible, and infected subpopulation.
Theorem 7. Bogdanov-Takens bifurcations occur at the following values and .
Proof. Based on the theory in [17] and the characteristic equation in (15), we obtain that Baogdanov-Takens bifurcation occurs when which is equivalent to . It is satisfied if . By solving with respect to , we obtain or . Then, by solving the equation formed after substituting to with respect to , we found and . Then, by solving the equation formed after substituting to with respect to , we found and .
4. Numerical Continuation
We do a numerical analysis to determine the other bifurcations that have not been identified analytically due to the complexity, such as the saddle-node bifurcation of cycle, branch point bifurcation of cycle, and homoclinic bifurcations. Analytical detection of saddle-node bifurcation of cycle and branch point bifurcation of cycle requires the determination of the cycle analytically which in general, in a three-dimensional system, the cycle is not known, so that it can only be identified numerically [18]. Analytical detection of homoclinic bifurcations needs the determination of homoclinic orbit analytically whereas it is usually hard to do because the occurrence of it is associated with a global bifurcation, so in general, we have to determine it numerically in order to analyze and continue its loci in some parameters [19]. Besides that, when we tried to apply the Melnikov method [20] to determine the homoclinic orbit; we had a hard time transforming the model such that it is reduced to the Melnikov condition. This kind of difficulty has also been presented in [21]. Hence, those discussions become some open problems in this research.
We use the bifurcations that have been identified analytically as the initial point to investigate the other bifurcations by conducting some continuations of the equilibrium points, bifurcation value, and limit cycle as the parameters vary using MATCONT. Furthermore, some phase portraits which describe the different phenomena that occurred when the parameters vary are illustrated using MATLAB. Some specific phenomena, such as homoclinic orbits and an indication of chaotic dynamic are also found and investigated in the simulation. We do not investigate the bifurcations with neutral saddle properties, because they have no dynamical interpretation [22].
4.1. Effect of Infected Subpopulation Recovery Rate Variation
Continuation of virus-absence endemic equilibria as the infected subpopulation recovery rate varies generates a branch point bifurcation of equilibria and two Hopf bifurcations with neutral saddle properties, see Figure 2.
(a)
(b)
(c)
Branch point bifurcation of equilibria will occur at by increasing the infected subpopulation recovery rate when saddle type of disease-free and virus-absence endemic equilibrium exists. The saddle type of sterile media endemic equilibria will not be biologically feasible, while the saddle type of disease-free equilibria is still biologically feasible if the recovery rate continues to increase.
From a biological point of view, increasing the infected subpopulation recovery rate implies the disappearance of endemics in the population with sterile media while maintaining the existence of disease-free condition. Branch point bifurcation of equilibrium value becomes the minimum limit of the recovery rate so that an endemic does not exist in the population when the media is sterile.
4.2. Effect of Difference of Susceptible-Coronavirus Interaction and Its Prevention Rate and Infected Subpopulation Recovery Rate Variation
In codimension-one bifurcation, continuation of nonsterile media endemic equilibria as a difference between susceptible-coronavirus interaction and its prevention rate varies triggers fold and Hopf bifurcation. In codimension-two bifurcation, a continuation of the Hopf point generates the Hopf curve, which triggers two Bogdanov-Takens when it intersects the fold curve generated by continuing the fold point. The fold curve also triggers zero-Hopf bifurcation that has neutral saddle properties as and varies, see Figure 3.
(a)
(b)
Increasing the preventing rate of human-nonhuman host interaction causes the collision of two endemic equilibria at the fold point, i.e., , and then they disappear. An unstable limit cycle is constructed at Hopf point, i.e., as the preventing rate decreases a bit. The first Bogdanov-Takens bifurcation is obtained at when the prevention rate increases and the infected subpopulation recovery rate decreases while the second Bogdanov-Takens bifurcation is obtained at when the prevention rate decreases and the infected subpopulation recovery rate increases.
The biological interpretation of the phenomena is described as follows. The fold bifurcation phenomenon represents the nonexistence of endemic conditions as the prevention rate increases. The unstable limit cycle is constructed when Hopf bifurcation is obtained as a disease cycle [23], i.e., a continuous fluctuation of susceptible and infected subpopulations in COVID-19 transmission. The Bogdanov-Takens bifurcations represent some complex phenomena indicated by homoclinic orbit and strange attractors that are constructed near them. Discussion of those phenomena will be presented in the related section.
4.3. Limit Cycle Continuation
The continuation of the unstable limit cycle formed at the Hopf point as the prevention rate varies generates two saddle-node bifurcations of cycles and one branch point bifurcation of cycle, see Figure 4.
(a)
(b)
(c)
A branch point bifurcation of cycle at is obtained by decreasing the prevention rate. Symmetric cycles are generated after the occurrence of this bifurcation [24]. After that, the first saddle-node bifurcation of cycle occurs at when the prevention rate keeps decreasing. At this stage, the unstable cycle that is generated at the Hopf point is annihilated, and stable cycles with increasing amplitude start to be gained until the second saddle-node bifurcation of cycle at . At the second saddle-node bifurcation of cycle, the stable cycle is annihilated and the unstable cycle is constructed.
The biological meaning of the mathematical results above is explained as follows. After the second saddle-node bifurcation of cycle, stable symmetric cycles with decreasing amplitude start to be formed until the first saddle-node bifurcation of the cycle constructs an unstable symmetric cycle as the rate of prevention increases. The stable symmetric cycles ensure that the fluctuation of susceptible and infected subpopulation are easier to control through their stability property and regular pattern. Increasing the prevention rate makes fluctuation even easier to control, because the amplitude of the stable symmetric cycles decreases. The two saddle-node bifurcations of the cycle become the indicator (checkpoint) [25] to generate cycles that play a role as the control of susceptible and infected subpopulation dynamics. The best value of the prevention rate to control the dynamic is at the first saddle-node bifurcation of cycle because the amplitude of the stable symmetric cycle formed at this state is the smallest one.
4.4. Homoclinic Orbit and Homoclinic Bifurcation Generated near the Bogdanov-Takens Bifurcation
The Bogdanov-Takens bifurcations trigger some homoclinic orbits through continuation. Continuation of these homoclinic orbit in codimension-two bifurcation generates some homoclinic bifurcations [26–29]. A family of homoclinic orbits near and some homoclinic bifurcations that were found by making a continuation of the homoclinic orbit as and vary are presented in Figure 5.
Orbit-flip with respect to the stable manifold (OFS), inclination-flip with respect to stable manifold (IFS), and neutrally-divergent saddle focus (stable) (NDS) bifurcation are obtained through a forward continuation of the homoclinic orbit generated near as and vary. Backward continuation of the homoclinic orbit generated near as and vary yields Neutral saddle (NS) and orbit-flip concerning the stable manifold (OFS) bifurcation.
A family of homoclinic orbits near and some homoclinic bifurcations that were found by making forward and backward continuation of a homoclinic orbit are shown in Figure 6.
Noncentral homoclinic to saddle-node (NCH), inclination-flip with respect to unstable manifold (IFU), inclination-flip with respect to stable manifold (IFS), and orbit-flip with respect to the stable manifold (OFS) bifurcation are found through a backward continuation of the homoclinic orbit constructed near .
The homoclinic orbits show a complex dynamic around the Bogdanov-Takens bifurcations. It joins a saddle equilibria to itself, i.e., they lay in the intersection of stable and unstable manifolds of equilibria.
The biological interpretation of the homoclinic orbits is described as follows. The homoclinic orbits represent a complex pattern of susceptible and infected subpopulation dynamics which occur when the prevention rate is equal to the infected subpopulation recovery rate or the prevention rate increases as the infected subpopulation recovery rate decreases.
4.5. Strange Attractor Generated near Homoclinic Bifurcation
Homoclinic orbit is possible to become an indication of a chaotic dynamic [30]. We especially investigate the attractors generated near an IFU bifurcation point , because the dynamic of the attractors is possible to be chaotic [31]. The possibility of the chaotic dynamic is investigated by calculating the Lyapunov exponent of coronavirus, susceptible, and infected subpopulation numerically because they are difficult or impossible to calculate analytically [32]. Justification of chaotic dynamic by analytical proof becomes an open problem in this research. The attractors and their Lyapunov exponent are illustrated in Figure 7.
(a)
(b)
(c)
(d)
Attractor generated at is possible to be chaotic because the Lyapunov exponent of the variable representing the coronavirus subpopulation has a positive value, i.e., . Dynamic difference of the phase portrait presented in Figures 7(a)–7(c) constructed by using the same parameter set and some slightly different initial values also show that the attractor has a sensitive response to the initial value, which is the characteristic of a chaotic attractor.
From a biological point of view, the phenomena represent an indication of the complexity that occurs in the dynamic. It illustrates that the fluctuation of susceptible and infected subpopulations have an irregular pattern that is difficult to predict when the prevention rate increases while the recovery rate decreases. This case represents an incident when the susceptible coronavirus prevention rate of humans increases, but on the other hand, the recovery rate of the infected subpopulation decreases.
5. Conclusions
Interaction between coronavirus attached in nonliving media, susceptible, and infected subpopulations in a population can be represented by a mathematical model. Numerical bifurcation analysis of the model gives some characteristics that represent the effect of infected subpopulation recovery rate and susceptible-coronavirus prevention rate variation on changing the dynamics of COVID-19 transmission.
Increasing infected subpopulation recovery rate implies the disappearance of endemic in the population while maintaining the biological feasibility of the disease-free condition. Branch point bifurcation of equilibrium value becomes the minimum limit of the recovery rate so that the endemic does not exist in the population when the media is sterile.
Increasing the rate of preventing interaction between viruses and susceptible subpopulations implies the collision of two endemic conditions at the fold point and then disappear. It represents the nonexistence of endemic conditions in the population when the media is not sterile as the prevention rate increases. At the condition represented by the Hopf point, there is a cycle that interprets a continuous fluctuation of the susceptible and infected subpopulations.
When the media is not sterile, increasing the susceptible-coronavirus prevention rate makes the fluctuation of infected and susceptible subpopulations easier to control by generating a stable cycle of the infected and susceptible subpopulations with regular patterns. The two saddle-node bifurcations of cycle become the sign to generate the cycles. The best value of the prevention rate to control the fluctuation of the infected subpopulation is at the first saddle-node bifurcation of cycle.
According to the indication of chaotic dynamic in the population when the media is not sterile, which represents an unpredictable fluctuation pattern of a susceptible and infected subpopulation, increasing the susceptible-coronavirus prevention rate and decreasing the recovery rate of the infected subpopulations should be avoided as a precaution against the phenomenon.
Data Availability
The data used to support the findings of this study are included in the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This research and publication article was funded by the Universitas Sebelas Maret through Hibah Penelitian dan Pengabdian kepada Masyarakat (P2M) Penerimaan Negara Bukan Pajak (PNBP) Universitas Sebelas Maret 2022 (254/UN.27.22/PT.01.03/2022).