Optimization of vaccination for COVID-19 in the midst of a pandemic

1.   Introduction

In the last years the study of collective behavior of multi-agent systems has attracted the interest of many researchers in different scientific fields, such as biology, physics, control theory, social sciences, economics. The celebrated Cucker-Smale model has been proposed and analyzed in [21,22] to describe situations in which different agents, e.g. animals groups, reach a consensus (flocking), namely they align and move as a flock, based on a simple rule: each individual adjusts its velocity taking into account other agents' velocities.

In the original papers a symmetric interaction potential is considered. Then, the case of non-symmetric interactions has been studied by Motsch and Tadmor [31]. Several generalizations and variants have been introduced to cover various applications' fields, e.g. more general interaction rates and singular potentials [8,10,19,27,30,32], cone-vision constraints [40], presence of leadership [17,37], noise terms [20,23,25], crowds dynamics [18,29], infinite-dimensional models [1,2,11,26,28,39], control problems [3,5,16,33]. We refer to [6,12] for recent surveys on the Cucker-Smale type flocking models and its variants.

It is natural to introduce a time delay in the model, as a reaction time or a time to receive environmental information. The presence of a time delay makes the problem more difficult to deal with. Indeed, the time delay destroys some symmetry features of the model which are crucial in the proof of convergence to consensus. For this reason, in spite of a great amount of literature on Cucker-Smale models, only a few papers are available concerning Cucker-Smale model with time delay [13,14,15,23,36]. Cucker-Smale models with delay effects are also studied in [34,35] when a hierarchical structure is present, namely the agents are ordered in a specific order depending on which other agents they are leaders of or led by.

Here we consider a distributed delay term, i.e. we assume that the agent $i,$ $i = 1, \dots, N,$ changes its velocity depending on the information received from other agents on a time interval $[t-\tau(t) , t].$ Moreover, we assume normalized communication weights (cf [13]). Let us consider a finite number $N\in\mathbb{N}$ of autonomous individuals located in $\mathbb{R}^d$, $d\geq 1$. Let $x_i(t)$ and $v_i(t)$ be the position and velocity of $i$th individual. Then, in the current work, we will deal with a Cucker-Smale model with distributed time delays. Namely, let $\tau :[0, +\infty)\rightarrow (0, +\infty),$ be the time delay function belonging to $W^{1, \infty}(0, +\infty)$. Throughout this paper, we assume that the time delay function $\tau$ satisfies

for some positive constant $\tau_*.$ Then, denoting $\tau_0: = \tau(0),$ we have

It is clear that the constant time delay $\tau(t) \equiv \bar\tau > 0$ satisfies the conditions above.

Our main system is given by

where $\phi(x_k(s), x_i(t))$ are the normalized communication weights given by

with the influence function $\psi: [0, \infty) \to (0, \infty).$ Throughout this paper, we assume that the influence function $\psi$ is bounded, positive, nonincreasing and Lipschitz continuous on $[0, \infty)$, with $\psi(0) = 1$. Moreover, $\alpha: [0, \tau_0] \to [0, \infty)$ is a weight function satisfying

and

We consider the system subject to the initial datum

i.e., we prescribe the initial position and velocity trajectories $x^0_i, v^0_i\in\mathcal C([-\tau_0, 0]; \mathbb{R}^d)$.

For the particle system (3), we will first discuss the asymptotic behavior of solutions in Section 2. Motivated from [13,26,31], we derive a system of dissipative differential inequalities, see Lemma 2.4, and construct a Lyapunov functional. This together with using Halanay inequality enables us to show the asymptotic velocity alignment behavior of solutions under suitable conditions on the initial data.

We next derive, analogously to [13] where the case of a single pointwise time delay is considered, a delayed Vlasov alignment equation from the particle system (3) by sending the number of particles $N$ to infinity:

where $f_t = f_t(x, v)$ is the one-particle distribution function on the phase space $\mathbb{R}^d \times \mathbb{R}^d$ and the velocity alignment force $F$ is given by

