text:
\section{The Cyclic Universe}

As a standard model of our universe, the big bang theory has
succeed in many aspects, especially in explaining the abundance of
primordial elements and prediction of the cosmic microwave
background radiation. However, its successes naturally lead the
physicists to another confused conception: the big bang, which
seems to be the beginning of our universe. Is it really the
beginning of our universe when everything suddenly emerged from
nothing or the bang? It is a radical notion for most physicists.
Actually, in the big bang theory, the big bang is nothing but a
moment when the temperature and density are extremely so high that
the physical theory, as far as we have known, is not available any
more. Then what is the universe before the big bang? The issue is
open.

Historically, an idea is that our universe is always existing with
no beginning and no ending, and the universe we live in today is
just one phase of infinite cyclic periods. Physicists revisited
this idea in the frame of the Einstein Gravitation, and developed a
kind of universe models which are called cyclic universe or
oscillating universe. In these models, the universe is homogeneous
with a periodical scalar factor and will orderly experience the
bounce, expansion, turnaround and contraction in each cycle.
Because the universe in cyclic models always rebirths in ``fire''
and will never end, it is also named as the Phoenix universe.

In different models, the evolution of cyclic universe is driven by
different physics. For example, in \cite{STS}, Steinhardt and Turok
introduce a cyclic model in high dimensional string theory with an
infinite and flat universe, while pointed by Brown\cite{Brown} and
Baum\cite{Baum} etc, the cyclic evolution of the universe can also
be realized with a modified gravitational equations. For other more
cyclic models one can see
\cite{BD,KSS,Piao04,Lidsey04,CB,Xiong,Xin,LS,Biswas,Biswas1,Cai0906},
and also see \cite{NB} for review. Although in some of models the
physics under such high energy is still uncertain, these studies
provide good insights and effective theories to describe the
dynamics of the cyclic universe.

As an alternate model of the early universe, the studies of cyclic
universe model not only bring a distinctive insight for the big
bang described in the standard universe model, but also promote the
study of and succeed in other problems. The explanation of the
primordial perturbations is just one of them.

\section{Cosmological Perturbations in Cycle}

Today, our universe has a well developed nonlinear structure which
takes the form of galaxies, clusters and super-clusters of
galaxies, and of filaments and voids on large scales. However, all
the inhomogeneities in the density distribution, when averaged over
hundreds mega-parsecs, are very small and can be seemed as the
perturbations of the homogeneous background. Moreover, if one
ascend to the time of recombination, it will be found that all of
these nonlinear inhomogeneous structures have not formatted and
only small perturbations, which are called primordial
perturbations, remain. By measuring the cosmic microwave
background, we can detect and even characterize these primordial
perturbations as gaussian, adiabatic fluctuations with a
scale-invariant power spectrum. There have been mature theories of
explaining how these primordial perturbations evolve and finally
lead the structure of our universe today. The details of this
process are so complicate that it can only be calculated
numerically, but fortunately the current supercomputers are strong
enough to deal with these calculations and have given very
beautiful results.

On the other hands, as the initial conditions of structure
formation, the primordial perturbations play an important role in
modern cosmology. However, what is the origin of primordial
perturbations, which is the seed of the structure formation, and
why these perturbations are gaussian, adiabatic and scale-invariant
are key challenges for the modern cosmology.

According to the present understanding, because the physical
wavelength of fixed co-moving scales is increasing less fast than
the Hubble radius in the late times (times later and including the
period of nucleosynthesis) when cosmological evolution is well
described in the big bang theory, we must look to the very early
universe to find an explanation for the observed structures.
Therefore, a successive model of the early universe must give an
appropriate explanation of origin of the primordial perturbations.

Cyclic universe models provide a good explanation about origin of
the primordial perturbations. In cyclic universe models, the prime
perturbations are generated in the previous cycle by quantum
fluctuations. During the contraction, since the Hubble radius
shrink faster than the wavelength of the perturbations with fixed
co-moving scales, the perturbations will leave out of the Hubble
radius. While after the bounce, which is seems like the big bang of
our cycle, the universe expands again, the perturbations leave our
the previous Hubble radius will reenter the horizon and then set as
the seeds of the structure formation.

\section{Proliferation of the Phoenix Universe}

