Mathematical Problems in Engineering
Volume 2015 (2015), Article ID 850124, 14 pages
http://dx.doi.org/10.1155/2015/850124

Research Article

Uncertainty Quantification in Control Problems for Flocking Models

Giacomo Albi,¹ Lorenzo Pareschi,² and Mattia Zanella²

¹Fakultät für Mathematik, Technische Universität München, Boltzmannstraße 3, Garching, 85748 Munich, Germany
²Department of Mathematics and Computer Science, University of Ferrara, Via N. Machiavelli 35, 44121 Ferrara, Italy

Received 19 January 2015; Revised 27 April 2015; Accepted 28 April 2015

Academic Editor: Angel Sánchez

Copyright © 2015 Giacomo Albi et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

The optimal control of flocking models with random inputs is investigated from a numerical point of view. The effect of uncertainty in the interaction parameters is studied for a Cucker-Smale type model using a generalized polynomial chaos (gPC) approach. Numerical evidence of threshold effects in the alignment dynamic due to the random parameters is given. The use of a selective model predictive control permits steering of the system towards the desired state even in unstable regimes.

1. Introduction

The aggregate motion of a multiagent system is frequently seen in the real world. Common examples are schools of fishes, swarms of bees and herds of sheep, natural phenomena that inspired important applications in many fields such as biology, engineering and economy [1]. As a consequence, the significance of new mathematical models, for understanding and predicting these complex dynamics, is widely recognized. Several heuristic rules of flocking have been introduced as alignment, separation, and cohesion [2, 3]. Nowadays these mathematical problems, and their constrained versions, are deeply studied from both the microscopic viewpoint [4–7] and their kinetic and mean-field approximations [8–13]. We refer to [1] for a recent introduction on the subject.

In an applicative framework a fundamental step for the study of such models is represented by the introduction of stochastic parameters reflecting the uncertainty due to wide range of phenomena, such as the weathers influence during an experiment, temperature variations, or even human errors. Therefore quantifying the influence of uncertainties on the main dynamics is of paramount importance to build more realistic models and to give better predictions of their behavior. In the modeling of self-organized system, different ways to include random sources have been studied and analyzed; see, for example, [3, 14–17]. In this paper we focus on the case where the uncertainty acts directly in the parameter characterizing the interaction dynamic between the agents.

We present a numerical approach having roots in the numerical techniques for uncertainty quantification (UQ) and model predictive control (MPC). Among the most popular methods for UQ, the generalized polynomial chaos (gPC) has recently received deepest attentions [18]. Jointly with Stochastic Galerkin (SG) this class of numerical methods is usually applied in physical and engineering problems, for which fast convergence is needed. Applications of gPC-Galerkin schemes to flocking dynamics, and their controlled versions, are almost unexplored in the actual state of art.

We give numerical evidence of threshold effects in the alignment dynamic due to the random parameters. In particular the presence of a negative tail in the distribution of the random inputs leads to the divergence of the expected values for the system velocities. The use of a selective model predictive control permits steering of the system towards the desired state even in such unstable regimes.

The rest of the paper is organized as follows. In Section 2 we introduce briefly a Cucker-Smale dynamic with interaction function depending on stochastic parameters and analyze the system behavior in the case of uniform interactions. The gPC approach is then summarized in Section 3. Subsequently, in Section 4 we consider the gPC scheme in a constrained setting and derive a selective model predictive approximation of the system. Next, in Section 5 we report several numerical experiments which illustrate the different features of the numerical method. Extensions of the present approach are finally discussed in Section 6.

2. Cucker-Smale Dynamic with Random Inputs

We introduce a Cucker-Smale type [10] differential system depending on a random variable with a given distribution . Let ,  , evolving as follows:

where is a time-dependent random function characterizing the uncertainty in the interaction rates and is a symmetric function describing the dependence of the alignment dynamic from the agents positions. A classical choice of space-dependent interaction function is related to the distance between two agents

where and are given parameters. Mathematical results concerning the system behavior in the deterministic case () can be found in [10]. In particular unconditional alignment emerges for . Let us observe that, even for the model with random inputs (1), the mean velocity of the system is conserved in time

