Ranking the effectiveness of non-pharmaceutical interventions to counter COVID-19 in UK universities with vaccinated population

Abstract

Several universities around the world have resumed in-person teaching after successful vaccination campaigns have covered 70/80% of the population. In this study, we combine a new compartmental model with an optimal control formulation to discover, among different non-pharmaceutical interventions, the best prevention strategy to maximize on-campus activities while keeping spread under control. Composed of two interconnected Susceptible-Exposed-Infected-Quarantined-Recovered (SEIQR) structures, the model enables staff-to-staff infections, student-to-staff cross infections, student-to-student infections, and environment-to-individual infections. Then, we model input variables representing the implementation of different non-pharmaceutical interventions and formulate and solve optimal control problems for four desired scenarios: minimum number of cases, minimum intervention, minimum non-quarantine intervention, and minimum quarantine intervention. Our results reveal the particular significance of mask wearing and social distancing in universities with vaccinated population (with proportions according to UK data). The study also reveals that quarantining infected students has a higher importance than quarantining staff. In contrast, other measures such as environmental disinfection seems to be less important.

The epidemic ascribed to the virus called Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) has greatly impacted society and the economy around the world. In December 2019, China reported cases of pneumonia of unknown cause in Wuhan. The World Health Organisation (WHO) declared the coronavirus disease 2019 (COVID-19) a pandemic in March 11, 2020¹. By 3 November 2021, WHO has reported over 246 million cases with more than 5 million deaths around the world¹.

The virus generally causes respiratory symptoms such as cough, sneezing, shortness of breath, along with other symptoms including fever, headache² and olfactory or gustatory dysfunctions³. Since direct contact and aerosol transmissions are two important ways of infection⁴, symptomatic individuals are high spreaders. Meanwhile, the virus can also survive on surfaces and then invade the human body through eyes, nose or mouth via touching^5,6. Some carriers are asymptomatic, but they can still infect other susceptible individuals^7,8.

Many countries have conducted successful vaccination campaigns. However, vaccination uptake in most of these countries have platooned at around 70/80% of their total population⁹. Moreover, SARS-CoV-2 has shown a relatively high ability to adapt¹⁰. Several variants have been reported¹¹ with four considered to be of main concern: Alpha (B.1.1.7), Beta (B.1.351), Gamma (P.1) and Delta (B.1.617.2). In particular, the Delta variant has shown to be more infectious and more severe, leading to second or third waves in the UK, India and South Africa¹², in addition to having partial resistance to vaccines^13,14. The fast emergence of viral mutations has also raised wide considerations on whether current vaccines will be effective on new lineages appearing in the future^15,16,17. Moreover, recent studies also show a decay of the protection that vaccines offer as time from vaccination increases^18,19, which prompted several counties to start a booster campaign. Therefore, in such a complex situation the implementation of non-pharmaceutical interventions such as mask wearing and social distancing have remained fundamental measures in containing the disease²⁰.

This work studies the situation in British universities in the Fall 2021 when the institutions re-opened (requiring that all non-self-isolating students attend lectures in person) while maintaining a combination of non-pharmaceutical measures in place. The objective of universities is to maximise on-campus activities while maintaining the spread of the disease under control. Universities are “small-environments” which have special features for which general purpose models may be inadequate. For instance, a university is composed of two fundamentally different populations, the students and the staff (e.g. professors, lecturers, technicians and administrators), which have different degrees of interaction, vaccination rates and serious symptoms. General purpose models focus mostly on modelling the spread of the disease, but here we are interested in maximising on-campus activity subject to limited spread. Our main contribution is two-fold: on one hand we provide a modelling framework to maximise safe on-campus activity. On the other hand, a ranking of non-pharmaceutical interventions and their fundamental importance in achieving the objective naturally emerges from our analysis. Moreover this emergent behaviour is shown to be relative robust to modelling parameters.

The first block of our framework is a compartmental model. Compartmental modelling is a popular choice for research on COVID-19. Cooper et al.²¹ used an SIR model with changing total population to estimate the growth of the epidemic in different nations. To explicitly model details such as incubation period, hospitalization, and quarantine, Leontitsis et al.²² proposed an SEAHIR model. Giordano et al.²³ created a far more complete SIDARTHE model, reflecting potential effects of non-pharmaceutical interventions implemented by the government. Various stochastic versions with higher complexity have also been designed to estimate the development of the pandemic under feasible countermeasures^24,25,26. However, most of these models consider the influence of prevention and control measures by tweaking the model’s parameter values. The use of control theory for analysis of COVID-19 countermeasures is then suggested by, e.g., Zhong et al.²⁷. Therein, a feedback control law is designed for a SIRD model to estimate the combined effects of three non-pharmaceutical interventions. Nevertheless, this study does not systematically evaluate and compare the importance of different interventions under different constraints. In this paper, we propose a new deterministic compartmental model for the spreading of COVID-19 in universities and investigate the importance of non-pharmaceutical countermeasures using optimal control techniques.

Consisting of a double SEIQR structure, our model distinguishes the epidemic evolution stages of students from those of staff. All individuals can potentially undergo five statuses: susceptible S (including the vaccinated), exposed E (asymptomatic), infected I (symptomatic), quarantined Q (hospitalized or isolated), and recovered R. Since all members are studying or working within the same confined space, the coronavirus can also transmit via the environment. An extra compartment C is used to represent the environmental virus concentration. The overall structure of the model is depicted in Fig. 1. This compartmental model automatically assumes that the total population is homogeneously mixed, which is reasonable because everyone belonging to the same department is following similar daily routines in the common confined spaces. The individuals in the compartments also have similar probability of infections due to the virus surviving in the environment. The model does not involve the isolation of asymptomatic subjects because it is less likely to occur in the context of a university population. Since the vaccinated and boosted individuals can still get infected, we regard them as susceptible but with lower infection rates (i.e. vaccine breakthrough infections). Because of the relatively short term of the time frame of this study, both the decrease of vaccine protection rate and the possibility of re-infection among recovered people are omitted. Death of vaccinated individuals is also negligible.

Figure 1

Model structure. Flow chart of five epidemic stages among students and staff in a university department: S, susceptible (including the vaccinated); E, exposed (asymptomatic); I, infected (symptomatic); Q, quarantined (hospitalized or mandatorily isolated); R, recovered. The subscript “y” stands for students while the subscript “s” denotes staff.

To assess the efficacy of the various intervention strategies under study we first evaluate a baseline scenario in which no countermeasures are implemented. The parameters for this scenario are taken from the literature. Some of the values of the parameters are obtained according to the proportion of vaccinated individuals and vaccine effectiveness reported in the UK^28,29,30. A detailed description of this procedure is reported in the “Methods”.

We consider five possible non-pharmaceutical interventions that can be implemented by the university after reopening: mask wearing, social distancing, environmental disinfection, quarantine on infected students and quarantine of infected staff. In the United Kingdom, surgical masks are the most commonly employed in the mask-wearing policy, while spraying of disinfectants in the common spaces and cleaning of frequently touched surfaces like desks and chairs in the lecture rooms are the most common environmental disinfection measures. To study the effectiveness of these measures in the model, we define five input control variables associated with each of these measures. We then use optimal control to study four desired scenarios: minimum number of cases, minimum intervention, minimum non-quarantine intervention and minimum quarantine intervention. In all scenarios, it is assumed that the University starts to control the epidemic 14 days after the asymptomatic carriers firstly appear. We also model the desire of keeping the infection under control as constraints in the optimisation. In particular, we impose constraints on the number of infected cases and on the number of days required to extinguish the epidemic. In the minimum-case scenario, the epidemic is controlled with no efforts spared, leading to the strongest minimization of COVID-19 cases. For minimum intervention, we study the possibility of minimizing the total effort of all control measures. In the minimum non-quarantine intervention, we minimise the use of mask wearing, social distancing and environmental disinfection. In the last scenario, minimum quarantine intervention, we minimize the use of quarantines. By comparing the obtained optimal trajectories in different scenarios, we identify an emergent behaviour that shows a ranking on the importance of the different non-pharmaceutical interventions.

Results

Baseline: no interventions

We first generate a baseline scenario in which no intervention is performed. We simulate the evolution of the epidemic in one department of the university. The results are shown in Fig. 2. As case study we consider a department where the total number of students is 1200 and the number of members of staff is 150 (i.e. similar to the EEE department of Imperial College London). We focus on the case study of students returning to in-person teaching after a period of lockdown or remote teaching. For this reason, the initial condition of the study corresponds to a small number of asymptomatic cases. Hence, we assume that 5 students and 2 staff members are initially asymptomatic (the exact numbers used for the initial condition do not matter as long as they are of the same order of magnitude). Based on our analysis of the model parameter values, the effective reproduction number $R_{0}$ on day 0 is 2.50 among non-immune population, but it reduces to 1.40 once we factor in the effect of vaccination. After two weeks, $1.8\%$ of students and $3.1\%$ of staff have caught the disease. Without countermeasures, increasingly more individuals will get infected during the coming 120 days. Simultaneously, $R_t$ keeps decreasing and becomes smaller than 1 after 76 days. At day 134, $56\%$ of students and $63\%$ of staff have been infected and $R_t = 0.718$. In the end (without considering further infections from outside the department), the epidemic last around 250 days, rendering $59\%$ of students and $66\%$ of staff infected. This prediction result reveals the necessity of imposing non-pharmaceutical measures at the current UK vaccination levels (see “Methods” for the exact percentages).

