{'fraction_non_alphanumeric': 0.043251276579529, 'fraction_numerical': 0.020272845057541344, 'mean_word_length': 4.3437181836693135, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
On the fractional SIRD mathematical model and control for the transmission of COVID-19: The first and the second waves of the disease in Iran and Japan

In this paper, a fractional-order SIRD mathematical model is presented with Caputo derivative for the transmission of COVID-19 between humans. We calculate the steady-states of the system and discuss their stability. We also discuss the existence and uniqueness of a non-negative solution for the system under study. Additionally, we obtain an approximate response by implementing the fractional Euler method. Next, we investigate the first and the second waves of the disease in Iran and Japan; then we give a prediction concerning the second wave of the disease. We display the numerical simulations for different derivative orders in order to evaluate the efficacy of the fractional concept on the system behaviors. We also calculate the optimal control of the system and display its numerical simulations.


Introduction
Coronaviruses, viruses that have a crown, have caused a variety of severe illnesses like Severe Acute Respiratory Syndrome (SARS) and Middle East Respiratory Syndrome (MERS). The disease COVID-19, which was caused by the virus SARS COV-2, was first identified in Wuhan City of China on December 2019. The COVID-19 disease appears in different people with different symptoms. In most people, fever, cough, and the shortness of breath have been the early symptoms; in some others, gastrointestinal symptoms such as anorexia, diarrhoea, and vomiting have been appeared. Additionally, other symptoms such as body aches, the loss of smell and taste senses, and some skin symptoms have been reported. The incidence of COVID-19 in children is lower than adults, and for people with some diseases such as diabetes, high blood pressure, asthma, and autoimmune diseases, it can lead to more serious conditions and sometimes death. Because of the new outbreak of COVID-19, there is no complete information on how it is transmitted, but according to the World Health Organization (WHO), coronavirus is spread through the sneezing and coughing of infected people, and it stays in the air for several hours and enters the body by breathing the polluted air. Also, it sits on the surface and is transmitted to a body by touching the infected environment and then the nose, eyes, or mouth. Due to these transmission kinds of this virus, the suggested ways of preventing the transmission include washing hands with soap, using masks outside the home, not presenting in the polluted environments, and reducing unnecessary presence in the community.
Many mathematicians and biologists have studied the dynamics of realistic systems by using the concept of mathematical modeling (see, for example, [1,2]). By extending the integerorder derivatives to the fractional-order ones, many researchers have studied the mathematical models of natural phenomena with fractional-order models [3,4]. Recently, the utilization of fractional operators has been expanded to model real-world systems, an efficient approach which investigates the transmission of diseases as well as their control [5,6]. With the expansion of COVID-19, researchers have also used mathematical models to simulate the spread of this disease in different environments (see, for example, [7][8][9][10]).
In this research, we study the transmission of COVID-19 by introducing a fractional SIRD model in which the Caputo derivative is used for the system formulation. With the onset of this disease, most countries made a series of measures to control its transmission, which led to a decline in the spread of the disease after the first peak. Iran was the first country that experienced the second peak of COVID-19. We provide a comparative study for the first and the second waves in Iran and the prediction of the disease transmission in the future. Japan was another country that experienced the second wave, while the mortality rate with COVID was low in that country. Thus, we examine the second wave of the disease in Japan and provide a forecast for it.
In the following, Section 2 provides some basic definitions with regard to the fractional calculus. The system of equations for the COVID-19 spread as well as the existence and uniqueness of the solution are presented in Section 3. The equilibrium points and the optimal control of the system are presented in Sections 4 and 5 , respectively. Numerical results and the simulations of the first and the second waves as well as the prediction of the second wave in Iran and Japan are discussed in Section 6. Finally, the paper is closed by some conclusions in Section 7.

Fractional calculus
In this section, we first recall some definitions of fractional operators; then we examine the important feature of fractional calculus called memory effect.

Basic definitions
in which υ in the fractional order. Also, for υ with Re(υ) > 0, the fractional integral of order υ is defined as Definition 2.2 ([11]). For the fractional-order Caputo operator of order υ, the Laplace transform is which can be written more simply as

Memory effect
Fractional calculus is one of the most powerful tools in the modeling of real-world phenomena with memory effect. The effect of memory in fractional derivatives has been studied by many researchers (see, for example, [13][14][15]). In the mathematical models with classical differential operators, the responses do not depend on the previous instance, so these models do not have memory. However, the fractional mathematical models include memory as their inherent property. Consider an integrable function w(t) in [0, T ] and time instances t 1 , t 2 ∈ [0, T ] such that 0 < t 1 < t 2 < T . If υ ̸ = 1, we have If υ = 1, then the first integral is canceled, and we obtain The above equations show that in the classical case υ = 1, the value of P depends only on [t 1 , t 2 ], but in the fractional case υ ̸ = 1, the value of P depends on [t 1 , t 2 ] and the previous what happened before t 1 ; this shows the effect of memory on the fractional process, while the integral in the classical case has nothing with the range before t 1 .
The nature of coronavirus transmission shows that the patient has no specific symptoms in the first few days of infection, but he/she is the carrier of the disease. Therefore, the population of infected people now affect the number infected individuals in the next few days. Hence, the population of infected people at time t is affected by the number of infected persons in the previous days, and the chain continues so that the number of infected people at time t is influenced by the number of infected people from the time 0 to t, indicating that the coronavirus transmission is a system with memory; therefore, it is more appropriate to model this phenomenon by fractional-order derivatives.

Model formulation
In this section, we provide a Caputo-type fractional SIRD model for the transmission of COVID-19. Generally, the SIRD models consider the population to be constant, but here we take into account the natural birth and death rates. We also consider that the total N(t) includes four sub-categories: susceptible persons S(t), infected people I(t), recovered individuals R(t), and died people with this disease D(t) so that N(t) = S(t)+I(t)+R(t)+D(t).
The diagram of the model is shown in Fig. 1. Also, the following compartmental model describes the spread of COVID-19 where (S(0), I(0), R(0), D(0)) = (S 0 , I 0 , R 0 , D 0 ) is a nonnegative initial state vector. The parameters of the model are also the per capita birth rate Ω, the recovery rate γ , the transmission rate β, the virus-induced average fatality rate µ, and the per capita natural death rate m.
Fractional-order systems are related to nonlocal effects, history, and memory, the features which exist in different biological systems. The integer-order system (1) suffers from the lack of internal memory effects regarding the COVID-19 transmission, so we are going to modify the aforesaid ordinary model to a fractional-order one. To this end, we replace the ordinary derivatives with the Caputo-type fractional operators. Also, an auxiliary parameter θ is used to solve the dimensional mismatching problem [16]. Thus, the fractional-order model for the spread of COVID-19 is given as follows where t ∈ [0, T ] and 0 < υ < 1.

Feasibility region
then we state that B is positively invariant by the next lemma. Proof. To prove this lemma, first we show that the solutions are non-negative by using the method presented in [17]. We have which are all nonnegative. Consider the auxiliary fractional differential system is the initial state vector. We will show that for all t ≥ 0, the solutions of the auxiliary system (S * e (t), I * e (t), R * e (t), D * e (t)) are nonnegative. If not, there exists a time instance t at which the solutions are negative. Consider Similarly, it can be shown that other components I * e , R * e , D * e are also nonnegative, which is a contradiction. Therefore, we obtain that (S * , for all t ≥ 0 as e → ∞; thus, the solutions are nonnegative. In the following, we first add the first three equations of the system (2). So Applying the Laplace transform results in where S(0), I(0), R(0) are the initial size of susceptible, infected and recovered population, respectively. With some calculations, we get Consequently, with respect to the system (2), the set B is positively invariant. □

Existence and uniqueness
To prove that the existence of a unique solution for the system (2), we first consider (2) in the form The integral operator is applied to the both sides, which results In the next theorem, we prove the Lipschitz and the contraction conditions for the kernels N 1 -N 4 .
Proof. According to the values of N 1 at S(t) and S 1 (t), one can write Therefore, the Lipschitz condition holds for N 1 , and if 0 ≤ βk 2 + m < 1, then N 1 is a contraction. □ Similarly, we can show that N 2 , N 3 , N 4 satisfy the Lipschitz It should be noted that if the system (2) has an equilibrium point at the origin, then one can use the method proposed in [18] to prove the existence and uniqueness analysis. Now suppose the following recursive formulas with respect to the system (2) . We take the norm of Φ 1,n as By the Lipschitz condition (4), we get Similarly, we obtain Then we can obtain In the next theorem, we prove the existence of the solution for the system (2).

Theorem 3.2. The fractional-order transmission model of COVID-19
(2) has a solution if there exists a time instant t e > 0 such that Proof. By considering a recursive technique and using Eqs. (5) and (7), we have This shows that the system has a continuous solution. We claim that a solution of the system (2) is constructed by the above functions. For this purpose, consider We repeat the method; then At t e , we get Then lim n→∞ ∥B 1,n (t)∥ → 0. Similarly, we have lim n→∞ ∥B j,n (t)∥ → 0, j = 2, 3, 4. □ In order to prove that the solution of the system (2) is unique, assume that another solutions S 1 (t), I 1 (t), R 1 (t), and D 1 (t) exist for the system (2). We can write Taking the norm, we get By the Lipschitz condition (4), we obtain Theorem 3.3. If the following inequality holds, then the solution of the spread of COVID-19 model (2) is unique Proof. From the above condition, we obtain With this inequality and (7), we conclude that ∥S(t) − S 1 (t)∥ = 0, so we obtain S(t) = S 1 (t). Following the same manner, we can derive

Reproduction number and equilibrium points
The We solve the above equations and obtain the equilibrium point E * = (S * , I * , R * ) such that Here we use the next generation method to find the basic reproduction number R 0 [19]. First, we consider the first three equations of the system (2) in a compact form as follows The values of the Jacobin matrices of F and V at the equilibrium point E 0 are Then the basic reproduction number is computed from R 0 = ρ(FV −1 ) where FV −1 is the next generation matrix. By performing simple algebraic calculations, the value of R 0 is obtained as follows The reproduction number is epidemiologically important and indicates the ability to transmit the infection and the continuation of the disease.

Stability analysis of equilibrium points
The Jacobian matrix of the fractional-order model (2) is computed by The Jacobian matrix J at the equilibrium point E 0 is also obtained as In the following, we first recall the definition of marginal stability.

Definition 4.1.
A system is marginally stable if all the poles on the imaginary axis are distinct, and all the remaining poles have negative real parts.
We state the necessary conditions for marginal stability by the following theorem.

Sensitivity analysis of R 0
To investigate the sensitivity of R 0 , we use the method used by Zafar et al. [3]. Since all parameters of the model (2) are positive, we have As can be seen, R 0 is sensitive to changes in the model parameters such that R 0 increases with Ω, β and decreases with m, µ, γ .

Optimal control
Currently, there is no vaccine to control the transmission of COVID-19, so the spread of the disease should be controlled by following health instructions and creating social distance. Meanwhile, creating social distance plays an important role in controlling the spread of COVID-19. We denote the social distance by 1 − u where u represents the social gathering. If we have social distance, then it is u = 0, and if there is not any social distance, then it will be u = 1. The aim here is to decrease the number of infected people by creating social distance while reducing the cost of implementing this control strategy. For this end, a fractional control system is considered as below If we remove the social distance, then it becomes u = 1; thus, the SIRD model (2) is recovered from the control system (8). In addition, 0 ≤ u(t) ≤ 1 represents the control space for t ∈ [0, T ].
For fractional-order systems, the optimal control theory is based on the Pontryagin's maximum principle [20]. Thus, we consider the performance index for the control system (8) as follows where 0 < b 1 , b 2 < ∞ are weighting coefficients corresponding to the number of infected individuals and the cost of control, respectively. We minimize the cost functional J(u) by the following theorem. Proof. We consider the Hamiltonian function as follows 2, 3, 4, are the co-state variables with λ i (T ) = 0, and they satisfy which is a fractional system of right Caputo derivative equations. By the Pontryagin principle, we obtain the following optimality condition Then according to the sign of ∂H ∂u , we consider the optimal control as follows: If ∂H ∂u < 0, then u * (t) = 0.
Therefore, by solving the boundary value problem (8) and (10) using the presented method in [21], the optimal control is derived. □

Simulation and numerical findings
Hereinafter, we employ the fractional Euler method to obtain the approximate solution of the transmission model (2) [22].

Numerical method
Let the model (2) be considered in the compact form where y = (S, I, R, D) ∈ (R + 0 ) 4 , g(t, y(t)) ∈ R is a vector function which is continuous and satisfies the Lipschitz condition, and y 0 = (S 0 , I 0 , R 0 , D 0 ) is the initial condition. We apply the fractional integral operator to Eq. (11); hence, we derive and t n = nq where n = 0, 1, 2, . . . , N and N is a natural number. Denote the approximation of y(t) at t = t n by y n ; then the fractional Euler method [22] provides where the coefficients are computed by the following equation Note that Theorem 3.1 in [22] proved the stability of the presented method. In addition, the solutions of the model (2) are calculated from where g 1 (t, y(t)) = Ω −mS(t)−βS(t)I(t), g 2 (t, y(t)) = −(m+µ+ γ )I(t) + βS(t)I(t), g 3 (t, y(t)) = −mR(t) + γ I(t), and g 4 (t, y(t)) = µI(t).

Case I: Iran
Iran was the first country that experienced the second wave among the others getting involved in the COVID-19 epidemic. To simulate the first and the second waves in Iran, we used the reported data for the infected cases from the world meter website [23]. The model (2) describes an endemic model (with vital dynamics) such that the births are added at the rate Ω into the class S and deaths reduce the classes at the rate m. Life expectancy in Iran is almost 73 years, so m −1 = 73 years and Ω = mN (the balance of births and natural deaths). Also, the average infectious period is γ −1 = 8 days, and the mortality rate with this disease is µ = 0.014. Additionally, the initial conditions for the first peak (started from February 22th) are S(0) = 83000000, I(0) = 34, R(0) = 6, and D(0) = 2. The auxiliary parameter is also considered as θ = 0.99 for all cases in the sequel. According to the reported cases by WHO for the first wave of infected people in Iran, we used the fmincon package of minimization and obtained β = 2.1 × 10 −9 for the model (2) with υ = 0.994. Since we used the fractional-order model for the fitting purpose, the obtained β cannot be the correct value for the integer-order case. Hence, we fitted the integerorder model separately and obtained β = 2.0973 × 10 −9 . In   Fig. 2, we plotted the results for the fractional transmission with υ = 0.994 and the integer transmission (υ = 1) in addition to the available real data. Also, Table 1 reported the errors related to each model. The results indicate that the fractional-order model has better performance and less errors. To show the influence of the fractional order, we plotted the state variables for the optimal order υ = 0.994 and the other orders υ = 0.98, 0.96, 0.94, 0.92 in Figs. 3 and 4. The figures indicate that the derivative order does not influence the behavior of the functions, although the obtained plots are different for various orders.
With the outbreak of COVID-19 in Iran, the government enforced quarantine rules until April 18th, so the number of infected people declined rapidly as shown in Fig. 2. This figure also depicts that with the continuation of this practice, the disease could be controlled in 100 days. However, the quarantine rules were lifted on April 18th; thus, by reopening of jobs, the second-wave of COVID-19 began on May 5th. In the second wave, according to the published information, the mortality rate was about µ = 0.034, and with the fitting technique, the amounts of β = 9.8 × 10 −9 and β = 9.702 × 10 −9 for the fractional-and the integer-order models were obtained, respectively. Fig. 5 shows the fractionalorder results with υ = 0.9952 and the classic integer responses with υ = 1 along with the reported cases of infected people from May 5th; Table 2 also indicates the errors for both cases.
The results indicate that the errors of the fractional-order responses are lower, a fact which confirms its superiority. Figs. 6 and 7 portray a forecasting, according to the second wave, for the spread of COVID-19 in Iran. According to these figures, it takes about eight months to pass the second wave, and about 70,000 people will die from this disease. To check the fractional order effect on the dynamics the model (2), the state variables were also plotted for the optimal order υ = 0.9952 and the other orders υ = 0.98, 0.96, 0.94, 0.92 in Figs. 6 and 7. The plots display that the derivative order does not influence the general functions' behavior, but the achieved values are different for various fractional orders, and the differences increase over the time.

Case II: Japan
We consider the second wave of COVID-19 transmission in Japan from June 30th. According to the available data, the mortality rate with respect to the COVID-19 transmission for the second wave in Japan was very low and about µ = 0.007. Life expectancy in Japan is almost 85 years, so m −1 = 85 years. Also, the other parameter is the average infectious period γ −1 = 8 days, and the initial conditions are S(0) = 19211669, I(0) = 7287, R(0) = 20534, and D(0) = 1799. According to the reported cases by WHO for the infected people in the second wave in Japan, we used the fmincon package of minimization and obtained β = 1.6 × 10 −9 for the fractional-order model with υ = 0.9961 and β = 1.573 × 10 −9 for the integer-order model with υ = 1. In Fig. 8, we plotted the results for the fractional and the integer derivatives along with the reported cases from June 30th to July 20th. As can be seen, both the integer-and fractional-order models with υ = 1 and υ = 0.9961, respectively, follow the real data very well, but Table 3 shows that the fractional-order model has less errors. Fig. 8 also forecasts the second wave of COVID-19 transmission in Japan, which will take about 5 months to subside. Figs. 9 and 10 show the results of the model in Japan, which indicates that the treatment process is well progressed, and eventually, by the end of the second wave, the number of deaths will have reached to nearly 1150 persons. Also, Figs. 9 and 10 display the results for different orders of fractional derivative; as can be seen, a small change in the derivative order makes a big difference in the resulting values, while the general functions' behaviors are the same.

Optimal control
Here we use some numerical simulations to investigate the effect of the suggested control strategy, social distance, on the outbreak of COVID-19. To investigate the effect of social distance on the second wave in Iran and Japan, we used the parameters presented in the previous section. Also, we consider n = 2, the weighting constants b 1 = 0.35, b 2 = 10, and the final intervention times T = 1000 and T = 100 days for Iran and Japan, respectively. To find the optimal control, we combined our suggested numerical method (Section 6.1) and an iterative algorithm proposed in the study [24]. Indeed, we implemented the proposed method in Section 6.1 forward and backward in time for the state and co-state equations, respectively. Note that the convergence and the stability of this technique were developed in [24]. In addition, the utilized scheme was employed in the previous literature for some practical cases such as the solution of general fractional optimal control problems [25], and the optimal control of a human respiratory syncytial virus infection [26]. Fig. 11 shows the population of susceptible and infected people in Iran in the presence and the absence of the aforesaid control strategy. It is apparent that the social distance reduces the number of infected people. Fig. 12 plots the variation in the number of susceptible and infected people in the presence and the absence of the control strategy in Japan. Our goal was to reduce the number of infected people; the results confirm that the number of infected people (peak of the disease) decreased, and since the initial number of infected people was small, this wave was ended faster, and the spread of the disease was controlled by the social distance strategy.

Conclusion
This research presented a fractional SIRD model in the Caputo sense for the spread of COVID-19 between humans. The equilibrium points, the basic reproduction number, and the feasibility region of the system (2) were determined. Also, we proved that a nonnegative solution exists for the model and checked the stability of the disease-free steady-state. An optimal control for the system under study was calculated considering the social  distancing strategy. In the numerical part, the fractional Euler method provided the approximate solutions of the system (2). Some numerical simulations were presented to investigate the first and the second waves of the disease in Iran and Japan. With the application of quarantine laws and the closure of jobs in Iran until April 18th, the first wave came to an end by May 2th.
However, by reopening the jobs on April 18th, the value of β was increased from 2.1×10 −9 to 9.8×10 −9 , and the second wave was started after two weeks. The second wave was continued for 2 months, and about 70,000 people were died at the end of the second wave. Also, comparing the results between the integer-order model, the fractional-order case, and the reported data showed that the fractional-order one follows the real data better. By plotting the results for different fractional orders, the effect of the fractional order on the results was investigated, which differed only in the resulting values and did not have any special effect on the general functions' behavior. Next, the second wave of the disease was simulated in Japan, but the amount of β produced was less than that in Iran, so we predicted that the second wave in Japan would take about 3 months. The number of infected cases was increased to 12,000, and the number of deaths from this disease was increased to 1200, which was much less than that in the second wave in Iran. This can be due to two reasons: first, because the value of β was low, the number of infected cases was slowly increased giving more time to the treatment department; second, the treatment system in Japan might be more successful. Additionally, some numerical simulations for the optimal control strategy in Iran and Japan were presented. The results showed that by creating social distance and observing basic health issues, the number of patients in Iran were decreased; also, if the desired measures continued, the number of patients would reach zero. In Japan, where the initial number of infected people was lower, social distance reduced both the peak and the duration of the disease, so the prevalence of the illness was controlled. The fractional differential system preserves the historical memory of the system, so in the second wave we should use the data from the beginning of the outbreak, but β is different in the two waves. Therefore, with the proposed model, we could not use the primary data for the second wave. Therefore, young researchers can examine how and by which method they can overcome this weakness.

Availability of data and materials
Data sharing not applicable to this article as no datasets were generated or analyzed during the current study.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.