Epidemic Synchronization in Robotic Swarms

Proceedings of the 2009 IFAC Workshop on Networked Robotics Golden, CO, USA, October 6-8, 2009

Epidemic Synchronization in Robotic Swarms Henrik Schioler ∗ Trung Dung Ngo ∗∗ Jens Dalsgaard ∗∗∗ ∗

Section for Automation and Control, Centre for Embedded SW Systems, Aalborg University, Denmark (e-mail: [email protected]). ∗∗ Section for Automation and Control, Dept. of Electronic Systems, Aalborg University, Denmark (e-mail: [email protected]) ∗∗∗ Section for Automation and Control, Dept. of Electronic Systems, Aalborg University, Denmark (e-mail: [email protected]) Abstract: Clock synchronization in swarms of networked mobile robots is studied in a probabilistic, epidemic framework. In this setting communication and synchonization is considered to be a randomized process, taking place at unplanned instants of geographical rendezvous between robots. In combination with a Markovian mobility model the synchronization process yields overall evolutionary dynamics for first and second conditional moments of synchronization error given geographical position. The established dynamics assume the shape of partial integrodifferential equations and the swarm is subsequently studied as an infinite-dimensional optimal control problem. Illustrative numerical examples are given and commented. Keywords: Mobile robots, Clock synchronization, Wireless communication, Random processes, Markov Models, Integro-Differential equations, Optimization 1. INTRODUCTION

synchronization in sparse ad-hoc networks and (Hu 2003) considers optimal multi-hop synchronization in dense sensor networks, where synchronization propagates in rings away from a reference clock and nodes in outer rings hear aggregated synchronization pulses from already synchronized nodes. This work is distinguished from those contributions, and to the authors knowledge, from any other works in this field by the inherent randomization of the synchronization dynamics, the inclusion of mobility modelling and finally the resulting overall probabilistic model obtained. In coherence with intuition, clock synchronization is improved along with the freqency of synchronization events. However excessive synchronization may lead to intolerable communication overhead as well as energy consumption associated to message transmission. Individual robots need some policy to determine their synchronization frequency. In this paper we study policies based on the current position of robots and seek to optimize such a policy w.r.t synchronization error and energy consumption. The paper is organized as follows; first a clock drift model including random as well as deterministic effects is presented and simplified to fit the purpose of this work, next the assumed randomized synchronization mechanism is explained and model dynamics for 1st. and 2nd. moment statistics, not including mobility effects, are elaborated. Mobility effects are introduced in the following section based in a suitable Markovian mobility model. Subsequently an optimization problem is introduced to express the trade-off between clock synchronization error and energy consumption. A numerical approximation scheme to handle the inherently infinitely dimensional optimization problem is outlined and a numerical results for a simplis-

Collaboration between individuals in groups is a frequently appearing phenomena in nature as well as among humans. It serves a variety of tasks for the benefit of the group as a whole as well as for the individuals. At a general level most studies in robotics are occupied with the transfer of biological/human skills and behaviour to the robotic counterpart; industrial robots replace human work force and newly developed social robots with imitated mimics substitute human company for elders or mentally disabled. In this perspective the transfer of biological group behaviour to groups/swarms of mobile robots seems to be a natural extension to this effort, which is not only motivated as a study object but indeed as a means of improving efficiency and increasing autonomy - not of the individual robots - but of the swarm as a whole. In this paper we study the process of randomized clock synchronization among networked mobile robots. Although clock synchronization does not seem to possess a natural biological origin, the underlying communication process facilitating synchronization is biologically inspired. In most biological groups communication is local and incidental, i.e. when members of the group come into proximity, they are likely to exchange information about the overall state of affairs such as location of food resources or dangers. Similar dynamics are found in broadcasting ad-hoc networks as well as rumour- and disease spread, hence the term epidemic. Clock synchronization is by now a well established body of theory and methodology widely applied in general purpose networks (Mills 2006) as well as in dedicated communication networks such as fieldbusses etc. (Kopetz 1997). For ad-hoc networks (Romer 2001) adresses the problem of

978-3-902661-53-1/09/$20.00 © 2009 IFAC