Figure 2

Prediction in the baseline scenario of no interventions. (a) Evolution of COVID-19 among students. (b) Evolution of COVID-19 among staff. Magnitudes are in proportion to the total number of students or staff.

Optimal interventions for four objectives

We now study the effect of non-pharmaceutical interventions. Which interventions to implement, when and how strongly are all decisions made by the optimisation method to optimise the objective. We have selected four different optimisation objectives. In all scenarios, interventions are introduced 14 days after the initial exposure of 5 students and 2 staff members. In all scenarios, we impose the constraints that the epidemic must end within 120 days and that at least $94\%$ of students and staff are not infected. Thus, the overall timeline is 134 days.

Minimum number of cases

In this scenario the optimisation objective is formulated to minimise infections, even though this may require that all non-pharmaceutical interventions are implemented at full strength. The optimal trajectories are depicted in Fig. 3. The epidemic is completely ended at around the 60th day when the individuals are only in two states: susceptible and recovered. More than $97\%$ of students and $96\%$ of staff do not get infected in this scenario. From the figure we see that this result is achieved by implementing all interventions unreservedly, from mask wearing to mandatory quarantine. This reduces $R_t$ to around 0.217. While initially there is a strong need for all countermeasures, after 30 days the optimal strategy relies mostly on distancing and masks. All cases have been quarantined by the 45th day. After approximately 60 days, the department reaches a steady state and the infection is stopped. Consequently, at this point the optimal strategy eases the interventions because the epidemic has been successfully contained.

Figure 3

Optimal trajectories for the minimum-case scenario. Optimal trajectories when the department spares no effort to contain the epidemic. (a) and (b) The epidemic evolution among students and staff, respectively. (c) The optimal strategies for mask wearing ($\kappa _{m}$), social distancing ($\kappa _{d}$), and environmental disinfection ($\kappa _{e}$). (d) The optimal strategies for mandatory quarantine on infected students ($\kappa _{qy}$) and infected staff ($\kappa _{qs}$).

Minimum intervention

In this case the objective function is the norm of all control variables. As a consequence, the aim here is to minimise the use of all interventions (including quarantine) while still satisfying the constraint that at least $94\%$ of the population is not infected. The resulting optimal trajectories are shown in Fig. 4. The epidemic is ended with $4.3\%$ of students and $6\%$ of staff having been infected. With respect to before we can see that there is a decrease in the strength of the interventions. We can also note that there is an emerging ranking between the interventions. Figure 4c shows that the environmental disinfection is far less important than mask wearing. For what concerns mandatory quarantine shown in Fig. 4d, we notice again that isolation of infected students plays a more significant role in controlling the spread of COVID-19 than the isolation of staff. All five control interventions are strongest at the beginning of the epidemic, and then their magnitude is attenuated gradually over time.

Figure 4

Optimal trajectories for the minimum intervention scenario. Optimal trajectories when the department would like to minimise the enforcement of control measures. (a) and (b) The epidemic evolution among students and staff, respectively. (c) The optimal strategies for mask wearing ($\kappa _{m}$), social distancing ($\kappa _{d}$), and environmental disinfection ($\kappa _{e}$). (d) The optimal strategies for mandatory quarantine on infected students ($\kappa _{qy}$) and infected staff ($\kappa _{qs}$).

Minimum use of non-quarantine interventions

In this case we want to minimize the use of masks, social distancing and environmental disinfection. As a result we expect an increase of the use of quarantines. The resulting optimal trajectories are shown in Fig. 5. As expected the figures show little use of non-quarantine interventions, and a strong use of quarantines. We stress that the primary objective of keeping $94\%$ of the susceptible population infection free is maintained. Again, Fig. 5d demonstrates the higher importance of quarantine among students with respect to staff. This shows the predominant role played by student quarantine in controlling the epidemic.

Figure 5

Optimal trajectories for minimum use of non-quarantine interventions. Optimal trajectories when the departments would like to minimise use of masks, social distancing and environmental disinfection. (a) and (b) The epidemic evolution among students and staff, respectively. (c) The optimal strategies for mask wearing ($\kappa _{m}$), social distancing ($\kappa _{d}$), and environmental disinfection ($\kappa _{e}$). (d) The optimal strategies for mandatory quarantine on infected students ($\kappa _{qy}$) and infected staff ($\kappa _{qs}$).

Minimum quarantine

In this scenario we want to minimise the use of quarantine, but we allow a strong use of mask wearing, social distancing and environmental disinfection. The resulting optimal trajectories are shown in Fig. 6. As a result, the quarantine control variables are zero and the other three interventions are major tools to resolve the epidemic in this scenario. Figure 6c clearly shows the primary role of mask wearing and social distancing in keeping the infection under control. The figure also shows that environmental disinfection in comparison play little role when strong mask wearing and social distancing are in place.

Figure 6

Optimal trajectories for the minimum quarantine scenario. Optimal trajectories when the department would like to minimise the enforcement of mandatory quarantines. (a) and (b) The epidemic evolution among students and staff, respectively. (c) The optimal strategies for mask wearing ($\kappa _{m}$), social distancing ($\kappa _{d}$), and environmental disinfection ($\kappa _{e}$). (d) The optimal strategies for mandatory quarantine on infected students ($\kappa _{qy}$) and infected staff ($\kappa _{qs}$).

Discussion

The figures show an emerging behaviour: non-pharmaceutical interventions have different importance and this importance arises mathematically from the evolution of the epidemic. In a typical university department composed of 1200 students and 150 staff, with a vaccination rate of $68\%$ for students and $78.8\%$ for staff (see “Methods”,²⁸) we see that the implementation of non-pharmaceutical interventions is still fundamental to reduce the number of infections to one tenth of the number of infections appearing in a completely uncontrolled scenario.

The ranking that arises from the study is as follows: wearing masks is the most effective measure among the considered interventions. Keeping social distance is ranked close second. This priority of mask wearing is reasonable in a university, where close contact is often unavoidable. Furthermore, we can also see that environmental disinfection seems to be far less necessary if both measures are strongly enforced. As for the enforcement of mandatory quarantines, the result yields that quarantine of symptomatic students is more significant than quarantine of staff. This ranking is robust with respect to model parameters. This is shown in the sensitivity study presented in the Supplementary Material.

Our ranking results between mask-wearing, social distancing, environmental disinfection and overall mandatory quarantine (or even lockdown) are similar to the results of other studies that consider the entire (country) population^27,31. The two main differences here are that we separate the small-environment population into two groups and that the policy is decided algorithmically in an optimal manner, rather than by a human policymaker. Our results also provide the optimal variations in implementation strength of different measures within different time periods.

From a practical perspective, the university should emphasize mask wearing and social distancing when on-campus teaching is resumed, especially among students. The study also suggests that the university should invest particular effort in identifying and quarantining infected students. Therefore, our findings demonstrably reflect the importance of different non-pharmaceutical interventions and help assess the trade-off between high-quality teaching and limiting COVID-19 infections. This is crucial at a time in which universities are under pressure to increase on-campus activities.

Our study indicates that optimal control theory can be used to determine the optimal combination of non-pharmaceutical interventions which do not just take into account the objective of no-infections, but also different goals, such as socioeconomic needs. In this regard, note that the optimal control idea can be applied to models that are completely different from the one we proposed, providing further flexibility and impact beyond this specific case study.

Our study also shows a gap in the epidemiological research regarding the evolution of the pandemic in small environments. There are plenty of small environments, such as primary schools, certain companies and nursing homes that have distinct populations (e.g. elderly and nurses) with homogeneous within-population characteristics but different cross-population characteristics (e.g. age, immunity) and require design tools to assess the best interventions to be implemented in order to maximise in-person activities while keeping the infections under control. This idea could even be used in studies of diseases different from COVID-19.

We also point out that some practical factors are not explicitly considered in our study. Firstly, the model does not consider the infections brought from outside the campus. We omitted this aspect because we wanted to focus on the study of the priority of different mitigation measures. The model could be modified to take into account the scenario of a sudden change in the external environment by either introducing another virtual compartment or using recently published results on the use of hybrid models to address exactly the problem of abrupt changes in the populations³². Another limitation is that the input variables (the interventions) in our optimal control problem are continuous in magnitude. It may be difficult to give practical significance to the numerical values representing the interventions. However, we stress that the optimal trajectories are used here only to compare the relative importance of different interventions. Further research can be done on discretizing the magnitude of the input variables into specific levels that correspond to scientifically-defined interpretable practical meanings. Finally, we omitted the effects of contact tracing measures, the analysis of which entails a network model that portrays detailed disease transmission paths among the population³³.

Methods