Before we discuss the proliferation of the Phoenix Universe, we
would like to address a little further about the behaviors of the
cosmological perturbations as preparation. Although there are many
cyclic universe models, for the generality of the results, we will
not involve the details of model building. Actually, it is not
difficult, since behaviors of the perturbations is mainly dependent
on the background evolution, namely the evolution of the scale
factor. To be more specifical, we would like to consider the
universe is dominated by matter with equation of state $\omega$,
and the scale factor can be described as $a \sim
\eta^{\frac{n}{n-1}}$ where $\eta$ is the conformal time. If we
solve the equation of the curvature perturbations $\zeta$ in
momentum space, we can find that the linear perturbation with each
$k$ evolves dependently. For the perturbations with wavelength
smaller than Hubble radius, the amplitude of $\zeta_k$ evolves in
inverse proportion to the scale factor. While for perturbations
with wavelength larger than Hubble radius, we can find that the
solution of $k$-mode curvature perturbation $\zeta_k$ can be
written into two terms which can be denote them as $D_1$ and $D_2$
term, and $D_1$ is constant term while $D_2$ is depended on the
background evolution.

Then we can discuss perturbations in the cyclic universe by
considering a cyclic universe model, starting with the turnaround.
The universe is homogeneous with small perturbations generated by
the quantum fluctuations inside the Hubble radius. During the
contraction with $n>{1\over 3}$, some perturbations leave out the
Hubble radius and the amplitude of $\zeta$ will be dominated by
$D_2$ term, which is increased during this time. After the bounce,
the universe comes into the expanding phase in which the $D_2$ term
is decreasing while the $D_1$ term dominates $\zeta_k$. Thus,
perturbations with wavelength large than Hubble radius keep
constant during expanding until they reenter the Hubble radius and
leading the formation of the large structure. In a word, for a
cycle of cyclic universe during the contraction $\zeta_k$ is
increased on super Hubble radius scale, until the bounce of
corresponding cycle, while during the expansion it becomes
constant. Therefore, the net result is that $\zeta_k$ on large
scale is inevitably amplified.

Since the net amplification is different for each $k$, we can find
that for curvature perturbations generated with quantum fluctuation
spectrum, this process will change the spectrum index as
\begin{eqnarray}
n_\zeta-1=3-\left|{3n-1\over n-1}\right|,
\end{eqnarray}
which is scale invariant for $n\simeq {2\over 3}$, i.e. the
contraction with $\omega \simeq 0$. The above estimation shows a
way of how the primordial perturbation with a scale invariant power
spectrum be gained in the cyclic models.

However, the picture may be not that simple as our usual
expectation. As we have concluded, the amplitude of the
perturbations may be amplified during the contraction, and in
\cite{Piao0901,Piao1001,Zhang:2010bb}we suggest that, this
amplification may change the global picture of our cyclic universe.

\begin{figure}
\includegraphics[width=8.5cm]{Graphf4.eps}
\caption{\label{fig:4} $\zeta(\vec{x})$ in position space, which
reflects the inhomogeneity of two different backgrounds. The length
scale is in unit of the Hubble radius at that time.}
\end{figure}

We would like to recheck the above process. One should notice that
during the hole process mentioned above, the universe is supposed
to be homogeneous, and evolves with the same background, from the
beginning of contraction to the time of turnaround, which denotes
the end of one cycle, and the net increasing of amplitude of the
perturbations in one cycle, say the $j^{th}$ cycle, can be
described as
\begin{eqnarray}
{\cal P}_{\zeta}^{j+1}(k,\eta_B^{j+1}) \sim k^{\Delta n_s}{\cal
P}_{\zeta}^{j}(k,\eta_B^{j}),
\end{eqnarray}
where $\Delta n_s$ is negative for most of cyclic models. We can
imagine that, due to the amplification, the modes on large scale
will have the amplitude be about order one after many cycles. This
will lead to ${\delta \rho\over \rho} \sim 1$ on corresponding
super horizon scale at this time. Here we use horizon scale, a
different conception with Hubble radius but with same scale in most
of cases, since region outside horizon is causally uncorrelated.

In this case, it is obviously impossible that the different regions
of global universe will evolve synchronously, even if it is
synchronous in previous cycle. Namely, the global universe at the
beginning time of this cycle will be separated into many different
parts, each of which will evolve independently of one another.
While inside any given part, all perturbation modes origin from the
interior of horizon, where is still causally correlative. Thus in
this sense, each of such parts actually corresponds to a new
universe. To illustrate this effect of the perturbation on
background, we show the distribution of curvature perturbation,
which is homogeneous at the very beginning, by transforming the
power spectrum into position space after many cycles in
Fig.\ref{fig:4}.

