SEIQRS model analysis and optimal control with two delays

Abstract

This paper addresses the limitations of assuming a bilinear infection rate in computer virus propagation models by proposing a more realistic nonlinear dual-delay SEIQRS model. The key focus of the research includes analyzing the local stability of the disease-free equilibrium point and the existence of positive equilibrium points with different delays. To further enhance control effectiveness, time-varying control terms are introduced, and a corresponding Hamiltonian function is constructed. The optimal control strategy is derived using the Pontryagin’s maximum principle. Numerical simulation experiments are conducted to validate the dynamic characteristics of the model, and the model’s properties are verified on real-time propagation networks. The experimental results demonstrate that appropriate control strategies can effectively suppress the spread of computer viruses.

1. Introduction

In the context of the current information age, the rapid advancements in modern science and technology, particularly in the realm of computers and electronic communications, have brought about significant transformations in both production processes and daily lifestyles. The Internet, in particular, has emerged as an indispensable and integral component of contemporary life. According to the 52nd Statistical Report on Internet Development in China [1] published by the China Internet Network Information Center (CNNIC), as of June 2023, the number of Internet users in China had reached 1.079 billion, reflecting an increase of 11.09 million since December 2022. Moreover, the Internet penetration rate had reached 76.4% in the country. The advent of the internet has brought about numerous conveniences for people, but it has also created favorable conditions for the rapid spread of computer viruses. Unlike biological viruses, computer viruses are artificially created program codes that search for other programs or media that meet their infection criteria when executed on a computer. Once a target is identified, the virus inserts its own code and proceeds to replicate itself. If an infected file is not promptly addressed, it can have a significant impact on human life. Therefore, researching the propagation patterns and control strategies of computer viruses holds great significance and profound implications. It is in this context that the field of computer virus epidemiology has emerged.

Computer virus epidemiology is an interdisciplinary field that investigates the propagation patterns and control strategies of computer viruses in networks. It draws upon methods and models from disciplines such as epidemiology, complex network theory, and dynamical systems. This field provides a theoretical foundation and technical means for preventing and mitigating network security threats. Scholars both domestically and internationally primarily focus on propagation models and mechanisms. They employ mathematical models and complex network theory to describe the spreading process of viruses in networks and explore the patterns and mechanisms of virus propagation. This research enables the development of corresponding defense strategies, which are of crucial significance in formulating virus defense strategies and enhancing network security.

1. (1). Research on computer virus propagation models

Research on the dynamics of computer virus propagation can be traced back to the late 1980s. In the late 1980s, F. Cohen [2] published a paper elucidating the existence of computer viruses and proposed the possibility of studying the propagation patterns of computer viruses using dynamic models borrowed from the field of infectious diseases in medicine. This work garnered significant attention from the academic community. In 1991, J.O. Kephart and S.R. White [3], drawing inspiration from a class of classical infectious disease models known as the SIS compartmental model, successfully established the first compartmental model for computer virus propagation. Subsequently, with the widespread adoption of the Internet, computer virus propagation entered a new stage. Researchers began collecting empirical data on actual virus propagation and conducting large-scale empirical analyses to understand the characteristics and patterns of virus propagation. During this stage, propagation models evolved beyond theoretical frameworks and became more practical, taking into account factors such as network topology and user behavior.

With the in-depth study of complex network theory, researchers have started applying it to computer virus propagation models. Complex network models can better describe the structure and characteristics of real networks, thereby improving the accuracy of virus propagation models. In 2020, Chenquan Gan et al. [4] studied the propagation behavior of computer viruses on complex networks under the joint influence of network topology structure and mobile storage media. They conducted a carefully designed set of experiments on the Oregon routing network and found that mobile storage media have a significant impact on virus propagation. Ronghua Shang et al. [5] proposed a dynamic node immunization model based on community structure and thresholds (NICT). This model considers nodes carrying harmful information as new nodes in the network and changes the way new nodes establish edges with original nodes based on the needs of different networks. They determine the propagation probability between nodes by using community structure information and a similarity function between nodes, and propose an improved immunization gain. Experimental results show that this method can provide better network immunity. In 2022, Zhang Ningbo et al. [6] proposed a hybrid resource allocation scheme and established a virus-resource asymmetric coupling propagation model to simulate the virus propagation process under the influence of resource recovery. They used a generalized discrete Markov chain approach to describe the dynamics of propagation. Juan Liu [7] studied the propagation patterns of malicious code in wireless sensor networks and established a prevalence model for wireless sensor networks with graded infection rates and two time delays. They discussed the impact of the two time delays on propagation, providing reference and effective suggestions for network administrators. In 2023, Rahat Zarin et al. [8] combined fractional calculus with the nonlinear interactions between nodes in wireless sensor networks and proposed a nonlinear fractal fractional mathematical model for analyzing the propagation of malicious code in wireless sensor networks. They elucidated the propagation patterns of viruses in vulnerable systems and potential countermeasures. In 2024, Kun Li et al. [9] combined evolutionary game theory with complex network theory to study the competition effect of different travel strategies on virus propagation. Simulation results show that self-isolation strategies can significantly inhibit the spread of infections in complex social networks, and the introduction of reward mechanisms further enhances this effect. Additionally, a higher network degree favors the prevalence of self-isolation, thereby impeding virus propagation.

In recent years, with the development of artificial intelligence and data science, researchers have started incorporating behavior modeling into computer virus propagation models. This type of model takes into account factors such as user behavior, social network relationships, etc., to more accurately predict the path and speed of virus propagation. In 2020, Dang Q.A. et al. [10] proposed a new method based on the Lyapunov stability theorem and its extension, combined with the global stability theorem of discrete-time nonlinear fractional-order systems, to theoretically prove that the proposed NSFD (Non-Standard Finite Difference) scheme preserves the global asymptotic stability of the original model. They also presented the positivity and global stability preservation of the NSFD scheme for computer virus mixed propagation models. In 2021, Shashank Awasthi et al. [11] considered the charging mechanism of low-energy sensor nodes and proposed the SILRD model, which describes the dynamic propagation of malicious software and the energy consumption of sensor nodes due to malicious software attacks. They also proposed defense techniques for preventing the propagation of malicious software in wireless sensor networks. In 2022, Nguyen Phuong Dong et al. [12] studied the uncertainty attack behavior of computer viruses in wireless sensor networks using fuzzy fractional derivatives with non-local Mittag-Leffler function kernels. They proposed a fuzzy fractional SIQR model to describe the dynamics of virus propagation in the network and discussed the concept of fuzzy fractional derivatives in the Caputo sense. They also discussed important properties of the fuzzy fractional SIQR model with fuzzy data. In 2023, Guowen Wu [13] proposed a new virus propagation model called STSIR based on the analysis framework of epidemiology and individual-group game theory. This model considers human behavior and introduces an individual-group game to establish an attack and defense model between infected SIoT (Social Internet of Things) nodes and susceptible SIoT nodes. By using the payoff matrix to obtain mixed Nash equilibrium solutions, the model provides a more reasonable description of virus propagation among devices.

