# Universal Critical Holography and Domain Wall Formation

Tian-Chi Ma<EMAIL_ADDRESS>Center for Gravitational Physics,
Department of Space Science, Beihang University, Beijing 100191, China  Han-
Qing Shi<EMAIL_ADDRESS>Center for Gravitational Physics, Department
of Space Science, Beihang University, Beijing 100191, China  Hai-Qing Zhang
<EMAIL_ADDRESS>Center for Gravitational Physics, Department of Space
Science, Beihang University, Beijing 100191, China  Peng Huanwu Collaborative
Center for Research and Education, Beihang University, Beijing 100191, China
Adolfo del Campo<EMAIL_ADDRESS>Department of Physics and Materials
Science, University of Luxembourg, L-1511 Luxembourg, Luxembourg  Donostia
International Physics Center, E-20018 San Sebastián, Spain

###### Abstract

Using holography, we study the universal scaling laws governing the coarsening
dynamics of strongly coupled domain walls. Specifically, we studied the
universal dependence of the length of the domain wall interfaces on the quench
rate. The relation satisfies the Kibble-Zurek scaling shortly after the
critical point. However, as time goes by, the coarsening dynamics suppresses
the Kibble-Zurek scaling in favor of a universal dynamical scaling of the
characteristic length and the adiabatic growth of the system. Theoretical
predictions of the universal scaling laws are consistent with numerical
findings in both regimes for both weak and strongly coupled systems.

A complete understanding of nonequilibrium phenomena involving many-body
systems is currently lacking. This is particularly the case in systems far
from equilibrium and at strong coupling. The Kibble-Zurek mechanism (KZM) is
one of the few universal paradigms in this context [1, 2, 3, 4, 5]. It uses
equilibrium scaling theory to describe the dynamics across a continuous or
quantum phase transition. In particular, it exploits the scaling laws for the
correlation length and the relaxation time,

$\xi(\varepsilon)=\frac{\xi_{0}}{|\varepsilon|^{\nu}},\quad\tau(\varepsilon)=\frac{\tau_{0}}{|\varepsilon|^{z\nu}}.$
(1)

Here, $\nu$ and $z$ are critical exponents and $\xi_{0}$ and $\tau_{0}$ are
system-dependent microscopic constants, while
$\varepsilon=(\lambda-\lambda_{c})/\lambda_{c}$ determines the proximity of
the control parameter $\lambda$ to its critical value $\lambda_{c}$. The
divergence of the relaxation time in the neighborhood of the critical point,
known as the critical slowing down, plays a key role in the KZM. Consider
driving the system by a linear modulation of the control parameter $\lambda$,
such that $\varepsilon=t/\tau_{Q}$, where $\tau_{Q}$ is the quench time, and
the critical point is reached at $t=0$. The timescale in which the
instantaneous equilibrium relaxation time matches the time elapsed after
crossing the critical point is known as the freeze-out time and plays a key
role in KZM. It is predicted to scale universally with the quench time
according to

$\hat{t}=\left(\tau_{0}\tau_{Q}^{\nu z}\right)^{\frac{1}{1+\nu z}}.$ (2)

Using it, the KZM predicts the nonequilibrium correlation length to be given
by

$\hat{\xi}=\xi[\varepsilon(\hat{t})]=\xi_{0}\left(\frac{\tau_{Q}}{\tau_{0}}\right)^{\frac{\nu}{1+z\nu}}.$
(3)

It follows that the density of topological defects of dimension $d$ in a
system in a $D$ dimensional space is given by.

$\rho=\frac{1}{\hat{\xi}^{D-d}}=\frac{1}{\xi_{0}^{D-d}}\left(\frac{\tau_{0}}{\tau_{Q}}\right)^{\frac{(D-d)\nu}{1+z\nu}}.$
(4)

A large body of accumulated evidence supports the KZM, relying on exact
solutions, numerical simulations, and experimental studies [6, 7]. In
particular, the last decade has witnessed remarkable experimental progress
thanks to advances in quantum simulation in platforms including trapped ions
[8, 9, 10, 11] and Rydberg gases [12], the development of programmable quantum
annealing devices [13, 14, 15], and the improved probing techniques in
condensed matter experiments [16, 17]. In addition, universal physics beyond
the scope of KZM has been discovered, governing the defect number statistics
[18, 11, 19, 14, 20, 21, 15] and the spatial distribution of defects [22, 23].

Despite this progress, KZM and its generalizations have been mostly studied in
scenarios characterized by point-like topological defects with $d=0$, such as
parity-breaking kinks [24, 25, 8, 9, 26] and solitons [27] in $D=1$, as well
as vortices in $D=2$ [28, 29, 16, 21, 30, 22, 31, 23]. Evidence supporting KZM
for extended topological defects with $d\geq 1$ remains limited. In such a
setting, deviations may arise due to the enhancement of coarsening [32, 33],
the role of the embedding geometry [34], and effects arising from the spatial
configuration of the extended defect, with no counterpart in the point-like
case. The latter can stem from conformational entropy, intra-defect
interactions, and excitations of modes propagating along the extended degrees
of freedom, such as Kelvin modes in vortex strings. Numerical simulations have
described 3D $U(1)$ vortex string formation following KZM scaling laws [35],
while recent simulations of skyrmions are better explained by coarsening [36].
Experiments on Rayleigh-Bénard convection with convective rolls and cells as
defects match the KZM prediction [37]. Domain formation on colloidal
monolayers, probing the Kosterlitz–Thouless–Halperin–Nelson–Young (KTHNY)
melting scenario, exhibits important deviations from the conventional KZM due
to the nature of the transition and suggests unusually large values of the
dynamic critical exponent [38, 39]. By contrast, recent experiments on domain
wall formation of 3D Ising systems are consistent with the KZM scaling [17];
see as well [40] for supporting numerics in 2D.

Coarsening is governed by a dynamic critical exponent $z_{d}$, which is
generally different from the equilibrium value $z$. According to Biroli et al.
[41], for $t\gg\hat{t}$, the nonequilibrium correlation length scales as

$R(t)=\xi[\varepsilon(t)]^{1-\frac{z}{z_{d}}}t^{\frac{1}{z_{d}}}.$ (5)

where $\xi[\varepsilon(t)]=\varepsilon(t)^{-\nu}=(t/\tau_{Q})^{-\nu}$.

In this study, we investigate the formation of domain walls in a holographic
setting [42]. Holography provides a framework to explore strongly-coupled
field theories, in and out of equilibrium, and has been used to explore the
validity of KZM [43, 44, 21, 45]. We focus on recent developments allowing for
the description parity symmetry breaking, which results in domain formation
analogous to that in a ferromagnet. In this setting, we establish the
universality of the critical dynamics by characterizing the domain wall
average length and its statistics beyond the scope of KZM, revealing a
crossover from KZM to coarsening-dominated scaling laws.

Boundary

Black Hole

AdS4 : ($\Psi$, $M_{\mu}$)

Figure 1: Schematic figure of the holographic mapping.

Background of gravity.  The schematic setting to probe the holographic domain
wall formation is shown in Fig. 1. The gravity background we choose is the
AdS4 black brane in Eddington-Finkelstein coordinates,

$\displaystyle
ds^{2}=\frac{1}{z^{2}}\left[-f(z)dt^{2}-2dtdz+dx^{2}+dy^{2}\right],$ (6)

in which $f(z)=1-(z/z_{h})^{3}$ with $z_{h}$ the horizon position. The AdS
infinite boundary is at $z=0$ where the field theory lives. The Lagrangian of
the model we adopt is the usual Abelian-Higgs model for holographic
superconductors,

$\displaystyle\mathcal{L}=-\frac{1}{4}F_{\mu\nu}F^{\mu\nu}-|D_{\mu}\tilde{\Psi}|^{2}-m^{2}|\tilde{\Psi}|^{2},$
(7)

where $F_{\mu\nu}=\partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}$ and $A_{\mu}$
is the $U(1)$ gauge field, $D_{\mu}=\nabla_{\mu}-iA_{\mu}$, and $\tilde{\Psi}$
is the complex scalar field. We work in the probe limit, in which the
equations of motion (EoMs) read

$\displaystyle D_{\mu}D^{\mu}\tilde{\Psi}-m^{2}\tilde{\Psi}=0,$
$\displaystyle\nabla_{\mu}F^{\mu\nu}=i\left(\tilde{\Psi}^{*}D^{\nu}\tilde{\Psi}-\tilde{\Psi}(D^{\nu}\tilde{\Psi})^{*}\right).$
(8)

To transform the Lagrangian with $U(1)$ symmetry to $Z_{2}$ symmetry, we must
transform the complex scalar fields into real ones. To this end, we make the
following transformations [46, 47]

$\displaystyle\tilde{\Psi}=\Psi e^{i\lambda},\qquad
A_{\mu}=M_{\mu}+\partial_{\mu}\lambda,$ (9)

where $\Psi,M_{\mu}$ and $\lambda$ are real functions. With these real
functions, the previous Eq. (Universal Critical Holography and Domain Wall
Formation) can be rewritten as

$\displaystyle\mathbb{D}_{\mu}\mathbb{D}^{\mu}\Psi-m^{2}\Psi=0,\qquad\nabla_{\mu}F^{\mu\nu}=2M^{\nu}\Psi^{2},$
(10)

where $\mathbb{D}_{\mu}=\nabla_{\mu}-iM_{\mu}$. We use the following ansätze
for the fields: $\Psi=\Psi(t,z,x,y)$, $M_{t}=M_{t}(t,z,x,y)$,
$M_{z}=M_{z}(t,z,x,y)$, $M_{x}=M_{x}(t,z,x,y)$, and $M_{y}=M_{y}(t,z,x,y)$.
The details of the EoMs can be found in [48].

