COVID-19 is a novel virus that has spread globally, and governments around the world often implement different strategies to prevent its spread. In the literature, several COVID-19 models have been studied with the bilinear incident rate. In this study, the *S*_{1}*V*_{1}*E*_{1}*I*_{1}*Q*_{1}*R*_{1} (susceptible-vaccinated-exposed-infective-quarantined-recovered) COVID-19 model is proposed. To investigate how the disease spreads in the population, an algorithm is used. The efficacy of the algorithm is used to calculate the disease-free equilibrium point. A next generation matrix technique is used to find *R*_{0}. Furthermore, to check the effect of parameters on the basic reproduction number (*R*_{0}), the sensitivity analysis is conducted. Numerical simulation displays that the disease spreads in the population by increasing the value of the contact rate *β* while the disease spread in the population reduces by increasing the value of the vaccination rate *θ*, quarantine rate *ϕ*, and recovery rate *γ*. Different optimal control strategies, such as social distance and quick isolation, are also implemented.

## I. INTRODUCTION

As of December 2019, COVID-19 had spread worldwide after first being reported in Wuhan, China. There have been many protective measures implemented against the transmission of this virus since it can be transmitted from person to person.^{1} These precautions include the use of masks, maintaining physical distance, and the administration of vaccination in recent times. Symptoms of COVID-19 include high fever, dyspnea, and invasive multilobed lesions, as seen in chest radiographs.^{2–9} According to reports, the symptoms of this virus appear for around five days after contracting COVID-19 infection.^{1} The COVID-19 virus has been linked to hundreds of hospitalizations due to respiratory issues, flu-like symptoms, and other severe complications.^{1} In some cases, these symptoms get worse with time, leading to death.^{2}

In order to minimize COVID-19 infections, researchers are continuously working on developing effective vaccines and remedies. There have been a number of research articles with varying viewpoints published in the literature regarding minimizing COVID-19 infection. A mathematical model is the most effective tool for analyzing infectious disease transmission dynamics and evaluating the efficacy of interventions in reducing the prevalence. Many mathematical methods are under investigation for reducing the COVID-19 disease load.^{10–14} Computational modeling was employed by Imai *et al.*^{15,16} to estimate COVID-19’s prevalence in Wuhan based on the person-to-person transmission of the virus. Approximately 60% of the disease transmission can be prevented through preventive measures, according to their findings. According to Ngonghala *et al.*,^{17} COVID-19 can be reduced by maintaining social distance. Ndaïrou *et al.*^{18} proposed a COVID-19 model. In their study, they take super-spreaders spreading the disease among the population into account. Mandal *et al.*^{19} established a mathematical model to minimize the outbreak of COVID-19 by incorporating quarantine compartments and government intervention measures. Prathumwan *et al.*^{20} investigated how confinement, hospitalization, and latent compartments affected COVID-19 transmission. Using a mask, social isolation, and early case detection, Okuonghae and Omame inspected the effectiveness of control measures in reducing COVID-19 dynamics.^{21} According to Khajanchi *et al.*,^{22} infection can be roughly divided into nine stages, and the most active method for stopping infection is the combination of pharmaceutical and non-pharmaceutical preventive measures. Rai *et al.*^{23} used non-pharmaceutical interventions to study the usefulness of social media promotion to stop COVID-19 outbreaks. In order to define vaccination schedules, Acuña-Zegarra *et al.*^{24} developed a mixed constraint optimal control problem. According to the findings, vaccines and reinfection periods are the main factors influencing COVID-19 mitigation. An investigation was conducted by Kurmi and Chouhan^{25} using vaccination as a control to determine the outcome of immunization on COVID-19 outbreaks. Venkatesh and Ankamma Rao^{26} investigated the impact of both effective and ineffective vaccination on the dynamics of COVID-19.

Numerous mathematical models have been examined based on the aforementioned literature in order to analyze the transmission dynamics of COVID-19 with the bilinear incidence rate. The transmission dynamics of COVID-19 with the effect of harmonic mean type incident for contact as well as for the vaccination rate have not been studied yet. As compared to other transmission rates, harmonic mean types are more realistic. In the present study, a harmonic mean type incident rate is considered for the transmission dynamics of COVID-19. In particular, the study uses harmonic rates to examine how vaccination affects COVID-19 transmission dynamics. Furthermore, different optimal control strategies, such as social distancing and quick isolation of infected individuals, are studied.