We show the global-in-time existence and stability of measure-valued solutions to (6) by employing the Monge-Kantorowich-Rubinstein distance. As a consequence of the stability estimate, we discuss a mean-field limit providing a quantitative error estimate between the empirical measure associated to the particle system (3) and the measure-valued solution to (6). We then extend the estimate of large behavior of solutions for the particle system (3) to the one for the delayed Vlasov alignment equation (6). For this, we use the fact that the estimate of large-time behavior of solutions to the particle system (3) is independent of the number of particles. By combining this and the mean-field limit estimate, we show that the diameter of velocity-support of solutions of (6) converges to zero as time goes to infinity. Those results will be proved in Section 3.

2.   Asymptotic flocking behavior of the particle system (3)

We start with presenting a notion of flocking behavior for the system (3), and for this we introduce the spatial and, respectively, velocity diameters as follows:

Definition 2.1. We say that the system with particle positions $x_i(t)$ and velocities $v_i(t)$, $i = 1, \dots, N$ and $t\geq 0$, exhibits asymptotic flocking if the spatial and velocity diameters satisfy

We then state our main result in this section on the asymptotic flocking behavior of the system (3).

Theorem 2.2. Assume $N>2.$ Suppose that the initial data $x_i^0, v_i^0$ are continuous on the time interval $[-\tau_0, 0]$ and denote

Moreover, denoted $\beta_N = \frac {N-2}{N-1}$, assume that

where $d_X$ and $d_V$ denote, respectively, the spatial and velocity diameters defined in (7). Then the solution of the system (3)–(5) is global in time and satisfies

and

for a suitable positive constant $\gamma$ independent of $t$. The constant $\gamma$ can be chosen also independent of $N.$

Remark 1. If the influence function $\psi$ is not integrable, i.e., it has a heavy tail, then the right hand side of (9) becomes infinite, and thus the assumption (9) is satisfied for all initial data and time delay $\tau(t)$ satisfying (1). This is reminiscent of the unconditional flocking condition for the Cucker-Smale type models, see [11,12,13]. On the other hand, if the influence function $\psi$ is integrable and the weight function $\alpha(s)$ is given by $\alpha(s) = \delta_{\tau}(s)$, then the system (3) with a constant time delay $\tau(t)\equiv \tau$ becomes the one with discrete time delay:

Note the above system is studied in [13]. For this system, the assumption (9) reduces to

Since $d_V(z) \geq 0$, the left hand side of the above inequality increases as the size of time delay $\tau>0$ increases. On the other hand, the integrand on the right hand side decreases as $\tau$ increases. This asserts that small size of the time delay provides a larger set of initial data.

Remark 2. Observe that our theorem above gives a flocking result when the number of agents $N$ is greater than two. The result for $N = 2$ is trivial for our normalized model.

Remark 3. Note that $0 < \beta_N \uparrow 1$ as $N \to \infty.$ Then, this implies that for any $\beta \in (0, 1)$ there exists $N_* \in \mathbb{N}$ such that $\beta_N \geq \beta$ for $N \geq N_*$.

For the proof of Theorem 2.2, we will need several auxiliary results. Inspired by [13], we first show the uniform-in-time bound estimate of the maximum speed of the system (3).

Lemma 2.3. Let $R_v> 0$ be given by (8) and let $(\mathbf{x}, \mathbf{v})$ be a local-in-time $\mathcal{C}^1$-solution of the system (3)–(5). Then the solution is global in time and satisfies

Proof. Let us fix $\epsilon>0$ and set

By continuity, $S^\epsilon\ne\emptyset.$ Let us denote $T^\epsilon : = \sup S^{\epsilon} > 0.$ We want to show that $T^\epsilon = +\infty.$ Arguing by contradiction, let us suppose $T^\epsilon$ finite. This gives, by continuity,

From (3)–(5), for $t<T^\epsilon$ and $i = 1, \dots, N,$ we have

Note that

for $s \in [t-\tau(t), t]$ with $t < T^\epsilon$. Thus we obtain

which gives

From (11) we obtain

This is in contradiction with (10). Therefore, $T^\epsilon = +\infty\, .$ Being $\epsilon>0$ arbitrary, the claim is proved.

In the lemma below, motivated from [13,38] we derive the differential inequalities for $d_X$ and $d_V$. We notice that the diameter functions $d_X$ and $d_V$ are not $\mathcal C^1$ in general. Thus we introduce the upper Dini derivative to consider the time derivative of these functions: for a given function $F = F(t)$, the upper Dini derivative of $F$ at $t$ is defined by

