text:
# GENERALIZED STRONG CURVATURE SINGULARITIES AND COSMIC CENSORSHIP

## 1 Introduction

The cosmic censorship hypothesis (CCH) of Penrose  says that in generic situations, all spacetime singularities arising from regular initial data are always hidden behind event horizons and hence invisible to outside observers (no naked singularities). This hypothesis plays a fundamental role in the theory of black holes and has been recognized as one of the most important open problems in classical general relativity. There exist many exact solutions of Einstein’s equations which admit naked singularities. However, Penrose  has stressed that the exact solutions with special symmetries have a rather limited value for verification of the CCH and what is required here is an understanding of the generic case.

One possible approach to this problem is to propose a class of generic singularities and then attempt to formulate and prove a censorship theorem which would constrain or prohibit the occurrence of naked singularities of the proposed class. Tipler et al.  and Królak  have argued that all singularities arising in generic situations should be of the strong curvature type. These singularities have the property that all objects approaching them are crushed to zero volume. The idea of a strong curvature singularity was introduced by Ellis and Schmidt  and defined in precise geometrical terms by Tipler  (a slightly different definition was given by Królak , see below). Królak  formulated and proved some censorship theorems that ruled out a class of naked singularities of the strong curvature type. Unfortunately, these results rely heavily on a further assumption (the so-called simplicity condition) which need not hold for generic spacetimes.

It should be stressed that there is no hope for finding a proof of the CCH using only the assumption that all singularities occurring in a given spacetime are of the strong curvature type. This follows just from the fact that naked singularities of this type do occur in certain exact solutions of Einstein’s equations. For instance, they occur in the Tolman-Bondi solution representing spherically symmetric inhomogeneous collapse of dust (see, e.g., Ref. ; see also Ref.  and references therein). Such naked singularities also occur in more general models of dust collapse — namely, in the Szekeres spacetimes which do not have any Killing vectors . Unnikrishnan  has argued that the existence of naked strong curvature singularities in the Tolman-Bondi solution can be ruled out by imposing certain reasonable constraints on the initial distribution of the energy density of dust. This argument, however, depends crucially on the spherical symmetry of the solution, and so cannot be applied to naked singularities occurring in more general cases — e.g. in the Szekeres spacetimes. It is worth recalling here that Bonnor  remarked, in another context, that “the Szekeres solution has a good deal of symmetry, even though it has in general no Killing vectors.” It is thus possible that the existence of naked singularities of the strong curvature type will always be accompanied by spacetime symmetries and/or instabilities of some sort, and so one can still hope to prove a formulation of the CCH involving — besides the assumption that all singularities are of strong curvature — a suitable criterion of genericity or stability.

To help identify such a criterion, it may be useful to establish and study various relations between strong curvature singularities and the causal structure in their neighborhood — such relations get to the heart of cosmic censorship. For this purpose, we need a new definition of a strong curvature singularity — one that not only describes the curvature strength of the singularity, but additionally enables one to relate the strength with properties of the causal structure in a neighborhood of the singularity. In Section II of this paper, we shall propose such a definition. Our definition is motivated by the definition of a strong curvature singularity given by Tipler and its modifications by Królak, but is significantly different and more general. The difference is in that our definition involves a certain focusing condition on solutions of the Raychaudhuri equation not only along a single causal geodesic, as in Tipler’s and Królak’s case, but rather along all causal geodesics in some small neighborhood about a given geodesic that reaches the singularity. All causal geodesics terminating at the singularities described by our definition, which we will refer to as generalized strong curvature singularities, will next be classified into three possible types; the classification is based on certain relations between the causal structure and the curvature strength of the singularities. Further on, in Sections III and IV, we shall formulate and prove a cosmic censorship theorem which shows that only one class of generalized strong curvature singularities, corresponding to a single type of geodesics according to our classification, can be naked.

This paper was inspired by ideas given in Ref. . Some of the results presented here are refinements of those announced without proofs in Ref. . The notation and fundamental definitions are as those in the monograph of Hawking and Ellis .

## 2 Generalized Strong Curvature Singularities

Before we give our definition, we need to recall one standard result on the behavior of geodesic congruences. Let $`\lambda (t)`$ be an affinely parametrized null (resp. timelike) geodesic. A congruence of null (timelike) geodesics infinitesimally neighboring $`\lambda (t)`$ and originating from a point on $`\lambda (t)`$ is characterized by two parameters: the expansion $`\theta `$ and the shear $`\sigma `$. The rate of change of the expansion $`\theta `$ along $`\lambda (t)`$ is given by the Raychaudhuri equation:

$$\frac{d\theta }{dt}=R_{ab}K^aK^b2\sigma ^2\frac{1}{n}\theta ^2,$$
(1)
where $`R_{ab}`$ is the Ricci tensor, $`K^a`$ is the tangent vector to $`\lambda (t)`$, $`n=2`$ if $`\lambda (t)`$ is a null geodesic and $`n=3`$ in the case when it is timelike (see pp. 78–88 of Ref. ).

Let us also recall the definition of a strong curvature singularity introduced by Królak .

Definition 1: Let $`\lambda `$ be a future-endless, future-incomplete null (timelike) geodesic. $`\lambda `$ is said to terminate in the future at a strong curvature singularity if, for each point $`p\lambda `$, the expansion $`\theta `$ of every future-directed congruence of null (timelike) geodesics emanating from $`p`$ and containing $`\lambda `$ becomes negative somewhere on $`\lambda `$.

This definition is equivalent to condition K given in Ref. . Królak’s definition generalizes Tipler’s definition of a strong curvature singularity because it implies weaker restrictions on the divergence of the curvature near the singularity .

We can now introduce our definition of a generalized strong curvature singularity.

Definition 2: Let $`\lambda `$ be a future-endless, future-incomplete null (timelike) geodesic. We say that $`\lambda `$ terminates in the future at a generalized strong curvature singularity if, for each sequence of endless null (timelike) geodesics, $`\{\lambda _n\}`$, for which $`\lambda `$ is a limit curve, at least one of the following conditions holds:

for each point $`p\lambda `$ and for each neighborhood $`N`$ of $`p`$ there exist a geodesic $`\stackrel{~}{\lambda }\{\lambda _n\}`$ and a point $`\stackrel{~}{q}\stackrel{~}{\lambda }N`$ such that the expansion $`\theta `$ of a future-directed congruence of null (timelike) geodesics emanating from $`\stackrel{~}{q}`$ and containing $`\stackrel{~}{\lambda }`$ becomes negative somewhere on $`\stackrel{~}{\lambda }`$;

$`I^{}(\lambda )=I^{}(\lambda _n)`$ for almost all geodesics belonging to $`\{\lambda _n\}`$.

Analogously, $`\lambda `$ is said to terminate in the past at a generalized strong curvature singularity if the time-reverse versions of the above conditions hold for $`\lambda `$.

Remark: By the term limit curve we mean in the above a curve that satisfies the definition given on p. 185 of Ref. .

It is easy to notice that every causal geodesic that terminates at a strong curvature singularity as defined by Królak also terminates at a generalized strong curvature singularity as defined above. To see this, consider, e.g., a null geodesic $`\lambda `$ that terminates in the future at Królak’s singularity. Next, take any sequence of endless null geodesics, $`\{\lambda _n\}`$, for which $`\lambda `$ is a limit curve. Let $`p`$ be an arbitrary point on $`\lambda `$, and let $`\{p_n\}`$ be an arbitrary sequence of points converging to $`p`$ such that, for each $`n`$, $`p_n\lambda _n`$. In accordance with Królak’s definition, the expansion $`\theta `$ of every future-directed congruence of null geodesics outgoing from $`p`$ and containing $`\lambda `$ must become negative. It follows by continuity that there must exist certain points $`p_k\{p_n\}`$ sufficiently close to $`p`$ such that, for each $`k`$, the expansion $`\theta `$ of a future-directed congruence of null geodesics outgoing from $`p_k`$ and containing $`\lambda _k\{\lambda _n\}`$ shall become negative as well. This clearly means that $`\{\lambda _n\}`$ shall obey condition (i) of Definition 2.

Let us now turn to condition (ii) of Definition 2. This condition has no direct equivalent in either Tipler’s or Królak’s definition; however, it appeals to the original idea behind the concept of a strong curvature singularity . To see this, observe first that the requirement $`I^{}(\lambda )=I^{}(\lambda _n)`$ found in condition (ii) means that both $`\lambda `$ and $`\lambda _n`$ must reach exactly the same point in the c-boundary of the spacetime (p. 218 of Ref. ). The geodesics $`\lambda _n`$ represent trajectories of point particles in motion towards the singularity within some small neighborhood about $`\lambda `$. Following the original idea of a strong curvature singularity, all small physical objects should be crushed to zero volume as they approach the singularity. Condition (ii) corresponds to the particular case where they are crushed to a single point in the $`c`$-boundary.

It is worth mentioning here that the $`c`$-boundary of cosmological models has been investigated by Tipler  in the context of limitations imposed on computation by general relativity. He has shown that a true universal Turing machine can be constructed only in a closed universe whose final singularity is a single point in the $`c`$-boundary topology. Thus all causal geodesics reaching this singularity will satisfy condition (ii) of Definition 2 and hence terminate in our generalized strong curvature singularity; moreover, they will be of type A according to the classification given below.

The following definition provides a complete classification of causal geodesics terminating at the singularities introduced in Definition 2.

Definition 3: Let $`\lambda `$ be a future-endless null (timelike) geodesic terminating in the future at a generalized strong curvature singularity, and let $`\{\lambda _n\}`$ be a sequence of endless null (timelike) geodesics for which $`\lambda `$ is a limit curve. $`\lambda `$ is said to be of type:

* $`A`$, if condition (2) of Definition 2 holds for each $`\{\lambda _n\}`$;
* $`B`$, if, for each $`\{\lambda _n\}`$ which does not satisfy condition (2) of Definition 2, there exist a geodesic $`\stackrel{~}{\lambda }\{\lambda _n\}`$ and a point $`\stackrel{~}{q}\stackrel{~}{\lambda }I^{}(\lambda )`$ such that the expansion $`\theta `$ of a future-directed congruence of null (timelike) geodesics emanating from $`\stackrel{~}{q}`$ and containing $`\stackrel{~}{\lambda }`$ becomes negative somewhere on $`\stackrel{~}{\lambda }`$;
* $`C`$, if $`\lambda `$ is neither of type $`A`$ nor $`B`$.