(a)(b)

Figure 2:  Time evolution of the order parameter and the birth of domain wall
in the AdS boundary. (a) Density plots of the evolving order parameter at
eight specific times with $\tau_{Q}$ = $e^{6.25}$. One can see how the order
parameters evolve from initial random configurations into domain structures.
(b) The average absolute value of the order parameter $\langle|O|\rangle$
during quenches. Each diamond corresponds to one of the snapshots in panel
(a).

Boundary conditions. The asymptotic behavior of fields near $z\rightarrow 0$
is described by $\Psi\sim{z}(\Psi_{0}(t,x,y)+\Psi_{1}(t,x,y)z+...)$,
$M_{t}\sim\mu(t,x,y)-\rho(t,x,y)z+...$, $M_{z}\sim
a_{z}(t,x,y)+b_{z}(t,x,y)z+...$, $M_{x}\sim a_{x}(t,x,y)+b_{x}(t,x,y)z+...$,
$M_{y}\sim a_{y}(t,x,y)+b_{y}(t,x,y)z+...$. We have set the scalar field mass
square as $m^{2}=-2$ and a unit AdS radius. By choosing the standard
quantization setting $\Psi_{1}\equiv 0$, it follows that $\Psi_{2}$ is related
to the condensate of the superconducting order parameter $O(t,x,y)$ in the
boundary field theory. The parameters $\mu$ and $\rho$ are interpreted as the
chemical potential and charge density, respectively, at the boundary. At the
horizon, we set $M_{t}(z_{h})=0$ and the regular finite boundary conditions
for other fields. In addition, we use the periodic boundary conditions for all
the fields along $(x,y)-$directions. The fourth-order Runge-Kutta method is
used to propagate the system with time step $\triangle t=0.01$. In the radial
direction $z$, the Chebyshev pseudo-spectral method with $21$ grid points is
used. In the $(x,y)$-direction, we set the length of the space as
$L_{x}=L_{y}=200$, and we use the Fourier decomposition with $201\times 201$
grid points. We add the random seeds of the fields in the bulk by satisfying
the distributions $\langle s(t,x_{i})\rangle=0$ and $\langle
s(t,x_{i})s(t^{\prime},x_{j})\rangle=h\delta(t-t^{\prime})\delta(x_{i}-x_{j})$,
with the amplitude $h=10^{-3}$.

Cooling the system. In the setup of a holographic superconductor, decreasing
$T$ is equivalent to increasing $\mu$. To linearly quench the temperature
across the critical point as $T(t)/T_{c}=1-t/\tau_{Q}$, we vary the chemical
potential as

$\displaystyle\mu(t)=\mu_{c}/(1-t/\tau_{Q}),$ (11)

where $\tau_{Q}$ is the quench rate and $\mu_{c}=4.06$ denotes the critical
potential in the static case. We drive the system from the initial temperature
$T_{i}=T_{c}$ to the superconducting state with the final equilibrium
temperature $T_{f}=0.8T_{c}$.

Following a thermal quench from high to low temperatures across the critical
point $T_{c}$, the system enters the superconducting phase, and domain walls
appear as predicted by the KZM. Figure 2 shows the evolution of the order
parameter in panel (a) and its average absolute value in panel (b) during the
transition from a homogeneous phase to the symmetry-breaking phase. Panel (a)
shows eight snapshots at different times of evolution for a fixed quench time
$\tau_{Q}=e^{6.25}$, each of which is marked by a diamond symbol in pane (b).
At the initial time $t/\tau_{Q}$ = 0, the order parameter exhibits a random
configuration similar to a mosaic picture and then evolves and grows rapidly
around $t/\tau_{Q}$ = 0.0859. At the time $t/\tau_{Q}$ = 0.1072, the system is
divided into domains, but the order parameter still exhibits smooth gradients.
From time $t/\tau_{Q}$ = 0.1294 to $t/\tau_{Q}$ = 1.1573, the average value of
$\langle|O|\rangle$ increases linearly and saturates at a plateau. As the
coarsening dynamics becomes more prominent, it reduces the small-scale
structures in the nonequilibrium configuration, leading to the coalescence of
domains into bigger ones.

Scaling of the average domain-wall length. In the domain formation associated
with the breaking of parity in two spatial dimensions, the number of defects
can be associated with the total domain wall length, which we denote by $L$,
averaged over an ensemble of realizations. Adapting the KZM scaling law in Eq.
(4), the generic universal dependence on the quench time $\tau_{Q}$ reduces to

$\displaystyle
L\propto\tau_{Q}^{-\left(D-d\right)\nu/\left(1+z\nu\right)}\propto\tau_{Q}^{-\nu/\left(1+z\nu\right)},$
(12)

in two spatial dimensions ($D=2$) for domain walls ($d=1$) in the AdS
boundary. From AdS/CFT, the boundary field theory is a mean-field theory. As a
result, $\nu=1/2$ and $z=2$, and thus $L\propto\tau_{Q}^{-1/4}$.

Figure 3(a) shows the relations between domain wall length $L$ and the quench
rate $\tau_{Q}$ for several values of the average absolute value of the order
parameter $\langle|O|\rangle$, using 2000 independent trajectories. Figure
3(a) indicates that for different $\langle|O|\rangle$, there exist different
scaling relations between $L$ and $\tau_{Q}$. Specifically, as
$\langle|O|\rangle=10^{-3}$, the fitted line has a power-law scaling
$L=28283\times\tau_{Q}^{-0.278}$, which is close to the KZM prediction
$-\left(D-d\right)\nu/\left(1+z\nu\right)=-1/4$ in Eq. (12). However, when
$\langle|O|\rangle=3.33$, the fitted line is $L=31571\times\tau_{Q}^{-0.488}$,
and the scaling deviates from the theoretical KZM prediction. This behavior is
attributed to coarsening, as it takes more time for a slowly quenched system
to achieve the same value of $\langle|O|\rangle$ than for a rapidly quenched
system. The coarsening process always exists after the formation of defects.
Therefore, to achieve the same $\langle|O|\rangle$ value, a slowly quenched
system must undergo a longer coarsening evolution than a rapidly quenched
system. As a result, coarsening competes with the KZM scaling.

(a)(b)

Figure 3: (a) Dependence of the domain wall length $L$ on the quench rate
$\tau_{Q}$ for the average absolute value of the order parameter
$\langle|O|\rangle$. (b) Lag time as a function of the quench rate for various
condensates. The circles, hexagons, stars, squares, and triangles are the
numerical data, while the straight lines are the best-fitted lines in the slow
quench regime.

According to KZM, the instant at which the system “freezes out” is set by
$\hat{t}$ in Eq. (2). In practice, one defines the lag-time $t_{L}$ by the
response of the system to represent the freeze-out time [49, 43, 44, 23]. In
theory, where $t_{L}\propto\tau_{Q}^{1/2}$ in the holographic setup. Figure
3(b) exhibits the relation between the lag time and the quench rate for
various condensates. When the condensate is small, $\langle|O|\rangle\approx
10^{-3}$, the scaling exponent is roughly $0.4375$, which is close to the KZM
prediction. By contrast, for large condensates (e.g., when
$\langle|O|\rangle\approx 3.33$), the scaling is roughly $0.98$, which rules
out the KZM prediction. This is an important advantage of our setting over
previous one-dimensional studies of parity-symmetry breaking in which the
power-law exponent associated with KZM [24, 50] and coarsening [26, 20] happen
to coincide, hindering the study of their interplay. We expect that the
deviations of the scaling in the $t_{L}\sim\tau_{Q}$ relation are also due to
the coarsening in the late stage of the evolution.

Domain coarsening. Coarsening involves gradient flow dynamics. The relaxation
of a non-equilibrium system follows the steepest descent in the energy
landscape. Adapting the prediction by Biroli et al. [41, 33], the domain wall
length scales as

$L\propto\frac{1}{R(t)^{D-d}}\propto\tau_{Q}^{\frac{(D-d)\nu}{z_{d}}(z-z_{d})}t^{-\frac{(D-d)}{z_{d}}[1+\nu(z-z_{d})]}.$
(13)

For $z=z_{d}=2$, $L\propto t^{-\frac{(D-d)}{z_{d}}}\propto t^{-1/2}$ and
$R(t)\propto t^{1/2}$, which is in agreement with the results in [51]. In Fig.
4(a), we present the length of domain walls $L$ as a function of $t$ for
several values of $\tau_{Q}$. Long after the quench, it can be found that
regardless of the value of $\tau_{Q}$, the time evolution of $L$ follows a
similar power law and is almost parallel, in log-log scale, to the reference
line $L\sim t^{-1/2}$, which agrees with the theoretical prediction in Eq.
(13). At this late stage, the system has reached a state of local equilibrium,
and coarsening dynamics completely govern the evolution of the system.

(a)(b)

Figure 4: (a) Time evolution of domain wall length for various quench rate
$\tau_{Q}$. The dashed line is the theoretical reference line with a slope of
$-1/2$; (b) Average condensate for slow quenches. The linear relation reflects
the adiabatic growth.