since the symmetry of implies

Therefore for each we have .

2.1. The Uniform Interaction Case

To better understand the leading dynamic let us consider the simpler uniform interaction case when , leading to the following equation for the velocities:

The differential equation (5) admits an exact solution depending on the random input . More precisely if the initial velocities are deterministically known we have that

where is the mean velocity of the system. In what follows we analyze the evolution of (6) for different choices of and of the distribution of the random variable .

Example 1. Let us consider a random scattering rate written in terms of the following decomposition:

where is a nonnegative function depending on . The expected velocity of the th agent is defined by

Then each agent evolves its expected velocity according to

For example, let us choose , where the random variable is normally distributed; that is, . Then, for each , we need to evaluate the following integral:

The explicit form is easily found through standard techniques and yields

From (11) we observe a threshold effect in the asymptotic convergence of the mean velocity of each agent toward . It is immediately seen that if

it follows that, for , the expected velocity diverges. In particular, if we have that the solution starts to diverge as soon as . Note that this threshold effect is essential due to the negative tail of the normal distribution. In fact, if we now consider a random variable taking only nonnegative values, for example, exponentially distributed for some positive parameter , from (9) we obtain

which corresponds to

and therefore as . Then independently from the choice of the rate we obtain for each agent convergences toward the average initial velocity of the system. Finally, in case of a uniform random variable we obtain

that is,

which implies the divergence of the system in time as soon as assumes negative values.

Example 2. Next we consider a random scattering rate with time-dependent distribution function; that is,

with . As an example we investigate the case of a normally distributed random parameter with given mean and time-dependent variance: . It is straightforward to rewrite as a translation of a standard normal-distributed variable ; that is,

where . The expected velocities read

which correspond to

Similarly to the case of a time-independent normal variable a threshold effect occurs for large times; that is, the following condition on the variance of the distribution

implies the divergence of system (5). As a consequence instability can be avoided by assuming a variance decreasing sufficiently fast in time. The simplest choice is represented by for some . Condition (21) becomes

For example, if the previous condition implies that the system diverges for each .

3. A gPC Based Numerical Approach

In this section we approximate the Cucker-Smale model with random inputs using a generalized polynomial chaos approach. For the sake of clarity we first recall some basic facts concerning gPC approximations.

3.1. Preliminaries on gPC Approximations

Let be a probability space, that is, an ordered triple with any set, a -algebra, and a probability measure on , where we define a random variable

with the Borel set of . Moreover let us consider ,  , and certain spatial and temporal subsets. For the sake of simplicity we focus on real-valued functions depending on a single random input

In any case it is possible to extend the setup of the problem to a -dimensional vector of random variables ; see [19]. Let us consider now the linear space of polynomials of of degree up to , namely, . From classical results in approximation theory it is possible to represent the distribution of random functions with orthogonal polynomials , that is, an orthogonal basis of :

with as the Kronecker delta function. Assuming that the probability law , involved in the definition of the introduced function , has finite second-order moment, then the complete polynomial chaos expansion of is given by

According to the Cameron-Martin theorem and to the Askey scheme, results that pave a connection between random variables and orthogonal polynomials [18, 20, 21], we chose a set of polynomials which constitutes the optimal basis with respect to the distribution of the introduced random variable in agreement with Table 1.

Table 1: The different gPC choices for the polynomial expansions.

Probability law of
Expansion polynomials
Support

Gaussian
Hermite

Uniform
Legendre

Beta
Jacobi

Gamma
Laguerre

Poisson
Charlier

Let us consider now a general formulation for a randomly perturbed problem

where we indicated with a differential operator. In general the randomness introduced in the problem by acts as a perturbation of , of the function , or occurs as uncertainty of the initial conditions. In this work we focus on the first two aspects assuming that initial positions and velocities are deterministically known.

The generalized polynomial chaos method approximates the solution of (27) with its th-order polynomial chaos expansion and considers the Galerkin projections of the introduced differential problem; that is,