Remark: If the $`\lambda `$ mentioned above is a null geodesic that admits a segment contained in the boundary, $`\dot{I}^{}(\lambda )`$, of its chronological past, one can always find a sequence of endless null geodesics, $`\{\lambda _n\}`$, for which $`\lambda `$ is a limit curve, such that none of the $`\lambda _n`$ will be contained in the closure of $`I^{}(\lambda )`$. Clearly, such a sequence cannot satisfy condition (ii) of Definition 2, and so $`\lambda `$ cannot be of type A. In general, however, there may exist a future-endless null geodesic $`\lambda `$ that never intersects the boundary of its chronological past. For example, if $`\lambda `$ terminates in the future at a curvature singularity and the curvature diverges fast enough along $`\lambda `$, every future-endless segment of $`\lambda `$ may admit a pair of conjugates points, which implies that no segment of $`\lambda `$ can be contained in the achronal boundary $`\dot{I}^{}(\lambda )`$. In this case $`\lambda `$ may be of type A.

## 3 Censorship Theorem

For the sake of simplicity, we shall restrict our considerations to weakly asymptotically simple and empty (WASE) spacetimes (p. 225 of Ref. ). Such a spacetime $`(M,g)`$ can be conformally imbedded in a larger spacetime $`(\stackrel{~}{M},\stackrel{~}{g})`$ as a manifold with boundary $`\overline{M}=MM`$, where the boundary $`M`$ consists of two null surfaces $`𝒥^+`$ and $`𝒥^{}`$ that represent future and past null infinity, respectively. Moreover, there exists an open neighborhood $`U`$ of $`M`$ in $`\stackrel{~}{M}`$ such that $`UM`$ coincides with part of an asymptotically simple and empty spacetime $`(M^{},g^{})`$, which means that all possible singularities of $`(M,g)`$ can only occur in the region $`MU`$. Since the CCH is concerned with singularities that develop from an initially non-singular state, we shall deal here only with such WASE spacetimes $`(M,g)`$ in which the region $`MU`$ does not extend arbitrarily far into the past and to a spatial infinity.aaaDue to its generality, the definition of a WASE spacetime given in Ref.  has the unwanted feature that the region $`MU`$ might extend to a spatial infinity. For a deeper discussion of this problem, see Ref. . To make this precise, we shall assume that $`(M,g)`$ admits a partial Cauchy surface $`S`$ for which the following two conditions hold:

$`S`$ has an asymptotically simple past (p. 316 of Ref. );

every null geodesic $`\mu `$ generating $`𝒥^+`$ admits a past-endless segment $`\mu _0`$ such that $`\mu _0I^+(S,\overline{M})`$ and $`I^{}(\mu _0,\overline{M})I^+(S)D^+(S)U`$.

Such a surface $`S`$ will be called a regular partial Cauchy surface. (This definition is a slight modification of that previously used in Ref. .)

Let us now recall that a WASE spacetime $`(M,g)`$ is said to be future asymptotically predictable from a partial Cauchy surface $`S`$ if the future null infinity $`𝒥^+`$ is contained in the closure of the future Cauchy development $`D^+(S,\overline{M})`$ (p. 310 of Ref. ). Future asymptotic predictability is a mathematically precise statement of cosmic censorship for $`(M,g)`$, since it ensures that there will be no singularities to the future of $`S`$ which are naked, i.e. which are visible from $`𝒥^+`$.

We are now in a position to state our cosmic censorship theorem. By this theorem, only singularities corresponding to null geodesics of type C can be naked.

Theorem 1: Let $`(M,g)`$ be a WASE spacetime admitting a regular partial Cauchy surface $`S`$. Suppose furthermore that the following conditions hold:

the null convergence condition, i.e. $`R_{ab}K^aK^b0`$ for every null vector $`K^a`$ of $`(M,g)`$;

the generic condition, i.e. every endless null geodesic of $`(M,g)`$ admits a point at which $`K_{[a}R_{b]cd[e}K_{f]}K^cK^d0`$, where $`K^a`$ is the tangent vector to the geodesic;

$`(M,g)`$ admits no naked points-at-infinity, i.e. for each point $`p𝒥^+`$, every future-endless, future-complete null geodesic of $`(M,g)`$ contained in $`\overline{I}^{}(p,\overline{M})D(S)`$ has a future endpoint on $`𝒥^+`$ in $`\overline{M}`$;

every incomplete null geodesic of $`(M,g)`$ terminates at a generalized strong curvature singularity.

If $`(M,g)`$ is not future asymptotically predictable from $`S`$, then there must exist a null geodesic $`\lambda `$ int$`D(S)`$ which terminates in the future at a generalized strong curvature singularity and is of type C.

Remarks: Conditions (i) and (ii) of Theorem 1 are reasonable requirements for any physically realistic model of a classical spacetime; they have been discussed extensively in the literature on the singularity theorems (see, e.g., Refs. ).

Condition (iii) is a slightly weakened version of that previously used by Newman and Joshi  in their censorship theorem. This condition ensures that any possible breakdown of future asymptotic predictability in $`(M,g)`$ must be associated with the occurrence of a naked singularity and not of a naked “point-at-infinity”. Known examples of WASE spacetimes with naked singularities contain no naked points-at-infinity in the sense of condition (iii). It should, however, be stressed that this condition restricts to some extent the generality of our result, since it is conceivable that some spacetimes might contain both naked singularities and naked points-at-infinity.

Note also that no causality conditions, other than those implied by the existence of $`S`$, are imposed on $`(M,g)`$ in Theorem 1. Thus this theorem may be applied to naked singularities associated with causality violation. The possible existence of such naked singularities results from Tipler’s singularity theorem .

## 4 Proof of the Theorem

Now we shall prove our theorem; the following two lemmas serve this purpose.

Lemma 1: Under the assumptions of Theorem 1, suppose that $`(M,g)`$ is not future asymptotically predictable from $`S`$. Then there must exist a past-incomplete null geodesic $`\eta H^+(S,\overline{M})`$ which has a future endpoint $`p𝒥^+`$. Moreover, the following conditions hold:

there exists a point $`r\eta M`$ such that the closure of $`I^{}(r)S`$ is compact;

for each point $`c_iI^{}(p,\overline{M})D^+(S)`$, there exists a null geodesic $`\lambda _i\dot{I}^+(c_i,\overline{M})`$ which extends from $`c_i`$ to some point on $`𝒥^+`$.

Proof: Let us first observe that from condition (ii) of the definition of $`S`$ it follows that every generator of $`𝒥^+`$ must intersect the closure of $`D^+(S,\overline{M})`$. In addition, as $`S`$ has an asymptotically simple past by condition (i), $`𝒥^{}`$ must be contained in the closure of $`D^{}(S,\overline{M})`$. This implies that $`(M,g)`$ is partially asymptotically predictable from $`S`$ as defined by Tipler . Thus, as $`(M,g)`$ is not future asymptotically predictable from $`S`$, by Proposition 2 of Ref.  there must exist a past-endless null geodesic generator $`\eta `$ of $`H^+(S)`$ with future endpoint $`p𝒥^+`$. Following steps as in the proof of Theorem 1 of Ref. , we can now easily show, by making use of conditions (i) and (ii) of our theorem, that $`\eta `$ must be incomplete in the past.

We shall now show that condition (a) holds. Let us choose an arbitrary point $`r\eta M`$; we shall demonstrate that the set $`R\overline{I^{}(r)S}`$ is compact. Since $`S`$ has an asymptotically simple past, the proof of this fact will be quite similar to the proof of the Theorem of Ref. .

Assume, to the contrary, that $`R`$ is not compact. Then, by Lemma 2 of Ref. , there must exist a sequence of null geodesics $`\eta _i`$ with future endpoints $`r_i`$ converging to $`r`$, such that each $`\eta _i`$ will be a generator of the achronal boundary $`\dot{I}^{}(r_i)`$ and will have a past endpoint on $`𝒥^{}`$; moreover, the geodesic $`\eta `$ must be a limit curve of the sequence $`\{\eta _i\}`$. As $`\eta `$ is past-incomplete, by condition (iv) of Theorem 1, it must terminate in the past at a generalized strong curvature singularity. Thus $`\{\eta _i\}`$ must fulfill at least one of the time-reverse versions of conditions (i) and (ii) of Definition 2. Since each of the $`\eta _i`$ has a past endpoint on $`𝒥^{}`$, for each $`\eta _i`$ one has $`I^+(\eta _i)I^+(\eta )`$, as $`\eta I^+(S)`$ while $`𝒥^{}I^{}(S,\overline{M})`$. It is thus clear that the time-reverse version of condition (ii) cannot be fulfilled, and so the time-reverse version of condition (i) must be. There must thus exist a geodesic $`\stackrel{~}{\eta }\{\eta _i\}`$ and a point $`\stackrel{~}{q}\stackrel{~}{\eta }`$ such that the expansion $`\theta `$ of a past-directed null geodesic congruence emanating from $`\stackrel{~}{q}`$ and containing $`\stackrel{~}{\eta }`$ must become negative somewhere on $`\stackrel{~}{\eta }`$. In this case, by Proposition 4.4.4 of Ref.  and condition (i) of our Theorem 1, there would exist a point conjugate to $`\stackrel{~}{q}`$ along $`\stackrel{~}{\eta }`$, since $`\stackrel{~}{\eta }`$ is past-complete as it has a past endpoint on $`𝒥^{}`$. However, by Proposition 4.5.12 of Ref. , the existence of a pair of conjugate points on $`\stackrel{~}{\eta }`$ would contradict the achronality of $`\stackrel{~}{\eta }`$. In view of this contradiction we must conclude that $`R`$ is compact.