Uncontrolled university model

One of the main premises of the work is that the model is constructed based on the mean-field assumptions on two groups of populations (students and staff) that have different characteristics. Students in the university are a strongly homogeneous population. They attend lectures in consistent groups and often live in the same accommodation (student halls). Most of their social interactions are with other students, through students societies. Apart for some exceptions, students also have the same age (18–22 years old), the same vaccination rates, the same response to the virus in terms of infection rates and morbidity. Moreover, the university is perfectly capable of imposing population-wide interventions (e.g. all students need to wear a mask). In contrast, staff are a bit less homogeneous population than students but share common characteristics, while being very different from the student population. For instance, as the staff group is older, it is more prone to infection and therefore they have higher infection rates. Another example is that the “exposed” staff have higher probability to present symptoms and so have enhanced infectiousness once infected. Thus, the mean-field assumption is an acceptable approximation in this case and all these characteristic differences have been factored in the various parameters of the model.

The proposed double SEIQR model is described by 11 differential equations, 5 associated to the epidemic evolution of students, 5 associated to the epidemic evolution of staff, and 1 representing environmental infection. We first introduce, describe and analyse the uncontrolled model. In “Control of the university model” section we modify the model by introducing the control variables that represent non-pharmaceutical interventions. The double uncontrolled SEIQR model is described by

$$\left\{ \begin{aligned} \frac{{\text{d}} S_{y}(t)}{{\text{d}} t}&= -\{\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)] + \beta _{cy}C(t)\}S_{y}(t) \\ \frac{{\text{d}} S_{s}(t)}{{\text{d}}t}&= -\left\{ \beta _{s}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] + \beta _{cs}C(t)\right\} S_{s}(t) \\ \frac{{\text{d}} E_{y}(t)}{{\text{d}} t}&= \{\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)] + \beta _{cy}C(t)\}S_{y}(t) \\&\quad - (\varepsilon _{y} + \xi _{y})E_{y}(t) \\ \frac{{\text{d}} E_{s}(t)}{{\text{d}} t}&= \left\{ \beta _{s}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] + \beta _{cs}C(t)\right\} S_{s}(t) \\&\quad - (\varepsilon _{s} + \xi _{s})E_{s}(t)\\ \frac{{\text{d}} I_{y}(t)}{{\text{d}} t}&= \varepsilon _{y}E_{y}(t) - (\eta _{y} + \rho _{y})I_{y}(t)\\ \frac{{\text{d}} I_{s}(t)}{{\text{d}} t}&= \varepsilon _{s}E_{s}(t) - (\eta _{s} + \rho _{s})I_{s}(t)\\ \frac{{\text{d}} Q_{y}(t)}{{\text{d}} t}&= \eta _{y}I_{y}(t) - \varphi _{y}Q_{y}(t)\\ \frac{{\text{d}} Q_{s}(t)}{{\text{d}} t}&= \eta _{s}I_{s}(t) - \varphi _{s}Q_{s}(t)\\ \frac{{\text{d}} R_y(t)}{{\text{d}} t}&= \xi _{y}E_{y}(t) + \rho _{y}I_{y} + \varphi _{y}Q_{y}(t)\\ \frac{{\text{d}} R_s(t)}{{\text{d}}t}&= \xi _{s}E_{s}(t) + \rho _{s}I_{s} + \varphi _{s}Q_{s}(t)\\ \frac{{\text{d}} C(t)}{{\text{d}} t}&= \mu _{y}(E_{y}(t) + kI_{y}(t)) + \mu _{s}(E_{s}(t) + kI_{s}(t)) - \delta C(t) \end{aligned}\right.$$

(1)

where all model parameters are denoted by Greek letters with specific biological meanings. We now provide a detailed explanation of each parameter. We stress that this model is uncontrolled, so the values discussed below are for the baseline scenario, i.e. no intervention is implemented. Also, the parameters are firstly introduced for unvaccinated population and then modified according to the UK vaccination proportions.

• Individual infection rates $\beta _{y}$, $\beta _{s}$ $k\beta _{y}$, $k\beta _{s}$

  The individual infection rates represent the average number of susceptible individuals who can be infected by a virus carrier via direct contacts in unit time. $\beta _y$ indicates student-to-student and staff-to-student infection rates. $\beta _s$ indicates staff-to-staff and student-to-staff infection rates. In other words, we assume that student-to-student and staff-to-student rates are the same and staff-to-staff and student-to-staff rates are the same. $\beta _s$ and $\beta _y$ are the infection rates of the asymptomatic compartments. Considering the age difference, it is reasonable to expect that $\beta _{y} < \beta _{s}$ because staff are more likely to get infected. $k\beta _s$ and $k\beta _y$ are the infection rates of the symptomatic group. Since sneezing and coughing play a major role in the direct transmission of the virus, the symptomatic carriers generally have larger infection rates, i.e. $k>1$. Here it is assumed that k has the same value among students and staffs. Referring to the scenario 3 of pandemic planning produced by the CDC³⁴, k is generally 4. However, in this confined environment case, k should be smaller because asymptomatic subjects can spread the virus more easily. We selected a value of $k = 1.5$. According to the study conducted by Leontitsis et al.²², the general infection rate $\beta$ is 0.1466. This value is expected to be larger in a confined space because of the higher number of direct contacts between people. In summary, putting together all these data and observations, the parameters have been selected as $\beta _{y} = 0.163$, $\beta _{s} = 0.225$, $k\beta _{y} = 0.2445$ and $k\beta _{s} = 0.3375$.

• Environmental infection rates $\beta _{cy}$, $\beta _{cs}$

  These parameters represent how many susceptible people are infected by the contaminated environment in unit time. They are properties of the virus in the environment. There is no clear value in the literature and we estimate the value of $\beta _{cy}$ to be 0.171 (based on the expected basic reproduction number). Moreover, it is reasonable to expect that the ratio $\beta _{s}/\beta _{y}$ equals the ratio $p = \beta _{cs}/\beta _{cy}$. Then $\beta _{cs} = p\beta _{cy} = 0.236$.

• Probability of becoming symptomatic $\varepsilon _{y}$, $\varepsilon _{s}$

  They are the inverse of the average incubation period. According to³⁵, this average period is 5 days. Considering that staff are of higher age, we set $\varepsilon _{y} = 1/5 = 0.2$ and $\varepsilon _{s} = 1/10 = 0.1$.

• Probability of recovery from an asymptomatic state $\xi _{y}$, $\xi _{s}$

  Similarly, the inverses of these quantities denote the average number of days spent by exposed/asymptomatic subjects to recover (i.e. they present no symptom during the whole period). Referring to³⁴, this portion accounts for $15\%$. Since $\varepsilon _{y} = 0.2$ and $\varepsilon _{s} = 0.1$ (corresponding to $85\%$ of the population). Then $\xi _{y} = 0.0353$ and $\xi _{s} = 0.0176$.

• Isolation rates $\eta _{y}$, $\eta _{s}$

  These parameters denote the proportion of symptomatic individuals who are isolated due to serious illness or mandatory quarantine. Since we initially model the unrestricted situation (i.e. no quarantines), infected individuals are at this point isolated mainly due to hospitalization. According to the survey conducted and reported by³⁶, the values are $\eta _{y} = 0.06$ and $\eta _{s} = 0.106$.

• Recovery rates of infected individuals $\rho _{y}$, $\rho _{s}$

  $\rho _{y}^{-1}$ and $\rho _{s}^{-1}$ indicate the average time for infected people which are not isolated to recover. If mandatory quarantine is not implemented, mild cases will not isolate. These two parameters mainly describe the recovery rate of this group. According to³⁷, the average length of recovery is approximately 14 days. In mild cases, we set this length at 10 days for students, and 18 days for staff. Thus, $\rho _y = 1/10$ and $\rho _s = 1/18$.

• Recovery rates of quarantined individuals $\varphi _{y}$, $\varphi _{s}$

  $\varphi _{y}^{-1}$ and $\varphi _{s}^{-1}$ indicate the average time for “quarantined” individuals to recover. In the baseline scenario this refers to hospitalised individuals. It may take six-nine weeks for severe cases to recover³⁸. We set $\varphi _{y} = 1/40$ while $\varphi _{s} = 1/55$.

• Virus shedding rates to the environment $\mu _{y}$, $\mu _{s}$, $k\mu _{y}$, $k\mu _{s}$

  These parameters measure the spread of the virus from asymptomatic/symptomatic individuals to the environment, with the effects brought by symptomatic subjects being higher. Similarly to the case of $\beta _{cy}$ and $\beta _{cs}$, there is no clear value in the literature for these parameters. In this “small environment” model we expect the values to be $\mu _{y} = \mu _{s} = \mu = 0.25$ (based on the expected basic reproduction number).

• Virus decaying rate in the environment $\delta$

  This measures the speed of decay of the virus in the small environment. Since the airborne virus could stay in aerosol for up to 1 day and survive on the surface for longer^39,40,41, this rate $\delta$ is set at 0.7.