43

10.3182/20091006-3-US-4006.0017

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

of Di corresponds to integral control action. Quantities Wi (t− ), Wj (t− ) account for white phase noise in clocks i and j respectively, which are only manifested when clocks are actually read. In between communication instants each clock value ci progresses according to (3) and

tic one-dimensional example are given and commented. Finally conclusions are made along with suggestions for future work. 2. CLOCK DRIFT MODEL

d/dt Di (t) = 0

Various models have been proposed to account for clock imprecision. We consider absolute synchronization, i.e. where nodes/robots seek a common understanding of time with an external reference authority. Thus leaving out relativistic effects, every clock may be related to some common reference defining the zero value. The clock ci for robot ”i” may be modelled by a random process, i.e. d/dt ci = Di + Ai t + Die + d/dt Ci (1)

We assume synchronization message transmission to be governed by a non-homogeneous Poisson process, i.e. in every infinesimal time interval of length h, transmission takes place with a probability proportional to h. Let dij be independent random booleans, then the overall synchronization dynamics for robot i is given by P ci (t + h) = (1 − j6=i dij ) ci (t) P + j6=i dij α (ci (t) + Wi (t)) + (1 − α)(cj (t) + Wj (t))) + h ∗ Di (t) (6)

where Di accounts for a constant drift caused by manufacturing, Ai accounts for ageing, Die is a random process accounting for environmental effects like temperature, pressure, etc., and Ci comprises the pure random noise driven effects. We shall assume Di to be tuneable and utilize this to accomodate for mean value clock frequency discrepancies as well as slow clock frequency variations. Also ageing and environmental terms are assumed to vary sufficiently slowly to be comprised in Di for the timespan considered. Traditionally (Rutman 1991) Ci has been specified in terms of its power spectral density Pi , i.e. Pi (f ) =

0 1 X hj f j 2π j=−4

(5)

and P Di (t + h) = (1 − j6=i dij ) Di (t) P + j6=i dij β Di (t) + (1 − β)(cj (t) + Wj (t) − ci (t) − Wi (t)) (7) Thus dij = 1 indicates that robot ”j” transmits a synchronization message, which is received by robot ”i” within the time interval [t, t+h]. Assuming robot positions xi and xj , we may state,

(2)

P (dij = 1) = hγ(xj )K(xi , xj )

where j = 0 corresponds to white phase noise and j = −2 to white frequency noise respectively. In general higher j-values account for faster clock variations. We assume variations for j ∈ {−3, −4} to be sufficiently slow to be comrised into Di . To maintain analytical tractability we disregard the so called phase flicker, i.e. j = −1 which, except for slow frequency variations, is similar to what is proposed in (Hurni 2007). Thus, expressed as a stochastic differential equation dci = dt Di + σi dBi (3)

(8)

where K is a kernel function modelling the locality of communication. K may be decreasing in the distance between robots, i.e. |xi − xj |. As an example K(xi , xj ) = I|xi −xj |
where Bi is a standard Brownian motion accounting for j = −2, i.e. white frequency variations. White phase variation (j = 0) only manifests itself when a particular clock is read, i.e. it appears as a measurement noise, modelled in this case as a zero mean normally distributed 2 random variable with variance σW .

3.1 First moment dynamics Taking conditional expectations given the position x of robot i and limits for h → 0 in equations (6) and (7), yields,

3. SYNCHRONIZATION MODEL

∂/∂t ci (x,P t) R = (1 − α) j6=i D K(x, y)γ(y)(cj (y, t) − ci (x, t))Lj (dy) + Di (x, t) (9)

In this work synchronization is assumed to take place upon message transmission. That is every message includes a timestamp indicating the clock of the transmitter ”j”, whereas clock synchronization is only performed at the recieving side ”i” according to equation (4),

∂/∂t Di (x,Pt) = (1 − β) j6=i R K(x, y)γ(y)(cj (y, t) − ci (x, t) − Di (x, t))Lj (dy) D