We shall now show that condition (b) holds. Let $`c_i`$ be an arbitrary point in $`I^{}(p,\overline{M})D^+(S)`$, and let $`\mu `$ be a generator of $`𝒥^+`$ that passes through $`p`$. We shall first show that $`\mu `$ must leave $`I^+(c_i,\overline{M})`$. To see this, suppose that $`\mu `$ were contained in $`I^+(c_i,\overline{M})`$. Then the past set $`X_{a\mu _0}[I^{}(a,\overline{M})M]`$, where $`\mu _0`$ denotes the past-endless segment of $`\mu `$ mentioned in condition (ii) of the definition of $`S`$, would be non-empty, as it would contain $`c_i`$. From the definition of $`X`$ it follows that any null geodesic generator of the boundary of $`X`$ cannot have its future endpoint on $`𝒥^+`$; otherwise, such a point would have to coincide with a past endpoint of $`\mu _0`$, which is impossible as $`\mu _0`$ is past-endless. But from the definition of the surface $`S`$ it follows that $`X`$ would have to be contained in an open neighborhood $`U`$ of $`𝒥^+𝒥^{}`$, such that $`UM`$ coincides with part of an asymptotically simple and empty spacetime, which implies that all null geodesics generating the boundary of $`X`$ would have to have their future endpoints on $`𝒥^+`$. This contradiction shows that $`\mu `$ must leave $`I^+(c_i,\overline{M})`$, so there must exist a point $`b\mu \dot{I}^+(c_i,\overline{M})`$. Since $`b\mu J^{}(p,\overline{M})`$, and $`pH^+(S,\overline{M})`$, it follows that $`b`$ must belong to the closure of $`D^+(S,\overline{M})`$. Therefore all the past-directed null geodesics outgoing from $`b`$, with the exception of $`\mu `$, must enter the interior of $`D^+(S)`$ immediately after leaving $`b`$; otherwise they would be generators of $`𝒥^+`$, which is impossible. Thus, as int$`D^+(S)`$ is a causally simple set, there must exist a null geodesic $`\lambda _i\dot{I}^+(c_i,\overline{M})`$ which extends from $`c_i`$ to $`b𝒥^+`$, as required in condition (b). $`\mathrm{}`$

Lemma 2: Let $`r`$ be a point on $`H^+(S)`$, where $`S`$ is a partial Cauchy surface. If the set $`\dot{I}^{}(r)J^+(S)`$ is compact, then every null geodesic generator of $`H^+(S)`$ through $`r`$ is geodesically complete in the past direction.

Proof: The proof of this lemma is identical to the proof of Lemma 8.5.5 of Ref. . In the course of the proof we only need to consider the set $`\dot{I}^{}(r)J^+(S)`$ instead of the Cauchy horizon $`H^+(S)`$. $`\mathrm{}`$

Proof of Theorem 1: Assume that $`(M,g)`$ is not future asymptotically predictable from $`S`$. Then, by Lemma 1, there exists a past-incomplete null geodesic $`\eta H^+(S)`$ which has a future endpoint $`p𝒥`$; moreover, there exists a point $`r\eta M`$ such that $`R\overline{I^{}(r)S}`$ is compact. As $`\eta `$ is past-incomplete, by Lemma 2 the set $`Q\dot{I}^{}(r)J^+(S)`$ cannot be compact. Let us put a timelike vector field on $`M`$ (such a field will always exist since $`M`$ admits a Lorentz metric $`g`$). If every integral curve of this field that intersects the set $`R`$ also intersects the set $`Q`$, we would have a continuous one-to-one mapping of $`R`$ onto $`Q`$, and hence $`Q`$ would have to be compact. Therefore there must exist a future-endless timelike curve $`\alpha I^{}(r)`$.

Let $`𝒫`$ be the family of all sets of the form $`I^{}(\alpha )`$, where $`\alpha `$ is a future-endless timelike curve contained in $`I^{}(r)`$; and let $`\widehat{P}`$ be a maximal chain determined in $`𝒫`$ by the relation of inclusion. Denote now by $`P_0`$ the set $`\{P:P\widehat{P}\}`$. This set is non-empty. To see this, note first that no member of $`\widehat{P}`$ can be contained in $`I^{}(S)`$; otherwise, there would exist a future-endless timelike curve $`\alpha I^{}(S)`$, which is impossible as $`I^{}(S)=D^{}(S)`$. Note now that, for each $`P\widehat{P}`$, the set $`\overline{PS}`$ must be compact, as it is a closed subset of the compact set $`R`$. This clearly implies, by the definition of $`P_0`$, that the set $`A\overline{P_0S}`$ is non-empty and compact. In addition, from the definition of $`P_0`$ it follows that there must exist a future-endless timelike curve $`\alpha _0`$ such that $`P_0=I^{}(\alpha _0)`$. Note also that $`P_0`$ must be a minimal element of $`\widehat{P}`$, i.e. we will have $`P_0=I^{}(\beta )`$ for every future-endless timelike curve $`\beta P_0`$.

Let us now fix some future-endless timelike curve $`\alpha _0`$ such that $`P_0=I^{}(\alpha _0)`$. Let $`\{c_i\}`$ be a sequence of points on $`\alpha _0D^+(S)`$ such that, for each $`i`$, $`c_{i+1}I^+(c_i)`$; assume also that $`\{c_i\}`$ has no accumulation point on $`\alpha _0`$ (such a sequence can always be found as $`\alpha _0`$ has no future endpoint). By condition (b) of Lemma 1, there must exist, for each $`i`$, a null geodesic $`\lambda _i`$ running from $`c_i`$ to some point on $`𝒥^+`$. Moreover, each $`\lambda _i`$ shall be a generator of the achronal boundary $`\dot{I}^+(c_i,\overline{M})`$ and shall be future-complete as it reaches $`𝒥^+`$.

Let us now extend each of the $`\lambda _i`$ maximally into the past. Since each $`\lambda _i`$ passes through $`c_i\alpha _0D^+(S)`$, each of the extended $`\lambda _i`$ must intersect $`S`$ at some point $`a_iA`$. As the set $`A`$ is compact, the sequence $`\{a_i\}`$ must have an accumulation point $`aA`$. Therefore, by Lemma 6.2.1 of Ref. , through $`a`$ there is a non-spacelike curve $`\lambda `$ which is future-inextendible in $`M`$ and which is a limit curve of the sequence $`\{\lambda _i\}`$. Since all the $`\lambda _i`$ are null geodesics, $`\lambda `$ must be a null geodesic as well. Note also that $`\lambda `$ must be contained in $`\overline{P_0}`$; otherwise $`\lambda `$ would have to intersect the curve $`\alpha _0`$ at some point that would be an accumulation point of the sequence $`\{c_i\}`$, which is impossible as $`\{c_i\}`$ has no accumulation point on $`\alpha _0`$. Moreover, as $`P_0`$ is a minimal element of $`\widehat{P}`$, we must have $`P_0=I^{}(\lambda )`$. Since $`\lambda `$ intersects the surface $`S`$, it cannot be a generator of $`H^+(S)`$. Thus, as $`\lambda \overline{P_0}\overline{D(S)}`$, we will have $`\lambda `$ int$`D(S)`$. Since $`\lambda \overline{P_0}\overline{I^{}(r)}`$, and $`rM`$, $`\lambda `$ cannot have a future endpoint on $`𝒥^+`$. In addition, as $`rJ^{}(p,\overline{M})`$, $`\lambda `$ must be contained in the closure of $`I^{}(p,\overline{M})`$. Therefore, as $`p𝒥^+`$, by condition (iii) of Theorem 1, $`\lambda `$ must be incomplete in the future.

By condition (iv) of Theorem 1, $`\lambda `$ must terminate in the future at a generalized strong curvature singularity. Thus the sequence $`\{\lambda _i\}`$ must satisfy at least one of conditions (i) and (ii) of Definition 2. Since $`\lambda \overline{P_0}`$, and each of the $`\lambda _i`$ must leave $`\overline{P_0}`$ as it has a future endpoint on $`𝒥^+`$, we will have $`I^{}(\lambda )I^{}(\lambda _i)`$ for each $`i`$. This means that condition (ii) of Definition 2 cannot hold for the $`\{\lambda _i\}`$, hence $`\lambda `$ cannot be of type A. Suppose that $`\lambda `$ were of type B. Then, according to Definition 3, there would exist a geodesic $`\stackrel{~}{\lambda }\{\lambda _i\}`$ and a point $`\stackrel{~}{q}\stackrel{~}{\lambda }I^{}(\lambda )`$ such that the expansion $`\theta `$ of a future-directed congruence of null geodesics outgoing from $`\stackrel{~}{q}`$ and containing $`\stackrel{~}{\lambda }`$ would be negative somewhere on $`\stackrel{~}{\lambda }`$. As $`\stackrel{~}{\lambda }`$ is future-complete, by condition (i) of our theorem and Proposition 4.4.4 of Ref. , there would thus exist some point $`xJ^+(\stackrel{~}{q})`$ conjugate to $`\stackrel{~}{q}`$ along $`\stackrel{~}{\lambda }`$. Consequently, by Proposition 4.5.12 of Ref. , there would also exist a timelike curve from $`\stackrel{~}{q}`$ to some point $`y\stackrel{~}{\lambda }J^+(x)`$. But this is impossible because from the above construction it follows that the points $`\stackrel{~}{q}`$ and $`y`$ would have to lie on the achronal segment of $`\stackrel{~}{\lambda }`$ contained in the boundary $`\dot{I}^+(\alpha _0\stackrel{~}{\lambda })`$. This contradiction shows that the geodesic $`\lambda `$ cannot be of type B. Thus, as $`\lambda `$ is neither of type A nor B, it must be of type C, which completes the proof. $`\mathrm{}`$

## 5 Concluding Remarks

In Theorem 1 we have demonstrated that the only possible naked singularities of strong curvature are those corresponding to null geodesics of type C, according to the classification proposed in this paper. As a consequence, one only needs to consider the stability of type C singularities as relevant to the cosmic censorship problem, which is a significant restriction. One argument showing that these singularities have a tendency to be unstable under generic perturbations was given in Ref. . We may conclude that a further study of the properties of type C singularities might bring us closer to formulating the genericity and stability criterium which is needed, as mentioned in the Introduction, to arrive at a satisfactory statement of the CCH.