Due to the Galerkin orthogonality between the linear space and the error produced in the representation of with a truncated series, from (28) we obtain a set of purely deterministic equations for the expansion coefficients . These subproblems can be solved through classical discretization techniques. From the numerical point of view through a gPC-type method it is possible to achieve an exponential order of convergence to the exact solution of the problem, unlike Monte Carlo techniques for which the order is where is the number of samples.

3.2. gPC Approximation of the Alignment Model

We apply the described gPC decomposition to the solution of the nonhomogeneous differential equation in (6) and to the stochastic scattering rate ; that is,

where

We obtain the following polynomial chaos expansion:

Multiplying the above expression by an orthogonal element of the basis and integrating with respect to the distribution of

we find an explicit system of ODEs

where and

We recall that the gPC numerical approach preserves the mean velocity of the alignment model (5). In fact, from (33) follows

thanks to the symmetry of . More generally it can be shown that if

where is time-independent, then its gPC decomposition is also mean-preserving since

Remark 3. The gPC approximation (33) can be derived equivalently without expanding the kernel function . In this way one obtains

where now

Note that since in general , the overall computational cost is .

4. Selective Control of the gPC Approximation

In order to stabilize the gPC approximation of the Cucker-Smale type model (1) with random inputs, we introduce an additional term which acts as control of the approximated dynamic. More specifically we modify the approximation of the alignment model (1) by introducing a control term to the th component of its gPC approximation

The control is a solution of

where we assume to have, for each , a bounded control having value in an admissible set ; for example, in the one-dimensional case we consider . Parameter is a regularization term and are the desired values for the gPC coefficients. For example,

where is a desired velocity.

Moreover the controller action is weighted by a bounded function,

Due to the dependence of the controller effect from the single agent velocity, we refer to this approach as selective control; see [22].

In order to tackle numerically the above problem, whose direct solution is prohibitively expansive for large numbers of individuals, we make use of model predictive control (MPC) techniques, also referred to as receding horizon strategy or instantaneous control [23]. These techniques have been used in [8, 9, 22] in the case of deterministic alignment systems.

4.1. Selective Model Predictive Control

Let us split the time interval in time intervals of length with with . The basic idea of the model predictive control approach is to consider a piecewise constant control

In this way it is possible to determine the value of the control , solving for a state the (reduced) optimization problem

Given the control on the time interval , we let evolve according to the dynamics

in order to obtain the new state . We again solve (45) to obtain with the modified initial data and we repeat this procedure until we reach .

The reduced optimization problem implies a reduction of the complexity of the initial problem since it depends only on the single real-valued variable . On the other hand the price to pay is that in general the solution of the problem is suboptimal with respect to the full one ((40)-(41)).

The quadratic cost and a suitable discretization of (46) allow an explicit representation of in terms of and , as a feedback controlled system, as follows:

where and . Note that since the feedback control in (47) depends on the velocities at time , the constrained interaction at time is implicitly defined. The feedback controlled system in the discretized form results in

Again the action of the control is substituted by an implicit term representing the relaxation toward the desired component of the velocity , and it can be inverted in a fully explicit system.

Considering the scaling for the regularization parameter , the previous scheme is a consistent discretization of the following continuous system:

Now the control is explicitly embedded in the dynamic of the th component of the gPC approximation as a feedback term, and the parameter determines its strength.

4.2. Choice of the Selective Control

For the specific choice of weight function we refer in general to nonselective control. Note that in this case the action of the control is not strong enough to control the velocity of each agent; indeed in this case we are able only to control the mean velocity of the system. In fact the control term is reduced to

where is the th coefficient of the expansion of ; that is,

Then, only the projections of the mean velocity are steered toward the respective components of the target velocity; that is, as soon as it follows that . Therefore, the choice of the selective function is of paramount importance to ensure the action of the control on the single agent.

In principle one can address directly the control problem on the original system (1) as

where the control term is solution of

Here is a target velocity, a regularization parameter, and the set of admissible control. Similarly to previous subsection, through the approach presented in [8, 9, 22], we can derive the time-continuous MPC formulation which explicitly embed the control term in the dynamic as follows:

Now the gPC approximation of (54) can be obtained as in Section 3 and leads to the set of ODEs

where