We denote $N_{t}$ to be the total population in the department. The total number of students is represented by $N_{y}$ and number of staff is labelled by $N_{s}$. A summary of the meaning the parameters is given in Table 1.

Table 1 Model parameter definitions.

These initial values refer to studies on the COVID-19 epidemic before vaccination. We now describe how the parameters are adapted to a vaccinated population. According to the UK data provided by²⁸, about $68\%$ of young adults between 18 and 24 years old have been vaccinated by the 1st November 2021. This quantity become $78.7\%$ among people aged between 25 and 64. Since vaccines utilized in the UK can reduce COVID-19 infections by around $65\%$^29,30, we can see that the $68\%$ vaccinated students will have $65\%$ less probability of getting infected. The same happens to the $78.7\%$ of staffs. Therefore, the average reductions in $\beta _{y}$ and $\beta _{cy}$ are $0.68 \times (1 - 0.65) + 0.32 = 0.558$. The average reductions in $\beta _{s}$ and $\beta _{cs}$ are $0.787 \times (1 - 0.65) + 0.213 = 0.488$. Consequently, due to vaccinations in the current situation, infection rates in this university model with mixed vaccinated/unvaccinated population becomes: $\beta _{y} = 0.0910$, $\beta _{cy} = 0.0954$, $\beta _{s} = 0.1098$, and $\beta _{cs} = 0.1152$. Since vaccines can also reduce the probability of symptomatic infections and of severe illness⁴², other parameter values are also tuned accordingly. This adjustment changed the $R_{0}$ from 2.50 (totally unvaccinated and no interventions) to 1.40 (mixed population but still no intervention). Since it is still greater than 1, the COVID-19 epidemic will still develop in the university model if no other control or prevention measures are introduced. Accordingly, the derived parameter values of the model before and after vaccinations are summarised in Table 2.

Table 2 Model parameter values.

Analysis of the uncontrolled university model

Equilibrium points

Denote $x = {[}S_{y},\;S_{s},\;E_{y},\;E_{s},\;I_{y},\;I_{s},\;Q_y,\;Q_s,\;R_y,\;R_s,\;C]^{\top }$. By equating all derivatives in (1) to zero, it is easy to determine the equilibria ${\bar{x}} = [{\bar{S}}_{y}, {\bar{S}}_{s}, 0, 0, 0, 0, 0, 0, {\bar{R}}_{y}, {\bar{R}}_{s}, 0]$, where

$${\bar{S}_{y}} + {\bar{S}}_{s} + {\bar{R}}_{y} + {\bar{R}}_{s} = 1, \quad {\bar{S}}_{y} \ge 0, \quad {\bar{S}}_{s} \ge 0, \quad {\bar{R}}_{y} \ge 0, \quad {\bar{R}}_{s} \ge 0$$

(2)

These equilibria imply that at the end of pandemic, the individuals are either susceptible to or recovered from the disease.

Basic reproduction number

The basic reproduction number is a crucial criterion to measure the average number of susceptible people that could potentially be infected by a primary case⁴³. This parameter is highly dependent on the fraction of the susceptible population and it provides information about the potential of the epidemic outbreak. If $R_0 < 1$ the disease will gradually disappear. If $R_0 > 1$ increasingly more people will be infected.

Derivation of the basic reproduction number for the uncontrolled university model is based on the next generation matrix method described by^44,45,46. The university model (1) has five infectious compartments: $E_{y}$, $E_{s}$, $I_{y}$, $I_{s}$ and C. We collect these in the infection state $x_{if} = [E_{y}, E_{s}, I_{y}, I_{s}, C]^{\top }$. Let ${\mathcal {F}}$ denote the rate of increase of secondary cases and ${\mathcal {V}}$ denote the progression rate. Accordingly, $x_{if}$ obeys the equation

$${\dot{x}}_{if} = {\mathcal {F} }- {\mathcal {V} }$$

(3)

where ${\mathcal {F}}$ and ${\mathcal {V}}$ given by

$$\begin{aligned} {\mathcal {F}} = \left[ \begin{array}{c} S_y \,{\left( C\,\beta _{cy} +\beta _{{y}} \,{\left( E_s +E_y +I_s \,k+I_y \,k\right) }\right) }\\ S_s \,{\left( C\,\beta _{cs} +\beta _{s} \,{\left( E_s +E_y +I_s \,k+I_y \,k\right) }\right) }\\ 0\\ 0\\ \mu _s \,{\left( E_s +I_s \,k\right) }+\mu _y \,{\left( E_y +I_y \,k\right) } \end{array}\right] \end{aligned}$$

(4)

and

$$\begin{aligned} {\mathcal {V}} = \left[ \begin{array}{c} E_y \,{\left( \varepsilon _y +\xi _y \right) }\\ E_s \,{\left( \varepsilon _s +\xi _s \right) }\\ I_y \,(\eta _y + \rho _{y}) -E_y \,\varepsilon _y \\ I_s \,(\eta _s + \rho _{s}) -E_s \,\varepsilon _s \\ C\,\delta \end{array}\right] . \end{aligned}$$

(5)

We linearise Eq. (3) around the equilibrium and we denote the Jacobians of ${\mathcal {F}}$ and ${\mathcal {V}}$ as F and V, respectively. Thus, we obtain

$${\dot{x}}_{if} = (F - V)x_{if}$$

(6)

where

$$\begin{aligned} F = \left[ \begin{array}{ccccc} {{\bar{S}} }_y \,\beta _{y} &{} {\bar{S} }_y \,\beta _{y} &{} {\bar{S} }_y \,\beta _{y} \,k &{} {\bar{S} }_y \,\beta _{y} \,k &{} {\bar{S} }_y \,\beta _{cy} \\ {\bar{S} }_s \,\beta _{y} \,p &{} {\bar{S} }_s \,\beta _{y} \,p &{} {\bar{S} }_s \,\beta _{y} \,k\,p &{} {\bar{S} }_s \,\beta _{y} \,k\,p &{} {\bar{S} }_s \,\beta _{cy} \,p\\ 0 &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} 0 &{} 0 &{} 0 &{} 0\\ \mu _y &{} \mu _s &{} k\,\mu _y &{} k\,\mu _s &{} 0 \end{array}\right] \end{aligned}$$

(7)

and

$$\begin{aligned} V = \left[ \begin{array}{ccccc} \varepsilon _y +\xi _y &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} \varepsilon _s +\xi _s &{} 0 &{} 0 &{} 0\\ -\varepsilon _y &{} 0 &{} \eta _y + \rho _{y} &{} 0 &{} 0\\ 0 &{} -\varepsilon _s &{} 0 &{} \eta _s + \rho _{s} &{} 0\\ 0 &{} 0 &{} 0 &{} 0 &{} \delta \end{array}\right] . \end{aligned}$$

(8)

Denoting $\sigma _y = \eta _y + \rho _y$ and $\sigma _s = \eta _s + \rho _s$, the next generation matrix K is therefore derived as

$$\begin{aligned} K&= FV^{-1}\\&= \left[ \begin{array}{ccccc} \frac{{\bar{S} }_y \,\beta _{y} \,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} &{} \frac{{\bar{S} }_y \,\beta _{y} \,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }} &{} \frac{{\bar{S} }_y \,\beta _{y} \,k}{\sigma _y } &{} \frac{{\bar{S} }_y \,\beta _{y} \,k}{\sigma _s } &{} \frac{{\bar{S} }_y \,\beta _{cy} }{\delta }\\ \frac{{\bar{S} }_s \,\beta _{y} \,p\,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} &{} \frac{{\bar{S} }_s \,\beta _{y} \,p\,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }} &{} \frac{{\bar{S} }_s \,\beta _{y} \,k\,p}{\sigma _y } &{} \frac{{\bar{S} }_s \,\beta _{y} \,k\,p}{\sigma _s } &{} \frac{{\bar{S} }_s \,\beta _{cy} \,p}{\delta }\\ 0 &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} 0 &{} 0 &{} 0 &{} 0\\ \frac{\mu _y \,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} &{} \frac{\mu _s \,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }} &{} \frac{k\,\mu _y }{\sigma _y } &{} \frac{k\,\mu _s }{\sigma _s } &{} 0 \end{array}\right] . \end{aligned}$$

(9)

The obtained matrix K is nonnegative and has rank 2. In particular, it has three zero eigenvalues and two positive eigenvalues. According to⁴⁴, $R_{0}$ is the spectral radius of K, i.e. its largest eigenvalue. By computing $\text {det}(\lambda I - K)$, the characteristic polynomial is

$$\begin{aligned} p_k(\lambda )&= \lambda ^3 \left[ \lambda ^2 - \left( \frac{{\bar{S} }_y \,\beta _{y} \,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} + \frac{{\bar{S} }_s \,\beta _{y} \,p\,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }}\right) \lambda \right. \\&\quad - \left. \frac{1}{\delta }\left( \frac{{\bar{S}}_y \,\beta _{cy}\,\mu _y{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} + \frac{{\bar{S} }_s \,\beta _{cy} \,p\,\mu _s{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }}\right) \right] . \end{aligned}$$