As we mentioned above, naked strong curvature singularities are present in certain exact solutions of the Einstein equations — solutions that represent gravitational collapse. In view of our result, it would therefore be of interest to verify whether in those solutions, the naked singularity is also associated in every case with the existence of type C geodesics. It should be stressed that this does not follow directly from Theorem 1, since it is not clear whether the solutions in question necessarily satisfy all the assumptions of that theorem.

## Acknowledgements

We wish to thank Andrzej Królak for many valuable discussions. This work was supported in part by the Polish Committee for Scientific Research (KBN) under Grant No. 2 P03B 130 16.


format: text

text:
# Anomaly in Spin Excitation Spectrum of Double-Exchange Systems with Randomness

One of the most important issues of colossal magnetoresistance (CMR) manganites AMnO3 is to understand effects of the A-site substitution. $`^{\text{?}\text{)}}`$ Generally, the A-site substitution has two different aspects. One is the carrier doping which is the substitution by A ions with different valences. The other is the ionic radius modification which is the A-site substitution with ions having the same valence but different ionic radius. In the latter, the averaged radius of A ions $`r_\mathrm{A}`$ affects length and angle of Mn-O-Mn bonds. This leads to a change of effective transfer energies between Mn sites. Therefore, it is widely accepted that the ionic radius modification gives a bandwidth control. Owing to a wide potential of these A-site substitutions, manganese oxides exhibit rich physics including the CMR phenomena.

Among various compounds AMnO3, recent developments in both experiments and theories have been revealed that La1-xSrxMnO3 (LSMO) at $`x1/3`$ is a canonical double-exchange (DE) systems. $`^{\text{?}\text{}\text{?}\text{)}}`$ The DE model $`^{\text{?}\text{)}}`$ can explain magnetic and transport properties in this compound quantitatively. For instance, the experimental values of Curie temperature $`T_\mathrm{C}`$ agree well with the theoretical estimate. $`^{\text{?)}}`$ Thus, the LSMO compound provides a well-established starting point to examine effects of the A-site substitutions.

Starting from LSMO as a reference, we focus on the ionic radius control by the A-site substitutions. This causes a systematic change in transport and magnetic properties. For instance, the value of $`T_\mathrm{C}`$ decreases systematically for compounds with smaller $`r_\mathrm{A}`$ and hence with ‘narrower-bandwidth’ such as Pr1-xSrxMnO3 and La1-xCaxMnO3 (LCMO). Hereafter, we call these A-site substituted compounds ‘low-$`T_\mathrm{C}`$ manganites’ whereas we call the compounds such as LSMO, which show relatively high $`T_\mathrm{C}`$, ‘high-$`T_\mathrm{C}`$ manganites’. An empirical relation between $`T_\mathrm{C}`$ and $`r_\mathrm{A}`$ has been proposed experimentally. $`^{\text{?}\text{)}}`$ The systematic change has been considered to be primarily due to the bandwidth change by the A-site substitution.

However, the picture of the bandwidth control has been doubted by some experimental results. One is a quantitative comparison between $`T_\mathrm{C}`$ and the bandwidth estimated from $`r_\mathrm{A}`$. The change of $`T_\mathrm{C}`$ is much larger than that of the bandwidth. For instance, from LSMO to LCMO, $`T_\mathrm{C}`$ decreases about 30% while the estimated bandwidth decreases less than 2%. $`^{\text{?}\text{)}}`$ This large decrease of $`T_\mathrm{C}`$ cannot be explained by the bandwidth control in the DE theory in which $`T_\mathrm{C}`$ should scale to the bandwidth. An extra degree of freedom as a hidden parameter seems to be necessary to explain such a change.

In this Letter, we examine a role of the randomness in the A-site substitution. In general, oxides are known to be far from perfect crystals. Moreover, in manganites, randomness is inevitably introduced by the A-site substitution since the compounds are solid solutions of different A ions. Random distribution of A ions induces a charge and structural disorder in the system, which introduces randomness in the electronic Hamiltonian.

Importance of the randomness effects has been discussed both experimentally and theoretically. One is that $`T_\mathrm{C}`$ is much affected not only by the averaged radius of A ions but also by the standard deviation of the ionic radii. A linear decrease of $`T_\mathrm{C}`$ as a function of the standard deviation has been reported. $`^{\text{?}\text{)}}`$ Another experimental fact is a correlation between the residual resistivity and $`T_\mathrm{C}`$. The A-site substitution from high-$`T_\mathrm{C}`$ to low-$`T_\mathrm{C}`$ manganites causes a systematic increase of the residual resistivity even in single crystals. $`^{\text{?}\text{}\text{?}\text{)}}`$ These results clearly indicate that it is necessary to take account of effects of the randomness in the A-site substitution. Recently, effects of the randomness have been discussed in the DE model. $`^{\text{?}\text{}\text{?}\text{}\text{?}\text{}\text{?}\text{)}}`$ A significant decrease of $`T_\mathrm{C}`$ by the randomness has been predicted. $`^{\text{?, ?, }\text{?}\text{}\text{?}\text{}\text{?}\text{)}}`$

Let us now consider another important result of the ionic radius modification which is the systematic change in the spin excitation spectrum. High-$`T_\mathrm{C}`$ manganites show a cosine-like spin-wave dispersion $`^{\text{?}\text{)}}`$ which is well described by the DE theory quantitatively. $`^{\text{?}\text{)}}`$ On the other hand, in low-$`T_\mathrm{C}`$ manganites, the spin excitation deviates from this behavior considerably, and shows some anomalies such as softening, broadening, anti-crossing and gap-opening. $`^{\text{?}\text{}\text{?}\text{}\text{?}\text{}\text{?}\text{)}}`$ The DE theory also fails to explain these anomalies in the spin dynamics.

In this Letter, we study effects of the randomness on the spin dynamics in the DE model. The spin excitation is calculated in the ground state by the spin wave approximation. We find clear deviation from the cosine-like spin-wave dispersion in the pure DE model and some anomalies in the spin excitation spectrum. Through comparisons with the experimental results in low-$`T_\mathrm{C}`$ manganites, our results strongly suggest that the randomness gives an important effect of the A-site substitution.

As a model for the A-site substituted manganites, we consider the DE model with randomness. The Hamiltonian is written by

$$=t\underset{ij,\sigma }{}(c_{i\sigma }^{}c_{j\sigma }+\mathrm{h}.\mathrm{c}.)\frac{J_\mathrm{H}}{S}\underset{i}{}\stackrel{}{\sigma }_i\stackrel{}{S}_i+\underset{i\sigma }{}\epsilon _ic_{i\sigma }^{}c_{i\sigma },$$
(1)
where the first two terms denote the DE Hamiltonian $`^{\text{?)}}`$ which consists of the nearest-neighbor hopping and the Hund’s-rule coupling between itinerant electrons and localized spins with a magnitude $`S`$. These two terms quantitatively describe physical properties in high-$`T_\mathrm{C}`$ manganites including the spin excitation. $`^{\text{?, ?, ?)}}`$ The last term denotes the on-site randomness. In this work, for simplicity, we incorporate the randomness by the A-site substitution as the on-site random potential although it may cause various other effects, for instance, a randomness in transfer integrals. Here we consider the binary-alloy type randomness, that is, $`\epsilon _i`$ takes $`\pm W_{\mathrm{imp}}/2`$ in equal probability in each site. We confirm that the results are independent of the type of randomness qualitatively. Although it is difficult to determine the actual magnitude of the randomness in manganites, the value of $`W_{\mathrm{imp}}`$ is roughly estimated to be the same order of magnitude as the bandwidth in low-$`T_\mathrm{C}`$ compounds. $`^{\text{?, }\text{?}\text{)}}`$ As an energy unit, we use the half bandwidth at the ground state for $`W_{\mathrm{imp}}=J_\mathrm{H}=0`$, i.e., $`W=6t=1`$ in three dimensions.

We study the spin excitation of the model (1) using the spin wave approximation. We consider a perfectly spin-polarized ferromagnetic ground state $`\stackrel{}{S}_i=(0,0,S)`$, and calculate the one-magnon excitation spectrum within the lowest order of $`1/S`$ expansion. In the absence of impurities, the formulation is given in Ref. ?.

Here we extend the method to the case with randomness. For a given configuration of the quenched randomness $`\{\epsilon \}`$ on a finite size cluster, we explicitly diagonalize the Hamiltonian matrix at the ground state. We have

$$\stackrel{~}{\phi }_{n\sigma }(i)=(\stackrel{~}{E}_n\sigma J_\mathrm{H})\stackrel{~}{\phi }_{n\sigma }(i),$$
(2)
where $`\stackrel{~}{\phi }_{n\sigma }`$ is the $`n`$-th wavefunction for the spin $`\sigma `$ electron. Since wavefunctions are given in real space, it is straightforward to formulate Green’s functions of electrons and spin wave using the real-space representation.

Electron Green’s function is given by

$$\stackrel{~}{G}_{ij,\sigma }(\omega )=\underset{n}{}\frac{\stackrel{~}{\phi }_{n\sigma }(i)\stackrel{~}{\phi }_{n\sigma }^{}(j)}{\omega (\stackrel{~}{E}_n\sigma J_\mathrm{H}\mu )+\mathrm{i}\eta \mathrm{sgn}\omega }.$$
(3)
Following Ref. ?, spin wave self-energy is calculated from electron’s spin polarization function depicted in Fig. 1. In the limit of $`J_\mathrm{H}/t\mathrm{}`$, which is realistic in manganites, we obtain