Fig. 4(b) shows that the condensates transit from exponential growth to a
common linear growth as a function of $t/\tau_{Q}$, for large $\tau_{Q}$. This
is an implication of the adiabatic evolution [43]. For large values of
$\tau_{Q}$, the condensate scales with the same function of $t/\tau_{Q}$. In
this regime, when $\langle|O|\rangle$ is large (e.g.,
$\langle|O|\rangle=3.33$) the lag-time is linearly proportional to $\tau_{Q}$.
This is consistent with the numerical relation in Fig. 3(b) that
$t_{L}\propto\tau_{Q}^{0.98}$ as $\tau_{Q}$ is large. From the relations
$L\sim t^{-1/2}$ and $t\sim\tau_{Q}$ after the coarsening dynamics, we deduce
that $L\sim\tau_{Q}^{-1/2}$, which is consistent with the relation between $L$
and $\tau_{Q}$ in Fig. 3(a) as $\langle|O|\rangle=3.33$.

(a)(b)

Figure 5: (a) Histogram of the probability of the domain wall length, for
different values of $\tau_{Q}$ at $\langle|O|\rangle=10^{-3}$, are well fitted
by a Gaussian distribution. (b) Histogram of the domain wall length
probability for different values of $\tau_{Q}$ at $\langle|O|\rangle=3.33$.
Coarsening-induced deviations from the Gaussian fit arise in the slow quench
limit, leading to a multimodal distribution.

Additional information is revealed by the domain wall statistics. Fig. 5(a)
shows the histogram of the probability of the domain wall lengths for three
kinds of quench rates, and a small condensate $\langle|O|\rangle=10^{-3}$. In
this regime, the nonequilibrium dynamics is consistent with the KZM. The
universality of average quantities with the quench rate carries over the full
distribution, according to the recent extension beyond KZM [18, 19]. Numerical
evidence to date supporting such extensions is limited to point-like defects
in one and two spatial dimensions [18, 11, 19, 14, 20, 21, 15]. For large
system sizes, the distribution approaches a normal distribution, justifying
the accuracy of the Gaussian fit in Fig. 5(a). This provides the first
numerical evidence regarding the universal character of the distribution of
domain walls as extended topological defects with $d=1$. Such universality
remains in larger condensates as shown in Fig. 5(b) for
$\langle|O|\rangle=3.33$ for quenches $\tau_{Q}=e^{6.5}$ and
$\tau_{Q}=e^{7.5}$. However, the distribution for slow quenches
$\tau_{Q}=e^{9.25}$ is dominated by coarsening; it becomes bimodal and no
longer universal, with increasing probability near the origin, with low defect
numbers.

Discussion and Conclusions. We have established the universal scaling laws
governing the critical dynamics leading to strongly coupled domain walls via
the AdS/CFT correspondence. The KZM describes the universal scaling of the
domain wall length with the quench time near the critical regime. Far away
from this regime, a different scaling holds that is no longer consistent with
KZM and reflects the coarsening dynamics of domain walls. A theoretical
analysis yields a universal prediction for the scaling dominated by coarsening
that is in agreement with numerical simulations. In [48], we further show the
similar phenomenon of domain wall in the Ginzburg-Landau model, which is a
weakly coupled theory. Therefore, we can conclude that our findings of the new
scalings in the far-from critical regime may apply universally in the
coarsening dynamics.

Note: Upon the completion of this work, Ref. [52] reported KZM deviations due
to coarsening in the digital-analog quantum simulation of a 2D XY ferromagnet.

Acknowledgements. This work was partially supported by the National Natural
Science Foundation of China (Grants No. 12175008). This research was funded in
part by the Luxembourg National Research Fund (FNR), grant No. 17132060. For
open access, the authors have applied a Creative Commons Attribution 4.0
International (CC BY 4.0) license to any Author Accepted Manuscript version
arising from this submission.

## References

  * Kibble [1976] T. W. B. Kibble, J. of Phys. A: Math. Gen. 9, 1387 (1976).
  * Kibble [1980] T. W. B. Kibble, Phys. Reports 67, 183 (1980).
  * Zurek [1985] W. H. Zurek, Nature 317, 505 (1985).
  * Zurek [1993a] W. H. Zurek, Acta Phys. Pol. B 24, 1301 (1993a).
  * Zurek [1993b] W. H. Zurek, Phys. Reports 276, 177 (1993b).
  * Dziarmaga [2010] J. Dziarmaga, Advances in Physics 59, 1063 (2010).
  * del Campo and Zurek [2014] A. del Campo and W. H. Zurek, Int. J. Mod. Phys. A 29, 1430018 (2014).
  * Ulm  _et al._ [2013] S. Ulm, J. Roßnagel, G. Jacob, C. Degünther, S. T. Dawkins, U. G. Poschinger, R. Nigmatullin, A. Retzker, M. B. Plenio, F. Schmidt-Kaler, and K. Singer, Nat. Comm. 4, 2290 (2013).
  * Pyka  _et al._ [2013] K. Pyka, J. Keller, H. L. Partner, R. Nigmatullin, T. Burgermeister, D. M. Meier, K. Kuhlmann, A. Retzker, M. B. Plenio, W. H. Zurek, A. del Campo, and T. E. Mehlstäubler, Nat. Comm. 4, 2291 (2013).
  * Cui  _et al._ [2016] J.-M. Cui, Y.-F. Huang, Z. Wang, D.-Y. Cao, J. Wang, W.-M. Lv, L. Luo, A. del Campo, Y.-J. Han, C.-F. Li, and G.-C. Guo, Sci. Rep. 6, 33381 (2016).
  * Cui  _et al._ [2020] J.-M. Cui, F. J. Gómez-Ruiz, Y.-F. Huang, C.-F. Li, G.-C. Guo, and A. del Campo, Communications Physics 3, 44 (2020).
  * Bernien  _et al._ [2017] H. Bernien, S. Schwartz, A. Keesling, H. Levine, A. Omran, H. Pichler, S. Choi, A. S. Zibrov, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Nature 551, 579 (2017).
  * Weinberg  _et al._ [2020] P. Weinberg, M. Tylutki, J. M. Rönkkö, J. Westerholm, J. A. Åström, P. Manninen, P. Törmä, and A. W. Sandvik, Phys. Rev. Lett. 124, 090502 (2020).
  * Bando  _et al._ [2020] Y. Bando, Y. Susa, H. Oshiyama, N. Shibata, M. Ohzeki, F. J. Gómez-Ruiz, D. A. Lidar, S. Suzuki, A. del Campo, and H. Nishimori, Phys. Rev. Research 2, 033369 (2020).
  * King  _et al._ [2022] A. D. King, S. Suzuki, J. Raymond, A. Zucca, T. Lanting, F. Altomare, A. J. Berkley, S. Ejtemaee, E. Hoskinson, S. Huang, E. Ladizinsky, A. J. R. MacDonald, G. Marsden, T. Oh, G. Poulin-Lamarre, M. Reis, C. Rich, Y. Sato, J. D. Whittaker, J. Yao, R. Harris, D. A. Lidar, H. Nishimori, and M. H. Amin, Nature Physics 18, 1324 (2022).
  * Lin  _et al._ [2014] S.-Z. Lin, X. Wang, Y. Kamiya, G.-W. Chern, F. Fan, D. Fan, B. Casas, Y. Liu, V. Kiryukhin, W. H. Zurek, C. D. Batista, and S.-W. Cheong, Nature Physics 10, 970 (2014).
  * Du  _et al._ [2023] K. Du, X. Fang, C. Won, C. De, F.-T. Huang, W. Xu, H. You, F. J. Gómez-Ruiz, A. del Campo, and S.-W. Cheong, Nature Physics 19, 1495 (2023).
  * del Campo [2018] A. del Campo, Phys. Rev. Lett. 121, 200601 (2018).
  * Gómez-Ruiz  _et al._ [2020] F. J. Gómez-Ruiz, J. J. Mayo, and A. del Campo, Phys. Rev. Lett. 124, 240602 (2020).
  * Mayo  _et al._ [2021] J. J. Mayo, Z. Fan, G.-W. Chern, and A. del Campo, Phys. Rev. Research 3, 033150 (2021).
  * del Campo  _et al._ [2021] A. del Campo, F. J. Gómez-Ruiz, Z.-H. Li, C.-Y. Xia, H.-B. Zeng, and H.-Q. Zhang, Journal of High Energy Physics 2021, 61 (2021).
  * del Campo  _et al._ [2022] A. del Campo, F. J. Gómez-Ruiz, and H.-Q. Zhang, Phys. Rev. B 106, L140101 (2022).
  * Thudiyangal and del Campo [2024] M. Thudiyangal and A. del Campo, Universal vortex statistics and stochastic geometry of bose-einstein condensation (2024), arXiv:2401.09525 [cond-mat.quant-gas] .
  * Laguna and Zurek [1998] P. Laguna and W. H. Zurek, Phys. Rev. D 58, 085021 (1998).
  * Sabbatini  _et al._ [2011] J. Sabbatini, W. H. Zurek, and M. J. Davis, Phys. Rev. Lett. 107, 230402 (2011).
  * Krapivsky [2010] P. L. Krapivsky, Journal of Statistical Mechanics: Theory and Experiment 2010, P02014 (2010).
  * Damski and Zurek [2010] B. Damski and W. H. Zurek, Phys. Rev. Lett. 104, 160404 (2010).
  * Weiler  _et al._ [2008] C. N. Weiler, T. W. Neely, D. R. Scherer, A. S. Bradley, M. J. Davis, and B. P. Anderson, Nature 455, 948 (2008).
  * Griffin  _et al._ [2012] S. M. Griffin, M. Lilienblum, K. T. Delaney, Y. Kumagai, M. Fiebig, and N. A. Spaldin, Phys. Rev. X 2, 041022 (2012).
  * Goo  _et al._ [2021] J. Goo, Y. Lim, and Y. Shin, Phys. Rev. Lett. 127, 115701 (2021).
  * Rabga  _et al._ [2023] T. Rabga, Y. Lee, D. Bae, M. Kim, and Y. Shin, Phys. Rev. A 108, 023315 (2023).
  * Libál  _et al._ [2020] A. Libál, A. del Campo, C. Nisoli, C. Reichhardt, and C. J. O. Reichhardt, Phys. Rev. Res. 2, 033433 (2020).
  * Samajdar and Huse [2024] R. Samajdar and D. A. Huse, Quantum and classical coarsening and their interplay with the kibble-zurek mechanism (2024), arXiv:2401.15144 [quant-ph] .
  * Stoop and Dunkel [2018] N. Stoop and J. Dunkel, Soft Matter 14, 2329 (2018).
  * Antunes  _et al._ [1999] N. D. Antunes, L. M. A. Bettencourt, and W. H. Zurek, Phys. Rev. Lett. 82, 2824 (1999).
  * Reichhardt and Reichhardt [2023] C. Reichhardt and C. J. O. Reichhardt, Phys. Rev. Res. 5, 033221 (2023).
  * Casado  _et al._ [2006] S. Casado, W. González-Viñas, and H. Mancini, Phys. Rev. E 74, 047101 (2006).
  * Deutschländer  _et al._ [2015] S. Deutschländer, P. Dillmann, G. Maret, and P. Keim, Proc. Nat. Acad. Sciences 112, 6925 (2015).
  * del Campo [2015] A. del Campo, Proceedings of the National Academy of Sciences 112, 6780 (2015).
  * Dziarmaga [2004] J. Dziarmaga, Acta Phys. Pol. B 35, 2205 (2004).
  * Biroli  _et al._ [2010] G. Biroli, L. F. Cugliandolo, and A. Sicilia, Phys. Rev. E 81, 050101 (2010).
  * Hartnoll  _et al._ [2018] S. Hartnoll, A. Lucas, and S. Sachdev,  _Holographic Quantum Matter_ , The MIT Press (MIT Press, 2018).
  * Chesler  _et al._ [2015] P. M. Chesler, A. M. Garcia-Garcia, and H. Liu, Phys. Rev. X 5, 021015 (2015).
  * Sonner  _et al._ [2015] J. Sonner, A. del Campo, and W. H. Zurek, Nat. Comm. 6, 7406 (2015).
  * Xia  _et al._ [2023] C.-Y. Xia, H.-B. Zeng, C.-M. Chen, and A. del Campo, Phys. Rev. D 108, 026017 (2023).
  * Horowitz  _et al._ [2011] G. T. Horowitz, J. E. Santos, and B. Way, Phys. Rev. Lett. 106, 221601 (2011), arXiv:1101.3326 [hep-th] .
  * Li  _et al._ [2023] Z.-H. Li, H.-Q. Shi, and H.-Q. Zhang, Phys. Rev. D 108, 106015 (2023), arXiv:2207.10995 [hep-th] .
  * [48] See the Supplemental Material.
  * Das  _et al._ [2012] A. Das, J. Sabbatini, and W. H. Zurek, Sci. Rep. 2, 352 (2012).
  * Chiara  _et al._ [2010] G. D. Chiara, A. del Campo, G. Morigi, M. B. Plenio, and A. Retzker, New Journal of Physics 12, 115003 (2010).
  * Bray [1994] A. J. Bray, Adv. Phys. 43, 357 (1994).
  * Andersen and et al. [2024] T. I. Andersen and et al., Thermalization and criticality on an analog-digital quantum simulator (2024), arXiv:2405.17385 [quant-ph] .