## II. MODEL FORMULATION

An infectious disease spreading in a population is described by a system of differential equations. Each parameter in the equations represents a specific function in the context of the disease and represents the rate of change of different groups of individuals in the population. The total population *N*(*t*_{1}) = *S*_{1}(*t*_{1}) + *V*(*t*_{1}) + *E*_{1}(*t*_{1}) + *I*_{1}(*t*_{1}) + *Q*_{1}(*t*_{1}) + *R*_{1}(*t*_{1}) is subdivided into six classes: *S*_{1}(*t*_{1}) susceptible, *V*_{1}(*t*_{1}) vaccinated, *E*_{1}(*t*_{1}) exposed, *I*_{1}(*t*_{1}) infectious, *Q*_{1}(*t*_{1}) quarantined, and *R*_{1}(*t*_{1}) recovered population. Vaccinated and quarantined people are included in the model of COVID-19 presented in this study. The parameter Λ represents the birth or immigration rate of susceptible individuals into the population, while *β* is the transmission rate of the disease, indicating how likely an individual is to become infected through contact with an infected individual, The term $2\beta S1I1S1+I1$ represents the rate at which susceptible individuals become infected. Individuals susceptible to getting vaccinated are represented by *λ*, where *μ* is the natural death rate in each class. $2\beta 1\u2212\theta V1I1V1+I1$ represents the rate at which vaccinated individuals become infected. *α*_{1} is the rate at which exposed individuals become infectious, while *γ* represents the infected people isolated from the society and become quarantined, and *η* is the rate at which infected people die due to the disease. Quarantined individuals get proper treatment during isolation, denoted by $\gamma 1\u2212\delta I1$, and transfer to the recovered class at a rate *γδ*. Infected individuals become healthy during quarantine while having proper treatment, denoted by *ϕ*, and transfer to the recovered class. The transmission dynamics of Eq. (1) are presented in Fig. 1.

*S*

_{1}

*V*

_{1}

*E*

_{1}

*I*

_{1}

*Q*

_{1}

*R*

_{1}model,

### A. Positive characterization

The solution under consideration—positive $S1,E1,V1,I1,Q1,Ro\u2208R6+$ of system (1)—is seen at any time *t* ≥ 0 when the rate of change of state variables is non-negative during the trivial stage.

#### 1. Invariant region

*R*

^{6}is

*S*

_{1}

*V*

_{1}

*E*

_{1}

*I*

_{1}

*Q*

_{1}

*R*

_{1}modal at time

*t*

_{1}is specified as

### B. Disease-free equilibria and reproduction number

*E*

_{0},

*R*_{0}< 1, it indicates that the disease is dying out from the population. This helps in preventing the outbreak of epidemics. On the other hand, when

*R*> 1, it indicates that the disease is spreading in the population. The next-generation matrix process $\rho F\u2217V\u22121$ is used to compute the value of

*R*

_{0}as follows:

## III. LOCAL STABILITY OF DISEASE-FREE EQUILIBRIA

When *R*_{0} < 1, the SIR model (1) is stable locally asymptotically at disease-free equilibrium *E*^{0}.

*J*(

*E*

_{0}) is set to zero, we obtain the Jacobian matrix of model (1), which gives

According to the Routh–Hurwitz criterion of order two, $\gamma +\eta +2\mu +\alpha 1>0$ and $\u22122\beta \alpha 12\u2212\theta +\mu +\alpha 1\gamma +\eta +\mu >0$.

*R*

_{o}< 1, so DFE is locally asymptotically stable.

## IV. SENSITIVITY ANALYSIS

*R*

_{0}. If the

*R*

_{0}value is huge, then it is challenging to control a disease. By targeting the parameters that produce a greater deviation in the basic reproduction number, we aim to reach a greater level of reproducibility. Variations in a parameter can cause the associated variance in a state variable. Based on the definition given in Ref. 27, these indices have been calculated. Given below is the partial derivative representation of the sensitivity index:

The sensitivity analysis shown in Table I demonstrates that certain parameters are more sensitive to reproduction numbers than others. *μ*, ** θ**,