$`\stackrel{~}{\mathrm{\Pi }}_{ij}(\omega )`$ $`=`$ $`{\displaystyle \frac{1}{2SN}}{\displaystyle \underset{mn}{}}f_n\stackrel{~}{\phi }_n(i)\stackrel{~}{\phi }_n(j)^{}\stackrel{~}{\phi }_m(j)\stackrel{~}{\phi }_m(i)^{}`$ (4)
$`\times (\stackrel{~}{E}_m\stackrel{~}{E}_n\omega ),`$
where $`N`$ is the number of lattice sites and $`f_n`$ is the fermi distribution function for up-spin states. Spin wave Green’s function is given by $`\stackrel{~}{D}_{ij}(\omega )=(\omega \stackrel{~}{\mathrm{\Pi }}_{ij}(\omega )+\mathrm{i}\eta )^1`$. Spin wave excitations are obtained from the poles of the Green’s function. Since $`\stackrel{~}{\mathrm{\Pi }}_{ij}1/S`$, positions of poles are at $`\omega O(1/S)`$. Then, as discussed in Ref. ?, spin wave excitations within the lowest order of $`1/S`$ expansion can be determined from the static part of the self-energy $`\stackrel{~}{\mathrm{\Pi }}_{ij}(\omega =0)`$.

The spectral function for spin wave is defined by

$$\stackrel{~}{A}(𝒒,\omega )=\frac{1}{N}\underset{ij}{}\frac{1}{\pi }\mathrm{𝐼𝑚}\stackrel{~}{D}_{ij}(\omega )\mathrm{exp}\left[\mathrm{i}𝒒(𝒓_i𝒓_j)\right].$$
(5)
Using the eigenvalues $`\stackrel{~}{\omega }_l`$ and eigenvectors $`\stackrel{~}{\psi }_l(j)`$ of $`\stackrel{~}{\mathrm{\Pi }}_{ij}`$,

$$\underset{j}{}\stackrel{~}{\mathrm{\Pi }}_{ij}(\omega =0)\stackrel{~}{\psi }_l(j)=\stackrel{~}{\omega }_l\stackrel{~}{\psi }_l(i),$$
(6)
the matrix inversion of $`(\omega \stackrel{~}{\mathrm{\Pi }}_{ij})`$ is obtained. We have

$`\stackrel{~}{A}(𝒒,\omega )`$ $`=`$ $`{\displaystyle \underset{l}{}}\stackrel{~}{A}_l(𝒒)\delta (\omega \stackrel{~}{\omega }_l),`$ (7)
$`\stackrel{~}{A}_l(𝒒)`$ $`=`$ $`{\displaystyle \frac{1}{N}}\left|{\displaystyle \underset{j}{}}\stackrel{~}{\psi }_l(j)\mathrm{exp}(\mathrm{i}𝒒𝒓_j)\right|^2.`$ (8)
Finally, we average $`\stackrel{~}{A}(𝒒,\omega )`$ for random configurations of the potential $`\{\epsilon \}`$ to obtain the spectral function $`A(𝒒,\omega )`$.

Figure 2 shows the spectral function for spin excitations at $`x=0.3`$ in three dimensions. Here, $`x`$ is the hole density in model (1). System size $`N`$ is $`16\times 16\times 16`$ with the periodic boundary conditions and the random average is taken for $`24`$ different configurations of the random potential for each value of $`W_{\mathrm{imp}}`$. For small values of $`W_{\mathrm{imp}}`$, the spectrum is similar to the cosine-like dispersion for $`W_{\mathrm{imp}}=0`$ (gray curve) although it shows some broadening. However, when the value of $`W_{\mathrm{imp}}`$ becomes moderate, some anomalies appear conspicuously in the spectrum. Apart from the large broadening, a remarkable feature is a branching which is clearly seen, for instance, along the $`\mathrm{\Gamma }`$-M and $`\mathrm{\Gamma }`$-R lines for $`W_{\mathrm{imp}}/W=0.5`$. Near the zone boundary, extra branches with considerable weight appear well below the cosine-like branch. If one follows only the extra lower branch, the spin excitation spectrum appears to show a significant softening.

In order to consider the origin of the anomalies, we perform detailed calculations in one dimension. We take $`N=256`$ and the random average for $`500`$ configurations. Figure 3 shows the results for different values of $`x`$ at $`W_{\mathrm{imp}}/W=0.5`$ ($`W=2t`$ in one dimension). The precise calculations in one dimension reveal the detailed structure of the anomalies; an anti-crossing structure with a gap opening, shadow bands and a significant broadening. These remarkable features become conspicuous for the values of $`x`$ far from the quarter filling $`x=0.5`$.

An important point in Fig. 3 is that the anomalies appear near the Fermi wave number $`k_\mathrm{F}`$ which is given by $`\pi (1x)`$ in one dimension. This indicates that the anomalies come from a singularity of $`2k_\mathrm{F}`$ momentum transfer between $`qk_\mathrm{F}`$ and $`qk_\mathrm{F}`$. There are also small anomalies at $`q\pi k_\mathrm{F}`$ and $`\pi +k_\mathrm{F}`$. We note that these anomalies can be seen even in the case with a single impurity. Namely, positions of anomalies in $`q`$-space scale with $`k_\mathrm{F}`$, while the intensities of the anomalous parts depend on impurity strength. From these observations, it is clear that the origin of the anomalies is the fermionic responses of the itinerant electrons to impurities.

We speculate that randomness induces modulation of effective exchange coupling due to the Friedel oscillation. The random potential causes the $`2k_\mathrm{F}`$ singularity in the charge channel, which appears in the effective ferromagnetic coupling through the magnon-electron interaction. Scattering of magnons by the oscillating components of the exchange couplings creates the anti-crossing structure. The Friedel oscillation occurs universally in itinerant electron systems with randomness although its effect might be weaker in higher dimensions. Thus, the anomalies in three-dimensional systems in Fig. 2 may be an analog of the anti-crossing due to the Friedel oscillation.

We discuss our findings in comparison with experimental results. Our results reproduce the anomalies of spin excitations reported in low-$`T_\mathrm{C}`$ compounds qualitatively: (i) The spectrum shows a significant broadening near the zone boundary. $`^{\text{?, ?, ?, ?)}}`$ (ii) The softening near the zone boundary $`^{\text{?, ?, ?)}}`$ can be interpreted as a consequence of the branching or the anti-crossing in our results. The lower-energy branch has a substantial weight and is easily identified in neutron scattering experiments, while the higher-energy one may be difficult to identify due to a deficiency of energy resolution in experiments. (iii) The anti-crossing or the gap opening $`^{\text{?, ?)}}`$ have been reproduced. Thus, the randomness is reasonably considered as a hidden parameter to give a qualitative description for the anomalous spin dynamics in low-$`T_\mathrm{C}`$ manganites.

Several other scenarios have also been proposed to explain the spin excitation spectrum in low-$`T_\mathrm{C}`$ manganites. As an origin of the softening, some additional elements to the DE mechanism have been examined such as magnon-orbital $`^{\text{?}\text{)}}`$ and magnon-phonon couplings. $`^{\text{?}\text{)}}`$ Magnetic origin has also been proposed. $`^{\text{?}\text{)}}`$ For the origin of the broadening, magnon-phonon $`^{\text{?)}}`$ or magnon-electron couplings $`^{\text{?}\text{)}}`$ have been discussed. These arguments attribute the origin of the spectrum change by the A-site substitution to the modification of the bandwidth through some additional couplings to the spin degrees of freedom. As discussed previously, however, the actual changes of the bandwidth in these compounds are small. Therefore, it seems difficult for these scenarios to explain the large changes of the spin excitations in a quantitative manner. Moreover, neither the anti-crossing nor the gap opening has been reproduced by them.

Our results give a comprehensive scenario to understand the systematic changes of spin excitations in the A-site substituted manganites, since the change of the randomness are expected to be large, as previously discussed in reproducing the systematic decrease of $`T_\mathrm{C}`$. We suggest that the randomness gives important effects on various thermodynamic properties of A-site substituted manganites.

We now propose some experimental tests for our scenario. One is the sample dependence. Our scenario predicts a correlation between the sample quality and the anomalies in the spin wave spectrum. Therefore, even in compounds with the same chemical formula, the extent of the anomalies would depend on the purity of samples which is measured independently, for instance, by the residual resistivity. The second test is to investigate the location of the anomalies in the whole Brillouin zone, and to make a comparison with information of the Fermi surface, which is obtained by independent experiments, for instance, by the angle-resolved photoemission. The anomalies appear at the Fermi wave number in our scenario. Doping dependence is also crucial in this test. The third is to observe the high energy excitation near the band top $`𝒒(\pi ,\pi ,\pi )`$. In our results, there remains an excitation with substantial weight and relatively sharp peak-width near the band top even under strong randomness. This may be observed by the pulse-neutron scattering. All these experimental tests are feasible and may help to discuss the origin of the anomalies in the spin excitation as well as the role of the A-site substitution in CMR manganites.

Finally, we comment on higher-order corrections in the $`1/S`$ expansion. The aim of this Letter is to study the systematic changes of the spectrum by the A-site substitution of manganites. In this case, the parameter $`1/S`$ should be kept fixed, and the $`1/S`$ expansion does not contribute to the spectrum changes. Furthermore, it has been shown that the contribution of the higher order terms are very small at $`x1/3`$ where we are interested in, at least in the pure case. $`^{\text{?}\text{}\text{?}\text{)}}`$ Therefore, we conclude that the lowest order terms in the $`1/S`$ expansion is sufficient to discuss the randomness effects as a first step, although quantitative estimates of the $`1/S`$ corrections in the presence of randomness are left for further studies.

To summarize, we have investigated effects of randomness on the spin excitation spectrum in the double-exchange model. The spin dynamics in the ground state is calculated by the spin wave approximation in the lowest order of the $`1/S`$ expansion. Randomness induces some anomalies in the spectrum near the Fermi wave number such as the broadening, the anti-crossing and the gap opening. We speculate that the origin of the anomalies is modulation of effective exchange coupling by the Friedel oscillation. The anomalous spin dynamics in CMR manganites is reproduced qualitatively. Since the anomalies in experiments become conspicuous in A-site substituted compounds which have low Curie temperature, our results suggest that the randomness is an important element in the A-site substitution which has previously been understood as the bandwidth control. We have discussed our results in comparison with other theories and proposed several crucial experimental tests.