1. (2). Research on computer virus spread control

In the context of biological virus models in the real world, the essential factors for infectious diseases are the infectious source, transmission routes, and susceptible individuals. In the biological realm, susceptible individuals refer to the organisms that can be infected. Similarly, for models of computer virus propagation, the elements of computer virus propagation in a network are infectious computer hosts, transmission routes, and susceptible computers. Therefore, most computer virus control models are studied from these aspects.

In 2015, Ahn et al. [14] proposed the C-SEIRA model for the propagation of malicious worms. They applied optimal control theory to minimize the number of infected hosts and the system processing cost of isolating infected computers from the network as objectives. This approach effectively controlled the spread of malicious worms. In 2018, Abu Arqub et al. [15] proposed a novel approach to simplify many problems encountered in various scientific and engineering fields by discretizing integral-differential algebraic equation systems or their corresponding system sequences. This approach provided a new perspective for establishing computer virus propagation models. In 2020, Upadhyay et al. [16] introduced two different frameworks, namely attack class and target class models, to study the mobility of viruses in target networks. They proposed a new model for controlling DDoS attacks and utilized optimal control theory to regulate virus propagation and reduce infection levels. In 2021, Xingguo Li et al. [17] presented an intermittent isolation-immune strategy for controlling virus propagation in the presence of information diffusion. In this model, information and viruses propagate and interact in different subnetworks. They developed a heterogeneous mean-field approach with time delay to study the dynamics of virus spread, and systematically investigated the propagation dynamics using Monte Carlo simulations. The results demonstrate that the optimal information transmission probability of the virus will be significantly suppressed during relatively short intermittent periods. Conversely, increasing the probability of information transmission can control virus propagation and inhibit the increase of intermittent periods when the periods are very long. Abu Arqub et al. [18] employed the reproducing kernel Hilbert space method with error estimation and convergence analysis to solve the Fredholm-constrained optimal control problem, providing strong theoretical support for controlling the spread of computer viruses. In 2022, Shahini et al. [19] proposed a method for solving the optimal control problem (OCP) of nonlinear time-delay differential equations at Legendre-Gauss-Radau points. They utilized operational matrices and collocation method, combined with a hybrid of Legendre polynomials and block pulse functions to discretize the model, transforming it into a large-scale finite-dimensional nonlinear programming problem. Experimental results indicate that strategies that increase failure rates and retrieval rates while reducing prevalence rates and the number of new computers have a significant impact on controlling and suppressing the spread of computer viruses. Hind Sweis et al. [20] transformed the fractional delay problem into a Volterra integral equation. They established a local existence theorem using the Arzelà-Ascoli theorem and Schauder fixed-point theorem and investigated ρ-order fractional delay differential equations with Atangana-Baleanu fractional derivatives in the Caputo sense. Their research provided strong support for stability studies of computer virus models. In 2023, Ruiling Wang et al. [21] introduced the exposed partitions of computers and USB devices based on the potential characteristics of worms. They proposed a series of computer worm propagation models with saturation occurrence rates and studied the qualitative behavior of the proposed models. Furthermore, they extended the worm propagation models to include three control strategies aimed at minimizing the number of computers and USB devices carrying computer worms. Qingyi Zhu et al. [22] presented a new mathematical model for the propagation of industrial viruses in SCADA systems and established an optimal control system for the model. Through optimal control, they achieved better control of virus propagation with limited resources. Haesung Ahn et al. [23] addressed the significant challenges posed by popular malicious software such as Stuxnet to autonomous vehicles. They developed a mathematical model for the temporal and spatial propagation of Stuxnet-style malware and discussed the future applications of the model in the temporal and geographical spread of infectious malware based on AVs, providing new insights for future research. Qianhui Zhuang et al. [24] proposed a malicious virus propagation model with diffusion and time delay for Cyber-Physical Systems (CPS) and introduced a proportional derivative (PD) feedback control scheme. Experimental results demonstrated that the PD control scheme successfully suppressed the instability of Turing machines and achieved the goal of regulating Hopf bifurcation, providing guidance for predicting and controlling the spread of malicious software. Avcı Derya [25] considered the possibility of failure in recycled computers due to mechanical reasons. They established a fractional-order SEIR computer virus propagation model and obtained numerical solutions for the optimal system using a fractional Euler method combined with forward-backward sweep algorithm. Their approach aimed to minimize the cost of antivirus software installation and eliminate the harm caused by virus propagation in computer networks. Pushpendra et al. [26] derived a mathematical model for the Flash Worm using an improved Caputo-type fractional-order derivative. They analyzed the stability of the model under different values of new computer arrival rate, damage rate, virus transmission rate, and natural clearance rate. Optimal control theory was employed to control the spread of the Flash Worm. Shen Shigen et al. [27] addressed the characteristics of simple edge devices and limited central computer resources in Industrial Internet of Things (IIoT). They proposed a two-layer malicious software propagation patch model called IIPV based on a hybrid patch allocation method. They solved the optimization problem with subjective effort parameters using optimization theory to obtain the optimal control strategy for devices against malicious software and patches. Additionally, they integrated deep reinforcement learning algorithms into the IIPV model and designed a new algorithm called DDQN-PV, which effectively suppressed the propagation of malicious software in IIoT during experimental processes, providing a reference for future research.

Currently, many scholars have conducted extensive research on the dynamic models of computer virus propagation based on the mechanisms of computer virus transmission. However, the majority of these studies assume a bilinear infection rate. In reality, this assumption is not accurate because the increase in susceptible computers and infected computers makes it impossible to achieve a bilinear infection rate. Nonlinear infection rates are more realistic. In light of this, this paper constructs an SEIQRS computer virus propagation model with preventive security controls and conducts research on the following aspects of virus transmission:

1. 1. Considering the influence of virus latency period and temporary immunity period on transmission: This paper investigates the impact of virus latency period and temporary immunity period on the speed and scale of virus transmission. By introducing these factors, the model more accurately reflects the reality of virus propagation.