In general systems (55) and (49), without further assumptions on the selective function , are not equivalent. In addition to the nonselective case, there exists at least one choice of selective control that makes the two approaches totally interchangeable. In fact, taking

and if , we have that is bounded and the control term in (55) takes the following form:

Similarly the control term in (49) reduces to

and therefore system (55) coincides with (49). Note that as both systems are driven towards the controlled state which implies a strong control over each single agent.

In Figure 1 we summarize the two equivalent approaches. In the case of nonselective control and of selective function given by (57) the constrained gPC system can be obtained from our initial unconstrained model (1) through two different but equivalent methods. The first approximates the solution of the Cucker-Smale type model via the gPC projection and then introduces a control on the coefficients of the decomposition through a MPC approach in order to steer each component to , whereas the second method considers a constrained Cucker-Smale problem (52), introduces its continuous MPC approximation, and then computes the gPC expansion of the resulting system of constrained differential equations.

Figure 1: Sketch of the two numerical approaches to solve the control problem with uncertainty, combining MPC and gPC. In both cases, of nonselective control, that is, , and of selective control with defined in (57), the two approaches are equivalent.

Remark 4. We remark that the choice of stated in (57), for which the two approaches sketched in Figure 1 are identical, is equivalent to consider the constrained dynamic (52), modified as follows:

where the control term, , for each agent , is given by the minimization of the following functional:

Since the functional is strictly convex, applying the (MPC) procedure on a single time interval for the discretized dynamic of (60)-(61), we obtain in terms of feedback control

Thus the same considerations on the equivalence of the approaches hold.

5. Numerical Tests

We present some numerical experiments of the behavior of the flocking model in the case of a Hermite polynomial chaos expansion. This choice corresponds to the assumption of a normal distribution for the stochastic parameter in the Cucker-Smale type equation (1) and in its constrained behavior (49). Numerical results show that the introduced selective control with the weight function (57) is capable of driving the velocity to a desired state even in case of a dynamic dependent on a normally distributed random input, with fixed or time-dependent variance. In the uniform interaction case, since the effect of agents’ positions does not influence the alignment we report only the results of the agents’ velocities.

5.1. Unconstrained Case

In Figures 2 and 3 we present numerical results for the convergence of the error using the gPC scheme described in (33) for and solved through a -order Runge-Kutta method. In particular Figure 2 shows the behavior of the error with respect to increasing terms of the gPC decomposition. Here we considered the average in time of the error for the mean and the variance at time in the norm

where

with and defined in (6) and (8). Observe that if the scattering rate is of the from described in (7) with and then, in addition to the explicit evolution for the expected velocity as in (9), we can obtain the exact version for the evolution of the variance of the th agent

(a)

(b)

Figure 2: Error convergence for increasing number of polynomials in the gPC decomposition approximation. (a) Convergence of the mean error at two fixed times and . (b) Convergence of the variance error. In both cases we considered a random time-independent scattering , where the random variable is normally distributed . The system of ODEs is solved through a -order Runge-Kutta with .

(a)

(b)

Figure 3: Evolution of the variance error defined in (63) for the gPC decomposition for the unconstrained model (5) with over the time interval with and time step .

In (63) we indicated with the approximated variance

It is easily seen how the error decays spectrally for increasing value of ; however the method is not capable of going above a certain accuracy and therefore for large a threshold effect is observed. This can be explained by the large integration interval we have considered in the numerical computation and by the well-known loss of accuracy of gPC for large times [19]. In the case of the error of the variance (Figure 3) the gPC approximation exhibits a slower convergence with respect to the convergence of the mean. Next in Figure 4 we see how for large times the solution of the differential equation (5) diverges and the numerical approximation is capable of describing accurately its behavior only through an increasing number of Hermite polynomials.

(a)

(b)

Figure 4: (a) Sixth-order Hermite gPC decomposition solved through a -order Runge-Kutta. (b) -order Hermite gPC decomposition solved through a -order Runge-Kutta. In both cases the final time considered is , with time step .

5.2. Constrained Uniform Interaction Case