This work is supported by “a Grant-in-Aid from the Ministry of Education, Culture, Sports, Science and Technology”.


text:
# Simulation of evacuation processes using a bionics-inspired cellular automaton model for pedestrian dynamics

## 1 Introduction

Methods from physics have been successfully used for the investigation of vehicular traffic for a long time . On the other hand, pedestrian dynamics has not been studied as extensively . Due to its generically two-dimensional nature, pedestrian motion is more difficult to describe in terms of simple models. However, many interesting collective effects and self-organisation phenomena have been observed (see  for an overview and a comprehensive list of references), e.g. jamming and clogging, lane formation and oscillations at bottlenecks in counterflow or collective patterns of motion at intersections. These phenomena will be discussed in sec. 2.4.

The model takes its inspiration111Such “learning from nature” is the central idea of a field called Bionics. from the process of chemotaxis (see  for a review). Some insects create a chemical trace to guide other individuals to food places. This is also the central idea of the active-walker models used for the simulation of trail formation. In the approach of  the pedestrians also create a trace. In contrast to trail formation and chemotaxis, however, this trace is only virtual although one could assume that it corresponds to some abstract representation of the path in the mind of the pedestrians. Its main purpose is to transform effects of long-ranged interactions (e.g. following people walking some distance ahead) into a local interaction (with the “trace”). This allows for a much more efficient simulation on a computer.

The basic idea of our approach might be used for studying a variety of problems, especially from biology . Here we want to apply this model to a simple evacuation process with people trying to escape from a large room. Such a situation can lead to a panic where individuals apparently act irrationally. A nice discussion of empirical results can be found in . Our motivation here is rather the determination and classification of the different types of behaviour exhibited by the model than a realistic application.

The phenomena observed during panics can be quite different from those found in “normal” situations. Nevertheless it is desirable to have a model which is able to describe the whole spectrum of possible pedestrian behaviour in a unified way. So far mainly the social-force model  has been used which allows to reproduce the observed behaviour  quite accurately. In this continuum model the pedestrians interact by a repulsive (social) force which decays exponentially with the distance between them. This means that in each step of a simulation of $`N`$ individuals $`O(N^2)`$ interaction terms have to be evaluated. Furthermore, in complex geometries it occurs quite frequently that two pedestrians are rather close to each other but do not interact since they are separated by a wall. Therefore in principle one has to check for all pairs of individuals whether an interaction is possible or not. For large crowds this becomes very time consuming. In contrast, in the model used here pedestrians only interact with the floor field in their immediate neighbourhood. Therefore one has only $`O(N)`$ interaction terms. A further advantage is the discreteness of the model which allows for a very efficient implementation for large-scale computer simulations.

We start with a short summary of the models basic concepts.

## 2 Basic principles of the model

In the model the space is discretised into small cells which can either be empty or occupied by exactly one pedestrian. Each of these pedestrians can move to one of its unoccupied neighbour cells at each discrete time step $`tt+1`$ according to certain transition probabilities (see fig. 1).

The probabilities are given by the interaction with two floor fields . These two fields $`S`$ and $`D`$ determine the transition probabilities in such a way that a particle movement is more likely in direction of higher fields. They will be defined in the following subsections where also the basic update rules of the model are summarised.

### 2.1 The static floor field $`S`$

The static floor field $`S`$ does not evolve with time and is not changed by the presence of the pedestrians. Such a field can be used to specify regions of space which are more attractive, e.g. an emergency exit or shop windows. In case of the evacuation processes considered here, the static floor field describes the shortest distance to a an exit door, lying at the middle of the top wall of the room. Fig. 2 shows graphical representations of $`S`$ for different geometries.

$`S`$ is calculated due to a certain distance metric for each lattice site so that the field values are increased in the direction to the door. The field values are highest for the door cells. The explicit construction of $`S`$ can be found in Appendix A.

### 2.2 The dynamic floor field $`D`$

The dynamic floor field $`D`$ is a virtual trace left by the pedestrians and has its own dynamics through diffusion and decay. It is used to model an attractive interaction between the particles. Fig. 3 shows 3-dimensional plots of $`D`$ for three stages during evacuation process.

At $`t=0`$ for all sites $`(i,j)`$ of the lattice the dynamic field is zero, i.e. $`D_{ij}=0`$. Whenever a particle jumps from site $`(i,j)`$ to one of the neighbouring cells, $`D`$ of the starting place is increased by one: $`D_{ij}D_{ij}+1`$.

Therefore $`D`$ has only non negative integer values and can be compared to a bosonic field222A simplified one-dimensional model with a fermionic field has been used recently for the description of chemotaxis ., i.e. the bosons dropped by the pedestrians during their movement create the virtual trace. Thus the field value $`D_{ij}`$ corresponds to $`D_{ij}`$ bosons. The dynamic floor field is time dependent, it has diffusion and decay controlled by two parameters $`\alpha `$ $`[0,1]`$ and $`\delta `$ $`[0,1]`$, which means broadening and dilution of the trace. In each time step of the simulation each single boson of the whole dynamic field $`D`$ decays with probability $`\delta `$ and diffuses with probability $`\alpha `$ to one of its neighbouring cells. Finally this yields $`D=D(t,\delta ,\alpha )`$.

Note that the dynamic floor field $`D`$ is only altered by moving particles and therefore it corresponds to a virtual velocity density field, rather then a particle density field. In sec. 3.2 we discuss an alternative definition of the dynamic floor field $`D`$ as a pure particle density field and its consequences. In the following we give a short summary of the specific update rules used in the simulations.

### 2.3 Update rules

The update rules of the full model including the interaction with the floor fields have the following structure :

1. The dynamic floor field $`D`$ is modified according to its diffusion and decay rules (for details, see sec. 2.2).
2. For each pedestrian, the transition probabilities $`p_{ij}`$ for a move to an unoccupied neighbour cell $`(i,j)`$ (see fig. 1) are determined by the local dynamics and the two floor fields. The values of the fields $`D`$ and $`S`$ are weighted with two sensitivity parameters $`k_S[0,\mathrm{}[`$ and $`k_D[0,\mathrm{}[`$. This yields

$$p_{ij}=N\mathrm{exp}\left(k_DD_{ij}\right)\mathrm{exp}\left(k_SS_{ij}\right)(1n_{ij})\xi _{ij},$$
(1)
with

occupation number: $`n_{ij}=0,1,`$
obstacle number: $`\xi _{ij}=\{\begin{array}{cc}0\hfill & \text{for forbidden cells, e.g. walls}\hfill \\ 1\hfill & \text{else}\hfill \end{array},`$
normalisation: $`N=\left[{\displaystyle \underset{(i,j)}{}}\mathrm{exp}\left(k_DD_{ij}\right)\mathrm{exp}\left(k_SS_{ij}\right)(1n_{ij})\xi _{ij}\right]^1.`$
3. Each pedestrian chooses a target cell based on the transition probabilities $`p_{ij}`$ determined in the previous step.
4. The conflicts arising by any two or more pedestrians attempting to move to the same target cell are resolved by a probabilistic method333One can as well solve the arising conflicts using the procedure described in . It can be shown  that conflicts are important for a correct description of the dynamics.. The pedestrians which are allowed to move execute their step.
5. $`D`$ is increased by all moving particles.

The above rules have to be applied to all pedestrians at the same time, i.e. we use parallel update.

Note that we do not use here the so called matrix of preference  because it encodes a direction of preferred motion of the pedestrians, which they usually not have at the beginning of an evacuation process. In the following all information about the desired direction of motion is obtained from the floor fields.

In  it has been shown that simple, local update rules for the pedestrians and the two floor fields are sufficient to yield a richness of complex phenomena (see sec. 2.4). Obviously this route is superior concerning the computational efficiency and even allows for faster-than-real-time simulations of large crowds , e.g. in evacuation processes in public buildings.

### 2.4 Collective Phenomena

Pedestrian dynamics exhibits are variety of fascinating collective effects :

* Jamming: At large densities various kinds of jamming phenomena occur, e.g. when many people try to leave a large room at the same time . This clogging effect is typical for a bottleneck situation. Other types of jamming occur in the case of counterflow where two groups of pedestrians mutually block each other . This happens typically at high densities and when it is not possible to turn around and move back, e.g. when the flow of people is large.
* Lane formation: In counterflow, i.e. two groups of people moving in opposite directions, a kind of spontaneous symmetry breaking occurs. The motion of the pedestrians can self-organise in such a way that (dynamically varying) lanes are formed where people move in just one direction . In this way, strong interactions with oncoming pedestrians are reduced and a higher walking speed is possible.
* Oscillations: In counterflow at bottlenecks, e.g. doors, one can observe oscillatory changes of the direction of motion . Once a pedestrian is able to pass the bottleneck it becomes easier for others to follow him in the same direction until somebody is able to pass (e.g. through a fluctuation) the bottleneck in the opposite direction.
* Patterns at intersections: At intersections various collective patterns of motion can be formed . A typical example are short-lived roundabouts  which make the motion more efficient. Even if these are connected with small detours the formation of these patterns can be favourable since they allow for a “smoother” motion.
* Trail formation: Although human and animal trails are formed for rather different purposes their structures have some similarities . Often human trails are formed as a short-cut which makes it attractive to leave a paved path. Animal trails usually are related to chemotaxis and mark the way to food places.
* Panics: In panic situations many counter-intuitive phenomena (e.g. “faster-is-slower” and “freezing-by-heating” effects ) can occur. For a thorough discussion we refer to  and references therein.

In  it has been shown that the new CA model described in the previous subsection is – despite its simplicity – able to reproduce these observed collective effects. This is essential if one intends to use the model for real applications, e.g. the optimisation of evacuation procedures.

## 3 Evacuation Simulations

In the following we describe results of simulations of a typical situation, i.e. the evacuation of a large room (e.g. in case of fire). At this we focus on the influence of the sensitivity parameters $`k_D`$ and $`k_S`$ on the evacuation times in order to identify the different classes of behaviour exhibited by the model. As we will see interesting collective phenomena between the pedestrians lead to a nontrivial dependence of the evacuation times on $`k_D`$ and $`k_S`$. In sec. 3.2 we investigate possible alternative definitions of the dynamical floor field. Sec. 3.3 is devoted to a simple optimisation problem, namely evacuation from a room with two doors .