In principle, each of these new universes will experience the
contraction, bounce and expansion, hereafter all or some of them
will enter into next cycle and proliferate again, and then the
above process may be repeated. In another word, the proliferation
will inevitably occur cycle by cycle, and we can have a cyclic
multiverse scenario. In this cyclic multiverse, the experience of
each universe after proliferation is generally not expected to be
synchronous. Thus when some universe are in a period of matter
domination, it is possible that there are many other universes
which are in the period of contraction or bounce or others. There
is also the proliferation of global universe in chaotic eternal
inflation, in which it is induced by the large quantum fluctuation
of inflation field in its horizon scale. Here, however, the
proliferation is induced by the cyclical amplification of
perturbation on super horizon scale, which is in classical sense,
thus it occurs cycle by cycle.

In conclusion, we have showed that the global configuration of
cyclic universe is more complex than people expected ever, which
actually shows itself a cyclic multiverse. It is found that if the
contracting phase with $w\simeq 0$ is included in each cycle of a
cycle universe, the increasing mode of metric perturbation is
inherited by the constant mode of curvature perturbation that lead
the amplification of metric perturbation on super horizon scale
cycle by cycle. Therefore after few cycles the universe will be
inevitably separated into many parts independent of one another,
each of which corresponds to a new universe and evolve up to next
cycle, and then is separated again. Therefore a cyclic multiverse
scenario, in which the universe proliferates cycle by cycle, is
actually presented. If the cyclic universe is like the Phoenix,
rebirth in the ``fire'' and will never ended, then we suggest that
this Phoenix will also proliferate eternally.

\textbf{Acknowledgments} The author would like to thank the editors
of the Journal of Cosmology for inviting this paper and Prof.
Yun-Song Piao for helpful discussions and comments.

\newpage

redpajama_set_name: RedPajamaArXiv


text:
\section{Introduction}

Over the past decades various models have been proposed for simulating traffic flow
(see e.g. \cite{Maerivoet2005,Chakroborty2006,Chandler1958,hoogendoorn2001,Brackstone1999,Chowdhury2000} and references therein),
which range on the  phenomenon treatment detail level (microscopic, mesoscopic and macroscopic models),
are based on individual events or global temporal behavior (discrete and/or stochastic),
implementing real time or off-line simulations.

Among these, the so called ``car-following'' models play an important role, as they capture
much of the realism of the phenomenon, with a high level of detail.
In particular, the
``Intelligent Driver Model'' (IDM) proposed by  Treiber and collaborators \cite{treiber1999,treiber2000}
has attracted much attention, as it is capable of reproducing diverse realistic effects.
This model is based on a system of coupled non-linear differential equations, where motion of each vehicle
depends on the position and speed of the adjacent vehicles. The solution of this system of equations
is only possible through numerical methods.

However, none of these models has up to date shown the ability to reproduce the
empirically observed effect of a maximum in the vehicle flux versus density curve,
which represents a fundamental phenomenon from the point of view of urban traffic planning.

In this work we implement a modified version of the IDM model, by assuming a truncated Gaussian
distribution for the ``desired" individual vehicle speed values, on top of the real street map of the suburb of Afogados
in the city of Recife, Brazil, where we also implement real traffic rules.
It is found that the composite effect of intermittent traffic lights, multiple merging lanes,
and random ``desired routes" of individual vehicles does produce a maximum (albeit weak) in the flux versus density curve
when a single avenue is considered as part of a larger composite system,
whereas this effect vanes off when the measurements are made on the district system as a whole.

To this end, we have developed a simulator with a ``Graphical User Interface'' (GUI), using animation to facilitate the correction of eventual problems, visualize the operation and the simultaneous interactions, and aid in a presentation which gives credibility to the model. Most importantly, this approach enables one to interactively adapt the model variables,
which greatly aids the realization of numerical experiments.

In the following section we give a brief overview of the necessary aspects of the traffic theory
and the IDM model, in the subsequent section we describe the numerical experiment and present the results,
and finally, we draw the conclusions.

\section{Modeling urban traffic}

\subsection{Basic concepts}