2. 2. Studying the stability of equilibrium points and the Hopf bifurcation phenomenon: This paper analyzes the local stability of the disease-free equilibrium point and the positive equilibrium point with delays in the model, and discusses the possible occurrence of Hopf bifurcation. These analyses contribute to understanding the dynamic characteristics of virus transmission.

3. 3. Proposing optimal control strategies: This paper formulates an optimal control problem for the model, aiming to determine the optimal security control strategies to effectively mitigate the speed and scale of virus transmission. By applying mathematical methods, the existence of optimal control is determined, and practical and effective transmission control measures are provided.

4. 4. Numerical simulation validation: To validate the effectiveness and accuracy of the model, numerical simulation experiments are conducted in this paper. By simulating virus transmission processes under different scenarios, the dynamic characteristics of the model and the effectiveness of control strategies are verified.

Through research in these areas, this paper aims to gain a better understanding of the propagation patterns of computer viruses and provide scientific basis for the development of prevention and control strategies. This paper is orgnized as follows. In Section 2, a two-delay SEIQRS with nonlinear incidence rates is used to represent the propagation of computer viruses. In Section 3, the local stability of the model’s disease-free equilibrium point and the positive equilibrium point with various time delays are examined, and its Hopf bifurcation is discussed. The optimal control problem is presented in Section 4 by adding control criteria, and a workable and efficient control strategy for computer virus propagation is created. By using numerical simulations, Section 5 verifies the model’s dynamic properties and demonstrates how, after enforcing an appropriate control approach, it is possible to effectively control the spread of computer viruses. The findings are summarized and discussed in Section 6.

2. Model establishment

V. MadhuSudanan et al. [28] developed a computer viral model with a nonlinear infection rate and infection cycle delay, taking into account that the network architecture in virus diffusion induces nonlinear infection.:

(1)

Fangfang Yang et al. [29] introduced the isolation strategy based on the work of V. MadhuSudanan and R. Geetha and established a double-delay SIQRS computer virus model:

(2)

The incubation period refers to the time interval between the invasion of a pathogen into the human body and the onset of clinical symptoms. It is influenced by various factors, including the quantity of the pathogen, its virulence, reproductive capacity, and the body’s resistance. Similarly, computer viruses also possess the significant characteristic of an incubation period. In real-life situations, computer viruses require the execution of specific programs to become activated, and if these programs are not executed, the computer virus can remain dormant for an extended period. Therefore, when studying the propagation patterns of computer viruses, it is necessary to consider the impact of the incubation period on virus transmission. Furthermore, similar to biological viruses, coming into contact with an infected individual does not necessarily result in immediate infection. It takes some time of interaction with the infected individual to enter the incubation period. Hence, there is a time delay effect from the susceptible state to the latent state that needs to be adequately considered. However, existing models overlook these issues and are not sufficiently aligned with the real-life situation.

This work was inspired by the research of V. MadhuSudanan [28] and Fangfang Yang [29] as well as the aforementioned model. According to the warehouse theory [30], in this paper, all state nodes are divided into five compartments: , , , and . Where, is the susceptible node, indicating the number of healthy nodes that do not have immunity and is vulnerable to virus infection at time . is a latent node, which means that after the node is infected with the virus at time , there is no immediate outbreak and the number of nodes that cannot infect. is the infected node, representing the number of nodes infected by the virus and having a specific infection capacity at time . is the isolation node and represents the number of nodes that disconnect from the network and enter the isolation state at time . is the immune node, representing the number of nodes not capable of infection but entering the immune state at time .

The state transition diagram of the computer virus spreading in the network is shown in Fig 1:

[Figure omitted. See PDF.]

The rules for node status transfer are as follows:

1. 1). Susceptible node Latent node : Some susceptible nodes enter the latent state to be activated due to virus infection, and temporarily have no attack capability;

2. 2). Latent node Susceptible node : Some latent nodes are still not activated after the end of the network operation and become susceptible nodes again.

3. 3). Latent node Immune node : Some latent nodes are transformed into immune state due to anti-virus software.

4. 4). Latent node Infected node : Viruses in some latent nodes are activated and become exposed nodes that can infect other network nodes.

5. 5). Infected node Immune node : Some infected nodes are transformed into immune nodes by anti-virus software.

6. 6). Infected node Isolated node : The infected node in virus defence disconnects and enters the isolation state.

7. 7). Isolated node Immune node : The node that enters the isolation state after being infected with the virus has the ability to antiviral attack infection after autonomic immune function, and enters the immune state;

8. 8). Immune node Susceptible node : Due to virus mutation and other reasons, the immune node cannot resist the new virus and is transformed into a susceptible node again.

The model description is shown in system (3), and the parameter description is shown in Table 1:

[Figure omitted. See PDF.]

(3)

The SEIQRS computer virus spreading model studied in this paper follows the dynamic equation of system (3) and assumes that the initial conditions that system (3) needs to satisfy are:  Among them:

Let . According to system (3) we get:

When , there is . Let , is the forward invariant set of system (3).

3. Study on stability of equilibrium point and Hopf bifurcation

Let , , , , then system (3) is simplified to system (4):

(4)

Let , , , , , The disease-free equilibrium point of system (4) can be obtained: . Let , , , , , the positive equilibrium point of system (4) can be obtained:

with. In the propagation dynamics, the primary reproduction number is the critical threshold to determine whether there is a virus in the system. According to the previous method [15], We can get the transfer-in matrix and transfer-out matrix of system (4) as follows:

According to , it is easy to obtain the primary reproduction number of the system (4).

3.1. Stability of disease-free equilibrium point

Theorem 1 When , the disease-free equilibrium of system (4) is locally asymptotically stable.

Proof: Jacobian matrix at the disease-free equilibrium point of system (4) is:

(5)

The characteristic equation of matrix (5) is:

(6)

So we can get , and it can be obtained by calculation: , , , , ; that is, the genuine part of all characteristic roots of equation (6) is negative, so when , the disease-free equilibrium point is locally asymptotically stable, and Theorem 1 is proved.

3.2. Stability of positive equilibrium point

When , the system (4) is linearized as:

(7)

where , , , , , , , , , , , , , , , , ; Then the Jacobian matrix at the positive equilibrium point of system (7) is:

(8)

The characteristic equation of matrix (8) is:

(9)

According to equation (9):

(10)

where

, ,

, ,

, ,

,