### 3.1 The impact of sensitivity parameters

The value of $`k_S`$, the coupling to the static field, can be viewed as a measure of the knowledge of the pedestrians about the location of the exit. A large $`k_S`$ implies a motion to the exit on the shortest possible path. For vanishing $`k_S`$, on the other hand, the people will perform a random walk and just find the door by chance. So the case $`k_S1`$ is relevant for processes in dark or smoke-filled rooms where people do not have full knowledge about the location of the exit.

The parameter $`k_D`$ for the coupling to the dynamic field controls the tendency to follow the lead of others. A large value of $`k_D`$ implies a strong herding behaviour which has been observed in the case of panics .

We consider a grid of size $`63\times 63`$ sites with an exit of one cell in the middle of one wall. The particles are initially distributed randomly and try to leave the room. The only information they get is through the floor fields. Fig. 4 shows typical stages of the dynamics for an initial particle density of $`\rho =0.3`$, which means 1116 particles.

In the middle picture of fig. 4 a half-circle jamming configuration in front of the door is easy to spot. A typical feature of the dynamics is a radial motion of ‘holes’ created by particles escaping through the door.

The evacuation times $`T`$ and their variances $`\sigma `$ are averaged over $`500`$ samples and strongly depend on the sensitivity parameters $`k_D`$ and $`k_S`$. Fig. 5 (a) shows the evacuation times for fixed sensitivity parameter $`k_D`$ of the dynamic field and variable sensitivity parameter $`k_S`$ of the static field and fig. 5 (b) shows the evacuation times for fixed $`k_D`$ and variable $`k_S`$. Fig. 6 shows all corresponding variances $`\sigma `$ and relative variances $`\sigma /T`$ of the evacuation times of fig. 5. The averaged evacuation times are always measured in update time steps. With the generic time-scale of the model  of about $`0.3`$ sec/timestep it easy to translate that into a real time value.

Let us first consider the case of $`k_D=0`$, i.e. no coupling to the dynamic field. In fig. 5(a) one can see the influence of $`k_S`$. For $`k_S0`$ the pedestrians do not sense the strength of the field. Therefore they do not have any guidance through the surroundings and perform a pure random walk which leads to a maximal value of the evacuation time for $`k_S=0`$. For $`k_S\mathrm{}`$ they have full information about the shortest distance to the door and the evacuation time converges towards a minimal value. The movement of the particles becomes almost deterministic. Therefore $`k_S`$ can be interpreted as some kind of inverse temperature for the degree of information about the inanimate surrounding of the pedestrians.

In the same way the sensitivity parameter $`k_D`$ of the dynamic field works as an inverse temperature for the information about the virtual trace. If $`k_S`$ is turned on from zero to infinity, a non-zero value of $`k_D`$ only means additional noise to the pedestrians and evacuation times increase for higher coupling strength to $`k_D`$ (see fig. 5(a)).

Much more interesting is the behaviour for fixed $`k_S`$ (fig. 5(b)). The evacuation times saturate at maximal values for growing sensitivity parameter $`k_D`$ of the dynamic field. The most interesting point is non-monotonic behaviour of $`T(k_D)`$ with the occurrence of minimal evacuation times for non-vanishing small values of the sensitivity parameter $`k_D`$ of the dynamic field. Therefore a small interaction with the dynamic field, which is proportional to the velocity density of the particles, represents some sort of minimal intelligence of the pedestrians. They are able to detect regions of higher local flow and minimise their waiting times. This effect is most pronounced for intermediate coupling ($`k_S0.4`$) to the static field $`S`$, but very weak for strong coupling to $`S`$ (e.g. $`k_S10`$), where the evacuation is dominated by the orientation to the inanimate surrounding. It vanishes again for very weak coupling to $`S`$ ($`k_S=0.1`$, not shown in fig. 5(b)) where the movement of the pedestrians is similar to a random walk.

If the coupling to the dynamic field is further increased, the evacuation times increase again and saturate at maximal values. The interaction with other pedestrians becomes more and more unfavourable similar to the arising of a panic situation. The weaker the coupling to the static field $`S`$ is, the higher are the evacuation times: the particles than have less information of the inanimate surrounding (for example they cannot find the way to the door because of smoke in a fire situation).

The variances $`\sigma `$ of the evacuation times (fig. 6) show qualitatively the same behaviour as the averaged time values. For strong coupling to $`S`$ ($`k_S10`$, $`k_D<2`$) $`\sigma `$ tends to zero, which means that the process is nearly deterministic. For weak coupling to $`S`$ the variance $`\sigma `$ increases strongly and also the relative variance $`\frac{\sigma }{T}`$ becomes rather large for $`k_S0`$. From a practical point of view, i.e. evacuations of real buildings, one has to ask whether it makes sense to specify safety just by the average evacuation time alone when (relative) variances can become large. For fixed $`k_S`$ and varying coupling $`k_D`$ to the dynamic field the variance $`\sigma `$ also behaves qualitatively as $`T`$, i.e. it increases with increasing $`k_D`$ (see fig. 6(b)).

In the following we discuss the influence of the decay and diffusion parameters $`\delta `$ and $`\alpha `$.

In fig. 7 (a) one can see that for $`k_S=0.4`$ and a high density of $`\rho =0.3`$ the effect of minimal evacuation times for non-vanishing small values of $`k_D`$ becomes most pronounced in the limit $`\alpha 0`$. Besides that all evacuation times are increased with increasing $`\alpha `$. That implies that diffusion of the field only means additional noise and no advantageous information for the particles due to the fast broadening and dilution of the trace for large $`\alpha `$.

For decay parameter $`\delta 0`$ the dynamic field persists for a long time and the regime of minimal evacuation times is shifted towards smaller $`k_D`$ values. Evacuation times for fixed $`k_D`$ increase monotonically with $`\delta `$ (see fig. 7(b)). This is a collective effect where the herding behaviour helps to overcome the insufficient knowledge about the location of the exit. In contrast, the evacuation times in the panic regime (i.e. high $`k_D`$ values) are increased with smaller $`\delta `$. With the field values of $`D`$ the memory of previous steps grows which increases the noise in the system for the particles. Therefore for large values of $`k_D`$ a strong herding behaviour would be unfavourable since it tends to ‘confuse’ the pedestrians which already have a good knowledge about the geometry. This reversed behaviour in comparison with the case $`k_D0`$ leads to a crossing of the curves in fig. 7(b). It is interesting to note that this crossing happens in a rather narrow interval around $`k_D1`$.

The influence of the diffusion parameter $`\alpha `$ of the dynamic field in the regime of minimal evacuation times strongly depends on the particle density $`\rho `$.

As in fig. 7(a) in fig. 8(a) one can see that in the high density regime an increasing diffusion parameter $`\alpha `$ always increases the evacuation times for arbitrary decay parameter $`\delta `$, since the diffusion of the field increases the noise for the particles. In the low density regime one finds the opposite behaviour. Here the diffusion of the dynamic field leads to a higher degree of information for the particles and to favourably long-ranged interactions between them (fig. 8(b)). For large $`\alpha `$ the evacuation time increases with decreasing $`\delta `$ whereas for small $`\alpha `$ one finds just the opposite behaviour. This leads to a crossing of the curves in fig. 8(a). The main reason is that for larger densities the evacuation time is mostly determined by the clogging which occurs at the exit whereas for small densities clogging is negligible. Furthermore for small $`\alpha `$ the traces are rather sharp and so a larger $`\delta `$ just reduces the field strength. In contrast, for large $`\alpha `$ the traces broaden quickly and are thus diluted such that they form an almost constant background.

Thus three main regimes for the behaviour of the particles can be distinguished. For strong coupling to $`k_S`$ and very small coupling to $`k_D`$ we find an ordered regime where particles only react to the static floor field and the behaviour than is in some sense deterministic. The disordered regime, characterised by strong coupling to $`k_D`$ and weak coupling to $`k_S`$, leads to a maximal value of the evacuation time. The behaviour here is typical for panic situations, e.g. the herding tendency dominates. Between these two regimes an optimal regime exists where the combination of interaction with the static and the dynamic floor fields minimises the evacuation time. Here the individuals have some knowledge about the location of the exits, but through some sort of cooperation the evacuation time $`T`$ can be optimised.

Fig. 9(a) shows the influence of an increasing particle density for the evacuation times in the three main regimes.

For very small densities one finds a very strong increase of the evacuation times with increasing density (fig. 9(b)). If there are only very few particles in the system, which are distributed randomly over the lattice at $`t=0`$, they reach the door nearly independently from each other during the evacuation. Than the evacuation times are nearly proportional to the largest initial distance to the exit. In this regime diffusion of the dynamic field $`D`$ provides advantageous information to the particles (see fig. 8 (b)). In the density regime $`\rho 0.005`$ a queue in front of the door begins to form that strongly controls the evacuation time which then increases nearly linearly with growing density due to clogging.

### 3.2 Alternative definitions of the dynamic floor field

In the following we want to give a brief discussion of the consequences of several variations and extensions of the definition of the dynamic field $`D`$ and the corresponding coupling parameter $`k_D`$. First we concentrate on the effects of a negative sensitivity parameter $`k_D<0`$, corresponding to repulsive interactions between the individuals, on the evacuation times. Fig. 10(a) shows the averaged evacuation times $`T`$ for a high density of $`\rho =0.3`$ and weak coupling $`k_S=0.4`$ to the static field for several diffusion parameter values $`\alpha `$ ($`\delta =0.3`$ fixed).

For all diffusion parameters $`\alpha `$ the averaged evacuation time is increasing for decreasing sensitivity parameter $`k_D<0`$. This is obvious since the dynamic field $`D`$ marks the regions of highest flow in the system (e.g. nearby the exit), because of its relation to the velocity density, and therefore an avoidance of this regions should lead to an increased evacuation time. An increasing $`\alpha `$ value weakens the sharpness of $`D`$ in the high flow regions and so the effect of the avoidance of this regions is moderated, i.e. the evacuation times decrease for increasing $`\alpha `$ and $`k_D<0`$.