There are many different aspects of vehicle properties, driver behavior, and the surrounding traffic infrastructure that may be
taken into account when modeling urban traffic.
For example, the vehicle weight has direct impact on the ability of changing the vehicle speed \cite{Drew1968}. The movement of a vehicle $\alpha$ in a given lane of an avenue may be described by considering the vehicle size $l_{\alpha}$, its position in relation to the lane beginning $x_{\alpha}$, speed $v_{\alpha}$ and the acceleration  $a_{\alpha}$.
In addition, individual driver characteristics influence the unique form in the vehicle conduction \cite{Drew1968, Boer1999}, where the driver perception and the reaction time are the key elements for the efficient vehicle operation \cite{Drew1968}.
Furthermore, individual driver factors in the physiologic context \cite{Groeger1998} may be taken into account.

Generally speaking, urban traffic represents a highly complex non-linear phenomenon,
directly related with the vehicle and driver properties, the interactions between the vehicles,
the properties and configuration of the traffic infrastructure (streets, lanes, traffic lights etc.)
and the quantity of vehicles that enter and exit the observed region per unit time \cite{manual2000}.

The traffic flow $q$ is defined
\cite{Maerivoet2005, manual2000, Ashton1966}
as the number of vehicles N which pass through a detector during a given time interval $\Delta t$
\begin{equation}
q = \frac{N}{\Delta }\ ,
\label{eq:fluxo}
\end{equation}
where in the case of multiple lanes, the composite vehicle flow on $L$ lanes is simply given by
\begin{equation}
q = \frac{1}{\Delta t}\sum^{L}_{i=1}N_{i}\ .
\label{eq:fluxonf}
\end{equation}
Depending on the vehicle density and the street capacity,
traffic flow may be classified as belonging to one of the three basic regimes: free flow, synchronized traffic flow (when the street capacity is approached),
and ``wide moving jam'' (when the street capacity is exhausted) \cite{Maerivoet2005, helbing2002}.

Vehicles in a given lane are ordered due to the impossibility of overtaking within the lane \cite{Maerivoet2005}, the first in will be the first out, unless lane switching is permitted. The distance $S_{\alpha}$ between the vehicle $\alpha$ and the vehicle $\alpha - 1$ immediately ahead, is given by the distance from the back bumper of the vehicle in front to the front bumper of the following vehicle, represented by
\begin{equation}
S_{\alpha}  = x_{\alpha-1} - x_{\alpha} - l_{\alpha-1}\ ,
\label{eq:distancia}
\end{equation}
where $x_{\alpha}$ represents the vehicle $\alpha$ relative position within the avenue, and $l_{\alpha}$ is the vehicle length,
as presented schematically in Fig.~\ref{fig:veiculos}.
\begin{figure}
       \includegraphics[scale=0.30]{carros.eps}
       \caption{\label{fig:veiculos}
       Consecutive vehicles in a lane.}
\end{figure}
Taking the instantaneous speed of all the vehicles in the region at a given instant,
the average speed $v_{s}$ is calculated \cite{Maerivoet2005}
as the arithmetic mean of these values
\begin{equation}
v_{s}  = \frac{1}{N} \sum^{N}_{\alpha=1}v_{\alpha}\ .
\label{eq:velmedia}
\end{equation}
Finally, the vehicle density $\rho$ is defined \cite{Maerivoet2005} as the fraction of the space occupied by all the vehicles
\begin{equation}
\rho = \frac{1}{K} \sum^{N}_{\alpha=1}l_{\alpha}\ ,
\label{eq:ocupacao}
\end{equation}
where $K$ is the lane length.

Empirical studies of traffic flow indicate a deterministic relationship between density $\rho$, flow $q$ and the average vehicle speed $v_s$ \cite{Drew1968, Maerivoet2005, manual2000, Ashton1966, Guan2008}. This relation is fundamental to the traffic flow theory, expressed by the equation $q=\rho v_{s}$ and graphically represented by the Fig.~\ref{fig:relacao}, where knowledge of two of these three variables permits determining the third.
However, the functional relationship $q(\rho)$ remains unknown.
\begin{figure}[htb]
\centering
\includegraphics[scale=0.30]{relacao.eps}
\caption{Fundamental diagram relating the density $\rho$ to the flow $q$. The flow capacity $q_{cap}$ is reached at the critical density $\rho_{c}$. The average vehicle speed $v_{s}$ is found as the inclination of the straight line that passes through the origin.}
\label{fig:relacao}
\end{figure}