*α*_{1},

**, and**

*η**γ*reduce the spread of infection in the population, while on the other hand,

**increases the spread of the infection when it gets high (Fig. 2).**

*β*No. . | Parameters . | Analytical values . | Numerical values . |
---|---|---|---|

1 | β | 1 | 1.000 00 |

2 | μ | $\u2212\mu \gamma +\eta +\mu \u2212\mu \mu +\alpha 1$ | −0.947 804 993 048 809 4 |

3 | θ | $\theta \u22122+\theta $ | −0.052 631 578 947 368 42 |

4 | α_{1} | $\mu \mu +\alpha 1$ | 0.460 000 115 000 028 8 |

5 | γ | $\u2212\gamma \gamma +\eta +\mu $ | −0.365 853 658 536 585 4 |

6 | η | $\u2212\eta \gamma +\eta +\mu $ | −0.146 341 463 414 634 17 |

No. . | Parameters . | Analytical values . | Numerical values . |
---|---|---|---|

1 | β | 1 | 1.000 00 |

2 | μ | $\u2212\mu \gamma +\eta +\mu \u2212\mu \mu +\alpha 1$ | −0.947 804 993 048 809 4 |

3 | θ | $\theta \u22122+\theta $ | −0.052 631 578 947 368 42 |

4 | α_{1} | $\mu \mu +\alpha 1$ | 0.460 000 115 000 028 8 |

5 | γ | $\u2212\gamma \gamma +\eta +\mu $ | −0.365 853 658 536 585 4 |

6 | η | $\u2212\eta \gamma +\eta +\mu $ | −0.146 341 463 414 634 17 |

## V. NUMERICAL SIMULATION OF SVEIQR MODEL

The following figures show the effects of primary and secondary parameters of the *S*_{1}*V*_{1}*E*_{1}*I*_{1}*Q*_{1}*R*_{1} model specially on the infectious class to discuss and observe the effects of different parameters. In addition, we find the range of the following parameters to find if the effects are directly or inversely proportional to the infectious class and if the disease spreads to or emanates from the population.

Figure 3(a) depicts the relationship with the contact rate *β*. As the contact rate beta increases, the susceptible population likely decreases, indicating a faster spread of disease due to more frequent contacts. In Fig. 3(b), a higher *β* shows an increase in exposed individuals, representing an increase in people who have contracted the disease but are not yet infectious. Both graphs illustrate the impact of the contact rate on disease dynamics within a population.

As shown in Fig. 4(a), as the vaccination rate *λ* increases, the proportion of the susceptible population is expected to decrease. This decline occurs because more individuals are receiving vaccinations, reducing the number of people vulnerable to the disease. In contrast, Fig. 4(b) displays a direct positive correlation: as *λ* increases, the number of vaccinated individuals in the population increases accordingly. These graphs together illustrate the critical role of vaccination in controlling disease spread by decreasing susceptibility and increasing immunity within a population.

Figure 5(a) shows that as *θ* increases, the vaccine effectiveness decreases, potentially leading to more breakthrough infections, but the total number of vaccinated individuals remains unaffected by *θ*. Meanwhile, Fig. 5(b) shows an upsurge in the number of exposed individuals within the vaccinated population as *θ* increases. This indicates that higher rates of infection among vaccinated individuals correspond to more of them being exposed to the disease. These graphs collectively highlight the importance of vaccine efficacy and the implications of breakthrough infections in understanding disease dynamics in vaccinated populations.

Figure 6(a) illustrates how changes in *α*_{1} affect exposed individuals. An increase in *α*_{1} might indicate a faster transition from being exposed to becoming infectious, which could initially lead to a decrease in the exposed population as they quickly move to the infected category. Figure 6(b) shows that a higher *α*_{1} would typically correspond to an increase in the infected population. This is because a greater *α*_{1} means exposed individuals are becoming infectious more rapidly, thereby swelling the ranks of the infected group. Together, these graphs demonstrate the dynamics of disease progression within a population, particularly how changes in the rate of progression from exposed to infectious impact both the exposed and infected populations. It underscores the importance of this transition rate, considering and managing the spread of infectious diseases.