, , ,

, ,

, , ,

, ,

, , ,

, ,

Next, we will discuss the stability and Hopf bifurcation of positive equilibrium points in different states.

Case 1 When , equation (10) can be reduced to:

(11)

Theorem 2 According to Routh-Hurwitz criterion, if the assumption condition (R) holds, when , the positive equilibrium point is locally asymptotically stable.

Case 2 When , equation (10) can be reduced to:

(12)

Suppose is a root of equation (12), then:

(13)

Here

[Figure omitted. See PDF.]

It can be obtained from (13):

According to the formula of the sum of squares of the trigonometric function: , we can obtain the following:

(14)

Assuming that is a positive root of equation (14), then:

(15)

For equation (12), taking the derivative of for gives:

(16)

Therefore:

(17)

With

[Figure omitted. See PDF.]

Theorem 3 If , according to Hopf bifurcation theory [31]: When , if is true, then the positive equilibrium point is locally asymptotically stable; when , the positive equilibrium point is unstable, and Hopf bifurcation is generated.

Case 3 When , equation (10) can be reduced to:

(18)

Suppose that is a root of equation (18), then:

(19)

With

[Figure omitted. See PDF.]

It can be obtained from (19):

According to the formula of the sum of squares of the trigonometric function: , we can obtain the following:

(20)

Assuming that is a positive root of equation (20), then:

(21)

For equation (18), taking the derivative of for gives:

(22)

Therefore:

(23)

With

[Figure omitted. See PDF.]

Theorem 4 If , according to Hopf bifurcation theory [31]: When , if is true, then the positive equilibrium point is locally asymptotically stable; when , the positive equilibrium point is unstable, and Hopf bifurcation is generated.

Case 4 When , equation (10) can be reduced to:

(24)

Multiply both sides of the equation by , and we get:

(25)

Suppose that is a root of equation (25), then:

(26)

With

[Figure omitted. See PDF.]

It can be obtained from (26):

According to the formula of the sum of squares of the trigonometric function: , we can obtain the following:

(27)

Assuming that is a positive root of equation (27), then:

(28)

For equation (25), taking the derivative of for gives:

(29)

Therefore:

(30)

With

,

Theorem 5 If , according to Hopf bifurcation theory [31]: When , if is true, then the positive equilibrium point is locally asymptotically stable; when , the positive equilibrium point is unstable, and Hopf bifurcation is generated.

Case 5 When , this situation is similar to case 2, equation (10) can be reduced to:

(31)

Suppose that is a root of equation (31), then:

(32)

With

[Figure omitted. See PDF.]

It can be obtained from (32):

According to the formula of the sum of squares of the trigonometric function: , we can obtain the following:

(33)

Assuming that is a positive root of equation (33), then:

(34)

For equation (31), taking the derivative of for gives:

(35)

Therefore:

(36)

With

[Figure omitted. See PDF.]

Theorem 6 If , according to Hopf bifurcation theory [31]: When , if is true, then the positive equilibrium point is locally asymptotically stable; when , the positive equilibrium point is unstable, and Hopf bifurcation is generated.

Case 6 When , this situation is similar to case 3, equation (10) can be reduced to:

(37)

Suppose that is a root of equation (37), then:

(38)

With

[Figure omitted. See PDF.]

It can be obtained from (38):

According to the formula of the sum of squares of the trigonometric function: , we can obtain the following:

(39)

Assuming that is a positive root of equation (39), then:

(40)

For equation (37), taking the derivative of for gives:

(41)

Therefore:

(42)

With

[Figure omitted. See PDF.]

Theorem 7 If , according to Hopf bifurcation theory [31]: When , if is true, then the positive equilibrium point is locally asymptotically stable; when , the positive equilibrium point is unstable, and Hopf bifurcation is generated.

4. Optimal control problem

Real-world effective control of the spread of computer viruses typically necessitates sizable control efforts, resulting in significant financial outlays. In this section, an optimal control system based on system (3) is established based on Pontryagin’s optimal control theory [32] to control the spread of the virus so that the number of infected nodes reaches the minimum within a predetermined time frame. As shown in the following formula, this paper introduces a time-dependent control variable based on formula (3):

(43)

Here, the variables , and are replaced by the control variable . By adjusting the control variable, such as the effectiveness parameters of antivirus software (including enhancing the detection and removal capabilities of antivirus software, reducing the growth rate of infected nodes, increasing the recovery rate, etc.), the dynamic behavior of the system can be modified. This allows the recovery rate of computers in the latent state e, infected state i, and isolated state q to become time-varying, represented by the control variable .

The goal of implementing optimal control is to minimize the number of infected individuals. However, there is a specific cost associated with putting any control strategy into action. As a result, the objective function is defined as follows when taking into account the treatment cost and control cost:

(44)

Here, is the weight coefficient of the infected node; is the weight coefficient of the control strategy; is the cost of the corresponding control strategy. The optimal control problem is to find the optimal control so that the objective function take the minimum value, that is, , where is the control set.

Theorem 8 For a state system with a known initial value (43) and a given objective function (44), there must exist the following optimal control , making .

Proof: According to the optimal existence theory [33], it is only necessary to prove:

1. a). The set of control variables and corresponding state variables is non-empty;

2. b). The control set is a convex closed set;

3. c). The right side of the state system is bounded for state variables and control variables;

4. d). The integrand in the objective function is a convex function for the permissible control;

5. e). There are constants , , and so that the integrand in the objective function has an upper bound .

Firstly, conditions (a) and (c) are satisfied according to the previous dynamics discussion of system (3). By definition, the control set is non-empty, and the integrand in the objective function is a convex function for the control set, so the conditions (b) and (d) are satisfied. In addition, because , , so take , and , condition (e) is established. In summary, Theorem 8 is proven.

According to the Pontryagin minimum principle [32], the necessary conditions for optimal control can be obtained as follows.

The Lagrange function is defined as: .

Define the Hamiltonian function as:

Theorem 9 For the given optimal control , let the corresponding optimal solution be , then there is a co-state variable that satisfies:

1. a). Costate equation:

1. b). Cross section conditions: .

2. c). According to the Pontryagin minimum principle, by solving the partial derivative of the Hamiltonian equation with respect to the control variable, the control point whose derivative is equal to zero can be obtained, and the optimal control solution can be obtained. That means . And the optimality condition can be obtained:

This way, a win-win optimal control can be achieved so that the number of infected computers is the least and the cost is the least.

5. Numerical simulation

