Selected results from the first five years of RHIC data taking are reviewed with emphasis on evidence for thermalization in central Au+Au collisions at $\sqrt{s_{NN}}=200$ GeV.

We present here a brief summary of new results on heavy quarks and heavy quarkonia from the PHENIX experiment.

Based on AMPT model, the elliptic flow $v_{2}$ of $\phi$ mesons which is reconstructed from $K^{+}K^{-}$ at the Relativistic Heavy Ion Collider (RHIC) energy has been studied. The results show that reconstructed $v_{2}$ of $\phi$ meson can keep the earlier information before $\phi$ decays and it seems to obey the number of constituent quark scaling as other mesons and baryons. This result indicates that the $\phi$ $v_2$ mostly reflects the parton level collectivity developed during the early stage of the collisions and the strange and light up/down quarks have similar azimuthal angular anisotropy properties at the hadronization.

We can establish a new picture, the perfect fluid sQGP core and the dissipative hadronic corona, of the space-time evolution of produced matter in relativistic heavy ion collisions at RHIC. It is also shown that the picture works well also in the forward rapidity region through an analysis based on a new class of the hydro-kinetic model and that this is a manifestation of rapid increase of entropy density in the vicinity of QCD critical temperature, namely deconfinement.

Elliptic flow at RHIC is computed event-by-event with NeXSPheRIO. Reasonable agreement with experimental results on $v_2(\eta)$ is obtained. Various effects are studied as well: reconstruction of impact parameter direction, freeze out temperature, equation of state (with or without crossover), emission mecanism.

The large elliptic flow observed in Au-Au collisions at RHIC is often put forward as a compelling evidence for the formation of a strongly-interacting quark-gluon plasma. The main argument is that the measured elliptic flow is as large as the value given by fluid-dynamics models that assume complete thermalization. It is argued that this claim may not be justified, since a detailed examination of experimental data rather suggests that the system created is not fully equilibrated at the time when anisotropic flow develops.

We investigate how thermalization of gluons depends on the initial conditions assumed in ultrarelativistic heavy ion collisions at RHIC. The study is based on simulations employing the pQCD inspired parton cascade solving the Boltzmann equation for gluons. We consider independently produced minijets with $p_T > p_0=1.3 \sim 2.0$ GeV and a color glass condensate as possible initial conditions for the freed gluons. It turns out that full kinetic equilibrium is achieved slightly sooner in denser system and its timescale tends to saturate. Compared with the kinetic equilibration we find a stronger dependence of chemical equilibration on the initial conditions.

Three-gluon to three-gluon scatterings lead to rapid thermalization of gluon matter created in central Au-Au collisions at RHIC energies. Thermalization of quark matter is studied from three-quark to three-quark scatterings.

We argue that isotropization and, consequently, thermalization of the system of gluons and quarks produced in an ultrarelativistic heavy ion collision does not follow from Feynman diagram analysis to any order in the coupling constant. We conclude that the apparent thermalization of quarks and gluons, leading to success of perfect fluid hydrodynamics in describing heavy ion collisions at RHIC, can only be attributed to the non-perturbative QCD effects not captured by Feynman diagrams.
We proceed by modeling these non-perturbative thermalization effects using viscous hydrodynamics. We point out that matching Color Glass Condensate inital conditions with viscous hydrodynamics leads to continuous evolution of all components of energy-momentum tensor and, unlike the case of ideal hydrodynamics, does not give a discontinuity in the longitudinal pressure. An important consequence of such a matching is a relationship between the thermalization time and shear viscosity: we observe that small viscosity leads to short thermalization time.

In the initial stage of the bottom-up picture of thermalization in heavy ion collisions, the gluon distribution is highly anisotropic which can give rise to plasma instability. This has not been taken account in the original paper. It is shown that in the presence of instability there are scaling solutions, which depend on one parameter, that match smoothly onto the late stage of bottom-up when thermalization takes place.

Because the initial shape of the QGP in a heavy ion collision is anisotropic, the momentum distribution becomes anisotropic after a short time. This leads to plasma instabilities, which may help explain how the plasma isotropizes. We explain the physics of instabilities and give the latest results of numerical simulations into their evolution. Nonabelian interactions cut off the size to which the soft unstable fields grow, and energy in the soft fields subsequently cascades towards more ultraviolet scales. We present first results for the power spectrum of this cascade.

I discuss recent advances in the understanding of non-equilibrium gauge field dynamics in plasmas which have particle distributions which are locally anisotropic in momentum space. In contrast to locally isotropic plasmas such anisotropic plasmas have a spectrum of soft unstable modes which are characterized by exponential growth of transverse (chromo)-magnetic fields at short times. The long-time behavior of such instabilities depends on whether or not the gauge group is abelian or non-abelian. I will report on recent numerical simulations which attempt to determine the long-time behavior of an anisotropic non-abelian plasma within hard-loop effective theory. For novelty I will present an interesting method for visualizing the time-dependence of SU(2) gauge field configurations produced during our numerical simulations.