ci (t+ ) = α (ci (t− ) + Wi (t− )) + (1 − α)(cj (t− ) + Wj (t− )) Di (t+ ) = β Di (t− ) + (1 − β)(cj (t− ) + Wj (t− ) − ci (t− ) − Wi (t− )) (4)

(10)

where D is the domain of operation and Lj is a probability measure on D capturing the stationary location distribution of robot j, which is elaborated in the next section. Assuming a homogeneous situation, where Lj = L, Di (x, t) = D(x, t) and cj (x, t) = c(x, t), ∀j gives

where α and β are positive reals slightly less than one. In feedback control terms we may think of the update of ci as a proportional feedback, whereas the update

44

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

to randomized movement. We adopt a so called velocity jump model (Keshet 2001), which as Brownian motion is a Markovian mobility model, where velocity remains constant between time instants {t n } of velocity change. Time instants {tn } are assumed to ∂/∂t D(x, t) R be a homogenous Poisson process with intensity parameter = (1−β)(N −1) D K(x, y)γ(y)(c(y, t)−c(x, t)−D(x, t))L(dy) λ. (12) At times {tn }, of velocity change, velocity is drawn independently from a distribution LQ (x(tn )) in this case where N is the number of robots in the swarm. A number assumed to be concentrated on a discrete set of velociof zero drift reference clocks located at points y1 , .., yM ties {v1 , .., vW }. As a consequence a time dependent joint may be included as in equation (13) and (14). Reference position/velocity distribution L (A, {v }, t) expressing the i j clocks are assumed to be transmitted at a constant Poisson kinematic state of robot i may be found, which is assumed rate γR . to be composed as follows, Z ∂/∂t c(x, t) R L (A, {v }, t) = Li ({vj }|x, t) · fLi (x, t)dx (19) i j = (1 − α)(N − 1) D K(x, y)γ(y)(c(y, t) − c(x, t))L(dy) PM A + (1 − α)γR k=1 K(x, yk )(0 − c(x, t)) + D(x, t) (13) i.e. the distribution of position possesses a density function fL . The velocity jump model constitutes a Markov process ∂/∂t D(x, t) R transforming kinematic state distribution as follows, = (1−β)(N −1) D K(x, y)γ(y)(c(y, t)−c(x, t)−D(x, t))L(dy) PM ∂/∂t(Li ({vj }|x, t) · fLi (x, t)) + (1 − β)γR k=1 K(x, yk )(0 − c(x, t) − D(x, t)) (14) = λ(x)(L ({v }, x) − L ({v }|x, t))f i (x, t) ∂/∂t c(x, t) R = (1 − α)(N − 1) D K(x, y)γ(y)(c(y, t) − c(x, t))L(dy) + D(x, t) (11)

Q

A stationary solution c(x, ∞) = c, exists, which is constant over D, for which (13) gives c(1 − α) · γR

M X

K(x, yk ) = D(x, t)

(15)

and from (14)

L

(20)

We give a 2-dimensional numerical example, where D = [−1, 1] with 2 discrete velocities {v1 = −1, v2 = 1}. We let LQ ({−1}, x) = LQ ({1}, x) = 1/2 for x ∈ D, i.e. robots tend to move backwards and forwards with equal probability within the domain of operation. Outside D we let velocities unanimously point towards D, i.e. LQ ({1}, x < −1) = LQ ({1}, x > 1) = 1, for x ∈ / D. Having λ(x) assume very high values outside D together with the definition of LQ , models reflection on the borders of D. Specifically fL (x) becomes uniformly concentrated on D for λ(x) approaching infinity outside D. Altogether in the limit Li ({vj }|x) = LQ ({vj }, x) = fL (x) = 1/2 for x ∈ D. Solutions for a variety of values for λ(x), x ∈ / D er shown in figure (1). The stationary solution constitutes the basis of the analysis clock synchronization dynamics in the following.

(16)

(17)