##  Appendix A Appendix A: Explicit forms of the equations of motions

The EoMs (10) in the main text can be decomposed into

$\displaystyle\nabla_{\mu}\nabla^{\mu}\Psi-M_{\mu}M^{\mu}\Psi-m^{2}\Psi$
$\displaystyle=$ $\displaystyle 0,$ (S1)
$\displaystyle\left(\nabla_{\mu}M^{\mu}\right)\Psi+2M^{\mu}\nabla_{\mu}\Psi$
$\displaystyle=$ $\displaystyle 0,$ (S2) $\displaystyle\nabla_{\mu}F^{\mu\nu}$
$\displaystyle=$ $\displaystyle 2M^{\nu}\Psi^{2}.$ (S3)

The Eqs. (S1) and (S2) follow respectively from the real and imaginary parts
of the scalar equations in Eq.(10). We note that Eqs. (S2) and (S3) are not
independent, since

$\displaystyle\nabla_{\nu}(\nabla_{\mu}F^{\mu\nu})=0=2(\nabla_{\nu}M^{\nu})\Psi^{2}+4M^{\nu}\Psi(\nabla_{\nu}\Psi).$
(S4)

Dividing $2\Psi$ on the right side of the Eq. (S4), we can readily get the
left side of the Eq. (S2). Therefore, there are five real functions with five
independent equations, which implies that our ansätze for the fields are
consistent. Using them in the frame of the line-element (6), the EoMs involve
the following contributions:
1) The gauge fields (S3) part,

$\displaystyle 0$ $\displaystyle=$
$\displaystyle-\frac{2\Psi^{2}M_{t}}{z^{2}}+\partial^{2}_{x}M_{t}+\partial^{2}_{y}M_{t}+f\partial^{2}_{z}M_{t}-\partial_{tx}M_{x}-\partial_{ty}M_{y}-\partial_{tz}M_{t}-f\partial_{tz}M_{z}+\partial^{2}_{t}M_{z},$
(S5) $\displaystyle 0$ $\displaystyle=$
$\displaystyle-\frac{2\Psi^{2}M_{z}}{z^{2}}+\partial^{2}_{x}M_{z}+\partial^{2}_{y}M_{z}-\partial_{zx}M_{x}-\partial_{zy}M_{y}+\partial^{2}_{z}M_{t}-\partial_{tz}M_{z},$
(S6) $\displaystyle 0$ $\displaystyle=$
$\displaystyle-\frac{2\Psi^{2}M_{x}}{z^{2}}+\partial^{2}_{y}M_{x}-f^{\prime}\partial_{x}M_{z}-\partial_{xy}M_{y}+f^{\prime}\partial_{z}M_{x}+\partial_{zx}M_{t}-f\partial_{zx}M_{z}+f\partial^{2}_{z}M_{x}+\partial_{tx}M_{z}-2\partial_{tz}M_{x},$
(S7) $\displaystyle 0$ $\displaystyle=$
$\displaystyle-\frac{2\Psi^{2}M_{y}}{z^{2}}+\partial^{2}_{x}M_{y}-f^{\prime}\partial_{y}M_{z}-\partial_{xy}M_{x}+f^{\prime}\partial_{z}M_{y}+\partial_{zy}M_{t}-f\partial_{zy}M_{z}+f\partial^{2}_{z}M_{y}+\partial_{ty}M_{z}-2\partial_{tz}M_{y}.$
(S8)

2) The real part of scalar fields (S1),

$\displaystyle 0$ $\displaystyle=$
$\displaystyle-\frac{m^{2}\Psi}{2z^{2}}-\frac{1}{2}\Psi
M_{x}^{2}-\frac{1}{2}\Psi M_{y}^{2}+\Psi M_{t}M_{z}-\frac{1}{2}\Psi
fM_{z}^{2}+\frac{1}{2}\partial^{2}_{x}\Psi+\frac{1}{2}\partial^{2}_{y}\Psi$
(S9)
$\displaystyle-\frac{f\partial_{z}\Psi}{z}+\frac{1}{2}\partial_{z}(f\partial_{z}\Psi)+\frac{\partial_{t}\Psi}{z}-\partial_{tz}\Psi.$

3) The imaginary part of scalar fields (S2),

$\displaystyle 0$ $\displaystyle=$
$\displaystyle\frac{2\Psi}{z}(fM_{z}-M_{t})-\Psi\left(\partial_{z}(fM_{z})+\partial_{x}M_{x}+\partial_{y}M_{y}-\partial_{t}M_{z}-\partial_{z}M_{t}\right)$
(S10) $\displaystyle+2\left(M_{t}\partial_{z}\Psi+M_{z}\partial_{t}\Psi-
M_{x}\partial_{x}\Psi-M_{y}\partial_{y}\Psi-fM_{z}\partial_{z}\Psi\right),$

in which $f^{\prime}=f^{\prime}(z)$.

##  Appendix B Appendix B: Domain wall formations in a time-dependent
Ginzburg-Landau model

We start with a time-dependent Ginzburg-Landau model (TDGL) of the real non-
conserved scalar field $\phi$ undergoing Langevin dynamics

$\ddot{\phi}+\eta\dot{\phi}-\partial_{xx}\phi-\partial_{yy}\phi+\partial_{\phi}V\left(\phi\right)=\xi\left(x,y,t\right),$
(S11)

where $\eta$ is a global damping constant, $V(\phi)$ is a generic Ginzburg-
Landau potential given by
$V\left(\phi\right)=\frac{1}{8}\left(\phi^{4}-2\epsilon\phi^{2}+1\right)$ and
$\xi\left(x,y,t\right)$ is a real white noise with zero mean. This potential
has two possible minima given by $\phi_{\rm min}=\pm\sqrt{\epsilon}$. The
system undergoes a second-order phase transition as $\epsilon$ changes sign
from negative to positive. For simplicity, we linearly quench $\epsilon$ as
$\epsilon(t)=t/\tau_{Q}$, where $\tau_{Q}$ is the quench rate. Specifically,
we quench the system from $t/\tau_{Q}=0$ to $t/\tau_{Q}=10$ and then maintain
it at $\epsilon_{f}=10$. When the system reaches a state of local equilibrium,
the order parameter saturates at the values
$\phi=\pm\sqrt{\epsilon_{f}}=\pm\sqrt{10}$. It is not possible to solve the
Eq. (S11) analytically, and we resort to a numerical approach. In the
numerics, we set $\eta=1$ and the amplitude of the noise equal to $10^{-3}$.
The length of the system is $L_{x}=L_{y}=200$, and we use the Fourier
decomposition with $201\times 201$ grid points for the periodic boundary
conditions. We use the fourth-order Runge-Kutta method with the time step
$\triangle t=0.01$ in the time direction.