(10)

We recall that we have assumed that $\mu _y = \mu _s = \mu$, which means that students and staff have equal rates of spreading the virus into the environment. The two non-zero eigenvalues can be derived by finding the roots of the polynomial

$$\begin{aligned} \tilde{p}_{k}(\lambda )&= \lambda ^2 - \beta _{y} \left( \frac{{\bar{S} }_y \,\left( {\sigma _y} +\varepsilon _y \,k\right) }{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} + \frac{{\bar{S} }_s p\,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }}\right) \lambda \\&\quad - \frac{\mu \,\beta _{cy}}{\delta }\left( \frac{{\bar{S} }_y \,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} + \frac{{\bar{S} }_s p\,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }}\right) . \end{aligned}$$

(11)

Denote

$$a = \beta _{y}\left( \frac{{\bar{S} }_y \,{\left( \sigma _y +\varepsilon _y \,k\right) }}{\sigma _y \,{\left( \varepsilon _y +\xi _y \right) }} + \frac{{\bar{S} }_s p\,{\left( \sigma _s +\varepsilon _s \,k\right) }}{\sigma _s \,{\left( \varepsilon _s +\xi _s \right) }}\right) ,\qquad c = \frac{\mu \beta _{cy}}{\delta \beta _{y}}.$$

(12)

Then

$$\tilde{p}_{k}(\lambda ) = \lambda ^2 - a\lambda - ac$$

(13)

Since the constants a and c are both positive, the quadratic equation has two real roots and $R_{0}$ will be the larger one. Therefore, we can express $R_{0}$ as

$$R_{0} = \frac{a(1+\sqrt{1+4c/a})}{2}.$$

(14)

In the next section we show how $R_0$ is related to the stability of the equilibrium point.

Stability analysis

Proposition 1

If $R_{0} < 1$, the equilibrium points ${\bar{x}}$ of the uncontrolled university model (1) is asymptotically stable.

Proof

Model (1) can be reformulated into a feedback interconnection. Compartments $E_y$, $E_s$, $I_y$, $I_s$, $Q_y$, $Q_s$, C form a positive linear subsystem with output feedback topology. Defining $x_{l} = [E_y, E_s, I_y, I_s, Q_y, Q_s, C]^\top$, $y_s = C_{s}x_l(t)$, $y_R = C_{R}x_l(t)$, the subsystem can be formulated as

$$\begin{aligned} {\dot{x}}_{l}(t) &= Ax_{l}(t) + Bu(t)\\ &= \left[ \begin{array}{ccccccc} -r_y &{} 0 &{} 0 &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} -r_s &{} 0 &{} 0 &{} 0 &{} 0 &{} 0\\ \varepsilon _y &{} 0 &{} -\sigma _y &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} \varepsilon _s &{} 0 &{} -\sigma _s &{} 0 &{} 0 &{} 0\\ 0 &{} 0 &{} \eta _y &{} 0 &{} -\varphi _y &{} 0 &{} 0\\ 0 &{} 0 &{} 0 &{} \eta _s &{} 0 &{} -\varphi _s &{} 0\\ \mu _y &{} \mu _s &{} k\mu _y &{} k\mu _s &{} 0 &{} 0 &{} -\delta \end{array}\right] x_{l}(t) + \left[ \begin{array}{ccc} 1 &{} 0 \\ 0 &{} 1 \\ 0 &{} 0 \\ 0 &{} 0 \\ 0 &{} 0 \\ 0 &{} 0 \\ 0 &{} 0 \end{array}\right] u(t) \end{aligned}$$

(15)

$$\begin{aligned} y_s(t) = C_{s}x_l(t) = \left[ \begin{array}{ccccccc} \beta _{y} &{} \beta _{y} &{} \beta _{y} \,k &{} \beta _{y} \,k &{} 0 &{} 0 &{} \beta _{cy} \\ \beta _{s} &{} \beta _{s} &{} \beta _{s} \,k &{} \beta _{s} \,k &{} 0 &{} 0 &{} \beta _{cs} \end{array}\right] x_l(t) \end{aligned}$$

(16)

$$\begin{aligned} y_R(t) = C_{R}x_l(t) = \left[ \begin{array}{ccccccc} \xi _y &{} 0 &{} \rho _y &{} 0 &{} \varphi _y &{} 0 &{} 0 \\ 0 &{} \xi _s &{} 0 &{} \rho _s &{} 0 &{} \varphi _s &{} 0 \\ \end{array}\right] x_l(t) \end{aligned}$$

(17)

$\square$

where $r_y = \varepsilon _y + \eta _y$, $r_s = \varepsilon _s + \eta _s$. Since output $y_{R}$ does not contribute to the variations of the state variable $x_{l}$, we can represent the overall system dynamics by the output feedback loop as

$$u(t) = K_s(t)y_s(t) = \left[ \begin{array}{ccc} S_y(t) &{} 0 \\ 0 &{} S_s(t) \\ \end{array}\right] y_s(t).$$

(18)

Derivatives of the remaining compartments can be calculated as

$$\left[ \begin{array}{c} {\dot{S}}_y(t)\\ {\dot{S}}_s(t)\\ \end{array}\right] = -\left[ \begin{array}{ccc} S_y(t) &{} 0 \\ 0 &{} S_s(t) \\ \end{array}\right] y_s(t)$$

(19)

$$\left[ \begin{array}{c} {\dot{R}}_y(t)\\ {\dot{R}}_s(t)\\ \end{array}\right] = y_R(t).$$

(20)

Note that the system has a time-varying feedback $K_s(t)$. Around the equilibrium ${\bar{x}}$, the system behaviour is determined by using the constant feedback term ${\bar{K}}_{s}$ as

$${\bar{K}}_{s} = \left[ \begin{array}{ccc} {\bar{S}}_y &{} 0 \\ 0 &{} {\bar{S}}_s \\ \end{array}\right] y_s(t).$$

(21)

Therefore, the dynamics of the original model is equivalent to that of this closed-loop system. To study its stability, we firstly derive the closed-loop system matrix $A_{cl}$ as

$$\begin{aligned} A_{cl}&= A + BK_{s}C_{s}\\&= \left[ \begin{array}{ccccccc} {\bar{S} }_y \,\beta _{y} -r_y &{} {\bar{S} }_y \,\beta _{y} &{} {\bar{S} }_y \,\beta _{y} \,k &{} {\bar{S} }_y \,\beta _{y} \,k &{} 0 &{} 0 &{} {\bar{S} }_y \,\beta _{cy} \\ {\bar{S} }_s \,\beta _{y} \,p &{} {\bar{S} }_s \,\beta _{y} \,p-r_s &{} {\bar{S} }_s \,\beta _{y} \,k\,p &{} {\bar{S} }_s \,\beta _{y} \,k\,p &{} 0 &{} 0 &{} {\bar{S} }_s \,\beta _{cy} \,p\\ \varepsilon _y &{} 0 &{} -\sigma _y &{} 0 &{} 0 &{} 0 &{} 0\\ 0 &{} \varepsilon _s &{} 0 &{} -\sigma _s &{} 0 &{} 0 &{} 0\\ 0 &{} 0 &{} \eta _y &{} 0 &{} -\varphi _y &{} 0 &{} 0\\ 0 &{} 0 &{} 0 &{} \eta _s &{} 0 &{} -\varphi _s &{} 0\\ \mu _y &{} \mu _s &{} k\mu _y &{} k\mu _s &{} 0 &{} 0 &{} -\delta \end{array}\right] . \end{aligned}$$

(22)

To determine its closed-loop poles, we compute its characteristic equation, i.e. $\text {det}(\lambda I - A_{cl})$. Since all parameters are positive, $A_{cl}$ must have two negative eigenvalues at $-\varphi _y$ and $-\varphi _s$. The remaining five eigenvalues are the roots of the polynomial $p_5(\lambda )$:

$$\begin{aligned} p_{5}(\lambda )&= \left\| \begin{array}{ccccc} \lambda +r_y -{\bar{S} }_y \,\beta _{y} &{} -{\bar{S} }_y \,\beta _{y} &{} -{\bar{S} }_y \,\beta _{y} \,k &{} -{\bar{S} }_y \,\beta _{y} \,k &{} -{\bar{S} }_y \,\beta _{cy} \\ -\bar{p} \,{\left( \lambda +r_y \right) } &{} \lambda +r_s &{} 0 &{} 0 &{} 0\\ -\varepsilon _y &{} 0 &{} \lambda +\sigma _y &{} 0 &{} 0\\ 0 &{} -\varepsilon _s &{} 0 &{} \lambda +\sigma _s &{} 0\\ -\mu _y &{} -\mu _s &{} -k\,\mu _y &{} -k\,\mu _s &{} \delta +\lambda \end{array}\right\| \\&= \left\| \begin{array}{ccccc} \lambda +r_y -{\bar{S} }_y \,\beta _{y} &{} -\lambda -r_y &{} 0 &{} 0 &{} -{\bar{S} }_y \,\beta _{cy} \\ -\bar{p} \,{\left( \lambda +r_y \right) } &{} \lambda +r_s +\bar{p} \,{\left( \lambda +r_y \right) } &{} -k\,{\left( \lambda +r_s \right) } &{} 0 &{} 0\\ -\varepsilon _y &{} \varepsilon _y &{} \lambda +\sigma _y &{} -\lambda -\sigma _y &{} 0\\ 0 &{} -\varepsilon _s &{} \varepsilon _s \,k &{} \lambda +\sigma _s &{} 0\\ -\mu &{} 0 &{} 0 &{} 0 &{} \delta +\lambda \end{array}\right\| \end{aligned}$$