To test the above conclusions, this section revolves around the basic reproduction number and conducts numerical simulations based on different parameters. The experimental environment used in this study is MATLAB R2021. All experiments and computational tasks were performed on a computer equipped with an Intel (R) Core (TM) i7-7500U CPU and 8GB RAM.

5.1. Stability of disease-free equilibrium point

Select parameters: , , , , , , , , , , , , , , , . Through MATLAB calculation, can be obtained. The disease-free equilibrium exists in system (4). And the disease-free equilibrium point is locally asymptotically stable. The simulation results are shown in Fig 2.

[Figure omitted. See PDF.]

w_. (a) The change of each node with respect to time when (b) SQR phase diagram for (c) EIQ phase diagram for (d) Number of hosts infected with different models when _.

The phase diagram can provide insights into the dynamics of virus spread, such as the existence of stable equilibrium points, limit cycles, or other patterns of behavior. It allows for visualizing the transitions between different compartments and understanding the long-term behavior of the system. The phase diagram takes different state variables as axes and shows how these variables change over time. By drawing the phase diagram of SQR model, we can observe the transition, equilibrium point and stability characteristics between different states, so as to better understand the dynamic behavior of computer virus propagation.

Fig 2a draws the waveform diagram of each node changing with time when , and Fig 2b and 2c illustrate the SQR phase diagram and EIQ phase diagram at this time. As can be seen from the above Fig: when , no matter how the initial value of hosts in each state changes, over time, the susceptible node and the immune node will eventually tend to be stable, while the latent node , the infected node , and the isolated node will ultimately converge to 0. Thus, when , the disease-free equilibrium point of system (4) tends to be stable. To better verify the effectiveness of model (4), keep the node state transition parameters of V. MadhuSudanan’s SIRS model [28] and Fangfang Yang’s SIQRS model [29] consistent with the SEIQRS model, let , , , , the changes in the number of infected nodes in each model are shown in Fig 2d. Through comparison, it is found that the result is consistent with Theorem 1, that is, when , the system is stable at the disease-free equilibrium point, the virus is eliminated in the computer network, and the infected host will eventually return to normal under anti-virus software. It is also apparent from Fig 2d that the infected node in the SEIQRS model tends to zero the fastest, which also means that the virus can be eliminated faster in the SEIQRS model.

5.2. Stability of positive equilibrium

Select parameters: , , , , , , , , , , , , , . Through MATLAB calculation, can be obtained. System (4) has a unique positive equilibrium point . Figs 4–9, respectively, describe the number of hosts of SEIQRS model with different time delays when . In addition, in order to better verify the validity of model (4), keep the node state transition parameters of V. MadhuSudanan’s SIRS model [28] and Fangfang Yang’s SIQRS model [29] consistent with the SEIQRS model, let , , , the changes in the number of infected nodes in each model are shown in Fig 3. Through comparison, it is found that when , the system is stable at the positive equilibrium point. It can also be clearly seen from Fig 3 that when the system finally stabilizes at the positive equilibrium point, the number of infected nodes in the SEIQRS model is the least, which also indicates that the computer virus is well controlled in the SEIQRS model.

[Figure omitted. See PDF.]

_.

[Figure omitted. See PDF.]

, the positive equilibrium is locally stable. (a) with respect to time (b) SQR phase diagram.

[Figure omitted. See PDF.]

, the positive equilibrium point becomes unstable. (a) with respect to time (b) SQR phase diagram.

[Figure omitted. See PDF.]

, the positive equilibrium is locally stable. (a) with respect to time (b) SQR phase diagram.

[Figure omitted. See PDF.]

, the positive equilibrium point becomes unstable. (a) with respect to time (b) SQR phase diagram.

[Figure omitted. See PDF.]

, the positive equilibrium is locally stable. (a) with respect to time (b) SQR phase diagram.

[Figure omitted. See PDF.]

, the positive equilibrium point becomes unstable. (a) with respect to time (b) SQR phase diagram.

1. (1). Fig 4 shows that, when , the system (4) is locally asymptotically stable at the positive equilibrium point; Fig 5 shows that when , the system is out of control and unstable and Hopf bifurcation occurs. When , case 5 is similar to case 2. When , the system is locally asymptotically stable. When , the system is out of control and produces Hopf bifurcation. The experimental simulation diagrams of those situations are shown in Figs 4 and 5.

1. (2). Fig 6 shows that, when , the system (4) is locally asymptotically stable at the positive equilibrium point; Fig 7 shows that when , the system is out of control and unstable and hopf bifurcation occurs. When , case 6 is similar to case 3. When , the system is locally asymptotically stable. When , the system is out of control and produces hopf bifurcation. The experimental simulation diagrams of those situations are shown in Figs 6 and 7.

1. (3). Fig 8 shows that, when , the system (4) is locally asymptotically stable at the positive equilibrium point; Fig 9 shows that when , the system is out of control and unstable and Hopf bifurcation occurs.

5.3. Systems with optimal control

This section performs a numerical simulation using the fourth-order Runge-Kutta method and the forward and backward scanning method to assess the efficacy of the suggested optimal control strategy. An appropriate control strategy can successfully stop the spread of computer viruses, as shown by the experimental simulation. Let , , , , , , , , , , , , , , , , , , , , , . The optimal control function diagram is obtained, as shown in Fig 10a. And Fig 10b shows the value of the target function at different . Fig 10c shows the changes of infected nodes over time under different control degrees. It can shows the change of infected nodes over time under different control levels. It can be seen from the Fig that when , that is, without any control, the spread of computer viruses is uncontrollable. Although the number of infected nodes seems to be lower than the optimal control when , , more importantly, the optimal control strategy can minimize the value of the objective function and achieve the balance between control objective and control cost. Table 2 shows the values of the objective function of V. MadhuSudanan’s SIRS model [28], Fangfang Yang’s SIQRS model [29] and SEIQRS model relative to the five constant control and optimal controls. Fig 10d reflects the change in the number of infected nodes in different models over time under optimal control. According to these data, we can see that the value of an objective function under optimal control is the smallest, and it can be seen that the time-varying optimal control strategy is more effective in controlling computer infection and virus spread. It can be seen from Table 2 and Fig 10d that after adopting the optimal control scheme, the three systems of SIRS [28], SIQRS [29] and SEIQRS can all achieve good control effect with minimal consumption.

[Figure omitted. See PDF.]

[Figure omitted. See PDF.]