Figure 7(a) exhibits the impact of varying recovery rates on the infected population over time. As the recovery rate *γ* increases from 0.05 to 0.21, the peak of the infected population occurs earlier and at a lower level, indicating that a higher recovery rate can significantly reduce the spread of an infection and shorten the duration of an epidemic. Figure 7(b) depicts the trajectory of a quarantined population over time, given different values of parameter *γ*. Each curve represents a scenario with a distinct *γ* value, showing how the size of the quarantined population changes over a period of 100 days. The peak of each curve suggests the maximum number of people quarantined at any given time, and this peak occurs at different times depending on the value of *γ*. Higher values of *γ* result in a quicker and lower peak, indicating a faster rate of quarantine and a smaller maximum quarantined population. Conversely, lower *γ* values lead to a slower response and a higher peak. After reaching the maximum, each curve shows a decline, suggesting that fewer people remain in quarantine over time.

Figure 8(a) exhibits an increase in *ϕ* (the rate of recovery of an infected individual becoming healthy during quarantine), which might initially lead to a larger quarantine population due to effective detection and isolation. However, as *ϕ* becomes higher, indicating more successful treatments, the duration of individuals in quarantine could shorten, affecting the overall quarantine population count. Conversely, Fig. 8(b) shows a clear positive trend, with higher *ϕ* values leading to a greater number of recoveries. This indicates that effective treatment and speedy recovery from infection significantly increase the recovered population, highlighting the crucial role of medical intervention in managing disease outbreaks.

Figure 9(a) shows that as *δ* (effectiveness of treatment for quarantined individuals) increases, the quarantine population might initially increase due to more effective isolation. However, with better treatment, the duration in quarantine may decrease, leading to quicker recoveries and a fluctuating quarantine population. Figure 9(b) shows a positive trend, with an increase in *δ* resulting in more individuals recovering from the disease. This suggests that efficient and effective treatment significantly boosts the recovery rate, highlighting the vital role of medical care in managing disease outbreaks and improving patient outcomes.

## VI. FORMULATION OF THE OPTIMAL CONTROL SVEIQR MODEL

*u*

_{1}(

*t*) include social distancing and the frequent use of masks. The aim is to prevent infection from spreading to healthy individuals. COVID-19 virus isolates people quickly by using control

*u*

_{2}(

*t*). A objective functional

*J*is developed to determine the best methods for controlling the disease. In addition to reducing infectious individuals,

*u*

_{1}(

*t*) and

*u*

_{2}(

*t*) are focused on minimizing the cost of implementation,

*M*

_{1},

*M*

_{2}, and

*M*

_{3}describe the positive weights. The objective function we have described aims to reduce the costs of the above-mentioned controls while reducing the number of infections. There are two controls, $u1*$ and $u2*$, which can be found as follows:

*U*, which is called the control set and is described by

^{28}

### A. An analysis of optimal control and its existence

With the aid of powerful classical procedures, it is possible to verify the existence and analysis of an optimal control. The following hypotheses should be examined according to Ref. 29:

H1: The control set and corresponding variables described in the hypothesis are nonempty.

H2: The convexity and closeness properties of U must be satisfied.

H