(a)(b)

Figure S1: Time evolution of the order parameter and the birth of domain walls
in the TDGL model. (a) Density plots of the evolving order parameter at eight
specific evolution times with $\tau_{Q}$ = $e^{4.5}$. One can see how the
order parameters evolve from initial random configurations into domain
structures. (b) The average absolute value of the order parameter
$\langle|O|\rangle$ during quenches. The eight diamond symbols correspond to
the eight snapshots in panel (a).

Figure S1(a) shows eight snapshots of the time evolution of the domain
profiles of the order parameter. Their corresponding average values
$\langle|O|\rangle$ are shown as eight diamonds in the right panel (b) in Fig.
S1. At $t/\tau_{Q}=0$, the scalar field takes very small random values that
serve as the inhomogeneous seeds for the time evolution of the system. At the
early time $t/\tau_{Q}=0.1734$ and $t/\tau_{Q}=0.2823$, domains formed because
of the KZM. From $t/\tau_{Q}=0.3756$ to $t/\tau_{Q}=3.6639$, due to the
influence of coarsening, it can be found that the domains merge and form a
larger domain.

(a)(b)

Figure S2: (a) Double logarithmic relation between the domain wall length $L$
and quench rate $\tau_{Q}$ for the average absolute value of the order
parameter $\langle|O|\rangle$. (b) The relation between the lag and quench
times for various condensates. The circles, hexagons, stars, squares, and
triangles are the numerical data, while the straight lines are the best-fitted
lines in the slow quench regime (larger $\tau_{Q}$).

Figure S2(a) displays the length of the domain walls as a function of the
quench rate for several values of the average absolute value of the order
parameter $\langle|O|\rangle$. As in the holographic case shown in Fig. 3 in
the main text, we see that when $\langle|O|\rangle$ is small (e.g., for
$\langle|O|\rangle=10^{-3}$), the power-law exponent for slow quenches is
$-0.246$, which is close to the KZM prediction $-1/4$ in Eq. (12) in the main
text. However, for larger condensates (say, when $\langle|O|\rangle=2.67$),
the scaling becomes $-0.499$, which rules out the KZM scaling. Fig. S2(b)
exhibits the dependence of the lag time on the quench rate for various
condensates. As in the case in holography (see Fig. 3(b)), for small
condensates and fast quenches, the scaling relation is
$t_{L}\sim\tau_{Q}^{0.5012}$, which agrees with theoretical predictions Eq.
(2) in the main text. However, as $\langle|O|\rangle$ is large, this relation
becomes $L\sim\tau_{Q}^{0.9923}$, which differs from the KZM prediction.
Therefore, as the system is close to the critical point, the scaling relations
are consistent with the KZM. However, as the system is driven away from the
critical point, the scaling relations deviate from the KZM predictions. We
explain this phenomenology in the following.

(a)(b)

Figure S3: (a) Time evolution of domain wall length for various quench rates
$\tau_{Q}=e^{2},e^{5}$ and $e^{8}$. The black dashed line is the theoretical
reference line with a slope of $-1/2$. (b) The average condensate as a
function of $t/\tau_{Q}$ for relatively slow quenches. The overlap at the late
stage indicates the adiabatic growth of $\langle|O|\rangle$.

Figure S3(a) shows the time evolution of domain wall length $L$ for three
different values of $\tau_{Q}$. In the later stage of the evolution, the
scaling of $L$ is consistent with the theoretical relation $L\sim t^{-1/2}$,
which is similar to the case in holography in the main text. Figure S3(b)
shows the time evolution of the average condensate $\langle|O|\rangle$ for
slow quenches. As $\langle|O|\rangle$ is relatively large, the condensates for
various quenches overlap due to the adiabatic growth for the slow quench.
Therefore, as $\langle|O|\rangle$ is relatively large, such as for
$\langle|O|\rangle=2.67$, the overlapping of the lines indicates that the lag-
time is linearly proportional to the quench time, i.e., $t_{L}\sim\tau_{Q}$,
for the slow quenches. This explains the linear relation between $t_{L}$ and
$\tau_{Q}$ in the Fig. S2(b). Furthermore, from the relations $L\sim t^{-1/2}$
and $t_{L}\sim\tau_{Q}$ in the large $\tau_{Q}$ and large $\langle|O|\rangle$,
the scaling $L\sim\tau_{Q}^{-1/2}$ explains the numerical results presented in
Fig. S2(a).

(a)(b)

Figure S4: (a) Histogram of the domain wall lengths with small condensates
$\langle|O|\rangle=10^{-3}$; (b) Histogram of the domain wall lengths with
large condensates $\langle|O|\rangle=2.67$. The data involves $2000$
independent simulations.

Figure S4(a) exhibits the histogram of the length of the domain wall for
various quench rates with small condensates $\langle|O|\rangle=10^{-3}$. For
the three values of $\tau_{Q}$ (i.e., $\tau_{Q}=e^{5},e^{6}$ and $e^{8}$), the
distributions of the domain wall length with large samplings satisfy the
Gaussian distribution. Panel (b) of Fig. S4 shows the histogram of the domain
wall lengths with large condensate $\langle|O|\rangle=2.67$. When the quench
is not too slow, the distributions satisfy the Gaussian distribution. In the
opposite limit (e.g., for $\tau_{Q}=e^{8}$, and for a large condensate), the
distribution strongly deviates from the Gaussian approximation if the quench
is very slow, as coarsening leads to the merging of domains. As a result,
events with no domain walls at all become likely, resulting in a multimodal
distribution.


# Modeling Interactions of Multimodal Road Users in Shared Spaces

Fatema T. Johora1 and Jörg P. Müller2 This research has been supported by the
German Research Foundation (DFG) through the Research Training Group
SocialCars: Cooperative (De-) centralized Traffic Management (GRK 1931). We
cordially thank Federico Pascucci and Bernhard Friedrich for providing road
user datasets.This paper has been published with copyright by IEEE, accessible
at https://ieeexplore.ieee.org/document/8569687.

###### Abstract

In shared spaces, motorized and non-motorized road users share the same space
with equal priority. Their movements are not regulated by traffic rules, hence
they interact more frequently to negotiate priority over the shared space. To
estimate the safeness and efficiency of shared spaces, reproducing the traffic
behavior in such traffic places is important. In this paper, we consider and
combine different levels of interaction between pedestrians and cars in shared
space environments. Our proposed model consists of three layers: a layer to
plan trajectories of road users; a force-based modeling layer to reproduce
free flow movement and simple interactions; and a game-theoretic decision
layer to handle complex situations where road users need to make a decision
over different alternatives. We validate our model by simulating scenarios
involving various interactions between pedestrians and cars and also car-to-
car interaction. The results indicate that simulated behaviors match observed
behaviors well.

##  I INTRODUCTION

In recent years, the shared space design paradigm introduced by Dutch traffic
engineer Hans Monderman in the 1970s [1], has received considerable attention
as an alternative approach to classic traffic design. Shared spaces are mixed
traffic environments, in which heterogeneous road users such as pedestrians,
cyclists, and motorized vehicles are not separated by temporal and spatial
division of road resources, but share the same space; their movements are
conducted based on social protocols (e.g. courtesy behavior) and informal
rules (e.g. colored road surface) instead of traffic regulations. They
interact more frequently to coordinate trajectories, by performing actions
they choose through negotiating with others.

Despite the ongoing debate on efficiency (average road user delays and road
capacity) and safeness of shared space, many cities in Europe already have
adopted the shared space design principles to reconstruct city centers [2].
This makes the investigation of traffic performance in shared spaces critical.
Modeling and analyzing the motion behavior of multimodal road users including
interactions among them can reproduce the operation of shared spaces, which
provides a basis for measuring aspects of their performance.

Modeling mixed traffic interactions is challenging as it is the result of
complex human decision-making processes. To the best of our knowledge, until
now, only little work has been done on modeling and simulating shared spaces.
Previous works on shared spaces modeling use a physics-based model, most
prominently the classical social force model (SFM) [3] to model pedestrian
movement behaviors and extend SFM using long-range conflict avoidance
mechanisms [4] [5], rule-based constraints [6], a multinomial logit model [7],
or a game theoretic approach[8] for modeling multimodal interaction. Cellular
automata are also used to model mixed traffic in presence of traffic
regulations like in [9], for modeling car-following and lane-changing
behaviors of cars and motorcycles or in [10], for modeling pedestrians-
vehicles interactions at crosswalks.

All these works cover interactions among multimodal road users at a certain
level of details. A restriction of the state of the art is that only single
bilateral conflicts are considered. This means that, at a time, each road user
can only handle a single conflict with another user, but not all conflicts she
is encountering at that time. We define a conflict in accordance with [11] as
“an observable situation in which two or more road users approach each other
in time and space to such an extent that there is a risk of collision, if
their movements remain unchanged”. These models also do not consider courtesy
behavior of car drivers for car-car interactions.

This paper describes a first step to address this gap, by proposing a
multiagent-based model for more realistically describing, recognizing and
dealing with real-world interactions of pedestrians and cars in shared spaces.
It includes interactions involving two or more road users and the courtesy
behavior of car drivers.