(a) Control variable as a function of time (b) The value of the target function at different (c) Change of infected node with time under different (d) Changes in the number of infected hosts of different models under optimal control.

At the same time, to verify the reliability of the proposed SEIQRS model and the effectiveness of optimal control, we conduct corresponding experiments on a real-time network. In the simulation experiments, we randomly set the initial positions and velocities of the nodes. We determine the initial infection rate by randomly selecting a subset of media and infecting them completely, while the remaining nodes are in a completely healthy state.

In this section, we create a network system with 100 nodes, where each node is connected to any 4 nodes in the network graph. We randomly select five nodes from the system to be infected with the virus. Each node in the network structure is assigned a certain healing rate and infection rate, with the healing rate set to 1 and the infection rate set to 0.1. We set the maximum number of time steps to 1000, and within a certain number of time steps, all nodes in the system eventually get healed. This is shown in Fig 11. We use Matlab to conduct simulation experiments on the computer virus model, demonstrating the dynamic process of node infection in the computer virus model. Since the experimental results are a dynamic changing process, in this section, we capture four stages from the dynamic change graph to demonstrate the process of node infection and healing. Fig 11 shows the four stages of node infection and healing in the dynamic network system. (a) represents the initial stage, (b)(c) represent the intermediate stages, and (d) represents the final stage.

[Figure omitted. See PDF.]

(a) The initial state (b) The intermediate stages (c) The intermediate stages (d) The final stage.

In Fig 11, we show the changes in susceptible and infected nodes for each stage. Susceptible nodes are shown in gray, latent nodes in yellow, infected nodes in red, isolated nodes in blue, and recovered nodes in green. In Fig 11, (a) shows the initial state of the simulation experiment. Four nodes are randomly selected from the network structure as already infected nodes, shown in red. The remaining uninfected nodes are shown in green. (b)(c) show the intermediate stages of the simulation experiment, illustrating the process of infection for uninfected nodes and the process of healing for already infected nodes. (d) shows the final stage of the simulation experiment, where the previously infected nodes have been healed and all nodes are shown in gray.

Through the simulation experiment with 100 nodes shown in Fig 11, we demonstrate the dynamic process of node infection and healing in the computer virus model, providing a theoretical visualization basis for computer virus control.

As the states of nodes in the dynamic network change, the number of infected and susceptible nodes in the network system also changes over time. Fig 12 shows the line graphs of the changing quantities of susceptible and infected nodes over time in the dynamic graph of the computer virus model simulation experiment. The green line represents the change in the quantity of susceptible nodes, while the red line represents the change in the quantity of infected nodes. The changing curves in Fig 12 illustrate the process of the changing quantities of susceptible and infected nodes in the dynamic graph shown in Fig 11.

[Figure omitted. See PDF.]

(a) The initial state (b) The intermediate stages (c) The intermediate stages (d) The final stage.

In particular, (a)(b)(c)(d) correspond to the node quantity change curves extracted from the four images in Fig 11. (a) represents the change in the quantity of infected and susceptible nodes in the initial state. (b)(c) represent the changing quantities of susceptible and infected nodes in the intermediate states of the system. (d) represents the changing quantities of susceptible and infected nodes in the final state of the system. From the changes in the red curve, it can be observed that under the initial conditions, when five nodes in the computer virus simulation system are infected and infect nearby nodes at a certain infection rate, the quantity of infected nodes increases while the quantity of susceptible nodes decreases. However, under certain healing rate conditions, when all the infected nodes are eventually healed, the system enters a stable state.

6. Conclusion

This paper presents a class of computer virus transmission models with nonlinear infection rates and double delay SEIQRS. The primary regeneration number is obtained through an in-depth analysis of the model, and the local stability of disease-free equilibrium is proved. The regional stability of positive balance under different time delays and Hopf bifurcation are discussed. And on this basis, the optimal control problem of this model is studied.

These results show that when the primary regeneration number is less than 1, the system will eventually stabilize at the disease-free equilibrium point. As time goes on, the infected node converges to 0, and the virus will eventually be eliminated. The system will eventually stabilize at the positive equilibrium point when the primary regeneration number is greater than 1 and the time delay is less than the critical value. At this time, the installation of anti-virus software and other measures can effectively inhibit the spread of computer viruses. Once the time delay exceeds the critical value, Hopf bifurcation appears in the system, and the spread of the virus will become uncontrollable. Under the optimal control condition, the trend of the number of infected computers is in line with the expectation of the objective function. The time-varying optimal control strategy can control the infection of computers and the spread of the virus more effectively. In the absence of control, the spread of computer viruses is uncontrollable. Meanwhile, applying different controls in SIRS [28], SIQRS [29], and SEIQRS proves that optimal control achieves the balance between control objectives and control costs. And at the same time, we carry out corresponding experiments on the real-time network to verify the reliability of the proposed SEIQRS model and the effectiveness of the optimal control.

References

1. 1. China Internet Network Information Center. The 52nd Statistical Report on China’s Internet Development. 2023.

2. 2. Cohen F. Computer viruses Theory and experiments. Computers & Security. 1987, 6(1): 22–35.

3. 3. Kephart JO, White SR, Chess DM. Computers and epidemiology. IEEE spectrum. 1993;30(5):20–6.

* View Article

* Google Scholar

4. 4. Gan C, Feng Q, Zhu Q, et al. Analysis of computer virus propagation behaviors over complex networks: a case study of Oregon routing network. China Engineering Res. 2020;100(2):1725–40.

* View Article

* Google Scholar

5. 5. Ronghua Shang WZ, Licheng J, Xiangrong Z, Rustam S. Dynamic Immunization Node Model for Complex Networks Based on Community Structure and Threshold. Joint International Research Laboratory of Intelligent. 2020;52(3):1–14.

* View Article

* Google Scholar

6. 6. Zhang N, Yang Q, Zhu X. The impact of social resource allocation on epidemic transmission in complex networks. Applied Mathematics and Computation. 2022;433:127405.

* View Article

* Google Scholar

7. 7. Liu J. Periodic solutions of a delayed model for propagation of malicious codes in WSNS with graded infection rate. Results in Physics. 2022;38:105631.

* View Article

* Google Scholar

8. 8. Zarin R, Ullah N, Khan A, et al. A numerical study of a new non-linear fractal fractional mathematical model of Malicious Codes Propagation in Wireless Sensor Network. Computers & Security. 2023; 135:103484.

* View Article

* Google Scholar