Note that the Dini derivative coincides with the usual derivative when the function is differentiable at $t$.

Lemma 2.4. Let $R_v> 0$ be given by (8) and let $(\mathbf{x}, \mathbf{v})$ be the global $\mathcal{C}^1$-solution of the system (3)–(5) constructed in Lemma 2.3. Then, for almost all $t>0$, we have

Proof. The first inequality is by now standard, Then, we concentrate on the second one. Due to the continuity of the velocity trajectories $v_i(t),$ there is an at most countable system of open, mutually disjoint intervals $\{\mathcal{I}_\sigma\}_{\sigma\in\mathbb{N }}$ such that

and thus for each ${\sigma\in\mathbb{N }}$ there exist indices $i(\sigma)$, $j(\sigma)$ such that

Then, by using the simplified notation $i: = i(\sigma)$, $j: = j(\sigma)$ (we may assume that $i\ne j$), we have for every $t\in \mathcal{I}_\sigma$,

Set

Note that, from the definition (4) of $\phi (x_k(s), x_i(t)),$ it results $\bar\phi_{ij}(s, t)<1$ for all $i\ne j,$ $1\le i, j\le N$ and for all $t>0, \ s\in [t-\tau(t), t].$ Then, we get

where

Observe that $a_{ij}^k(s, t)\ge 0$ and $\sum_{k = 1}^N a_{ij}^k(s, t) = 1.$ Thus, if we consider the convex hull of a finite velocity point set $\{ v_1(t), \dots, v_N(t)\}$ and denote it by $\Omega (t),$ then

This gives

which, used in (13), implies

Now, by using the first equation in (3), we estimate for any $1 \leq i, k \leq N$ and $s \in [t-\tau (t), t]$,

Then, Lemma 2.3 gives

and due to the monotonicity property of the influence function $\psi,$ for $k\ne i,$ we deduce

On the other hand, we find

Then, from (15), we obtain

Using the last estimate in (14), we have

that, used in (12), concludes the proof.

Lemma 2.5. Let $u$ be a nonnegative, continuous and piecewise $\mathcal C^1$-function satisfying, for some constant $0 < a < 1$, the differential inequality

Then we have

with $\gamma = 1 - \tau_0^{-1}H(a\tau_0 \exp(\tau_0))$, where $H$ is the product logarithm function, i.e., $H$ satisfies $z = H(z) \exp(H(z))$ for any $z \geq 0$.

Proof. Note that the differential inequality (16) implies

Then, the result follows from Halanay inequality (see e.g. [24,p. 378]).

We are now ready to proceed with the proof of Theorem 2.2.

Proof of Theorem 2.2. For $t > 0$, we introduce the following Lyapunov functional for the system (3)–(5):

where $R_v$ is given by (8) and the diameters $d_X(t)$, $d_V(t)$ are defined in (7). We first estimate $\mathcal{L}_1$ as

for almost all $t \geq 0$, due to (1) and the first inequality in Lemma 2.4. Here $\mbox{sgn}\, (\cdot)$ is the signum function defined by

Analogously, we also get

for almost all $t \geq 0$. Then, from Lemma 2.4, (17) and (18), we have, for almost all $t>0$,

On the other hand, since $h^\prime(t) = \alpha(\tau(t))\tau^\prime(t) \leq 0$, we have $D^+\left (h(t)\mathcal{L}(t)\right )\le 0,$ namely

Moreover, it follows from the assumption (9) that there exists a positive constant $d_* > d_X(-s) + R_v\tau_0$ such that

This, together with (19) and (2), implies

Now, observe that, if $d_X(t-s)\geq d_X(-s),$ then

Similarly, when $d_X(t-s)<d_X(-s),$ then

Thus, from (20), (21) and (22) we deduce that

Note that, for $s \in [0, \tau(t)]$,

due to Lemma 2.3 and the first inequality of Lemma 2.4. Analogously, we also find for $s \in [0, \tau(t)]$ that

This gives

Thus we get

Combining this and (20), we obtain

Since the left hand side of the above inequality is positive, we have

We then again use (23) to find

for $s \in [0, \tau(t)]$ and $t \geq 0$. Hence, by Lemma 2.4 together with the monotonicity of $\psi$ we have