3.2 Second moment dynamics Squaring equation (6), taking conditional expectations given the position of robot ”i” and taking limits for h → 0 yields the following dynamics for the second conditional moment c2(x, t), ∂/∂t c2(x, t) R = (1 − α)2 (N − 1) D K(x, y) γ(y) c2(y, t) L(dy) R 2 + (1 − α)2 σW ((2(N − 1) D K(x, y) γ(y) L(dy) PM + γR k=1 K(x, yk )) R − (2 − α) c2(x, t) ((N − 1) D K(x, y) γ(y) L(dy) PM + γR k=1 K(x, yk )) + σ 2

j

4.1 Numerical example

which altogether gives

c = D(x, t) = 0

i

Transformation (20) under mild assumptions defines a unique stationary solution Li ({vj }|x) · fLi (x), which, by suitable selection of λ(x) and LQ , may be set to model some observed kinematic distribution of robots, e.g. to be uniform within D.

k=1

0 = (1 − β)D(x, t) · R ((N − 1) D K(x, y)γ(y)L(dy) + PM γR k=1 K(x, yk ) + 1/(1 − α))

j

− < vj , ∇x [Li ({vj }|·, ·)fLi ](x, t) >

4.2 Mobility Effect on Syncronization Even though a stationary location distribution is assumed mobility has an effect on the clock statistics c and c2 since robots over time carry their individual clocks between different sub areas within the domain of operation D. Thus robots, which for some time have been in the vicinity of the reference clock and thus being well synchronized, may move to areas distant to the reference clock and effectively spread this information by movement. To combine the mobility dynamics with the synchronization dynamics in equations (9) to (18) we need to embrace the combined

(18)

4. MOBILITY MODEL To facilitate probabilistic modelling we adopt a randomized view on robot motion. Although robot motion is likely to be deterministically governed by overall objectives and individual robot policies, we assume that the complexity of the swarm as a whole lead to overall statistics similar

45

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

to the selection of model parameters, which raises the question of how to select these in a reasonable way. At least parameters α, β, γ, λ seem to be candidates for design choice, since they may be programmed into the behaviour of individual robots. In this case, for brevity, we shall exclusively focus on the parameter γ(x) which as indicated allows position dependence. γ scales the rate of synchronization events depending on the current position of a particular robot. Large values of γ gives a higher rate of synchronization events and in turn a higher energy consumption rate, whereas the converse is rather obvious. Thus γ should be set sufficiently high to assure a sufficient synchronization rate and on the other hand sufficiently low not to exceed tolerable energy consumption. Altogether it seems appropriate to consider the problem in an optimization setting. Several alternative strategies seem natural;

0.5 0.45 0.4 0.35

Density

0.3 0.25

lambda=2 lambda=4

0.2

lambda=8 lambda=16

0.15 0.1 0.05 0 −2

−1.5

−1

−0.5

0 Position

0.5

1

1.5

2

Fig. 1. Stationary densities for a variety of λ-values Markov state of the two dynamics, i.e. synchronization dynamics should additionally comprise velocity information. We define ci (x, v, t) to be the conditional expected clock value in robot ”i”, when found in position x and velocity v at time t. For the homgenous case mobility dynamics transform equation (13) into (21), i.e.

• a weigthed performance function J(γ, ρg , ρc ) may be defined by Z J(γ, ρc , ρg ) = ρc c2(x, γ) + ρg γ(x) dx (23)

∂ fL (x)L({vj }|x) ∂t c(x, vj , t) = fL (x)L({v }|x) ((1 − α) · j R ((N − 1) D K(x, y)γ(y) (c(y, t) − c(x, vj , t))L(dy) PM + γR k=1 K(x, yk )(0 − c(x, vj , t))) + Di (x, vj , t)) − < vj , ∇x [ci fL L](x, vj , t) > + λ (LQ ({vj }, x) c(x, t) − L({vj }|x) c(x, vj , t))fL (x)

D

where c2(·, γ) denotes the stationary conditional second moment of syncronization error for a particular choice of γ. Since γ scales transmission rate it constitutes a reasonable measure for energy consumption related to communication/synchronization. Thus J aggregates energy consumption and synchronization error into one overall performance measure, where the tradeoff between energy and error is determined by the tuning parameters ρc , ρg .

(21) A similar transformation valid for equation (14) is not shown due to space limitations. Its is readily shown that the inclusion of mobility preserves unbiasedness, i.e. c(x, vj , ∞) = 0, when clock frequency tuning is assumed. For conditional 2nd. moment further arguments are needed. The synchronization and mobility effects on clock synchronization dynamics are assumed to be probabilistically independent processes. Additionally 2nd. moments of changes for infinitesimally small time steps h are assumed to exhibit certain smoothness properties for h = 0. Under these assumptions independent effects may be combined by summation of terms into overall dynamics for the 2nd. moment, i.e.

• As one alternative to the above, an average transmission rate J(γ, 0, 1) across the entire population of robots may be held below an acceptable level during optimization of J(γ, 1, 0) w.r.t. γ. • Finally may J(γ, 0, 1) be mimimized for J(γ, 1, 0) held below a fixed level.

In the sequel we pursue to latter two options since they remain independent on the choice of tuning parameters ρc and ρg . Alternating those options yield a sequence of non increasing values of J(γ, 0, 1) and J(γ, 1, 0) converging to ∂ fL (x)L({vj }|x) ∂t c2(x, vj , t) = a local Pareto optimum. R 2 fL (x)L({vj }|x) ((1−α) (N −1) D K(x, y) γ(y)c2(y, t) L(dy) The optimization problem so established constitutes an R 2 infinite dimensional constrained optimization problem, + (1 − α)2 σW (2 (N − 1) D K(x, y)γ(y) L(dy) PM where constraints are defined by imposing the stationarity + γR k=1 K(x, yk )) R condition to (22) as well as bounding J(γ, 0, 1) or J(γ, 1, 0) − (2 − α) c2(x, vj , t) ((N − 1) D K(x, y)γ(y) L(dy) i.e. PM + γR k=1 K(x, yk )) + σ 2 ) − < vj , ∇x [c2 fL L](x, vj , t) > ∂/∂t c2(x, vj , t) = 0 for ∀vj ∈ {v1 , .., vW } + λ (LQ ({vj }, x) c2(x, t) − L({vj }|x) c2(x, vj , t))fL (x) J(γ, 0, 1) ≤ J¯ or (22) J(γ, 1, 0) ≤ J¯ (24) 5. OPTIMAL CONTROL 5.1 Numerical approximation procedure Under stationary conditions we find clocks to be mutually unbiased as well as unbiased w.r.t. a reference clock as In this work we do not aspire for an analytical solution expressed in (17). This result is valid for any meaningful γ ∗ to the presented optimization problem. Instead we choice of model parameters α, β, γ, λ, K, N etc. Second consider a sequence of approximate solutions {γn∗ ∈ Γ} moment stationary statistics however remain sensitive such that

46

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

lim J(γn∗ ) = inf J(γ)

n→∞

γ∈Γ

(25)

0.35

where Γ denotes the set of admissible γ-functions/controls. The sequence {γn∗ } is found by defining the approximation mappings An : Γ → Γn ⊆ Γ and Bn : Γn → Cn ⊆ C, where C denotes an appropriate solution space to (22) and (5). As an example Γn may be the space of piecewise constant functions defined on a partition Dn of the domain of operation D, where Dn become increasingly fine for increasing n. The mapping Bn may be defined by Bn : γ 7→ c2n for some discrete approximate solution c2n to (22) and (5). Next define Jn : Γn → R by Z Jn (γ, ρc , ρg ) = ρc c2n (x, γ) + ρg γ(x) dx (26)

0.3 v=−1 v=1

Finally we

min Jn (γ)

γ∈Γn

Conditional 2nd. moment

0.2

0.15

0.1

0.05

0 −1

−0.8

−0.6

−0.4

−0.2

0 Position

0.2

0.4

0.6

0.8

1

Fig. 2. Solution for conditional second R R moments constant value of γ = 0.001, γ = 0.002, c2 = 0.36

D

define γn∗ by γn∗ = arg

gamma

0.25

is significant however not excessive. Despite a relatively small gain in performance γ seems to undergo a significant transformation. The optimal γ effectively confines message transmission to a regular subgrid of D with a spacing of 0.2, which matches transmission range d = 0.1 so that every position within D is within range of transmission except a small interval next to the reference clock location x = −1. The edge effect at the opposite interval end x = 1 has been covered through a high transmission level close to x = 0.9. Fixing J(γ, 1, 0) below 0.3 during minimization of J(γ, 0, 1) yields a result shown in figure (4). In this case no significant improvement in J(γ, 0, 1) is obtained. Therefore the result depicted in figures (3) and (4) is assumed to be a local Pareto optimum. With an initial solution as shown in figure (2) J(γ, 1, 0) is kept below 0.36 during minimization of J(γ, 0, 1) from 0.002 to 0.0017. The result is shown in figure (5). In this case transmission is eliminated in proximity of the reference clock and intensified in the opposite end of D to neutralize edge effects. From the solution shown in (5) J(γ, 0, 1) is kept below 0.0017 during mimization of J(γ, 1, 0), which yields no significant performance improvement. However the shown solution (figure (6)) has undergone significant transformation to show close resemblence to results shown in figures (2) and (4). Thus the comments stated above pertain to this result as well. Since no significant performance change is observed between figures (5) and (6) both results are assumed to be close to Pareto optimal.

(27)

If a sequence of mappings {An : Γ → Γn } exists such that lim Jn (An (γ)) = J(γ) (28) n→∞

uniformly in Γ, it is readily shown that (25) is fulfilled. In the sequel it is without proof assumed that the chosen approximation mappings An and Bn exist and fulfill (28). 6. NUMERICAL EXAMPLE We extend the example introduced above for the mobility model. Additionally we set N = 50, α = 0.95, λ(x) = 0.5 for x ∈ D, γR = 3, σ = 1, σW = 0.1, K(x, y) = I|x−y|<0.1 , i.e. transmission range is limited to 0.1. One reference clock is assumed to be located at the leftmost interval end, i.e. x = −1. In this case (22) reduces to ∂/∂t c2(x, vj , t) = R (1-α)2 (N −1)/4 D K(x, y)γ(y) (c2(y, v1 , t)+c2(y, v2 , t)) dy  R 2 (N − 1)/2 D K(x, y)γ(y) dy +(1-α)2 σW + γR K(x, −1) R -(2-α) c2(x, vj , t) ((N − 1)/2 D K(x, y)γ(y) dy + γR K(x, −1)) + σ 2 -sign(vj ) ∂/∂x c2(x, vj , t) +λ ((c2(x, v1 , t) + c2(x, v2 , t))/2 − c2(x, vj , t)) (29) 6.1 Results

7. CONCLUSION

Solving (29) for a constant γ(x) = 0.001 yields results as shown in figure (2). Conditional clock drift variances for both predefined velocities {−1, 1} are significantly low in proximity of the reference clock located at x = −1. Both conditional variance functions increase with the distance to the reference clock. However since robots moving towards the left (v1 = −1) carry clock variance from a region further away from the reference clock a higher variance for this velocity is observed, whereas the converse argument holds forR(v2 = 1). In figure (3) J(γ, 0, 1) = D γ(x)dx is kept below the level associated to (2) through optimization of J(γ, 1, 0) = R c2(x, {−1}) + c2(x, {1})dx w.r.t. γ itself. Between figD ures (3) and (2) J(γ, 1, 0) is reduced from 0.36 to 0.3, which

A probabilistic model for epidemic clock synchronization, among mobile robots in a swarm, is presented including random clock drift, random mobility and Poisson generated clock messages. Slow variations are captured by tuning individual clock frequency, whereas rapid variations such as white phase and white frequency effects are handled through individual clock phase tuning. It is shown how frequency tuning in a static clock drift scenario produces zero mean clock drift. This result is carried along to the analysis of fast variations, where zero mean is assumed to give an overall model for rapid clock drift variance, which is adopted as an overall performance measure for the swarm synchonization. A velocity jump model is adopted

47

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

0.35

0.35 v=−1

v=−1

0.3

0.3

v=1

v=1

gamma

0.25

Conditional 2nd. moment

Conditional 2nd. moment

gamma

0.2

0.15

0.1

0.05

0 −1

−0.8

−0.6

−0.4

−0.2

0 0.2 Position

0.4

0.6

0.8

0.1

v=1 gamma

0.2

0.15

0.1

0.05

−0.6

−0.4

−0.2

0 0.2 Position

0.4

0.6

0.8

1

Fig. R4. Optimal Rconditional second moments for constant c2 = 0.3, γ = 0.002 0.35 v=−1

0.3

−0.6

−0.4

−0.2

0 0.2 Position

0.4

0.6

0.8

1

A numerical approximation procedure is presented, where optimization over a sequence of discrete approximations is used to obtain an asymptotically optimal solution for an aggregate performance measure of variance and energy. Variance and energy may be aggregated through a weighted sum. Alternatively the one may be minimized for a fixed value of the other or vice-versa. Numerical results are shown in a simplistic one-dimensional example. Although performance is not improved to an excessive level results reveal a surprising and yet intuitively coherent shape for message transmission frequency, where transmission is effectively confined to a regular grid with a granularity matching transmission range. Various extensions to the presented work seem natural, such as extension to 2- and 3-dimensional examples as well as carrying out large scale statistical experiments and simulation to validate analysis. Letting syncronization parameters such as α and β depend on estimated clock drift variance of the transmitting robot represents one possible extension out of a variety to the synchronization policy presented along with methods for estimating the quality of recieved synchronization messages.

0.25

−0.8

−0.8

Fig. R6. Optimal conditional second moments for constant R γ = 0.0017, c2 = 0.36

v=−1

0.3

Conditional 2nd. moment

0.15

0 −1

1

0.35

v=1 gamma

Conditional 2nd. moment

0.2

0.05

Fig. R3. Optimal conditional second moments for constant R γ = 0.002, c2 = 0.3

0 −1

0.25

0.25

0.2

0.15

REFERENCES

0.1

[Mills 2006] D. L. Mills, Computer Network Time Synchronization: the Network Time Protocol, CRC Press 2006, ISBN 978-0-8493-5805-0 [Kopetz 1997] H. Kopetz, Real-Time Systems Design Principles for Distributed Embedded Applications The Springer International Series in Engineering and Computer Science , Vol. 395, 1997,ISBN: 978-0-79239894-3 [Romer 2001] K. Romer, Time synchronization in ad hoc networks, Proceedings of the 2nd ACM international symposium on Mobile ad hoc networking, 2001, ISBN:1-58113-428-2 [Hu 2003] A. Hu and S. D. Servetto, Asymptotically optimal time synchronization in dense sensor networks Proceedings of the 2nd ACM international conference on Wireless sensor networks and applications [Rutman 1991] J. Rutman and F. L. Walls, Characterization of frequency stability in precision frequency sources, Proc. IEEE, vol. 79, June 1991. [Hurni 2007] P. Hurni Unsynchronized Energy-Efficient Medium Access Control and Routing in Wireless Sensor Networks Masterarbeit der Philosophisch-

0.05

0 −1

−0.8

−0.6

−0.4

−0.2

0 0.2 Position

0.4

0.6

0.8

1

Fig. R5. Optimal conditional second moments for constant R c2 = 0.36, γ = 0.0017

to account for robot mobility, assumed to take place within a predefined set of dicrete velocities. Conditional clock drift expectation and variance are associated to each predefined velocity given robot position and velocity at time t. Stationary solutions are obtained by solving associated partial integro-differential equations for infinite time. Since clock drift variance is adopted as an overall performance measure, expected conditional variance functions are chosen as objective functions for optimization along with energy consumed for message transmission. Transmission energy is controlled through transmission frequency, i.e. the average number of message transmission per time, which is allowed to depend on current robot position.

48

IFAC NetRob 2009 Golden, CO, USA, October 6-8, 2009

Naturwissenschaftlichen Fakultat der Universitat Bern, 2007 [Keshet 2001] L. Edelstein-Keshet Mathematical models of swarming and social aggregation. Paper, International Symposium on Nonlinear Theory and Its Applications (NOLTA 2001), Miyagi, Japan, 2001

49