9. 9. Li K, Chen Z, Cong R, et al. Simulated dynamics of virus spreading on social networks with various topologies. Applied Mathematics and Computation. 2024;470:128580.

* View Article

* Google Scholar

10. 10. Dang QA, Hoang MT. Positivity and global stability preserving NSFD schemes for a mixing propagation model of computer viruses. Journal of Computational and Applied Mathematics. 2020;374:112753.

* View Article

* Google Scholar

11. 11. Awasthi S, Kumar N, Srivastava PK. An epidemic model to analyze the dynamics of malware propagation in rechargeable wireless sensor network. Journal of Discrete Mathematical Sciences and Cryptography. 2021;24(5):1529–43.

* View Article

* Google Scholar

12. 12. Dong NP, Long HV, Giang NL. The fuzzy fractional SIQR model of computer virus propagation in wireless sensor network using Caputo Atangana–Baleanu derivatives. Fuzzy Sets and Systems. 2022;429:28–59.

13. 13. Wu G. STSIR: An individual-group game-based model for disclosing virus spread in social internet of things. Journal of Network and Computer Applications. 2023;214:103608.

* View Article

* Google Scholar

14. 14. Ahn I, Oh H, Park J. Investigation of the C-SEIRA model for controlling malicious code infection in computer networks. Applied Mathematical Modelling. 2015;39(14):4121–33.

* View Article

* Google Scholar

15. 15. Abu Arqub O, Rashaideh H. The RKHS method for numerical treatment for integrodifferential algebraic systems of temporal two-point BVPs. Neural Computing and Applications. 2018;30(8):2595–606.

* View Article

* Google Scholar

16. 16. Upadhyay R, Singh P. Modeling and control of computer virus attack on a targeted network. Physica A: Statistical Mechanics and its Applications. 2020;538:122617.

* View Article

* Google Scholar

17. 17. Li X, Luo X, Wang Y. Controlling Virus Spread Using an Intermittent Quarantine Strategy on Multiplex Networks. Security and Communication Networks. 2021;2021:1–9.

* View Article

* Google Scholar

18. 18. Abu Arqub O, Shawagfeh N. Solving optimal control problems of Fredholm constraint optimality via the reproducing kernel Hilbert space method with error estimates and convergence analysis. Mathematical Methods in the Applied Sciences. 2021;44(10):7915–32.

* View Article

* Google Scholar

19. 19. Shahini M, Ebrahimzadeh A, Khanduzi R. A Spectral Collocation Method for Computer Virus Spread Case of Delayed Optimal Control Problem. Bulletin of the Iranian Mathematical Society. 2022;48(2):507–35.

* View Article

* Google Scholar

20. 20. Sweis H, Shawagfeh N, Arqub OA. Fractional crossover delay differential equations of Mittag-Leffler kernel: Existence, uniqueness, and numerical solutions using the Galerkin algorithm based on shifted Legendre polynomials. Results in Physics. 2022;41:105891.

* View Article

* Google Scholar

21. 21. Wang R, Xue Y, Xue K. Dynamic analysis and optimal control of worm propagation model with saturated incidence rate. Frontiers in Physics. 2023;10:1098040.

* View Article

* Google Scholar

22. 22. Zhu Q, Zhang G, Luo X, Gan C. An industrial virus propagation model based on SCADA system. Information Sciences. 2023;630:546–66.

* View Article

* Google Scholar

23. 23. Ahn H, Choi J, Kim YH. A Mathematical Modeling of Stuxnet-Style Autonomous Vehicle Malware. IEEE Transactions on Intelligent Transportation Systems. 2023;24(1):673–83.

* View Article

* Google Scholar

24. 24. Zhuang Q, Xiao M, Ding J, et al. Spatiotemporal Evolution Control of Malicious Virus Propagation in Cyber Physical Systems Via PD Feedback Control. IEEE Transactions on Control of Network Systems. 2023:1–12.

* View Article

* Google Scholar

25. 25. Avcı D, Soytürk F. Optimal control strategies for a computer network under virus threat. Journal of Computational and Applied Mathematics. 2023;419:114740.

* View Article

* Google Scholar

26. 26. Kumar P. Fractional mathematical modeling of the Stuxnet virus along with an optimal control problem. Ain Shams Engineering Journal. 2023;14(7):102004.

* View Article

* Google Scholar

27. 27. Shen S, Xie L, Zhang Y, et al. Joint Differential Game and Double Deep Q-Networks for Suppressing Malware Spread in Industrial Internet of Things. IEEE Transactions on Information Forensics and Security. 2023;18:5302–15.

* View Article

* Google Scholar

28. 28. Madhusudanan V, Geetha R. Dynamics of Epidemic Computer Virus Spreading Model with Delays. Wireless Personal Communications. 2020;115(3):2047–61.

* View Article

* Google Scholar

29. 29. Yang F, Zhang Z. Dynamics of a nonlinear SIQRS computer virus spreading model with two delays. AIMS Mathematics. 2021;6(4):4083–104.

* View Article

* Google Scholar

30. 30. Mishra BK, Saini D. Mathematical models on computer viruses. Applied Mathematics and Computation. 2007;187(2):929–36.

* View Article

* Google Scholar

31. 31. Hopf E. Abzweigung einer periodischen Lösung von einer stationären eines Differentialsystems. 1943;95(1):3–22.

32. 32. Pontryagin LS, Boltyanskii VG, Gamkrelidze RV, et al. The mathematical theory of optimal processes. Interscience Publishers: John Wiley & Sons, Inc., New York, London; 1962.

33. 33. Kennedy J, Eberhart R. Particle Swarm Optimization; Proceedings of the ICNN95 - International Conference on Neural Networks, F; 1995.

Citation: Wang J, Zhong L, Chang X (2025) SEIQRS model analysis and optimal control with two delays. PLoS One 20(6): e0319417. https://doi.org/10.1371/journal.pone.0319417

About the Authors:

Junling Wang

Contributed equally to this work with: Junling Wang

Roles: Writing – review & editing

Affiliation: School of Information Engineering, Jiangxi University of Science and Technology, Ganzhou, Jiangxi, China

Lei Zhong

Roles: Writing – original draft

E-mail: [email protected]

¶‡ LZ and XC also contributed equally to this work.

Affiliation: School of Information Engineering, Jiangxi University of Science and Technology, Ganzhou, Jiangxi, China

ORICD: https://orcid.org/0000-0002-2883-7430

Xinxin Chang

Roles: Data curation, Visualization

¶‡ LZ and XC also contributed equally to this work.