(23)

where we have defined ${\bar{p}} = \frac{{\bar{S}}_{s}p}{{\bar{S}}_y}$ and used that $\mu _y = \mu _s = \mu$ according to the previous analysis. The polynomial is finally derived as

$$p_5(\lambda ) = \lambda ^5 + \alpha _4\lambda ^4 + \alpha _3\lambda ^3 + + \alpha _2\lambda ^2 + \alpha _1\lambda + \alpha _0,$$

(24)

where

$$\begin{aligned} \alpha _4&= \delta +r_s +r_y +\sigma _s +\sigma _y -{\bar{S} }_y \,\beta _{y} -{\bar{S} }_s \,\beta _{y} \,p \\ \alpha _3&= \delta \,r_s +\delta \,r_y +\delta \,\sigma _s +\delta \,\sigma _y +r_s \,r_y +r_s \,\sigma _s +r_s \,\sigma _y +r_y \,\sigma _s +r_y \,\sigma _y +\sigma _s \,\sigma _y\\&\quad -{\bar{S} }_y \,\beta _{y} \,\delta -{\bar{S} }_y \,\beta _{cy} \,\mu -{\bar{S} }_y \,\beta _{y} \,r_s -{\bar{S} }_y \,\beta _{y} \,\sigma _s -{\bar{S} }_y \,\beta _{y} \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\varepsilon _y \,k\\&\quad -{\bar{S} }_s \,\beta _{y} \,\delta \,p-{\bar{S} }_s \,\beta _{cy} \,\mu \,p-{\bar{S} }_s \,\beta _{y} \,p\,r_y -{\bar{S} }_s \,\beta _{y} \,p\,\sigma _s -{\bar{S} }_s \,\beta _{y} \,p\,\sigma _y\\&\quad -{\bar{S} }_s \,\beta _{y} \,\varepsilon _s \,k\,p \\ \alpha _2&= \delta \,r_s \,r_y +\delta \,r_s \,\sigma _s +\delta \,r_s \,\sigma _y +\delta \,r_y \,\sigma _s +\delta \,r_y \,\sigma _y +\delta \,\sigma _s \,\sigma _y +r_s \,r_y \,\sigma _s +r_s \,r_y \,\sigma _y \\&\quad +r_s \,\sigma _s \,\sigma _y +r_y \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,r_s -{\bar{S} }_y \,\beta _{y} \,\delta \,\sigma _s -{\bar{S} }_y \,\beta _{y} \,\delta \,\sigma _y -{\bar{S} }_y \,\beta _{cy} \,\mu \,r_s\\&\quad -{\bar{S} }_y \,\beta _{cy} \,\mu \,\sigma _s -{\bar{S} }_y \,\beta _{cy} \,\mu \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,r_s \,\sigma _s -{\bar{S} }_y \,\beta _{y} \,r_s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\sigma _s \,\sigma _y\\&\quad -{\bar{S} }_y \,\beta _{y} \,\delta \,\varepsilon _y \,k-{\bar{S} }_y \,\beta _{cy} \,\varepsilon _y \,k\,\mu -{\bar{S} }_y \,\beta _{y} \,\varepsilon _y \,k\,r_s -{\bar{S} }_y \,\beta _{y} \,\varepsilon _y \,k\,\sigma _s\\&\quad -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,r_y -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,\sigma _s -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,\sigma _y -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,r_y -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,\sigma _s\\&\quad -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,\sigma _y -{\bar{S} }_s \,\beta _{y} \,p\,r_y \,\sigma _s -{\bar{S} }_s \,\beta _{y} \,p\,r_y \,\sigma _y -{\bar{S} }_s \,\beta _{y} \,p\,\sigma _s \,\sigma _y\\&\quad -{\bar{S} }_s \,\beta _{y} \,\delta \,\varepsilon _s \,k\,p-{\bar{S} }_s \,\beta _{cy} \,\varepsilon _s \,k\,\mu \,p-{\bar{S} }_s \,\beta _{y} \,\varepsilon _s \,k\,p\,r_y -{\bar{S} }_s \,\beta _{y} \,\varepsilon _s \,k\,p\,\sigma _y \\ \alpha _1&= \delta \,r_s \,r_y \,\sigma _s +\delta \,r_s \,r_y \,\sigma _y +\delta \,r_s \,\sigma _s \,\sigma _y +\delta \,r_y \,\sigma _s \,\sigma _y +r_s \,r_y \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,r_s \,\sigma _s\\&\quad -{\bar{S} }_y \,\beta _{y} \,\delta \,r_s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{cy} \,\mu \,r_s \,\sigma _s -{\bar{S} }_y \,\beta _{cy} \,\mu \,r_s \,\sigma _y \\&\quad -{\bar{S} }_y \,\beta _{cy} \,\mu \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,r_s \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,\varepsilon _y \,k\,r_s -{\bar{S} }_y \,\beta _{y} \,\delta \,\varepsilon _y \,k\,\sigma _s\\&\quad -{\bar{S} }_y \,\beta _{cy} \,\varepsilon _y \,k\,\mu \,r_s -{\bar{S} }_y \,\beta _{cy} \,\varepsilon _y \,k\,\mu \,\sigma _s -{\bar{S} }_y \,\beta _{y} \,\varepsilon _y \,k\,r_s \,\sigma _s -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,r_y \,\sigma _s\\&\quad -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,r_y \,\sigma _y -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,\sigma _s \,\sigma _y -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,r_y \,\sigma _s -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,r_y \,\sigma _y\\&\quad -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,\sigma _s \,\sigma _y -{\bar{S} }_s \,\beta _{y} \,p\,r_y \,\sigma _s \,\sigma _y -{\bar{S} }_s \,\beta _{y} \,\delta \,\varepsilon _s \,k\,p\,r_y -{\bar{S} }_s \,\beta _{y} \,\delta \,\varepsilon _s \,k\,p\,\sigma _y\\&\quad -{\bar{S} }_s \,\beta _{cy} \,\varepsilon _s \,k\,\mu \,p\,r_y -{\bar{S} }_s \,\beta _{cy} \,\varepsilon _s \,k\,\mu \,p\,\sigma _y -{\bar{S} }_s \,\beta _{y} \,\varepsilon _s \,k\,p\,r_y \,\sigma _y \\ \alpha _0&= \delta \,r_s \,r_y \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,r_s \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{cy} \,\mu \,r_s \,\sigma _s \,\sigma _y -{\bar{S} }_y \,\beta _{y} \,\delta \,\varepsilon _y \,k\,r_s \,\sigma _s\\&\quad -{\bar{S} }_y \,\beta _{cy} \,\varepsilon _y \,k\,\mu \,r_s \,\sigma _s -{\bar{S} }_s \,\beta _{y} \,\delta \,p\,r_y \,\sigma _s \,\sigma _y -{\bar{S} }_s \,\beta _{cy} \,\mu \,p\,r_y \,\sigma _s \,\sigma _y\\&\quad -{\bar{S} }_s \,\beta _{y} \,\delta \,\varepsilon _s \,k\,p\,r_y \,\sigma _y -{\bar{S} }_s \,\beta _{cy} \,\varepsilon _s \,k\,\mu \,p\,r_y \,\sigma _y \end{aligned}$$

(25)

The conditions to obtain a stable equilibrium point can be determined using the Routh-Hurwitz stability criterion. The Routh table is given in Table 3.

Table 3 Routh table.

From the Routh-Hurwitz stability criterion follows that the equilibrium point is asymptotically stable if and only if $\alpha _4 >0$, $b_{31}>0$, $b_{21}>0$, $b_{11}>0$ and $\alpha _0>0$. We evaluate these coefficients numerically for different values of $R_0$. The results are shown in Table 4 for $R_0<1$ and in Table 5 for $R_0>1$. Since $\alpha _0 > 0$ in Table 4 and $\alpha _0<0$ in Table 5, the Routh-Hurwitz stability criterion confirms that the $R_0$ found in (14) is consistent with the expected epidemic dynamics.

Table 4 Values of the first column of the Routh table when $R_{0} < 1$.

Table 5 Values of the first column of the Routh table when $R_{0} > 1$.

Control of the university model

Formulation of the controlled university model

In this study we consider five non-pharmaceutical interventions that the university can implement.

• 1. Compulsory mask wearing: how strongly this measure is implemented is represented by the normalised variable $0\le \kappa _{m} \le 1$.

• 2. Keep safe social distance: how strongly this measure is implemented is represented by the normalised variable $0\le \kappa _{d} \le 1$.