_{3}: The described system is bounded in the state of control by a linear function on the right-hand side.- H
_{4}: A convex integrand on*U*bounds the objective functional*J*,According to Lukes,(33)$d1\u2211i=13ui22\u2212d2,whered1,d2>0,andt1>1.$^{30}the solution in system (29) exists if the result is valid. Thus, the hypothesis given above can be verified. Bounded coefficients and solutions are in accordance with H1, and boundedness demonstrates that*U*fulfills the second hypothesis. The right-hand side of (29) fulfills the third hypothesis H3 because system (39) is bilinear in*u*_{1}and*u*_{2}and because the solutions are bounded. We can choose constants*M*_{1},*M*_{2},*M*_{3},*C*_{1},*C*_{2},*d*_{1}, and*d*_{2}that are all positive and*t*_{1}> 1. The inequality holds since control functions*u*_{1}and*u*_{2}are convex and the number of infectious hosts and vectors is also convex,Thus, we can also verify the last condition. As a result, we arrive at the following theorem:(34)$M1E1+M2I1+M3Q1+C1u12+C2u22\u2265d1\u2211i=13ui2t12\u2212d2.$Based on the above-mentioned formula, one must construct the least possible value. Defining the Hamiltonian this way gets us to objective results: Let us examine $X=S1,V1,E1,I1,Q1,R1,U=u1,u2$ and $\psi =\psi 1(t1),$ $\psi 2(t1),\psi 3(t1),\psi 4(t1),\psi 5(t1),\psi 6(t1)$, with the intention of obtaining(35)$LE1,I,1Q1,u1,u2=M1E1+M2I1+M3Q1+C1u12+C2u22\u2265d1u12+u2\u221222.$For determining the essential optimal control conditions, Pontryagin’s maximum principle is used as reported in Ref. 31. Using the following formula, this can be calculated: For optimized solution $u1*,u2*$ of optimal control problem (30), a vector function $\psi =\psi 1(t1),\psi 2(t1),$ $\psi 3(t1),\psi 4(t1),\psi 5(t1),\psi 6(t1)$, which is nontrivial, prevails, which fulfills the ensuing conditions. In the statement, it is stated that(36)$HX,U,\psi =M1E1+M2I1+M3Q1+C1u12+C2u22+\psi 1(t1)\Lambda \u22122\beta 1\u2212u1S1I1S1+I1\u2212\lambda S1\u2212\mu S1+\psi 2(t1)\lambda S1\u22122\beta 1\u2212\theta V1I1V1+I1\u2212\mu V1+\psi 3(t1)2\beta 1\u2212u1S1I1S1+I1+2\beta 1\u2212\theta V1I1V1+I1\u2212\alpha 1E1\u2212\mu u2E1\u2212\omega u2E1+\psi 4(t1)\alpha 1E1\u2212\gamma I1\u2212\eta I1\u2212\mu I1+\psi 5(t1)\gamma 1\u2212\delta I1\u2212\varphi Q1\u2212\mu Q1+\omega u2E1+\psi 6(t1)\varphi Q1+\gamma \delta I1\u2212\mu R1,.$Currently,(37)$\u2202xdt1=\u2202\u2202\psi Ht1,U*,\psi (t1),0=\u2202\u2202uHt1,U*,\psi (t1),\u2202\psi dt1=\u2212\u2202\u2202xHt1,U*,\psi (t1).$*H*is subject to mandatory conditions.

For $U=u1*,u2*$ and solution $X=S1,V1,$ $E1,I1,Q1,R1$ of system (1), variables $\psi =\psi 1(t1),\psi 2(t1),\psi 3(t1),$ $\psi 4(t1),\psi 5(t1),\psi 6(t1)$.

*H*with describe to variables $X=S1,V1,E1,I1,Q1,R1$, Hence, we have

*ψ*

_{1}(

*T*) =

*ψ*

_{2}(

*T*) =

*ψ*

_{3}(

*T*) =

*ψ*

_{4}(

*T*) =

*ψ*

_{5}(

*T*) =

*ψ*

_{6}(

*T*) = 0. As a result, the property set

*U*of controls and the optimality conditions give the result

### B. Numerical simulation of optimal control

Parameters used for simulated optimal control in this subsection are mentioned in Table II. Moreover, we compare the results of the optimal control applied to the time intervals without and with control. Indicated by the weight of the treatment. control *u*_{1} and weight accompanying the infected population *u*_{2} are taken as unity. In this second attempt at solving the optimal system (10), we will use iterative Runge–Kutta procedures of fourth order. A Runge–Kutta fourth order procedure is applied to solve state (10), based on our initial guess for the state variables. Our primary objective is to address the deputy variables through backward Runge–Kutta fourth-order procedures while simultaneously including state variables and transverse conditions. In addition to providing a comprehensive understanding of system dynamics, this approach enhances the accuracy of our analysis.

Parameters . | Estimated values . | References . |
---|---|---|

Λ | 0.464 | 32 |

β | 0.503 | 33 |

λ | 0.4 | 32 |

μ | 0.04 | Assumed |

θ | 0.1 | Assumed |

α_{1} | 0.046 956 5 | Assumed |

γ | 0.030 | Assumed |

η | 0.012 | 33 |

δ | 0.2 | Assumed |

ϕ | 0.082 | 34 |

ω | 0.3 | Assumed |