##  II Problem Statement and Requirements

As there are no explicit traffic regulations in shared spaces, road users
interact with each other for negotiating the priority of the integrated space.
To understand their interactions, we analyze the data of a shared space area
in Hamburg [7]. Based on our observation of the video data and on the
classification of road users’ behaviors given by Helbing and Molnàr [3], we
classify typical interactions into the following categories:

  * •

Simple Interaction: direct mapping of percepts to actions.

    * –

Reactive Interaction: When road users do not have time to think or plan their
actions, they behave re-actively. As an example, if a pedestrian suddenly
jumps in front of the vehicle, to avoid a serious collision, the car driver
will either break strongly or steer away immediately.

    * –

Car Following (car only): Even-though in shared space, there is no defined
lane for cars, but empirical observation indicates that cars merge into
assumed lanes and follow the car in front [6].

  * •

Complex Interaction: road users need to decide an action among different
alternatives.

    * –

Implicit Interaction between Two Road Users: Each road user chooses her best
action to avoid conflicts with another user by predicting their action. E.g.,
in reality, when a car comes with a moderate speed, a pedestrian may predict
that if she continues to walk, the car driver will decelerate and let her go.

    * –

Implicit Interaction among Multiple Road Users: Sometimes road user needs to
implicitly interact with multiple users at the same time as in pedestrian
crossing, where pedestrian needs to interact with cars coming from both
directions to make a precise decision to cross or not to cross the road.

    * –

Explicit Interaction: Road users also interact with other road users through
communicating.

Road users interact by performing some actions. Courtesy behavior is one of
them, where a road user decelerates and eventually stops to give priority to
others. We model this behavior both as reactive and complex interaction.

The simple interactions can be represented using forced-based [3, 6] or
cellular automata models [12, 13]. However, for complex interactions,
decision-theoretic models such as games are useful[3]. In a noncooperative
game-theoretic model, each decision maker makes a decision by predicting
others’ decision. We argue that such models can more precisely model implicit
interaction between road users. Indeed, they have been used to describe
merging-give way interaction [14], lane-changing interaction of cars [15], or
cyclist-car driver interaction at zebra crossings [16]. We note that modeling
explicit interaction is also missing in previous works. However, we do not
consider it in this paper.

##  III Multiagent-Based Simulation Model

We model the motion behaviors of multimodal road users, including the
interactions among them (discussed in Section II) in three layers: trajectory
planning, force-based modeling, and game-theoretic decision-making. Figure 1
summarizes the roles of each layer.

Figure 1: The conceptual structure of motion behaviors of multimodal road
users. Here, AF means Additional Force and SFM means Social Force Model
(adapted from [6],p.88).

###  III-A Trajectory Planning Layer

We plan the trajectories of road users by considering static obstacles like
walls or trees in the shared space environment. To perform a trajectory
planning algorithm, the environment needs to be represented as a graph. We
transform the environment into a visibility graph: where each edge connects
two outline vertices of any obstacle, which are visible to each other and
these edges do not collide with any line segment that outlines any
obstacle[17]. The origin and destination points of each road user are added to
the graph. We perform the A* algorithm[18] on the graph to plan free flow
trajectory for each road user. However, humans keep some distance from
obstacles. Hence, to model this behavior we fine-tune the placement of the
inner path vertices regarding their distance to obstacles.

###  III-B Force-Based Modeling Layer

We use the classical social force model (see also [19]) to model the driving
force of road users towards their destination, the interaction of road users
with static obstacles, the interaction between pedestrians, and also the
overtaking behavior of cars. In this model, the motion of a road user is
conducted by a set of simple forces that represent various influences she
faces in simple situations:

$\frac{d\vec{v}_{\alpha}(t)}{dt}=\vec{f}_{\alpha}^{o}+\sum\vec{f}_{\alpha
B}^{soc}+\sum_{\beta\neq V}\vec{f}_{\alpha
V}^{soc}+\color[rgb]{0,0,1}\vec{I}_{\alpha\beta}^{stopping}+\color[rgb]{0,0,1}\vec{I}_{\alpha\beta}^{courtesy}+\vec{I}_{\alpha\beta}^{following}+\vec{I}_{\alpha\beta}^{game}$

Here, $\vec{f}_{\alpha}^{o}$ represents the driving force of road user
$\alpha$ and can be computed by the formula:
$\vec{f}_{\alpha}^{o}(\vec{v^{*}}_{\alpha}(t),\vec{v}_{\alpha}(t))=\frac{\vec{v^{*}}_{\alpha}(t)-\vec{v}_{\alpha}(t)}{\tau}$
for a relaxation time $\tau$ and $\vec{v^{*}}_{\alpha}(t)$ and
$\vec{v}_{\alpha}(t)$ are the desired and current velocities of $\alpha$
respectively. Forces $\vec{f}_{\alpha\beta}^{soc}$ and $\vec{f}_{\alpha
B}^{soc}$ represent the human behavior of maintaining some distance from other
road users $\beta$ and static obstacles B. If $\alpha$ is a pedestrian then
$\vec{f}_{\alpha\beta}^{soc}$ applies between $\alpha$ and every other road
users, but if $\alpha$ is a car then $\vec{f}_{\alpha\beta}^{soc}$ only
applies between $\alpha$ and every other cars. These behaviors are formulated
by the following equations:
$f_{\alpha\beta}^{soc}(\vec{d}_{\alpha\beta}(t))=V_{\alpha\beta}^{o}\exp\bigg{[}\frac{-\vec{d}_{\alpha\beta}(t)}{\sigma}\bigg{]}\hat{n}_{\alpha\beta}F_{\alpha\beta}$
and $\vec{f}_{\alpha B}^{soc}(\vec{d}_{\alpha B}(t))=U_{\alpha
B}^{o}\exp\bigg{[}\frac{-\vec{d}_{\alpha B}(t)}{R}\bigg{]}\hat{n}_{\alpha B}$,
where $V_{\alpha\beta}^{o}$ and $U_{\alpha B}^{o}$ denote the interaction
strengths, and $\sigma$ and $R$ indicate the range of these repulsive
interactions. $\vec{d}_{\alpha\beta}(t)$ and $\vec{d}_{\alpha B}(t)$ are the
distances from $\alpha$ to $\beta$, or $\alpha$ to $B$ at a specific time.
$\hat{n}_{\alpha\beta}$ and $\hat{n}_{\alpha B}$ denote the normalized
vectors. $F_{\alpha\beta}$ represents the fact that human are mostly
influenced by the objects which can be captured within their field of view.

We extend the social force model to consider reactive interactions of cars
towards pedestrians ($\color[rgb]{0,0,1}\vec{I}_{\alpha\beta}^{stopping}$),
car following interactions
($\color[rgb]{0,0,1}\vec{I}_{\alpha\beta}^{following}$) and courtesy behavior
of car drivers ($\color[rgb]{0,0,1}\vec{I}_{\alpha\beta}^{courtesy}$):

##### Car to Pedestrian Reactive Interaction:

If the nearest conflict of a car driver $\alpha$ is with pedestrian(s) $\beta$
and $\beta$ has already started moving in front of $\alpha$ or is very close
to $\alpha$, then $\alpha$ decelerates to let $\beta$ pass. This type of
interaction only happens if $\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$
$\leq$ 15 $\parallel$ $\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\geq$
345.

##### Car Following

If $d_{\alpha\beta}$ $\geq$ $d_{\alpha}$, the car following interaction of
follower $\alpha$ with leader $\beta$ is $\vec{I}_{\alpha\beta}^{following}$ =
$\hat{n}_{p_{\alpha}x_{\alpha}}$, $\alpha$ maintains the direction of movement
towards $p_{\alpha}=\vec{x}_{\alpha}(t)+\hat{v}_{\beta}(t)\cdot d_{\alpha}$,
where $d_{\alpha}$ is the minimum vehicle distance, $d_{\alpha\beta}$ is the
distance between $\alpha$ and $\beta$, and $\vec{x}_{\alpha}(t)$ is the
current position of $\alpha$. Otherwise, $\vec{I}_{\alpha\beta}^{following}$ =
$\frac{s_{current}}{2}$, $\alpha$ slows down to maintain safe distance from
$\beta$. This interaction type only occurs if the leader $\beta$ is in the
movement direction of the follower $\alpha$ ($\vec{e}_{\alpha}$):
$\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\leq$ 10 $\parallel$
$\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\geq$ 350 and the movement
direction of $\beta$ ($\vec{e}_{\beta}$) is $\sim$ equal to
$\vec{e}_{\alpha}$: $\theta_{\vec{e}_{\beta}\hat{n}_{\beta\alpha}}$ $\leq$ 5
$\parallel$ $\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\geq$ 355\.
Here, $\hat{n}_{\beta\alpha}$ is the normalized vector from $\beta$ to
$\alpha$ and $\theta_{\vec{e}_{\beta}\hat{n}_{\beta\alpha}}$ is the angle
between vectors $\vec{e}_{\alpha}$ and $\hat{n}_{\beta\alpha}$. The force
$f_{\alpha\beta}^{soc}$ is set to zero when car following is active.

##### Courtesy Behavior