In Figure 5 we show different scenarios for the uniform interaction dynamic with constraints. In the first row we represent the solution for agents, whose dynamic is described by (49) with ; different values of originate different controls on the average of the system, which however do not prevent the system from diverging. In the second row we show the action of selective control (57). It is evident that, with this choice, we are able to control the system also in the case with higher variance.

(a)

(b)

(c)

(d)

Figure 5: Evolution of the uniform interaction alignment model (49) with agents, at distributed around with unitary variance, depending on a normal random parameter. Left column: . Right column: . The control term shows its ability to steer the system towards desired velocity , with different intensities and ; when the control has no influence. First row shows the action of the control acting just on the average velocity, . Second row shows the action of selective control with as in (57).

Observe that the numerical results are coherent with the explicit solution of the controlled equation. Let us consider the time-independent scattering rate ; then from the equation

we can compute the exact solution given :

The asymptotic behavior of the expected value of (68) can be studied similarly to what we did in Section 2.1. In other words in order to prevent the divergence of the leading term of the controlled expected exact solution we might study

which diverge if

Then for each fixed time we could select a regularization parameter so as to avoid the divergence of (68). Moreover we can observe that in the limit the introduced selective control is capable of correctly driving the system (67) for each .

In Figure 6 we consider the system with random time-dependent scattering rate . The dynamic shows how, for the choice of time-dependent variance described within Example 2 in Section 2.1, that is, with , the convergence depends on the mean value of the random input. In particular numerical experiments highlight the threshold effect for which we derived in Section 2. In the second figure we show that the action of the selective control (57), with desired velocity , is capable of stabilizing the system and driving the velocities towards the desired state.

(a)

(b)

Figure 6: Solution of the uniform interaction case with time-dependent random parameter distributed accordingly to a normal distribution , with a time-dependent standard deviation , and . (a) We see the threshold for different values of ; that is, for the system diverges. (b) Solution of the constrained model with ; observe that we are able to steer the system to the desired velocity , that is, the initial mean velocity of the system, using the selective control described in (57).

5.3. Constrained Space-Dependent Case

Next let us consider the full space nonhomogeneous constrained problem (1) with the interaction function defined in (2). In this case we assume that with . In Figures 7 and 8 we consider a system of agents with Gaussian initial position with zero mean and with variance and Gaussian initial velocities clustered around with variance . The numerical results for (33) have been performed through a -order gPC expansion. The dynamic has been observed in the time interval with . In Figure 8 we see how the selective control is capable of driving the velocity of each agent to the desired state . In fact in case of no control (see Figure 7) we have that the velocities of the system naturally diverge.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 7: Numerical solution of (49), with , , through a -order gPC Hermite decomposition (33) with with time step . The random input is normally distributed .

(a)

(b)

(c)

(d)

(e)

(f)

Figure 8: Numerical solution of (49), with , , through a 6th-order gPC Hermite decomposition for the selective control (33) with time step . Here we considered a normally distributed random input ; the desired velocity is and the control parameter is .

6. Conclusions

We proposed a general approach for the numerical approximation of flocking models with random inputs through gPC. The method is constructed in two steps. First, the random Cucker-Smale system is solved by gPC. The presence of uncertainty in the interaction terms, which is a natural assumption in this kind of problems, leads to threshold effects in the asymptotic behavior of the system. Next a constrained gPC approximation is introduced and approximated though a selective model predictive control strategy. Relations under which the introduction of the gPC approximation and the model predictive control commute are also derived. The numerical examples illustrate that the assumption of positivity of the mean value of the random input is not sufficient for the alignment of the system but that a suitable choice of the selective control is capable of stabilizing the system towards the desired state. Extension of this technique to the case of a large number of interacting agents through mean-field and Boltzmann approximations is actually under study.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

Giacomo Albi acknowledges the support of the ERC-Starting Grant High-Dimensional Sparse Optimal Control (HDSPCONTR-306274).

References

L. Pareschi and G. Toscani, Interacting Multiagent Systems: Kinetic Equations & Monte Carlo Methods, Oxford University Press, New York, NY, USA, 2013.

C. W. Reynolds, “Flocks, herds and schools: a distributed behavioral mode,” Computer Graphics, vol. 21, no. 4, pp. 25–34, 1987.