It is seen from Fig.~\ref{fig:relacao} that in the free traffic region, at low density values, the flow increases in a linear fashion with the  increase of density. The congested traffic begins when critical density is reached, which defines the flow capacity. From this point on, the flow begins to decrease with the density increase, corresponding to the regime of congested traffic.

At present there is no consensus in the literature as to the functional form of traffic flow versus density curve. Several approaches have been proposed up to date which end up extracting different functional forms from the same data sets. Also, there seems to be no general consent up to date as to the relevance of individual components of various existing models for the description of diverse effects of interest for the understanding of the urban vehicular traffic phenomenon.

\subsection{Intelligent Driver Model}

The ``Intelligent Driver Model'' (IDM) \cite{treiber1999,treiber2000} is governed by a system of non-linear
differential equations for the individual vehicles indexed by $\alpha$, given by
\begin{equation}
\label{eq:aceleracao}
\frac{\partial v_{\alpha}}{\partial t} = a \left[ 1 - \left(  \frac{v_{\alpha}}{v_{0}} \right) ^{\delta} - \left( \frac{S^{*}(v_{\alpha}, \Delta v_{\alpha}) } {S_{\alpha}} \right)^{2} \right],
\end{equation}
where $a$ is maximum (or intrinsic) acceleration, $\delta$ is the parameter controlling the rate at which the desired speed is approached,
 $v_{0}$ is the ``desired" vehicle speed, $S_{\alpha}$ is the distance to the preceding vehicle $\alpha-1$
\begin{equation}
\label{eq:follow}
S^{*}(v_{\alpha}, \Delta v_{\alpha}) = S_{0} + max \left(Tv_{\alpha} + \frac{v_{\alpha} \Delta v_{\alpha}}{2\sqrt{a b}} , 0\right) ,
\end{equation}
is the ``safety distance" between the vehicles $\alpha$ and $\alpha-1$, $S_{0}$ is the minimum security distance, $T$ is the breaking reaction time \cite{Leutzbach1988}, $\Delta v_{\alpha}$ is the difference between the velocities of the two vehicles, and $b$ is the deceleration parameter governing the influence of the velocity difference on the safety distance.

Beside an abundance of parameters (perhaps, over-abundance for the current purpose), this model contains intuitively attractive components for the description of vehicular behavior, and has been shown to display rather
realistic behavior in diverse situations (see \cite{treiber1999,treiber2000} and citing references).
The solution
of the model does not present physically unacceptable situations (such as negative velocities, vehicle overlapping etc.), and stable solutions
are attainable already with the Euler's method, while application of the Runge-Kutta method yields accurate results for describing effects such as e.g. traffic density waves propagation.

On the other hand, this model yields realistic effects only in the presence of obstacles or conflicting situations such as confluent lanes,
whereas in the case of free lane traffic it demonstrates trivial behavior (all vehicles transit in parallel with the same speed), without limit for the street capacity. It turns out that using a truncated Gaussian distribution \cite{Walck2007, boxmuller1958} for maximum desired speeds of individual vehicles, is sufficient to create a traffic congestion in an avenue without obstacles  \cite{Santos2008}, where street capacity becomes apparent through saturation of the flow curve.
However, the maximum of the flow curve at the critical density was not observed with this generalization of the original IDM model \cite{Santos2008}.

In order to investigate whether this phenomenon may be brought about by the additional confusion stemming from the complexity of the environment (multiple lanes, vehicle types, traffic rules etc.), in what follows we implement a numerical simulation implementing such a situation.

\section{Numerical experiment}

The area of study is located in the Recife Metropolitan Region (RMR), where a network of urban highways connects municipal districts and suburbs. According to Pernambuco State Traffic Department (DETRAN) data, in the year of $1990$ the RMR was traversed by some $250000$ vehicles, which has grown to over $700000$ in  $2007$, where over $400000$ of these are found in the city of Recife alone. This violent vehicular density increase
has been causing heavy traffic jams in several parts of the city over the past years, augmenting the attention of the authorities and researchers.
The study region itself represents one of the two principal transport corridors between the south metropolitan region and Recife downtown.