The courtesy behavior of driver $\alpha$ to another road user $\beta$ is
$\vec{I}_{\alpha\beta}^{courtesy}$ = $\frac{s_{current}}{2}$, here $\alpha$
decelerates and eventually stops when $\beta$ is passing in front of $\alpha$.
This type of interaction only happens if
$\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\leq$ 20 $\parallel$
$\theta_{\vec{e}_{\alpha}\hat{n}_{\beta\alpha}}$ $\geq$ 340 and 70 $<$
$\theta_{\vec{e}_{\beta}\hat{n}_{\beta\alpha}}$ $<$ 100 $\parallel$ 260 $<$
$\theta_{\vec{e}_{\beta}\hat{n}_{\beta\alpha}}$ $<$ 290\. The force
$f_{\alpha\beta}^{soc}$ is set to zero when courtesy behavior of $\alpha$ is
active.

This layer is also responsible for executing the game layer decisions. As an
example, if the result of a game-theoretic interaction between a car driver
and a pedestrian is that the pedestrian will continue to walk and the vehicle
driver will decelerate (see also Subsection III-C), the respective actions of
car and pedestrian will be executed in the Force-Based Modeling Layer. In this
paper, we consider the following possible actions for road users:

  * •

Continue: Car drivers continue the free-flow movement with their previous
speed. Pedestrians also continue free-flow movement but with maximum speed.

  * •

Decelerate: Road users decelerate and eventually stop.

  * •

Deviate (pedestrian only): We model this action as that any pedestrian
$\alpha$ first passes the car $\beta$ from behind from a position $p_{\alpha}$
and then continues moving towards the original destination. Here,
$p_{\alpha}=\vec{x}_{\beta}(t)$ \- 0.5$\vec{v}_{\beta}(t)$, if the car
continues moving and $p_{\alpha}=\vec{x}_{\beta}(t)$ \- $\vec{v}_{\beta}(t)$,
if the car stops. The pedestrian will move towards $p_{\alpha}$ until the car
stays within the field of view (FOV) of the pedestrian.

###  III-C Game-Theoretic Decision Layer

This layer handles complex interactions among road users. We use one-shot
Stackelberg games, which are sequential leader-follower games, where the
leader player commits to an action first; then followers choose an action to
optimize their utility considering the action chosen by the leader [8]. For a
single game, restrict the number of leaders to one, while there may be one or
more followers. For each implicit interaction only one game is played and the
games are independent of each other. In multiple players game, the leader
player makes a decision $s_{l}$ by considering the nearest follower and the
total number of followers. All followers follow the same action $s_{f}$.

Solving the Stackelberg game is done by finding the subgame perfect Nash
equilibrium (SPNE), which can be inferred by backward induction. Let $S_{l}$
and $S_{f}$ are the sets of all possible actions of the leader and follower
players respectively, and $u_{l}(s_{l},s_{f})$ and $u_{f}(s_{f}|s_{l})$ are
the utility of the leader and followers, when the leader choose the action
$s_{l}\in S_{l}$ and the followers react by the move $s_{f}\in S_{f}$.
$\forall$$s_{l}\in S_{l}$, $B_{f}(s_{l})$ = {$s_{f}\in
S_{f}|max(u_{f}(s_{f}|s_{l}))$} is the best answer of the follower. Then the
SPNE = $<${$s_{l}\in S_{l}|max(u_{l}(s_{l},B_{f}(s_{l})))$},
$B_{f}(s_{l})$$>$.

####  III-C1 Payoff Estimation

In the game, all actions of the players are ordinally valued from 2 to -100 as
shown in Figure 2a, with the assumption that every person prefers to reach
their destination safely and quickly. To model courtesy behavior and situation
dynamics, we consider six observable factors (in the following denoted by
Boolean variables $x_{1}\ldots x_{6}$) relevant for payoff estimation. They
are defined in the following and are used to calculate a set of parameters
$F_{1}\ldots F_{12}$, which are impacts of these factors. Let $\alpha$ is a
road user who needs to interact with another user $\beta$:

  * •

$x_{1}$: has value $1$, if current speed of $\beta$, $S_{current}$ $<$
$S_{normal}$. This factor determines the value of $F_{1}$, $F_{6}$ and
$F_{11}$.

  * •

$x_{2}$: has value $1$ if number of active interactions $<$ N (=2). This
factor determines the value of $F_{2}$ and $F_{7}$.

  * •

$x_{3}$: has value $1$ if the number of give ways already allowed by $\alpha$,
$G_{\alpha}$ $<$ M (=3) or not. This determines the value of $F_{3}$ and
$F_{8}$. We generate the value of $G_{\alpha}$ randomly.

  * •

$x_{4}$: has value $1$ if $\alpha$ is a car driver following another car. This
factor determines the value of $F_{4}$ and $F_{9}$.

  * •

$x_{5}$: has value $1$ if $\alpha$ is a pedestrian in a group. This factor
determines the value of $F_{5}$ and $F_{10}$.

  * •

$x_{6}$: has value $1$ if $\alpha$ is a car driver followed by another car.
This factor determines the value of $F_{12}$.

(a) Ordinal valuations of outcomes.

(b) The complete payoff matrix

Figure 2: The row player is a pedestrian with three possible actions; the
column player is a car driver with two possible actions.

For the scope of this paper, the parameters N and M were chosen in a “trial-
and-error” manner through simulation experiments.
The payoff of a car driver’s actions is calculated as follows:

  * •

Continue: $C_{c}$ = $F_{1}$ \+ $F_{2}$ \+ $F_{3}$ \+ $F_{4}$ \+ $F_{5}$

  * •

Decelerate: $C_{d}$ = $F_{6}$ \+ $F_{7}$ \+ $F_{8}$ \+ $F_{9}$ \+ $F_{10}$

The payoff value of the actions of a pedestrian is calculated as follows:

  * •

Continue: $P_{c}$ = $F_{1}$ \+ $F_{9}$ \+ $F_{10}$

  * •

Decelerate: $P_{d}$ = $F_{6}$ \+ $F_{4}$ \+ $F_{5}$

  * •

Deviate: $P_{dev}$ = $F_{11}$ \+ $F_{12}$

The utility calculation of the action pairs is shown in Figure 2b. As we were
not able to determine the parameters $F_{1}$, $F_{2}$, $F_{3}$, $F_{4}$,
$F_{5}$, $F_{6}$, $F_{7}$, $F_{8}$, $F_{9}$, $F_{10}$, $F_{11}$, and $F_{12}$
based on a substantial amount of real interaction scenarios, there remains a
considerable doubt about providing correct values for other and more complex
scenarios. Hence, while we did perform a sensitivity analysis on a certain
amount of interaction scenarios and analyzed the behavior of our model, we do
not present the values of these parameters in this paper. Automated
calibration and validation of these parameters (including N and M) is an
important task in our future work, see Section V.

###  III-D Recognition and Handling of Interactions

In our model, simple interactions are handled locally by the road users.
Remember that simple interactions are dealt with by local action and does not
require synchronization or communication. If a simple interaction situation is
detected by a road user within its FOV (set to 180° in this paper), it is
classified based on the interaction angle as given in Section II, and the
corresponding action is performed. Complex interactions are handled by a
(virtual) central mediation entity, further referred to as host agent. At
first, the host agent collects the nearest conflict (one conflict at a time)
and the total number of conflicts the user agent is facing currently for each
car agent in the environment by considering their FOV. The host agent checks
if this conflict can be solved without game playing, if yes then the conflict
will be handled by the respective car agent itself. Otherwise, the host
selects the car agent (faster agent) as the leader and organizes a game by
considering only the selected leader and its followers, not all other user
agent involved in this situation (if any), as then the host may need to
consider all possible subsets of conflicts as one game (worst-case scenario)
which will make the situation more complex. Hence, the impact of follower-to-
follower conflict and follower-to-other-agent conflict is not considered in
the scope of a single game. The workflow of a user agent and the host agent
for a single cycle are shown in Figure 7. The agent cycle is executed once per
simulation time unit including time to perform all interactions.

(a)  A Road User Agent

(b)  The Host Agent

Figure 3: The inner steps of the user and host agents with their plans and
goals.

##  IV Evaluation

To evaluate our simulation model, we simulate 117 conflict scenarios from the
shared space data and visualize the real and simulated trajectories and speed
difference of all road users involved in these scenarios. We select two
example scenarios among all those scenarios that involve multiple user
interaction, the courtesy behavior of car drivers towards pedestrians and also
the car following behavior to illustrate elaborately how our model works. We
also simulate another scenario revealing courtesy behavior of car drivers
towards other car drivers. However, for the last case we could not compare the
performance as real data is currently missing. We consider all these three
scenarios as complex situations because here the road users can choose their
actions among different alternatives like the pedestrians can deviate,
continue, or stop moving and the car driver can give way or continue driving.
We did not consider deviation as an option for car drivers because of the non-
holonomic constraints of the vehicular movement. To simulate different shared
space scenarios, we use Lightjason, a Java-based BDI multi-agent framework
[20]. All simulation runs were performed on an Intel$\circledR$ Core™i5
processor with 16 GB RAM.

##### Scenario 1

Here, an interaction between three pedestrians and one car driver is captured.
In this scenario, the car driver decelerates to let the pedestrians cross the
street. In our model, car handles this interaction re-actively as the angle
between the car and (at least) one of the pedestrians is $\leq$ 15 $\parallel$
$\geq$ 345 and reproduces the same solution.

(a) observation

(b) Comparison of the intensity of interaction

(c) Comparison of the trajectories. The arrows indicate the direction of
movement of the road users.

Figure 4: Comparison of a car-to-pedestrian interaction

We observe that in both simulation and reality, the car driver decelerates;
yet, the timing and the way of acting differ in both cases: The acceleration
and deceleration of simulated car driver is sharp compare to real car driver.
The mechanism of speed handing of cars in our model need to be improved. The
speed difference between the real and simulated cars are illustrated in Figure
4b. The comparison between the observed and simulated trajectories is given in
Figure 4c. In Figure 4c, the trajectories of both simulated and real car are
very similar, while the trajectories of pedestrians in simulation slightly
deviated from the real trajectories because of the repulsive force between
pedestrians when they come close to each other.