Figure 10(a) shows the number of people who are susceptible. Without interventions, the number of susceptible people decreases more slowly over time (red curve), indicating a slower disease spread. With interventions (blue curve), the number of susceptible people decreases more rapidly, suggesting that the interventions are operative in dropping the disease outbreak.

Figure 10(b) shows the number of vaccinated individuals. Without optimal control, the number of vaccinated people remains at zero (red curve), as expected. With optimal control, the number of vaccinated people increases and then decreases over time (blue curve), reflecting the implementation of a vaccination campaign. After reaching a peak, it declines as fewer susceptible people are left to vaccinate.

Figure 11(a) illustrates the population that has encountered the infection but has not yet reached the stage of being contagious. The red curve demonstrates a rapid increase in exposed individuals when no control measures are in place, reaching a peak before declining. The blue curve, where control measures are *u*_{1} ≠ 0 and *u*_{2} ≠ 0, shows a slower increase in exposed individuals, indicating that the measures are delaying the onset of the disease and reducing the number of exposed individuals.

Figure 11(b)’s graph represents the number of people who are infectious and can spread the disease. Similar to the left graph, the red curve without control measures shows a sharp increase to a peak, reflecting a rapid spread of the disease. The blue curve with control measures in place increases more slowly and peaks at a lower level, suggesting that the interventions are operative in dropping the spread of the disease and the peak infectious population.

Together, these graphs highlight the influence of interventions in controlling the spread of an infectious disease, showing that appropriate measures can significantly reduce the number of exposed and infectious individuals in a population over time.

Figure 12(a) demonstrates the number of people in quarantine over time. The red curve indicates the situation without any control measures *u*_{1} = 0, *u*_{2} = 0, showing a slow increase in quarantined individuals. The blue curve represents the scenario with control measures *u*_{1} ≠ 0, *u*_{2} ≠ 0, which leads to a much more pronounced upsurge in the number of quarantined individuals. This implies that with the implementation of control measures, more people are being identified and quarantined to stop further spread of the disease.

Figure 12(b) displays the number of people who have recovered from the disease over time. Again, the red curve shows the recovery trajectory without control measures, which results in a slower and lower number of recoveries. The blue curve, with control measures in place, shows a more rapid and higher number of recoveries, indicating that the control measures are effective in both preventing the spread of the disease and helping more people recover.

These graphs suggest that the interventions, which could include quarantine, isolation, vaccination, or treatment strategies, not only reduce the spread of the infection but also contribute to a faster and more effective recovery of the population.

## VII. CONCLUSION

Our study utilized a deterministic model for assessing the emergence and spread of COVID-19 cases, demonstrating its effectiveness. The positivity and boundedness of model (1) are determined. The reproduction number of model (1) is computed. Local and global stability of model (1) is also analyzed. By identifying the factors with the most significant impact on the basic reproduction number, such as the contact rate between susceptible and infectious individuals, our study highlights the crucial elements that need monitoring in pandemic response efforts. The sensitivity analysis conducted as part of our study is particularly valuable for understanding how different parameters interact and influence the spread of the virus. This can inform global public health policies, especially in regions with varying demographic and social characteristics. Our findings regarding the effectiveness of control measures such as social distancing, mask-wearing, and rapid quarantine are particularly relevant in the current global scenario. The use of control functions in our model to simulate these measures provides insights into their potential impact on dropping the spread of the virus. This aspect of our study can guide policymakers in implementing targeted interventions.

## ACKNOWLEDGMENTS

Princess Nourah bint Abdulrahman University Researchers Supporting Project number (PNURSP2024R404), Princess Nourah bint Abdulrahman University, Riyadh, Saudi Arabia.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Kamil Shah**: Conceptualization (equal); Data curation (equal); Investigation (equal); Methodology (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Jamal Shah**: Conceptualization (equal); Data curation (equal); Investigation (equal); Methodology (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Ebenezer Bonyah**: Conceptualization (equal); Data curation (equal); Investigation (equal); Methodology (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Tmader Alballa**: Investigation (equal); Methodology (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). **Hamiden Abd El-Wahed Khalifa**: Data curation (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). **Usman Khan**: Methodology (equal); Writing – original draft (equal). **Hameed Khan**: Investigation (equal); Writing – original draft (equal).

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

## REFERENCES

*Differential Equations: Classical to Controlled*