As already mentioned, the IDM model by itself does not yield results demonstrating critical behavior on the flux versus density scatter plot,
and the generalization of IDM which considers a distribution of desired velocities leads only to saturation of the flux versus density curve,
with no sign of diminishing flux when the density is further increased. In order to verify whether the complexity of an urban traffic system is able
to reproduce such an effect, in the current numerical experiment we have implemented various components encountered in the real life situation.
In particular, we consider the following basic elements:
\begin{itemize}
	\item The avenues are considered  in terms of length, number of lanes, crossings, and traffic lights;
  \item	Vehicles of different size and maximum acceleration are considered (trucks, busses and cars), with individual maximum and minimum ``desired" speed. Upon entry, each vehicle is (randomly) assigned a ``desired" itinerary.
  \item	Traffic rules are implemented in each avenue and lane, maneuvering is implemented in accord with the existing signalization, considering the possibility of changing lanes and overtaking slower vehicles, while avoiding collisions.
\end{itemize}
\begin{figure}[htb]
       \centering
       \includegraphics[width=3.4in]{afogados.eps}
       \caption{Street map of the Afogados suburb with indicated traffic flow directions.}
       \label{fig:afogados}
\end{figure}
We present in Fig.~\ref{fig:afogados} the street map of the Afogados suburb with six entries and seven exits, where the directions of the traffic flow are indicated.
All the avenue crossings were visited and examined physically, in order to take into account all the possible maneuvering options of a vehicle
in a given lane of a given avenue, which attempts to follow a given route, while honoring the existing traffic rules.

The traffic lights were modeled as elements that belong to the avenues, such that the vehicles detect the traffic light presence, and behave in accord with the stop or go state (in fact, this was implemented by inserting a virtual still vehicle at the crossroad when the traffic light is red, such that approaching vehicles sense its presence, and break according to the IDM rules).

As already mentioned, the vehicles entering the region are attributed (random) desired routes. In the attempt to follow this route, a vehicle chooses permitted maneuvers, attempting to choose a lane that optimizes access to a determined crossroad, and realizes lanes changing and overtaking in the attempt to reach the destination. To this end we implement a origin $\rightarrow$ destination route table taking into account all the possible  combinations of the entry and exit avenues that comply with the traffic rules.
While attempting to follow the desired route, the vehicle may also change lanes due to advantageous conditions (e.g. when directly in front of it there is a slower vehicle, and the adjacent lane is empty), while honoring traffic rules.
This basic setup creates a rather realistic behavior of individual vehicles, and it seems to capture the collective behavior
of the vehicle fleet, which is observed on the graphical interface, as shown in Fig.~\ref{fig:simulacao} .
\begin{figure}[htb]
       \centering
       \includegraphics[width=3.4in]{simulacao.eps}
       \caption{A snapshot of the graphic interface of the simulator during the simulation.}
       \label{fig:simulacao}
\end{figure}

In order to achieve a real time simulation on the current scale that permits realistic traffic behavior, we have opted for the Euler algorithm
with a $0.05s$ time increment, and a full scene refresh at the rate of two frames per second. While a higher order Runge-Kutta method yields higher
precision results, the computational demand on this simulation scale takes it out of the real time observation range on the current hardware, and we have verified that the overall behavior of the traffic system is not affected by the precision gain.

For the model parameters we adopt the values $a=1.5 m/s^2$, $b=2.0 m/s^2$, $T=1.2 s$, $S_{0}=2.0 m$, and $\delta =4$, while the values of the desired maximum speed for each vehicle are drawn from a normal distribution truncated at $V_{min}$ and $V_{max}$,
depending on the vehicle type, and in accordance with the local traffic rules,
with parameters given in Tab.~\ref{tab:velocidade}.
\begin{table}[h]
	\centering
		\caption{
		Parameters
		of the normal distribution truncated at $V_{min}$ and $V_{max}$,
		used to draw the desired maximum speed for different vehicle types
		(the values are given in $km/h$).
		}
		\begin{tabular}{l c c c c c}
		\hline
		Vehicle type & $V_{max}$ & $V_{min}$ & $\mu$ & $\sigma$ \\
		\hline
		cars              & 80 & 40 & 60  & 20 \\
		buses  & 50 & 30 & 40  & 5  \\
		trucks    & 40 & 30 & 35  & 2.5  \\
		\hline
		\end{tabular}
	\label{tab:velocidade}
\end{table}

The simulation is performed by inserting vehicles into all the incoming avenues (and lanes) at the edge of the simulation zone.
The insertion rate is gradually increased from $0.125$ cars per second per lane ($8 s$ delay between incoming vehicles), to continuous influx (no delay between inserting the vehicles), with successive decrements of $0.1 s$ for the delay between incoming vehicles. At each rate (altogether eighty values), ample equilibration time of $1800 s$ ($30 min$) is used to attain equilibrium, after which the vehicle count and individual speeds are recorded in each avenue segment bounded by street crossings and confluence points. This process was repeated $100$ times in order to establish the average values.