(a)

(b) The game matrices, left: Car2 to Ped4 and Ped5, right: Car4 to Ped6

(c) The real trajectories

(d) The simulated trajectories

(e) The simulated trajectories without game decision

Figure 5: Comparison of a car to pedestrian and a car to car interaction

##### Scenario 2

It includes two types of interactions, namely car-to-pedestrian and car-to-
car. In this scenario, Car2 has interaction with Ped4 and Ped5; Car3 and Car4
have car following interaction with Car2; and Car4 has interaction with Ped6.
Here, both Car2 and Car4 decelerate to let the pedestrians cross the street.
We handle both car-to-pedestrian interactions using game-theoretic model and
depict the game matrices in Figure 5b. Here, both Car2 and Car4 act as leader
and decide to decelerate; their payoff function reflects courtesy of the
drivers towards the pedestrians, as decelerating will maximize its value. We
visualize the real and simulated trajectories of these interactions at two
specific time steps in Figure 5c and Figure 5d respectively. In both Figure 5c
and Figure 5d, Ped indicates Pedestrian. Car3 and Car4 in Figure 5d decide to
stay behind Car2, instead of overtaking as in Figure 5c. We model this
behavior as reactive car following. Our simulation model reproduces real
interactions, but a large safety distance is maintained between pedestrian to
pedestrian, car to car and also car to pedestrian, compared to real situation.

We simulate this same scenario again but without game-theoretic layer. Figure
5e visualizes the simulated trajectories at a specific time step to show how
the classical social force model behaves in this situation. In Figure 5e all
cars deviate from their trajectory instead of stopping and pedestrians
sometimes go back and forth and sometimes stop.

##### Scenario 3

The third example studies courtesy behavior of a car driver towards another.
This is shown in Figure 6. In the example, when the driver of the blue car
perceives the green car, she slows down and lets the green car pass first.

Figure 6: Car-to-car interaction at different time steps. The blue car is
showing courtesy to the green car.

To sum up, in Scenario 1, our simulation model reproduces a multiple user
interaction among four road users realistically. Scenario 2 shows that our
model can handle three interdependent situations realistically. In the third
scenario, our model simulates the courtesy behavior between car drivers re-
actively. In all these example scenarios, the road users suitably detect and
classify the interaction and behave accordingly during simulation.

(a)  The real trajectories

(b)  The simulated trajectories

Figure 7: The trajectories of pedestrians and cars are indicated by blue lines
and red lines respectively.

Figure 7a and Figure 7b visualize the real and simulated trajectories of all
road users involved in 117 selected real scenarios, respectively. We start
simulating each road user five time steps before and end simulating five time
steps after than real data, because in some situations road users cross each
other without having any conflict, which is not the case in real scenario. The
reason behind that is the speed of road users varies in real and simulated
scenarios. The trajectories of real and simulated road users closely match.
However, the time step wise position of each road user differs in real
scenario and simulation because of speed difference between real and simulated
road users and also the fact that a large safety distance is maintained
between road users in simulation. The mean and standard deviation of speed of
the simulated and real road users are: 2.779 and 2.966 (simulated cars);
0.8469 and 0.7813 (simulated pedestrians); 1.27 and 1.007 (real cars), 0.502
and 0.2526 (real pedestrians). The calibration of the acceleration and
deceleration parameters is a part of our future work. Figure 8a represents
average distance between simulated and real trajectories of every car and
pedestrian involved in each scenario of our selected 117 scenarios. For more
than one car or pedestrian in a single scenario we take the average among them
to calculate the average distance. Figure 8b visualizes the individual average
speed difference of real and simulated road users, which is calculated in the
same way as trajectory distance.

(a) Distance between real and simulated trajectories

(b) Speed differences of real and simulated road users

Figure 8: The difference in trajectory and speed of real and simulated
pedestrians and cars.

##  V Conclusion and Outlook

This work-in-progress paper proposes a multiagent-based model for modeling
interactions of pedestrians and cars in shared spaces. We analyze mixed
traffic interactions, classify the interactions in terms of complexity and
then model different levels of interactions using the social force model and
Stackelberg games. Our model facilitates interaction recognition depending on
the situation context. First results support the hypothesis that this model
has good potential to simulate both complex and simple shared space scenarios,
as modeling multiple user interactions and interdependent interactions can
make the simulation model more realistic.

Future research will focus on calibrating and validating the model parameters
using genetic algorithms. Next, we will try to find out other factors which
might influence road user’s decision making. The explicit interaction and
pedestrian group dynamic will be integrated into our simulation model. Most
importantly, we shall study larger scenarios with larger numbers interactions.
Doing so, we shall explore the scalability of various interaction types for
different scenarios varying in number and heterogeneity of users as well as in
topology. In particular, more efficient algorithms for recognizing, framing
and resolving conflicts interactions will be required for realistically
simulating larger shared space scenarios.

## References

  * [1] E. Clarke, “Shared space-the alternative approach to calming traffic,” _Traffic engineering & control_, vol. 47, no. 8, 2006.
  * [2] R. Schönauer, M. Stubenschrott, H. Schrom-Feiertag, and K. Menšik, _Social and Spatial Behaviour in Shared Spaces_ , 2012.
  * [3] D. Helbing and P. Molnar, “Social force model for pedestrian dynamics,” _Physical review E_ , vol. 51, no. 5, p. 4282, 1995.
  * [4] F. Pascucci, N. Rinke, C. Schiermeyer, B. Friedrich, and V. Berkhahn, “Modeling of shared space with multi-modal traffic using a multi-layer social force approach,” _Transportation Research Procedia_ , vol. 10, pp. 316–326, 2015.
  * [5] N. Rinke, C. Schiermeyer, F. Pascucci, V. Berkhahn, and B. Friedrich, “A multi-layer social force approach to model interactions in shared spaces using collision prediction,” _Transportation Research Procedia_ , vol. 25, pp. 1249–1267, 2017.
  * [6] B. Anvari, M. G. Bell, A. Sivakumar, and W. Y. Ochieng, “Modelling shared space users via rule-based social force model,” _Transportation Research Part C: Emerging Technologies_ , vol. 51, pp. 83–103, 2015.
  * [7] F. Pascucci, N. Rinke, C. Schiermeyer, V. Berkhahn, and B. Friedrich, “Should i stay or should i go? a discrete choice model for pedestrian–vehicle conflicts in shared space,” Tech. Rep., 2018.
  * [8] R. Schönauer, M. Stubenschrott, W. Huang, C. Rudloff, and M. Fellendorf, “Modeling concepts for mixed traffic: Steps toward a microscopic simulation tool for shared space zones,” _Transportation Research Record: Journal of the Transportation Research Board_ , no. 2316, pp. 114–121, 2012.
  * [9] L. W. Lan and C.-W. Chang, “Inhomogeneous cellular automata modeling for mixed traffic with cars and motorcycles,” _Journal of Advanced Transportation_ , vol. 39, no. 3, pp. 323–349, 2005.
  * [10] Y. Zhang and H. Duan, “Modeling mixed traffic flow at crosswalks in micro-simulations using cellular automata,” _Tsinghua Science & Technology_, vol. 12, no. 2, pp. 214–222, 2007.
  * [11] D. Gettman and L. Head, “Surrogate safety measures from traffic simulation models,” _Transportation Research Record: Journal of the Transportation Research Board_ , no. 1840, pp. 104–115, 2003.
  * [12] C. Burstedde, K. Klauck, A. Schadschneider, and J. Zittartz, “Simulation of pedestrian dynamics using a two-dimensional cellular automaton,” _Physica A: Statistical Mechanics and its Applications_ , vol. 295, no. 3-4, pp. 507–525, 2001.
  * [13] K. Nagel and M. Schreckenberg, “A cellular automaton model for freeway traffic,” _Journal de physique I_ , vol. 2, no. 12, pp. 2221–2229, 1992\.
  * [14] H. Kita, “A merging–giveway interaction model of cars in a merging section: a game theoretic analysis,” _Transportation Research Part A: Policy and Practice_ , vol. 33, no. 3-4, pp. 305–312, 1999.
  * [15] N. Lütteken, M. Zimmermann, and K. J. Bengler, “Using gamification to motivate human cooperation in a lane-change scenario,” in _Intelligent Transportation Systems (ITSC), 2016 IEEE 19th International Conference on_. IEEE, 2016, pp. 899–906.
  * [16] T. Bjørnskau, “The zebra crossing game–using game theory to explain a discrepancy between road user behaviour and traffic rules,” _Safety science_ , vol. 92, pp. 298–301, 2017.
  * [17] A. Koefoed-Hansen and G. S. Brodal, “Representations for path finding in planar environments,” Ph.D. dissertation, 2012.
  * [18] I. Millington and J. Funge, _Artificial intelligence for games_. CRC Press, 2009.
  * [19] F. T. Johora, P. Kraus, and J. P. Müller, “Dynamic path planning and movement control in pedestrian simulation,” Sao Paulo, Brazil, May 2017. [Online]. Available: http://arxiv.org/abs/1709.08235
  * [20] M. Aschermann, P. Kraus, and J. P. Müller, “LightJason: A BDI Framework inspired by Jason,” in _Multi-Agent Systems and Agreement Technologies: 14th Europ. Conf., EUMAS 2016_ , ser. LNCS, vol. 10207. Springer, 2017, pp. 58–66. [Online]. Available: https://lightjason.github.io