Affiliation: School of Information Engineering, Jiangxi University of Science and Technology, Ganzhou, Jiangxi, China

[/RAW_REF_TEXT]

[/RAW_REF_TEXT]

[/RAW_REF_TEXT]

[/RAW_REF_TEXT]

[/RAW_REF_TEXT]

[/RAW_REF_TEXT]

References

1. China Internet Network Information Center. The 52nd Statistical Report on China’s Internet Development. 2023.

2. Cohen F. Computer viruses Theory and experiments. Computers & Security. 1987, 6(1): 22–35.

3. Kephart JO, White SR, Chess DM. Computers and epidemiology. IEEE spectrum. 1993;30(5):20–6.

4. Gan C, Feng Q, Zhu Q, et al. Analysis of computer virus propagation behaviors over complex networks: a case study of Oregon routing network. China Engineering Res. 2020;100(2):1725–40.

5. Ronghua Shang WZ, Licheng J, Xiangrong Z, Rustam S. Dynamic Immunization Node Model for Complex Networks Based on Community Structure and Threshold. Joint International Research Laboratory of Intelligent. 2020;52(3):1–14.

6. Zhang N, Yang Q, Zhu X. The impact of social resource allocation on epidemic transmission in complex networks. Applied Mathematics and Computation. 2022;433:127405.

7. Liu J. Periodic solutions of a delayed model for propagation of malicious codes in WSNS with graded infection rate. Results in Physics. 2022;38:105631.

8. Zarin R, Ullah N, Khan A, et al. A numerical study of a new non-linear fractal fractional mathematical model of Malicious Codes Propagation in Wireless Sensor Network. Computers & Security. 2023; 135:103484.

9. Li K, Chen Z, Cong R, et al. Simulated dynamics of virus spreading on social networks with various topologies. Applied Mathematics and Computation. 2024;470:128580.

10. Dang QA, Hoang MT. Positivity and global stability preserving NSFD schemes for a mixing propagation model of computer viruses. Journal of Computational and Applied Mathematics. 2020;374:112753.

11. Awasthi S, Kumar N, Srivastava PK. An epidemic model to analyze the dynamics of malware propagation in rechargeable wireless sensor network. Journal of Discrete Mathematical Sciences and Cryptography. 2021;24(5):1529–43.

12. Dong NP, Long HV, Giang NL. The fuzzy fractional SIQR model of computer virus propagation in wireless sensor network using Caputo Atangana–Baleanu derivatives. Fuzzy Sets and Systems. 2022;429:28–59.

13. Wu G. STSIR: An individual-group game-based model for disclosing virus spread in social internet of things. Journal of Network and Computer Applications. 2023;214:103608.

14. Ahn I, Oh H, Park J. Investigation of the C-SEIRA model for controlling malicious code infection in computer networks. Applied Mathematical Modelling. 2015;39(14):4121–33.

15. Abu Arqub O, Rashaideh H. The RKHS method for numerical treatment for integrodifferential algebraic systems of temporal two-point BVPs. Neural Computing and Applications. 2018;30(8):2595–606.

16. Upadhyay R, Singh P. Modeling and control of computer virus attack on a targeted network. Physica A: Statistical Mechanics and its Applications. 2020;538:122617.

17. Li X, Luo X, Wang Y. Controlling Virus Spread Using an Intermittent Quarantine Strategy on Multiplex Networks. Security and Communication Networks. 2021;2021:1–9.

18. Abu Arqub O, Shawagfeh N. Solving optimal control problems of Fredholm constraint optimality via the reproducing kernel Hilbert space method with error estimates and convergence analysis. Mathematical Methods in the Applied Sciences. 2021;44(10):7915–32.

19. Shahini M, Ebrahimzadeh A, Khanduzi R. A Spectral Collocation Method for Computer Virus Spread Case of Delayed Optimal Control Problem. Bulletin of the Iranian Mathematical Society. 2022;48(2):507–35.

20. Sweis H, Shawagfeh N, Arqub OA. Fractional crossover delay differential equations of Mittag-Leffler kernel: Existence, uniqueness, and numerical solutions using the Galerkin algorithm based on shifted Legendre polynomials. Results in Physics. 2022;41:105891.

21. Wang R, Xue Y, Xue K. Dynamic analysis and optimal control of worm propagation model with saturated incidence rate. Frontiers in Physics. 2023;10:1098040.

22. Zhu Q, Zhang G, Luo X, Gan C. An industrial virus propagation model based on SCADA system. Information Sciences. 2023;630:546–66.

23. Ahn H, Choi J, Kim YH. A Mathematical Modeling of Stuxnet-Style Autonomous Vehicle Malware. IEEE Transactions on Intelligent Transportation Systems. 2023;24(1):673–83.

24. Zhuang Q, Xiao M, Ding J, et al. Spatiotemporal Evolution Control of Malicious Virus Propagation in Cyber Physical Systems Via PD Feedback Control. IEEE Transactions on Control of Network Systems. 2023:1–12.

25. Avcı D, Soytürk F. Optimal control strategies for a computer network under virus threat. Journal of Computational and Applied Mathematics. 2023;419:114740.

26. Kumar P. Fractional mathematical modeling of the Stuxnet virus along with an optimal control problem. Ain Shams Engineering Journal. 2023;14(7):102004.

27. Shen S, Xie L, Zhang Y, et al. Joint Differential Game and Double Deep Q-Networks for Suppressing Malware Spread in Industrial Internet of Things. IEEE Transactions on Information Forensics and Security. 2023;18:5302–15.

28. Madhusudanan V, Geetha R. Dynamics of Epidemic Computer Virus Spreading Model with Delays. Wireless Personal Communications. 2020;115(3):2047–61.

29. Yang F, Zhang Z. Dynamics of a nonlinear SIQRS computer virus spreading model with two delays. AIMS Mathematics. 2021;6(4):4083–104.

30. Mishra BK, Saini D. Mathematical models on computer viruses. Applied Mathematics and Computation. 2007;187(2):929–36.

31. Hopf E. Abzweigung einer periodischen Lösung von einer stationären eines Differentialsystems. 1943;95(1):3–22.

32. Pontryagin LS, Boltyanskii VG, Gamkrelidze RV, et al. The mathematical theory of optimal processes. Interscience Publishers: John Wiley & Sons, Inc., New York, London; 1962.

33. Kennedy J, Eberhart R. Particle Swarm Optimization; Proceedings of the ICNN95 - International Conference on Neural Networks, F; 1995.