\section{Results}

Results of the simulation corresponding to the avenue labeled ``$00$" in Fig.~\ref{fig:afogados} are presented in Fig.~\ref{fig:results00},
where it is observed that the flux versus occupation density scatter plot does display signs of a maximum (and thus existence of critical density),
while the average speed versus occupation demonstrates large fluctuations in the low density regime, which diminish with density increase.

\begin{figure}
\centering
{
  \includegraphics[width=0.4\textwidth]{entrada00geral.eps}}
  \includegraphics[width=0.4\textwidth]{Graph23.eps}
  \caption{Simulation results for the avenue $00$.}
  \label{fig:results00}
\end{figure}

In order to alleviate the effect of large fluctuations and pinpoint the critical occupation density, in Fig.~\ref{fig:regressao} we display the same data as in Fig.~\ref{fig:results00}, together with the results of averaging the data in $2\%$ bins, and polynomial regression which indicates a critical occupation of $41.63\%$, with a maximum flow of $0.91998\; cars/s$ (or $3312\; cars/h$).

\begin{figure}[h]
       \centering
       \includegraphics[scale=0.23]{Graph21.eps}
       \caption{Polynomial regression of the third order for simulation data, avenue $00$.}
       \label{fig:regressao}
\end{figure}

Results integrated over the entire simulation region are displayed in Fig.~\ref{fig:geral},
which does not reveal a maximum in the flux versus occupation density scatter plot. Also, variation of
the average velocity is pronounced over the whole observed range of density occupation values.
This findings may be attributed to the fact (observed on the graphical interface during the simulation)
that the entire region never becomes entirely congested, rather, different regions take turns in
demonstrating free flow, and traffic jam states.

\begin{figure}[h]
  \centering
  \includegraphics[width=0.45\textwidth]{Graph3.eps}\label{fig:ocupacaogeral}

  \includegraphics[width=0.45\textwidth]{Graph6.eps}\label{fig:velocgeral}
  \caption{Collected gross data for the entire region.}
  \label{fig:geral}
\end{figure}

\section{Conclusions}

In this work we implement a modified version  of the Intelligent Driver Model on top of a real street map of a suburb in Recife, Brazil,
in order to verify whether the additional complexity brought about by the synchronized collective behavior of vehicles that attempt to follow
individual routes, is capable of producing the critical occupation density effect with a maximum of the flux versus density curve. This
behavior is empirically observed, and represents a fundamental effect from the point of view of urban traffic planning, but has not been reproduced
up to date by existing mathematical and numerical models.

The modification of the IDM model implemented in this work consists of several ingredients. In particular, we consider three types of vehicles
with different characteristics (length, maximum acceleration, and individual desired speed), multiple lanes with traffic rules (merging traffic, semaphores etc.), and predefined individual routes for every vehicle.

It turns out that different parts of the observed region go intermittently through phases of free flow and congested traffic behavior. If only a single avenue is observed, the composite effect of interaction with the neighboring regions is reflected in behavior that may be identified as critical
density (albeit weak in the current simulation). However, the phase offset between the regions (while a given region is congested, a neighboring region displays free flow, and vice versa) leads to canceling out of individual congestion effects, and the maximum in the flux versus density curve is not observed. Also, the average velocity as a function of occupation displays heteroscedastic behavior, with large fluctuations at low densities that diminish as the density is increased, whereas large average velocity fluctuations are observed for the entire occupation density range for the suburb as a whole.

We may conclude that the critical density phenomenon should be regarded as a local effect, brought about by the interplay of a given region with
the neighboring regions, as an exit point for one region represents an entry point for another. To the best of our knowledge, the current work represents the first report in the literature regarding a mathematical/numerical model capable of exhibiting this effect.

Further studies should be made as to the contribution of each of the implemented model components to the observed critical density effect, and whether it may be enhanced by inclusion of some others.
It should also be investigated whether other microscopic models are capable of demonstrating this phenomenon in a similar setup.

\section{Acknowledgements}
The authors acknowledge financial support of Brazilian federal agencies CNPq (project CASADINHO-620113/2008) and
CAPES (projects PROCAD-1396/2007 and PROCAD-2273/2008).