• 3. Environment disinfection: how strongly this measure is implemented is represented by the normalised variable $0\le \kappa _{e} \le 1$.

• 4. Mandatory quarantines: how strongly these measures are implemented is represented by the normalised variables $0\le \kappa _{qy} \le 1$ (for students) and $0\le \kappa _{qs} \le 1$ (for staff).

A group of five variables is initially defined to represent the reduction factors in the infection rates, shedding rates, environmental decaying rate, and isolation rates. Denoting these factors as $u = [u_{p},\;u_{m},\;u_{e},\;u_{qy},\;u_{qs}]^{\top }$, the university model at this stage is expressed as

$$\left\{ \begin{aligned} \frac{{\text{d}} S_{y}(t)}{{\text{d}} t}&= -\{(1-u_{p})\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)]\\&\qquad - (1-u_{m})\beta _{cy}C(t)\}S_{y}(t) \\ \frac{{\text{d}} S_{s}(t)}{{\text{d}} t}&= -\{(1-u_{p})\,p\,\beta _{y}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] \\&\qquad - (1-u_{m})\,p\,\beta _{cy}C(t)\}S_{s}(t) \\ \frac{{\text{d}} E_{y}(t)}{{\text{d}} t}&= \{(1-u_{p})\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)]\\&\qquad + (1-u_{m})\beta _{cy}C(t)\}S_{y}(t) - (\varepsilon _{y} + \xi _{y})E_{y}(t)\\ \frac{{\text{d}} E_{s}(t)}{{\text{d}} t}&= \{(1-u_{p})\,p\,\beta _{y}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] \\&\qquad + (1-u_{m})\,p\,\beta _{cy}C(t)\}S_{s}(t) - (\varepsilon _{s} + \xi _{s})E_{s}(t)\\ \frac{{\text{d}} I_{y}(t)}{{\text{d}} t}&= \varepsilon _{y}E_{y}(t) - u_{qy}\eta _{y}I_{y}(t)\\ \frac{{\text{d}} I_{s}(t)}{{\text{d}} t}&= \varepsilon _{s}E_{s}(t) - u_{qs}\eta _{s}I_{s}(t)\\ \frac{{\text{d}} Q_{y}(t)}{{\text{d}} t}&= u_{qy}\eta _{y}I_{y}(t) - \varphi _{y}Q_{y}(t)\\ \frac{{\text{d}} Q_{s}(t)}{{\text{d}} t}&= u_{qs}\eta _{s}I_{s}(t) - \varphi _{s}Q_{s}(t)\\ \frac{{\text{d}} R_y(t)}{{\text{d}} t}&= \xi _{y}E_{y}(t) + \varphi _{y}Q_{y}(t) \\ \frac{{\text{d}} R_s(t)}{{\text{d}} t}&= \xi _{s}E_{s}(t) + \varphi _{s}Q_{s}(t) \\ \frac{{\text{d}} C(t)}{{\text{d}} t}&= (1-u_{m})\,\mu \,(E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)) - u_{e}\,\delta \,C(t) \end{aligned}\right.$$

(26)

Note that $u_{p}$, $u_{m}$, $u_{e}$, $u_{qy}\eta _{y}$, and $u_{qs}\eta _{s}$ should vary within [0, 1].

Since a reduction factor can be influenced by multiple interventions, we need to identify the relationship between the reduction factors (u’s) and intervention variables ($\kappa$’s).

• 1. Reduction of interpersonal infection rates $\beta _{y}$ and $\beta _{s}$ The person-to-person infection rates are directly influenced by two measures: wearing masks and keeping social distance. Wearing masks could reduce the probability of infections during contacts while social distancing could reduce the number of direct contacts between individuals. Study conducted by Karaivanov et al.⁴⁷ argued that the mandatory mask-wearing policy in confined spaces could reduce the number of infected cases by up to $40\%$ weekly. Furthermore, Jarvis et al.⁴⁸ conducted a survey which showed that the physical distancing could reduce the number of direct contacts by $74\%$. However, since this survey might have selection and recall bias, the actual result should be lower than $74\%$ to obtain a conservative estimate. In this case, we set the maximum reduction at $65\%$. Hence, $1 - u_{p} = (1 - 0.4\kappa _{m})(1 - 0.65\kappa _{m})$.

• 2. Reduction of shedding rates and environment-to-person infection rates: $\mu _{y}$, $\mu _{s}$, $\beta _{cy}$, $\beta _{cs}$ Shedding rates as well as infections due to unclean environment can be reduced by the use of masks. According to laboratory-based investigations by^49,50, masks could block approximately $50\%$ to $70\%$ droplets and aerosol, greatly reducing the transmission of virus. In this case, we set the maximum reduction of shedding rate to be $60\%$ for a conservative estimate. Therefore, $1 - u_{m} = 1 - 0.6\kappa _{m}$.

• 3. Enhancement of virus environmental decaying rate $\delta$ The environmental decaying rate could be magnified by disinfection of the confined space. Recall that the uncontrolled decaying rate was $\delta =0.7$. It is hard to determine how much this rate is increased. A conservative estimate is that the rate is increased at the maximum by $30\%$, yielding a maximum new rate $u_{e}\delta = 0.91$. The linear relationship can be finally defined as $u_{e} = 1 + 0.3\kappa _{e}$.

While the analysis above is derived from considerations extracted from the literature, these still assumed an ideal enforcement of the interventions. Since the university cannot guarantee full compliance, we limit $\kappa _{m}$, $\kappa _{d}$, and $\kappa _{e}$ to $70\%$ of their values.

• 4. Quarantine enhancement Recall that in model (1) the quarantined populations were simply equivalent to the populations who developed serious symptoms and their rate were $\eta _{y}$ and $\eta _{s}$. To consider the effects of mandatory quarantines we replace $\eta _{y}$ and $\eta _{s}$ by $u_{qy}\eta _{y}$ and $u_{qs}\eta _{s}$, respectively. $(u_{qy}\eta _{y})^{-1}$ and $(u_{qs}\eta _{s})^{-1}$ now indicate the average time that unisolated symptomatic students/staff stay in campus before being detected by the university. We set these values to 2, meaning that the university takes two days in average to detect and isolate infected individuals after their symptoms develop. Thus, the maximum values of $(u_{qs}\eta _{s})^{-1}$ and $(u_{qs}\eta _{s})^{-1}$ are both 0.5. These values correspond to the situation where the enforcement of mandatory quarantines reaches the strongest degree, which means that $\kappa _{qy}$ and $\kappa _{qs}$ are 1. On the contrary, when mandatory quarantines are not implemented, i.e. $\kappa _{qy} = \kappa _{qs} = 0$, we want that the resulting quarantine rates $u_{qy}\eta _{y}$ and $u_{qs}\eta _{s}$ are still $\eta _{y}$ and $\eta _{s}$, respectively. These relations are formulated as $u_{qy}\eta _{y} = (0.5 - \eta _{y}) \kappa _{qy} + \eta _{y}$ and $u_{qs}\eta _{s} = (0.5 - \eta _{s}) \kappa _{qs} + \eta _{s}$.

In summary, the following relations hold

$$\begin{aligned} 1-u_{p}&= (1 - 0.4 \cdot 0.7\kappa _{m})(1 - 0.65 \cdot 0.7 \kappa _{d}) \\ 1-u_{m}&= 1 - 0.6 \cdot 0.7 \kappa _{m} \\ u_{e}&= 1 + 0.3 \cdot 0.7 \kappa _{e} \\ u_{qy}\eta _{y}&= (0.5 - \eta _{y}) \kappa _{qy} + \eta _{y} \\ u_{qs}\eta _{s}&= (0.5 - \eta _{s}) \kappa _{qs} + \eta _{s}. \end{aligned}$$

(27)

Substituting these relationships into Eq. (26), we obtain the final controlled university model