T. Vicsek, A. Czirók, E. Ben-Jacob, I. Cohen, and O. Shochet, “Novel type of phase transition in a system of self-driven particles,” Physical Review Letters, vol. 75, no. 6, pp. 1226–1229, 1995.

G. Albi, D. Balagué, J. A. Carrillo, and J. von Brecht, “Stability analysis of flock and mill rings for second order models in swarming,” SIAM Journal on Applied Mathematics, vol. 74, no. 3, pp. 794–818, 2014.

J. A. Carrillo, M. Fornasier, J. Rosado, and G. Toscani, “Asymptotic flocking dynamics for the kinetic Cucker-Smale model,” SIAM Journal on Mathematical Analysis, vol. 42, no. 1, pp. 218–236, 2010.

S. Motsch and E. Tadmor, “Heterophilious dynamics enhances consensus,” SIAM Review, vol. 56, no. 4, pp. 577–621, 2014.

H. Zhou, W. Zhou, and W. Zeng, “Flocking control of multiple mobile agents with the rules of avoiding collision,” Mathematical Problems in Engineering. In press

G. Albi, M. Herty, and L. Pareschi, “Kinetic description of optimal control problems in consensus modeling,” Communications in Mathematical Sciences, vol. 13, no. 6, pp. 1407–1429, 2015.

G. Albi, L. Pareschi, and M. Zanella, “Boltzmann-type control of opinion consensus through leaders,” Philosophical Transactions of the Royal Society A, vol. 372, no. 2028, 2014.

F. Cucker and S. Smale, “Emergent behavior in flocks,” IEEE Transactions on Automatic Control, vol. 52, no. 5, pp. 852–862, 2007.

P. Degond, J.-G. Liu, and C. Ringhofer, “Large-scale dynamics of mean-field games driven by local Nash equilibria,” Journal of Nonlinear Science, vol. 24, no. 1, pp. 93–115, 2014.

P. Degond, M. Herty, and J. G. Liu, “Mean-field games and model predictive control,” http://arxiv.org/abs/1412.7517

S. Motsch and E. Tadmor, “A new model for self-organized dynamics and its flocking behavior,” Journal of Statistical Physics, vol. 144, no. 5, pp. 923–947, 2011.

S. M. Ahn and S.-Y. Ha, “Stochastic flocking dynamics of the Cucker-Smale model with multiplicative white noises,” Journal of Mathematical Physics, vol. 51, no. 10, Article ID 103301, 2010.

S.-Y. Ha, K. Lee, and D. Levy, “Emergence of time-asymptotic flocking in a stochastic Cucker-Smale system,” Communications in Mathematical Sciences, vol. 7, no. 2, pp. 453–469, 2009.

Y. Sun, “Mean square consensus for uncertain multiagent systems with noises and delays,” Abstract and Applied Analysis, vol. 2012, Article ID 621060, 18 pages, 2012.

C. A. Yates, R. Erban, C. Escudero et al., “Inherent noise can facilitate coherence in collective swarm motion,” Proceedings of the National Academy of Sciences of the United States of America, vol. 106, no. 14, pp. 5464–5469, 2009.

D. Xiu, Numerical Methods for Stochastic Computations, Princeton University Press, Princeton, NJ, USA, 2010.

M. Gerritsma, J.-B. van der Steen, P. Vos, and G. Karniadakis, “Time-dependent generalized polynomial chaos,” Journal of Computational Physics, vol. 229, no. 22, pp. 8333–8363, 2010.

R. H. Cameron and W. T. Martin, “Transformations of Wiener integrals under translations,” Annals of Mathematics, vol. 45, no. 2, pp. 386–396, 1944.

D. Xiu and G. E. Karniadakis, “The Wiener-Askey polynomial chaos for stochastic differential equations,” SIAM Journal on Scientific Computing, vol. 24, no. 2, pp. 619–644, 2002.

G. Albi and L. Pareschi, “Mean-field selective control of flocking dynamics,” In preparation

H. Michalska and D. Q. Mayne, “Robust receding horizon control of constrained nonlinear systems,” IEEE Transactions on Automatic Control, vol. 38, no. 11, pp. 1623–1633, 1993.