for almost all $t>0$, where $\psi_* = \psi(d_* + 4R_v\tau_0)$. We finally apply Lemma 2.5 to complete the proof. Note that, in order to have an exponential decay rate of $d_V$ independent of $N,$ it is sufficient to observe that $\beta_N\ge 1/2$ for $N\ge 3.$

3.   A delayed Vlasov alignment equation

In this section, we are interested in the behavior of solutions to the particle system (3) as the number of particles $N$ goes to infinity. At the formal level, we can derive the following delayed Vlasov alignment equation from (3) as $N \to \infty$:

for $(x, v) \in \mathbb{R}^d \times \mathbb{R}^d, \ t > 0,$ with the initial data:

where $f_t = f_t(x, v)$ is the one-particle distribution function on the phase space $\mathbb{R}^d \times \mathbb{R}^d$. Here the interaction term $F[f_s]$ is given by

For the equation (24), we provide the global-in-time existence and uniqueness of measure-valued solutions and mean-field limits from (3) based on the stability estimate. We also establish the large-time behavior of measure-valued solutions showing the velocity alignment.

3.1.   Existence of measure-valued solutions

In this part, we discuss the global existence and uniqueness of measure-valued solutions to the equation (24). For this, we first define a notion of weak solutions in the definition below.

Definition 3.1. For a given $T >0$, we call $f_t \in \mathcal C([0, T); \mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$ a measure-valued solution of the equation (24) on the time interval $[0, T)$, subject to the initial datum $g_s\in \mathcal C([-\tau_0, 0]; \mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$, if for all compactly supported test functions $\xi \in \mathcal C_c^\infty( \mathbb{R}^d \times \mathbb{R}^d \times [0, T))$,

where $\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d)$ denotes the set of probability measures on the phase space $\mathbb{R}^d \times \mathbb{R}^d$ with bounded first-order moment and we adopt the notation $f_{t-\tau_0}:\equiv g_{t-\tau_0}$ for $t \in [0, \tau_0]$.

We next introduce the $1$-Wasserstein distance.

Definition 3.2. Let $\rho^1, \rho^2\in \mathcal{P}_1( \mathbb{R}^{d})$ be two probability measures on $\mathbb{R}^{d}$. Then $1$-Waserstein distance between $\rho^1$ and $\rho^2$ is defined as

where $\Pi(\rho^1, \rho^2)$ represents the collection of all probability measures on $\mathbb{R}^d \times \mathbb{R}^d$ with marginals $\rho^1$ and $\rho^2$ on the first and second factors, respectively.

Theorem 3.3. Let the initial datum $g_t \in \mathcal C([-\tau_0, 0]; \mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$ and assume that there exists a constant $R > 0$ such that

where $B^{2d}(0, R)$ denotes the ball of radius $R$ in $\mathbb{R}^d \times \mathbb{R}^d$, centered at the origin.

Then for any $T > 0$, the delayed Vlasov alignment equation (24) admits a unique measure-valued solution $f_t\in \mathcal C([0, T); \mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$ in the sense of Definition 3.1.

Proof. The proof can be done by using a similar argument as in [13,Theorem 3.1], thus we shall give it rather concisely. Let $f_t \in \mathcal C([0, T];\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$ be such that

for some positive constant $R > 0$. Then we use the similar estimate as in [13,Lemma 3.1] to have that there exists a constant $C > 0$ such that

for $(x, v) \in B^{2d}(0, R), \; t\in [0, T]$. Subsequently, this gives

and

for $(x, v) \in B^{2d}(0, R), \; t\in [0, T]$. This, together with an application of [4,Theorem 3.10], provides the local-in-time existence and uniqueness of measure-valued solutions to the equation (24) in the sense of Definition 3.1. On the other hand, once we obtain the growth estimates of support of $f_t$ in both position and velocity, these local-in-time solutions can be extended to the global-in-time ones. Thus in the rest of this proof, we provide the support estimate of $f_t$ in both position and velocity. We set

for $t \in [0, T]$, where supp$_x f_t$ and supp$_v f_t$ represent $x$- and $v$-projections of supp$f_t$, respectively. We also set

We first construct the system of characteristics

associated with (24),

where we again adopt the notation $f_{t-\tau_0}:\equiv g_{t-\tau_0}$ for $t \in [0, \tau_0]$. The system (26) is considered subject to the initial conditions

for all $(x, v)\in \mathbb{R}^d \times \mathbb{R}^d$. Note that the characteristic system (26)–(27) is well-defined since the $F[f_s]$ is locally both bounded and Lipschitz. We can rewrite the equation for $V$ as

Then, arguing as in the proof of Lemma 2.3, we get

due to (25). Using again a similar argument as in the proof of Lemma 2.3 and the comparison lemma, we obtain

which further implies $R^t_X \leq R^0_X + t R^0_V$ for $t \geq 0$. This completes the proof.

3.2.   Stability & mean-field limit

In this subsection, we discuss the rigorous derivation of the delayed Vlasov alignment equation (24) from the particle system (3) as $N \to \infty$. For this, we first provide the stability of measure-valued solutions of (24).

Theorem 3.4. Let $f^i_t \in \mathcal C([0, T);\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$, $i = 1, 2$, be two measure-valued solutions of (24) on the time interval $[0, T]$, subject to the compactly supported initial data $g^i_s\in \mathcal C([-\tau_0, 0];\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$. Then there exists a constant $C = C(T)$ such that

Proof. Again, the proof is very similar to [13,Theorem 3.2], see also [7,9]. Indeed, we can obtain

Then we have

for $t \in [0, T)$.

Remark 4. Since the empirical measure

associated to the $N$-particle system (3) with the initial data

for $s\in[-\tau_0, 0],$ is a solution to the equation (24), we can use Theorem 3.4 to have the following mean-field limit estimate:

where $C$ is a positive constant independent of $N$.

3.3.   Asymptotic behavior of the delayed Vlasov alignment equation (24)

In this part, we provide the asymptotic behavior of solutions to the equation (24) showing the velocity alignment under suitable assumptions on the initial data. For this, we first define the position- and velocity-diameters for a compactly supported measure $g\in\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d)$,

where supp$_x f$ denotes the $x$-projection of supp$f$ and similarly for supp$_v f$.

Theorem 3.5. Let $T>0$ and $f_t \in \mathcal C([0, T);\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$ be a weak solution of (24) on the time interval $[0, T)$ with compactly supported initial data $g_s\in \mathcal C([-\tau_0, 0];\mathcal{P}_1( \mathbb{R}^d \times \mathbb{R}^d))$. Furthermore, we assume that

Then the weak solution $f_t$ satisfies

where $C$ is a positive constant independent of $t$.

Let us point out that the flocking estimate at the particle level, see Section 2 and Remark 3, is independent of the number of particles, thus we can directly use the same strategy as in [11,13,26]. However, we provide the details of the proof for the completeness.

Proof of Theorem 3.5. We consider an empirical measure $\{g^N_s\}_{N\in\mathbb{N }}$, which is a family of $N$-particle approximations of $g_s$, i.e.,

where the $x_i^0, v_i^0 \in \mathcal C([-\tau_0, 0]; \mathbb{R}^d)$ are chosen such that

Note that we can choose $(x_i^0, v_i^0)_{i = 1, \dots, N}$ such that the condition (9) is satisfied uniformly in $N\in\mathbb{N}$ due to the assumption (28). Let us denote by $(x^N_i, v^N_i)_{i = 1, \dots, N}$ the solution of the $N$-particle system (3)–(5) subject to the initial data $(x_i^0, v_i^0)_{i = 1, \dots, N}$ constructed above. Then it follows from Theorem 2.2 that there exists a positive constant $C_1>0$ such that

with the diameters $d_V$, $d_X$ defined in (7), where $C_1>0$ is independent of $t$ and $N$. Note that the empirical measure

is a measure-valued solution of the delayed Vlasov alignment equation (24) in the sense of Definition 3.1. On the other hand, by Theorem 3.4, for any fixed $T>0$, we have the following stability estimate

where the constant $C_2>0$ is independent of $N$. This yields that sending $N\to\infty$ implies $d_V[f_t] = d_V(t)$ on $[0, T)$ for any fixed $T>0$. Thus we have

Since the uniform-in-$t$ boundedness of $d_X[f_t]$ just follows from the above exponential decay estimate, we conclude the desired result.