More subtle consequences arise for the dynamics of the model if the definition of the dynamic field $`D`$ is changed. Up to now we have considered a dynamic field which is related to the velocity density in the sense that it is only altered by moving particles (see sec. 2.2). As an alternative we investigate a field $`\stackrel{~}{D}`$ which is altered by all particles: in each time step each particle increases the dynamic field value of its site by one, i.e. $`\stackrel{~}{D}`$ is proportional to the particle density. Fig. 10(b) shows the evacuation times for fixed coupling strength $`k_S=0.4`$ to the static field, variable $`k_{\stackrel{~}{D}}[10,10]`$ and the two diffusion parameters $`\alpha =0.0`$ and $`\alpha =0.1`$ for this density field $`\stackrel{~}{D}`$ ($`\delta =0.3`$ fixed). As a comparison fig. 10(b) shows again the corresponding times for the velocity density field $`D`$ (see fig. 10(a)) with $`\alpha =0.0`$. For positive coupling parameter $`k_D>0`$, i.e. attractive interaction, one finds monotonously increasing evacuation times for all diffusion parameters $`\alpha `$. The time values for different $`\alpha `$ begin to branch for $`k_D>1`$ such that for $`k_D>5`$ the evacuation time for a diffusion parameter $`\alpha =0.0`$ is three times higher than for $`\alpha =0.1`$. The reason for that becomes obvious from fig. 11(c) where a two-dimensional visualisation of the density-dependent dynamic field $`\stackrel{~}{D}`$ is shown.

The darkest shaded sites correspond to regions of highest $`\stackrel{~}{D}`$ values. Since $`\stackrel{~}{D}`$ is related to the particle density, regions of jamming due to high particle concentration are signified by high $`\stackrel{~}{D}`$ values (e.g. a half-circle jamming in front of a door, see fig. 11(a)). This effect is strongest for $`\alpha =0.0`$, because a non-vanishing diffusion parameter only smoothes the field gradients. Therefore an attractive interaction to such a particle density-dependent field leads to very high evacuation times due to clustering and herding effects known from panic situations, especially for vanishing diffusion of $`\stackrel{~}{D}`$.

For negative sensitivity parameter $`k_{\stackrel{~}{D}}<0`$, i.e. repulsive interaction, one finds a different behaviour. In the region of $`k_{\stackrel{~}{D}}]1,0[`$ the evacuation time for the density-dependent field case is minimized. For smaller $`k_{\stackrel{~}{D}}`$ values the evacuation times increase again monotonically. This is similar to the effect of minimized evacuation times for the velocity density related field $`D`$ and $`k_D]0,1[`$ (see fig. 5(b) and inset of fig. 10(b)). A slight repulsion from regions of high density combined with a weak guidance through the static ground field $`S`$ ($`k_S=0.4`$) should lead to a larger flow and to smaller evacuation times, as well as a slight attraction to regions of high flow (inset of fig. 10(b)). A strong repulsive behaviour supersedes the directed walk through $`S`$ and therefore is counterproductive. It leads again to increasing evacuation times. For $`k_D<0`$ the evacuation times for the particle density dependent field $`\stackrel{~}{D}`$ are always smaller than for the velocity density dependent field $`D`$, since a repulsion from regions of higher particle density should always be more favourable than repulsion from regions of higher flow. However, fig. 10 indicates that a density dependent field with repulsive interactions behaves in some aspects as a velocity density dependent field with attractive interactions.

### 3.3 Room with two doors

As a simple example for safety estimations in architectural planing we investigate how evacuation times change if a gap between two doors is increased from zero to a maximal value. We consider again an ordinary room with no internal structure of grid size $`102\times 102`$. We start with one door of width two cells in the middle of one wall. Then the door is split into two doors of one cell each separated by gaps ranging from $`2`$ to $`98`$ (fig. 12(b) and (c)). Fig. 12 shows typical stages of the dynamics for different gap sizes.

Averaged evacuation times are measured for all three main regimes introduced in sec. 3.1. Depending on the regime, fig. 13(a) shows a strong influence of the size of the gap for small and for large gaps.

In the ordered regime, the evacuation time is almost uninfluenced by the gap size. The reason is the almost deterministic motion of the particles which do not interact with each other through the dynamic floor field $`D`$. In the other two regimes, however, the value of the gap can have a large influence due to the strong interaction effects through $`D`$. For both the optimal and the disordered regime one finds minimized evacuation times for a wide area of the gap length ranging from about $`20`$ to $`80`$ cells. For small gaps the evacuation times increase due to the interactions between the pedestrians. For large gaps the presence of the side walls has a negative effect (fig. 12(c)) and so again the evacuation times increase. For intermediate gaps the crowd of pedestrians will be subdivided into smaller groups (fig. 12(b)), leading to more favorable interactions between them. This picture is confirmed by the averaged flow values $`J`$ in fig. 13(b). A gap size of $`40`$ cells between the doors leads to maximised flow values in the optimal and disordered regime.

## 4 Conclusions

We have studied simple evacuation processes using a recently introduced stochastic cellular automaton for pedestrian dynamics which implements interactions between individuals using an idea similar to chemotaxis. Due to its simplicity the model allows very high simulation speeds and is very well suited for the optimization of evacuation procedures even in complex situations.

We have focussed on studying evacuation times in a very simple evacuation scenario. The main purpose was to elucidate the influence of the various model parameters, especially of the coupling strengths $`k_S`$ and $`k_D`$ to the static and dynamic floor fields, to facilitate their interpretation.

Most important are the coupling parameters $`k_D`$ and $`k_S`$ to the dynamic and static floor fields, respectively. $`k_S`$ is a measure for the knowledge of the individuals about the geometry, especially the location of the exits. For dominating coupling to the static field the pedestrians will choose the nearest exit without much detour. For $`k_S0`$, on the other hand, they have no knowledge at all and will find the doors just by chance.

The dynamic floor field $`D`$ is a measure for the velocity density of the pedestrians. The parameter $`k_D`$ controls the herding behaviour. For dominating coupling to the dynamic field the pedestrians have a strong tendency to follow in the footsteps of others, e.g. because they hope that others have more knowledge about the location of exits. Such a behaviour is relevant for panic situations where this herding tendency becomes important and has been observed empirically . In fact, the ratio $`k_D/k_S`$ has very similar effects as the panic parameter introduced in .

The investigation of alternative definitions of the dynamic floor field has shown that the most realistic results are obtained for the velocity density dependent field with repulsive interactions as used in the original definition of the model. However, a field related to the particle density with repulsive interactions behaves in some aspects very similarly.

An important result of our investigations is that for achieving optimal evacuation times a proper combination of herding behaviour and use of knowledge about the surrounding is necessary. Then through cooperative behaviour of the individuals evacuation times can be minimized.

We want to emphasize that here only a very simple scenario has been studied in order to identify the different regimes of the model. For realistic applications a procedure to determine the coupling parameters is needed. This work is in progress and results will be published elsewhere . In the Appendix one possible idea is discussed (see also ). Furthermore the effects of disorder (e.g. different individuals $`j`$ having different couplings $`k_S^{(j)}`$ and $`k_D^{(j)}`$) are important. Another interesting open question concerns the details of the dynamics of evacuation processes in panic situations. We have assumed, as is also done for the panic parameter in , that the coupling constants do not change during the process. However, in many realistic situations it might be more appropriate to use a panic parameter which increases with time. We leave this important problem for future study.

It is suprising that the properties of the model investigated here are in many respects very similar to the social-force model  although the interactions are very different. In our approach pedestrians interact with the velocity-density through an attractive coupling whereas the social-force model uses a repulsive density-density interaction. It would be interesting to get a deeper understanding of these similarities.

## Acknowledgment

We like to thank K. Nishinari, C. Burstedde, F. Zielen and D. Helbing for useful discussions.

## Appendix A Construction of the static floor field $`S`$

The values $`S_{ij}`$ for the static floor field for sites $`(i,j)`$ for the geometries used in this paper can be calculated in the following way: The lattice shall be surrounded by obstacle cells, which are not traversable, except for a number of door cells at $`\{(i_{T_1},j_{T_1}),\mathrm{},(i_{T_k},j_{T_k})\}`$. The pedestrians can only leave the room through these door cells. The explicit values of $`S`$ in the examples studied here are than calculated with a distance metric:

$$S_{ij}=\underset{(i_{T_s},j_{T_s})}{\mathrm{min}}\left\{\underset{(i_l,j_l)}{\mathrm{max}}\left\{\sqrt{(i_{T_s}i_l)^2+(j_{T_s}j_l)^2}\right\}\sqrt{(i_{T_s}i)^2+(j_{T_s}j)^2}\right\}.$$
(2)
This means that the strength of the static floor field depends on the shortest distance to an exit. $`\mathrm{max}_{(i_l,j_l)}\left\{\sqrt{(i_{T_s}i_l)^2+(j_{T_s}j_l)^2}\right\}`$, where $`(i_l,j_l)`$ runs over all cells of the lattice, is the largest distance of any cell to the exit at $`(i_{T_s},j_{T_s})`$. This is just a normalisation so that the field values increase with decreasing distance $`\sqrt{(i_{T_s}i)^2+(j_{T_s}j)^2}`$ to an exit and is zero for the cell farthest away from the door.

This is only one possible representation of the explicit calculation of the static floor field $`S`$. Alternative constructions of $`S`$ do not change the qualitative results of the investigations, provided that the strength of the static floor field is increased in the direction to the exit . The use of a Manhattan metric is one possible example for a successful construction of $`S`$ for more complex geometries .

Once the floor field is specified from the geometry of the problem the coupling constant $`k_S`$ has to be determined such that a quantitative agreement with experimental data is obtained. Currently we are investigating different strategies . A simple method would be to determine $`k_S`$ such that the observed velocity and velocity fluctuations of the motion of a single pedestrian are reproduced (see ).