Numerical 1D-3V solutions of the Wong-Yang-Mills equations with anisotropic particle momentum distributions are presented. They confirm the existence of plasma instabilities for weak initial fields and their saturation at a level where the particle motion is affected, similar to Abelian plasmas. The isotropization of the particle momenta by strong random fields is shown explicitly, as well as their nearly exponential distribution up to a typical hard scale, which arises from scattering off field fluctuations. By variation of the lattice spacing we show that the effects described here are independent of the UV field modes near the end of the Brioullin zone.

Based on hep-ph/0510121, we discuss further the numerical study of classical SU(2) 3+1-D Yang-Mills equations for matter produced in a high energy heavy ion collision. The growth of the amplitude of fluctuations as $\exp{(\Gamma \sqrt{g^2\mu \tau})}$ (where $g^2\mu$ is a scale arising from the saturation of gluons in the nuclear wavefunction) is shown to be robust over a wide range of initial amplitudes that violate boost invariance. We argue that this growth is due to a non-Abelian Weibel instability, the scale of which is set by a dynamically generated plasmon mass. We find good agreement when we relate $\Gamma$ to the prediction from kinetic theory.

The transport properties of certain strongly coupled thermal gauge theories can be determined from their effective description in terms of gravity or superstring theory duals. Here we provide a short summary of the results for the shear and bulk viscosity, charge diffusion constant, and the speed of sound in supersymmetric strongly interacting plasmas. We also outline a general algorithm for computing transport coefficients in any gravity dual. The algorithm relates the transport coefficients to the coefficients in the quasinormal spectrum of five-dimensional black holes in asymptotically anti de Sitter space.

We compute by numerical integration of the Dirac equation the number of quark-antiquark pairs produced in the classical color fields of colliding ultrarelativistic nuclei. Results for the dependence of the number of quarks on the strength of the background field, the quark mass and time are presented. We also perform several tests of our numerical method. While the number of qqbar pairs is parametrically suppressed in the coupling constant, we find that in this classical field model it could even be compatible with the thermal ratio to the number of gluons.

Calculations of nonequilibrium processes become increasingly feasable in quantum field theory from first principles. There has been important progress in our analytical understanding based on 2PI generating functionals. In addition, for the first time direct lattice simulations based on stochastic quantization techniques have been achieved. The quantitative descriptions of characteristic far-from-equilibrium time scales and thermal equilibration in quantum field theory point out new phenomena such as prethermalization. They determine the range of validity of standard transport or semi-classical approaches, on which most of our ideas about nonequilibrium dynamics were based so far. These are crucial ingredients to understand important topical phenomena in high-energy physics related to collision experiments of heavy nuclei, early universe cosmology and complex many-body systems.

Approximations based on the 2PI effective action are used to investigate the process of equilibration in phi^4 theory in 3+1 dimensions, both in the symmetric and broken phase. A special emphasis is put on the study of the kinetic and chemical equilibration.

We discuss how the dynamics of an exploding hot fireball of quark--gluon matter impacts the actual phase transition conditions between the deconfined and confined state of matter. We survey the chemical conditions prevailing at hadronization.

In a Relativistic Diffusion Model (RDM), the evolution of net-proton rapidity spectra with sqrt(s_NN) in heavy systems is proposed as an indicator for local equilibration and longitudinal expansion. The broad midrapidity valley recently discovered at RHIC in central Au + Au collisions at sqrt(s_NN)= 200 GeV suggests rapid local equilibration which is most likely due to deconfinement, and fast longitudinal expansion. Rapidity spectra of produced charged hadrons in d + Au and Au + Au systems at RHIC energies and their centrality dependence are well described in a three-sources RDM. In central collisions, about 19 % of the produced particles are in the equilibrated midrapidity region for d + Au.

The present status of the thermal model is reviewed and the recently discovered sharp peak in the $K^+/\pi^+$ ratio is discussed in this framework. It is shown that the rapid change is related to a transition from a baryon dominated hadronic gas to a meson dominated one. Further experimental tests to clarify the nature of the transition are discussed. In the thermal model the corresponding maxima in the Xi/pi and Omega/pi ratios occur at slightly different beam energies.

We discuss the role of noise and dissipation in the explosive spinodal decomposition scenario of hadron production during the chiral transition after a high-energy heavy ion collision. We use a Langevin description inspired by nonequilibrium field theory to perform real-time lattice simulations of the behavior of the chiral fields. Preliminary results for the interplay between additive and multiplicative noise terms, as well as for non-Markovian corrections, are also presented.

The transverse energy and the charged particle multiplicity at midrapidity are evaluated in a single-freeze-out model for different centrality bins at RHIC at $\sqrt{s_{NN}}=130$ and 200 GeV. The predictions of the model are done at the freeze-out parameters determined earlier from measured particle yields and $p_{T}$ spectra. The results agree qualitatively well with the experimental data.

We provide a renormalization procedure for Phi-derivable approximations in theories coupling different types of fields. We illustrate our approach on a scalar phi^4 theory coupled to fermions via a Yukawa-like interaction. The non-perturbative renormalization amounts to fixing the scalar coupling via a set of nested Bethe-Salpeter equations coupling fermions to scalars.