$$\left\{ \begin{aligned} \frac{{\text{d}} S_{y}(t)}{{\text{d}} t}&= -\{(1 - 0.28\kappa _{m})(1 - 0.455\kappa _{d})\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)] \\&\quad - (1 - 0.42\kappa _{m})\beta _{cy}C(t)\}S_{y}(t) \\ \frac{{\text{d}} S_{s}(t)}{{\text{d}} t}&= -\{(1 - 0.28\kappa _{m})(1 - 0.455\kappa _{d})\,p\,\beta _{y}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] \\&\quad - (1 - 0.42\kappa _{m})\,p\,\beta _{cy}C(t)\}S_{s}(t) \\ \frac{{\text{d}} E_{y}(t)}{{\text{d}} t}&= \{(1 - 0.28\kappa _{m})(1 - 0.455\kappa _{d})\beta _{y}[E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t)] \\&\quad + (1 - 0.42\kappa _{m})\beta _{cy}C(t)\}S_{y}(t) - (\varepsilon _{y} + \xi _{y})E_{y}(t) \\ \frac{{\text{d}} E_{s}(t)}{{\text{d}} t}&= \{(1 - 0.28\kappa _{m})(1 - 0.455\kappa _{d})\,p\,\beta _{y}\left[ E_{s}(t) + kI_{s}(t) + E_{y}(t) + kI_{y}(t)\right] \\&\quad + (1 - 0.42\kappa _{m})\,p\,\beta _{cy}C(t)\}S_{s}(t) - (\varepsilon _{s} + \xi _{s})E_{s}(t) \\ \frac{{\text{d}} I_{y}(t)}{{\text{d}} t}&= \varepsilon _{y}E_{y}(t) - [(0.5 - \eta _{y}) \kappa _{qy} + \eta _{y} + \rho _{y}]I_{y}(t) \\ \frac{{\text{d}} I_{s}(t)}{{\text{d}} t}&= \varepsilon _{s}E_{s}(t) - [(0.5 - \eta _{s}) \kappa _{qs} + \eta _{s} + \rho _{s}]I_{s}(t) \\ \frac{{\text{d}} Q_{y}(t)}{{\text{d}} t}&= [(0.5 - \eta _{y}) \kappa _{qy} + \eta _{y}]I_{y}(t) - \varphi _{y}Q_{y}(t) \\ \frac{{\text{d}} Q_{s}(t)}{{\text{d}} t}&= [(0.5 - \eta _{s}) \kappa _{qs} + \eta _{s}]I_{s}(t) - \varphi _{s}Q_{s}(t) \\ \frac{{\text{d}} R_y(t)}{{\text{d}} t}&= \xi _{y}E_{y}(t) + \rho _{y}I_{y}(t) + \varphi _{y}Q_{y}(t) \\ \frac{{\text{d}} R_s(t)}{{\text{d}} t}&= \xi _{s}E_{s}(t) + \rho _{s}I_{s}(t) + \varphi _{s}Q_{s}(t) \\ \frac{{\text{d}} C(t)}{{\text{d}} t}&= (1 - 0.42\kappa _{m})\,\mu \,(E_{y}(t) + kI_{y}(t) + E_{s}(t) + kI_{s}(t))\\&\quad - (1 + 0.21\kappa _{e})\,\delta \, C(t). \end{aligned}\right.$$

(28)

In this model, the variables $\kappa = [\kappa _{m},\; \kappa _{d},\; \kappa _{e},\; \kappa _{qy},\; \kappa _{qs}]^{\top }$ are the control variables that need to be optimised.

Formulation of the optimal control problem

Now we formulate the optimal control problem that, once solved, provides the optimal trajectories that give the best combination of interventions to contain the epidemic within the university. The optimal control problem is defined as

$$\begin{aligned} J^{*} = \min _{x(\cdot p ), \kappa (\cdot p )} \; J = \Vert x(t_{f})&-x_{tar}\Vert _{H}^{2} + \int ^{t_{f}}_{t_{0}} \left\| x(\tau )-x_{tar}\right\| _{Q}^{2} + \left\| \kappa _{k}(\tau )\right\| _{R}^{2} \;d\tau \\ {\text{s.t.}}\qquad \qquad \qquad \qquad {\dot{x}}(t)&= f(x(t), \kappa (t), t), \\ x_{\min }&\le x(t) \le x_{\max }, \\ x_{f,\min }&\le x(t_{f}) \le x_{f,\max }, \\ \kappa _{\min }&\le \kappa (t) \le \kappa _{\max }, \\ x(0)&= x_{14}, \end{aligned}$$

(29)

where $x_{tar}$ represents the target state vector, $\Vert \cdot \Vert _{H/Q/R}$ indicates the Euclidean matrix norm weighted by H/Q/R and the square matrices H, Q and R contain the weights for the final states, running states and running control variables, respectively. The values of $x_{tar}$ and the weights are changed depending on the scenario that we need to solve (see next section). The problem has also constraints on both states and control variables. In fact, we need that the states lie between 0 and 1. Thus, even without further requirement the state constraints are at least $x_{\min } = x_{f,\min } = {[}0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }$ and $x_{\max } = x_{f,\max } = [1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1]^{\top }$. Similarly, the lower bounds of the control variables should be $\kappa _{\min } = [0,\;0,\;0,\;0,\;0]^{\top }$ while their upper constrains are $\kappa _{\max } = [1,\;1,\;1,\;1,\;1]^{\top }$. The initial condition sets the starting point of the system. We assume that the university starts to react two weeks after the first exposed subjects appear among students or staff. Thus, we first run the simulation of the original un-controlled model for 14 days and get the resulting state values on day 14, namely $x_{14}$. Then this state is used as the initial condition for the optimisation problem, namely $x(0) = x_{14}$. We finally require that the epidemic is eliminated within 120 days, so $t_{0} = 0$ and $t_{f} = 120$.

Implementation of four different scenarios

Balancing the trade-off between controlling the spread of COVID-19 and resuming the normal campus activity is the main question considered in this study. According to different trade-off’s between these two objectives, four scenarios are studied: minimum number of cases, minimum intervention, minimum non-quarantine interventions, and minimum quarantine interventions. Each scenario corresponds to different weight matrices and different path constraints.

For the first scenario (minimum number of cases), the university does not impose any limitations on the strength of the interventions to lead to the fastest mitigation of the epidemic. Mathematically, the controller is designed to maximise the number of susceptible people during the whole period, i.e., the proportions of the susceptible students/staff are maintained as close to one as possible. In contrast, the values of the other compartments need to be as close as possible to zero in this scenario. Therefore, we set the target state as

$$x_{tar} = [1,\;1,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }.$$

(30)

Meanwhile, the weight matrices are set as

$$\begin{aligned} H&= diag([100;\;100;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0])\\ Q&= diag([100;\;100;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0;\; 0]) \\ R&= diag([10^{-2};\; 10^{-2};\; 10^{-2};\; 10^{-2};\; 10^{-2}]). \end{aligned}$$

(31)

where in the matrix R we use $10^{-2}$ instead of 0 to improve numerical stability of the solver. The constraints on both states and control variables remain the same as the basic requirements discussed before, namely

$$\begin{aligned} x_{\min } =&[0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }\\ x_{\max } =&[1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1]^{\top }\\ x_{f,\min } =&[0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }\\ x_{f,\max } =&[1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1]^{\top }\\ \kappa _{\min } =&[0,\;0,\;0,\;0,\;0]^{\top }\\ \kappa _{\max } =&[1,\;1,\;1,\;1,\;1]^{\top }. \end{aligned}$$

(32)

For the other three scenarios, we require that the number of susceptible students and staff cannot go below $94\%$ of the total population (i.e. the number of infected individuals is equal or below $6\%$). We can easily achieve this by setting more restrictive state constraints. Consequently, we do not need to use the weights H and Q (also because it is difficult to intuitively understand the meaning of the weights on the states in these scenarios). Hence, we select $H = 0$ and $Q = 0$ and $x_{tar}$ is not used. Additionally, to ensure that the epidemic is completely concluded after 134 days, the number of exposed and infected individuals should become zero at the final state. Recall that $N_{y}$, $N_{s}$, and $N_{t}$ indicate the total number of student, staff, and the total population, respectively. Then the constraints on states should be

$$\begin{aligned} x_{\min }&= {[}0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }\\ x_{\max }&= {[}1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1,\;1]^{\top }\\ x_{f,\min }&= {[}\frac{0.94N_y}{N_t},\;\frac{0.94N_s}{N_t},\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0,\;0]^{\top }\\ x_{f,\max }&= {[}1,\;1,\;10^{-4},\;10^{-4},\;10^{-4},\;10^{-4},\;1,\;1,\;1,\;1,\;1]^{\top }, \end{aligned}$$

(33)

while the input constraints remain unchanged

$$\begin{aligned} \kappa _{min}&= {[}0,\;0,\;0,\;0,\;0]^{\top }\\ \kappa _{max}&= {[}1,\;1,\;1,\;1,\;1]^{\top }. \end{aligned}$$

(34)

The only difference in the formulation of the other three scenarios is the value of the weight matrix R. In the minimum intervention scenario, the solver aims to minimize the total norm of the control signal while satisfying the $94\%$ constraint. Thus, all control variables are heavily weighted and the matrix R is selected as

$$R = diag({[}100;\; 100;\; 100;\; 100;\; 100]).$$

(35)

In the scenario of minimisation of non-quarantine interventions, the strength of first three interventions (wearing masks $\kappa _{m}$, keeping social distance $\kappa _{d}$ and disinfecting the environment $\kappa _{e}$) are minimised, while the quarantines are not. In this situation, the matrix R becomes

$$R = diag({[}100;\; 100;\; 100; \; 10^{-2}; \;10^{-2}]).$$

(36)

In the minimum quarantine scenario, the norm of quarantine rates for the infected students ($\kappa _{qy}$) and staff ($\kappa _{qs}$) are minimised. Since mandatory quarantine has the largest impact on campus activities, this scenario aims to find how the spread can be minimised while avoiding quarantine enforcement. Hence, the matrix R is

$$R = diag({[}10^{-2};\; 10^{-2};\; 10^{-2}; \;100; \;100]).$$

(37)