Thermalization of holographic Wilson loops in spacetimes with spatial anisotropy

In this paper, we study behaviour of Wilson loops in the boost-invariant nonequilibrium anisotropic quark-gluon plasma produced in heavy-ion collisions within the holographic approach. We describe the thermalization studying the evolution of the Vaidya metric in the boost-invariant and spatially anisotropic background. To probe the system during this process we calculate rectangular Wilson loops oriented in different spatial directions. We find that anisotropic effects are more visible for the Wilson loops lying in the transversal plane unlike the Wilson loops with partially longitudinal orientation. In particular, we observe that the Wilson loops can thermalizes first unlike to the order of the isotropic model. We see that Wilson loops on transversal contours have the shortest thermalization time. We also calculate the string tension and the pseudopotential at different temperatures for the static quark-gluon plasma. We show that the pseudopotential related to the configuration on the transversal plane has the screened Cornell form. We also show that the jet-quenching parameter related with the average of the light-like Wilson loop exhibits the dependence on orientations.


Introduction
Wilson loops are known to play a key role as fundamental probes of gauge theories, in particular QCD. Owing to Wilson loops one can define many important quantities, for instance, we can derive the potential of a quark-antiquark interaction from the expectation value of the space-time rectangular Wilson loop. In the lattice QCD the Wilson loops are the prime observables and their expectation values are defined non-perturbatively [1]. One can also determine expectation values of Wilson loops in the framework of perturbative QCD after suitable renormalization [2].
At the same time, expectation values of Wilson loops are used to characterize properties of the quark-gluon plasma produced in heavy ion collisions (HIC). With the help of Wilson loops one can perform the analysis of radiative parton energy loss, quarkonium suppression, jet quenching, etc. [3].
In this paper, we investigate the behaviour of Wilson loops in the nonequilibrium spatially anisotropic quark-gluon plasma (QGP) produced in heavy-ion collisions in the framework of the holographic correspondence. It is widely appreciated that the holography, or the gauge/gravity duality, provides a powerful tool for studying dynamics of the strong coupling system, in particular, the QGP formed in heavy-ion collisions [3,4,5]. Main idea of this approach is using natural prescriptions of the generalized AdS/CFT correspondence to recover non-perturbative QCD phenomena, particularly, non-perturbative vacuum phenomena, finite temperature, high-dense and non-zero chemical potential phenomena. In this strategy, fitting parameters are ones specifying the form of the 5-dimensional metric. The 5-dimensional background is supposed to be a solution of Einstein equations with a suitable matter content, not necessary related with string theory.
According to the holographic approach, creation of the quark-gluon plasma in HIC can be represented here as a collision of shock waves in the 5-dimensional bulk, in which a black hole is formed [6]- [15]. After the collision the shocks slowly decay, leaving the plasma described by hydrodynamics. The formation of the black hole can be also described by an infalling shell [16] propagating in the 5-dimensional bulk. A gravitational collapse of the thin shell to the black hole provides also a gravitational dual description of a more general class of thermalization processes [17]- [19]. Note, that the holographic approach is convenient to incorporate anisotropic properties of the QGP created in heavyion collisions [20,21].
The holographic approach has been widely used to study Wilson loops in different settings. Expectation values of Wilson loops within the gauge/gravity duality have been calculated for the strongly coupled N = 4 super Yang-Mills theory [22,23,24,25]. The string dual description of the real QCD is unknown in spite of a lot of performed effort to find it [26,27,28]. However, suitable "bottom-up" holographic QCD models matching with experimental and lattice results have been proposed in [29,30,31,32,33,34,35,36,37]. Various physical quantities, in particular, expectation values of Wilson loops, have been calculated holographically. Wilson loops in static anisotropic backgrounds, static non-relativistic background, Lifshitz backgrounds and backgrounds with hyperscaling violation have been examined in [21], [38]- [44] and refs therein. As it is known following the holographic dictionary, gravity duals of light-like Wilson loops can be used for calculation of the jet quenching parameter [45,46] controlling the description of medium-induced energy loss for partons in QCD. Holographic evaluation of the jet-quenching parameter for anisotropic QGP was considered in [41,42,21].
The special feature of this paper is that we consider spatially anisotropic backgrounds, which have also boost invariance These backgrounds are characterized by the scale parameter R and anisotropic parameters a and ν. This spacetime is in accordance with the geometry of HIC, where one has a selected direction -the axis of ions collisions. As has been shown before [47], to fit the experimental form of the dependence of total multiplicity on energy, obtained at LHC, one can assume that the holographic model has boost invariance and a spatial anisotropy controlled by the so-called critical exponent. More precisely, it has been shown in [47] that the critical exponent should be taken ν = 4.45. The metrics of the form (1) are dual to models with so called Lifshitz-like fixed points [48,49]. These metrics are also occurs as the IR limit for the anisotropic background suggested for studies of QGP in [50]. It's also worth to be noted that the anizotropic background (1) differs from the Lifshitz metrics [51] by the anisotropic scaling of spatial coordinates. The Vaidya shell in the background (1) has been found in [53] and used to estimate the thermalization time of the 2-point correlators and the holographic entanglement entropy in HIC [52,53]. The holographic model with the shell in (1) was used to explore the behaviour of the quark-antiquark potential in [54] with the method developed previously for in the AdS-Vaidya background [55].
As already mentioned, the information about processes during HIC can be read off from the expectation values of Wilson loops. Thus, it is natural to study the behaviour of Wilson loops during the HIC within the same holographic model we used to fit the energy dependence of total multiplicity. In particular, it is reasonable to investigate thermalization of Wilson loops and their behaviour in the end of the thermalization process. By thermalization of Wilson loop we mean the thermalization of its expectation value.
In the present work we calculate holographically Wilson loops in the backgrounds (1) with the Vaidya shell and black brane. We consider spatial and lightlike configurations, which represent rectangles with two infinitely long sides and two sides of finite lengths. By virtue of the metric (1) possessing a spatial anisotropy, the expectation values of Wilson loops depend on the orientation of the corresponding configuration. Further, we will study potentials and its evolution during the thermalization process. We show that the order of thermalization in our background is the following: first thermalizes two-point correlator, then Wilson loops and the last thermalization occurs for the entanglement entropy. Our anisotropy reduces the thermalization time for non-local observables as compared to the isotropic case. We also find an analytical representation for the holographic light-like Wilson loops in the background (1) and derive the jet quenching parameter.
The paper is organized as the following. In Sec. 2 we briefly remind the holographic description of Wilson loops, gravitational backgrounds and the notations. In Sec. 3 we analyze the static Wilson loops as well as calculate pseudopotentials, string tensions for different orientations of Wilson loop and the jet quenching parameter. In Sec. 4 we study the nonequilibrium dynamics of the same oriented Wilson loops and present the results. In Sec. 5 the thermalization time of Wilson loops is estimated and we compare it with the thermalization times of the entanglement entropy and two-point correlation functions. We conclude with a discussion of our results and further directions.

Wilson loops
In this work we consider rectangular Wilson loops. As already noticed, Wilson loops contain the information about the force between quarks. Following the holographic approach [22,23] the expectation value of the Wilson loop in the fundamental representation calculated on the gravity side reads as: where C in a contour on the boundary. More precisely, we mean that the contour is at the "regularized boundary" z = z 0 . On the gauge side we suppose to deal with the Wilson loop in the anisotropic gauge theory, i.e. theory on the 4-dimensional flat anisotropic background, The regularized version of this theory corresponds to the gauge theory on the anisotropic lattice [56], where spacing in the longitudinal and transversal directions, a and a ⊥ , are different so that a /a ⊥ = (aR/z 0 ) 1−1/ν . One can develop renormalizations in this theory and also the renormalization of the Wilson loops in the analogy with the isotropic theory. F means the fundamental representation (we will omit this symbol in what follows), S string is the minimal action of the string hanging from the contour C in the bulk. The Nambu-Goto action can be represented as with h αβ = g M N ∂ α X M ∂ β X N . In (4) σ 1 , σ 2 parametrize the worldsheet, g M N is the background metric, M, N = 1, . . . , 5 and X M = X M (σ 1 , σ 2 ) specify the string worldsheet. The pseudopotential of the interquark interaction can be extracted from the rectangular spatial Wilson loop of size 1 × L, for large L and its defines the so called pseudopotential V. We note, that for large L the behavior of W ( , L) is different from the behaviour of the time-like one. Then the pseudopotential can be straightforwardly extracted from the string action (4) as follows As it is known from the QCD lattice calculations the spatial Wilson loops obey the area law at all temperature, i.e. the spatial string tension σ s is given by The quantity σ s differs from the usual string tension which is defined from time-like Wilson-loops. By virtue to the non-Abelian Stokes formula equal time spatial Wilson loops [57] are related with the spatial components of the energymomentum tensor and by this reason σ s is also called the magnetic string tension. Spatial Wilson loops have been studied on the lattice [58,59], analytically [60], and also within the gauge/gravity duality [62,61].

The boost invariant anisotropic metrics
We will study spatial Wilson loops in gravity backgrounds with spatial anisotropy, given by (1), where the critical exponent ν controls the deviation from isotropic case. In [53] we have called this metric as an Lifshitz-like metrics. In this paper to avoid a misleading with the Lifshitz metrics we mainly call the metric (1) as the boost invariant spatial anisotropic metic. Note, that in (1) we used the standard dimensional coordinates. In what follows we take R 2 = 2πα , a = 1 and use the dimensionless coordinatest = t/(2πα ) 1/2 ,z = z/(2πα ) 1/2 etc., and we remove tilde on the top of coordinates. So, the metric take the form One can see that the background (8) with ν = 1 comes to be the 5-dimensional AdS spacetime. As already mentioned, the choice of this metric is motivated by the fact that holographic estimations of the total multiplicity performed in this background reproduce the experimental dependence of the multiplicity on the energy [47].
This background (9)-(10) describes holographically the anisotropic media on the boundary with the temperature corresponding to the Hawking temperature of the black brane: To study the thermalization process, corresponding to the black brane formation in the dual language, we will use the Vaidya generalization of solution (1) [53]: The metric (12) has been written in ingoing Eddington-Finkelstein coordinates (v, r). The function m(v) in (13) defines the thickness of the shell smoothly interpolating between the zero-temperature (1) at v = −∞ and black brane backgrounds (9) at v = ∞. We choose the following form for the function f where α is a parameter. For the calculations in this paper we keep α = 0.2.

Spatial Wilson loops in a time-independent background
In this work we consider rectangular Wilson loops in the static background (9) located in the spatial planes xy 1 (or xy 2 ) and y 1 y 2 . One can delineate the following possible configurations: • a rectangular loop in the xy 1 (or xy 2 ) plane with a short side of the length x in the longitudinal x direction and a long side of the length L y1 along the transversal y 1 direction, so that • a rectangular loop in the xy 1 plane with a short side of the length y1 in the transversal y 1 direction and a long side of the length L x along the longitudinal x direction, • a rectangular loop in the transversal y 1 y 2 plane with a short side of the length y1 in one of transversal directions (say y 1 ) and a long side of the length L y2 along the other transversal direction y 2 , namely In this section we perform all calculations in the static spacetime (9) with (10) (using Eddington-Finkelstein coordinates). We start from the rectangular Wilson loop in the xy 1 -plane assuming that the large extent is oriented in the y 1 -direction (see (15)). We parametrize the world-sheet of the string in the following way Taking into account (4) and (12) with the stationary f given by (10), the Nambu-Goto action can be presented as where it is supposed ≡ d dx . The subscript x, y 1(∞) in the LHS of (18) indicates the orientation of the loop contour.
The action (18) on the time-independent string configurations after division on the length of the Wilson loop in the y 1 -direction can be rewritten in the form The RHS of (19) defines a dynamical system that has the first integral Note, that the same is true for the action in the form (18). Indeed, it defines the dynamical system with two degrees of freedom, z = z(x) and v = v(x), that for the stationary case has two integral of motions and excluding v due to the conservation law one comes back to (20). Since we have the symmetric boundary conditions, the top of the configuration z = z(x) is at x = 0, i.e. z * = z(0) and z (0) = 0, and z * is related with the first integral, Equation (21) and the boundary condition z(± x /2) = 0 give us the relation between the top point z * and the length So the Nambu-Goto action (18) can be rewritten in the form Here we use the regularized boundary conditions z(± x /2) = z 0 . The action (23) has a divergent when z 0 → 0 and one can present it in the form Subtracting the divergent part ν z 1/ν 0 one gets the renormalized Nambu-Goto The length (22) between the ends of the string admits the removing of regularization z 0 → 0 without renormalization Then pseudopotential V x,y 1(∞) is given by: In Fig.1 we present the dependence of the pseudopotential V x,y 1(∞) (27) on the length (26). We see that for small x the pseudopotential has the Coulomb part deformed by the critical exponent, thus where C 1 is some constant dependent on ν. Putting m = 0 in (25) and (26), we get the leading term (28), while taking corrections on m one has for ν = 4 The asymptotics for arbitrary ν are presented in Appendix A. One should note that for small enough x the behaviour of V x,y 1(∞) extracted from (25) for all ν has a form of the deformed Coulomb law with the power equal to 1/ν reproducing the Coulomb behaviour of the pseudopotential in the case ν = 1 (the AdS case). For large distances x the pseudopotential V x,y 1(∞) behaves as a linearly increasing function To see this behaviour let us note that to get the large x we have to take z * near z h , or z * /z h = w 0 → 1, and in this case the denominator in the integrands in the RHS of (26) and (25) behaves as and we get the log-behaviour for x and V x,y 1(∞) which lead us to the asymptotic where From Fig.1.b we see that the asymptotics of a V x,y 1(∞) at large x is linear with the slope given by formula (35). The formula (35) exhibits the dependence on ν shown in the zoomed inset of Fig.1 The string tension (35) can be also seen by dimensional analysis keeping the anisotropic parameter a in (1). Indeed, since we are working with metric (9) and we can assign a factor A to t, x and z and factor A 1/ν to y1 and y 2 . Since the action S is invariant under this rescaling we have that σ 1,s ∼ S

Rectangular strip infinite along the x-direction
Another possible configuration in the xy 1 plane is the rectangular Wilson loop whose contour is infinite along the x-direction and has a finite size along the y 1 -direction (see (16)). We specify this type of the configuration by the subscript y 1 , x (∞) . For the parameterization we take v = v(y 1 ), z = z(y 1 ) with boundary conditions z(± y1 /2) = 0. By virtue to this assumption the Nambu-Goto action (4) reads where it is supposed ≡ d dy1 . Similar to the previous case this action can be rewritten as and the first integral is related with the top point as Due to this relation we get the expression for the action and the length where z 0 is the regularization. Similar to the previous case one can remove regularization in (39) and in (40) after renormalization The pseudopotential V y1,x (∞) is related to (41) as In Fig. 2 we show the dependence of the pseudopotential extracted from the action (41) on the length y1 for different values of the temperature and the dynamical exponent. Now the pseudopotential has a power-law dependence on ν for small y1 , so that with some constant C 2 dependent on ν. One gets asymptotics (43) with m = 0 in (41) and (40), and taking corrections on m we obtain for ν = 4 However, for large distances the pseudopotential represents a linear function and using estimation, similar to (32) and (33) we get that is also in agreement with the dimensional analysis. From Fig. 2.b we also see that the dependence of V y1, x (∞) ( y1 , ν) on y1 at large length is linear for all the dynamical exponent ν with the slops slightly deviating from the AdS case (ν = 1).

Wilson loop on the y 1 y 2 -plane
Now we come to the spatial rectangular Wilson loop located on the y 1 y 2plane. Let us assume that the loop contour is infinite along the y 2 -direction and has the finite extent of the length y1 in the y 1 -direction (see (17)). We specify this type of the orientation by the subscript y 1 , y 2(∞) , choosing only transversal coordinates for the parameterization of the worldsheet σ 1 = y 1 , σ 2 = y 2 .
Taking into account, that z(y 1 ) satisfying z(± y1 /2) = 0 one can represent the string action (4) in the following form where it is supposed ≡ d dy1 . Similar to the previous cases this action can be rewritten as and the first integral is related with z * point as Due to this relation we get the expression for the action and the length where z 0 is the regularization. Similar to the previous case one can remove regularization in (51) directly and in (50) after renormalization. The renormalized action (50) in terms of the w-variable takes the form The relation for the length is given by Finally, the pseudopotential V y1, y 2(∞) extracted from (52) reads as: In Fig.3 we display the behaviour of the pseudopotential (54) on the length (53).
It is easy to see that the behavior of V y1, y 2(∞) in Fig.3 is rather different from two previous cases. From Fig.3 (b) we observe, that now the dependence on ν is driven by some constant C 3 relying on ν. It should be noted that the pseudopotentials strongly deviate from the AdS case (ν = 1) both in the UV and the IR regions of y1 . Thus, one can write for small y1 where C 3 is some constant dependent on ν. For ν = 4 one can write down (55) with the corrections on m in the following form .
For large y1 we have As in the previous cases the large y1 behaviour is provided by the pole near to w ∼ 1 in (53) and we get, compare with (32) and (33), Looking at formulas (55)-(56) one can conclude that the pseudopotential corresponding to the configuration on the transversal plane reproduces the form of the Cornell potential.

Spatial string tension dependence on the orientation
It is interesting to analyze the behavior of the spatial string tension σ s (7) for different orientations of the Wilson loop, i.e. the behaviour of σ s,i , i = 1, 2, 3 given by (30), (45) and (56). The temperature dependences of the spatial string tension in the confining background, which reproduces the Cornell potential [31], and the deconfining one have been studied in [62] and [63], respectively. In [64] the universal behaviour of the spatial string tension for multiquark configurations was found. In the AdS/QCD model [31] the string tension dependence matching lattice data was found in [61].
In Fig.4 (the left panel) we present the dependence of the spatial string tension √ σ s as a function of T for all cases of the orientation and for ν = 1, 2, 3, 4.
We see, that for the configurations located on the xy 1 -plane (partially longitudinal orientations shown by solid and dashed lines) the temperature dependence of the string tension for different ν are rather similar. The deviations of solid lines from dashed ones increase with increasing T . We also see that the string tension corresponding to the Wilson loop in the y 1 y 2 -plane (the totally transversal orientation shown by the dotted lines) differs from the behaviour of the Wilson loop including the longitudinal direction, showing less dependence on the temperature with increasing ν. All these plots indicate that the structure of chromomagnetic fields in our holographic model has strong dependence on its orientation. In the right panel of Fig.4 the spatial string tension √ σ s for different orientations for ν = 4 is presented. We note that the behaviour of magnetic Wilson loops in heavy ion collisions was worked out in [65,66] and is put in the context of the color glass condensate model. The universal behavior of a large magnetic Wilson loop was found to have a nontrivial power-law dependence on the loop area. They have also argued that in contrast to usual Coulomb phase behaviour, magnetic flux does not propagates uniformly in the transverse plane, but instead, it is concentrated in small domains. In our work the Coulomb phases of pseudopotentials are modified for orientations different from the transversal one.

The holographic light-like Wilson loops and jet-quenching
One of the important characterizations of heavy-ions collisions is the jet quenching. The jet-quenching parameterq introduced in [67] is related with the average of the light-like Wilson loop in the adjoint representation [68] W where C is a rectangular contour with large extension L − in a light-like direction and small extension y1 in a transversal one. In the holographic approach W A (C) is equal to the classical string action S of a string worldsheet configuration stretched on the contour C on the boundary of the holographic background [45] W A (C) = e 2iS .
In this section we focus on holographic light-like Wilson loops in the black brane background (9). Choosing different light-like directions we obtain the dependence of the jet quenching parameterq on orientations.
Let us calculate the holographic light-like Wilson loop for the light-like coordinates related with the longitudinal direction, x ± = t±x √ 2 . The metric (9) in these light-cone coordinates takes the form where After introducing the parametrization for the string worldsheet by coordinates τ and σ such that the string action takes the form where i in front of the integral comes from the √ − det γ, since we deal with the Lorentz signature (here γ is the induced metric on the string). Note, that G −− is positive definite function when z < z h . For small y1 from this expression one gets [40] Evaluating the integral in (64) explicitly we get For ν = 1 formula (65) reproduces the T 3 dependence of the jet quenching parameter in the isotropic quark-gluon plasma [45], while in the case of ν > 1 the dependence on the temperature of the jet quenching parameter is caused by the anisotropic parameter.
Let us consider now the contour C with the light-like direction y ± 1 = t±y1 √ 2 . The metric in these coordinates has the form (60) with slightly different coefficients as compare to (61) .
The small side of the contour C can be oriented in x or y 2 directions. In both cases the corresponding string actions have the form (63) and contain the metric coefficient G −− as in (66). This G −− is non-positive definite below horizon z h , 0 < z < z h , for ν > 1. This makes the string action (63) complex that one can interpret as a suppression of the jet quenching parameter.

Spatial Wilson loops in a time-dependent background
Now we move to consider the thermalization of rectangular Wilson loops in the Vaidya background (12)- (14), which describes collapsing geometry in the special anizotropic spacetime (1). We proceed in a similar manner as in the static case studying three possible configurations of spatial Wilson loops.

4.1.
Wilson loops on the xy 1 -plane 4.1.1. Rectangular strip infinite along the y 1 -direction As in Sec. 3 we start from the spatial rectangular Wilson loop on the xy 1plane with the assumption that one side of the loop is infinite along the y 1direction and the other has finite size along the x-direction (see (15)). Here we suppose the dependence v = v(x), z = z(x). The Nambu-Goto action takes the form similar to (18) which for ν = 1 coincide with the Vaidya-AdS equations [17]. We have to consider eqs.(68) with the following boundary conditions z(± x /2) = 0, v(± x /2) = t, where x is the length of the Wilson loop along the x-direction. To solve numerically the equations of motion (68), it is convenient to impose the initial conditions For a given solution (v(x), z(x)) to eqs. (68) we can compute the functional for the Nambu-Goto action (67). We note that the dynamical system governed by (67) has the following integral of motion where we denote Taking into account (69)-(70) one can represent (67) in the following form where z * is the turning point defined from the requirements z = v = 0 and related with J as z 1/ν+1 * = J −1 . Coming to integration with respect to the z-variable the expression (71) can be represented as One can observe that b(z) tends to be −1 for z → 0 and divergencies at z = 0 are similar to the static configuration and by this reason we put the regularization z 0 . Making a substraction we come to the renormalized version of S x,y 1(∞) , in which we can remove the regularization The pseudopotential is expressed as : In Fig. 6 we present the behaviour of the renormalized pseudopotential V x,y 1(∞) derived from the action (74) as a function of x at fixed time moments for different values of ν. We see that for small distances the pseudopotential behaves similarly for different values of t. This dependence strengthens with increasing ν. For large times we see that the pseudopotential equilibrates to its thermal value. The pseudopotential reaches saturation for enough large size of the strip. The value of the thermalization time grows with increasing x for all values of ν. At the same time, we observe that the saturation is reached faster for large ν.

Rectangular strip infinite along the x-direction
Now we consider the rectangular Wilson loop on the xy 1 -plane with the assumption that its contour is infinite along the x-direction while it has finite stretch along the y 1 -direction (see (16)). As in the previous section, we specify [fm] Figure 7: The pseudopotential Vy 1 ,x (∞) as a function of y 1 at fixed values of t for ν = 2, 4 ((a),(b) respectively). The plot for ν = 1 is the same as in Fig.6.a. Different curves correspond to t = 0.1, 0.5, 0.9, 1.4, 2 from down to top. In (14) we take m = 1.
this type of the strip by the subscript y 1 , x (∞) . Thus, the corresponding Nambu-Goto action can be represented by with the notation ≡ d dy1 . The time-independent analogue of (76) is given by (36).
The equations of motion following from (76) are It is worth to be noted that eqs.(77) match with (68) taken with ν = 1 and also reproduce those for the AdS-case. The boundary conditions for eqs.(77) read z(± y1 /2) = 0, v(± y1 /2) = t, where y1 is the length of the Wilson loop along the y 1 -direction.
As in the previous case the action (76) can be simplified on equations of motions. For this purpose we note that the dynamical system governed by action (76) has the integral of motion where Taking into account (78)-(79) the action (76) is represented in the form with the turning point z * related with J as z 1/ν+1 * = J −1 . As in the previous case we present S y1, x (∞) as where b is defined by Performing renormalization we come to The renormalized pseudopotential derived from the action (83) as a function of y1 for different values of t and ν is demonstrated in Fig. 7. Here we again observe that for small y1 the behavior of the pseudopotential is similar for different values of t. As in the previous case the behavior intensifies with the increasing value of the dynamical exponent. However, in contrast to the previous case there is no substantial dependence, on a given scale, of the thermalization time on the dynamical exponent. Comparing Fig.6 and Fig.7 one can see that for the same scale y1 the thermalization of the Wilson loop occurs faster for the configuration with a long extent in the y 1 -direction. The dependence on the dynamical exponent ν is also stronger for the latter case of the orientation.

Wilson loop on the y 1 y 2 -plane
Finally, we come to the configuration located on the y 1 y 2 -plane. We assume that this infinite rectangular strip is invariant along the y 2 -direction, see (17). As in Sect. 3.2 we use y 1 , y 2,(∞) for the subscript of the action where we define ≡ d dy1 . The equations of motion corresponding to (84) can be written down in the following form One can check that eqs. (85) coincide with (68) and (77) for ν = 1 as well as come to be the equations for the AdS case. The boundary conditions to be satisfied by eqs.(85) have the standard form z(± y1 /2) = 0, v(± y1 /2) = t, where y1 is the length of the Wilson loop along the y 1 -direction. Following our strategy we compute the functional (84) on a given solution to (85). We note that the integral of motion for the configuration governed by the action (84) reads where we define Plugging (86) into (84) we come to the following form for the functional of the Nambu-Goto action S y1, y 2,(∞) = L y2 where the turning point z * is related with J as z 2/ν * = J −1 . Changing the variable of integration from y 1 to z we get The renormalized action in terms of the z-variable reads Here the pseudopotential is expressed from (91) as: The dependence of the pseudopotential V y1, y 2(∞) on the length y1 is shown in Fig.8. As for the previous configurations of Wilson loop located on the xy 1plane, the pseudopotential V y1, y 2(∞) tends to its thermal value for large t. We note that the influence of the critical exponent on the rate of the thermalization process for V y1, y 2(∞) is even higher than for V x, y 1(∞) and V y1, x (∞) . In all these cases the saturation time increases with y1 and ν.

Thermalization times of spatial Wilson loops
In this section, we compare the thermalization time for spatial Wilson loops with different orientations and its dependence on the value of the dynamical exponent ν. To simplify these estimations we consider the thin shell limit. We are interested in the value of the boundary time t therm when the string profile is totally covered by the thin shell, i.e.: where is the length between the string endpoints on the boundary given by (22), (40) or (51) for different orientations. In Fig.9 we plot the dependence on of the thermalization time for two configurations in the xy 1 -plane. In Fig.10(a) the behavior of the thermalization time as a function of for the configuration in the transverse y 1 y 2 -plane is presented. One can see that the thermalization time decreases with increasing ν for all cases plotted in Fig. 9 and Fig. 10. The dependence on the length for the loop in the xy 1 -plane with the short extent in the x-direction is linear. At the same time, the dependence for the loop in the same plane, but with the short extent in the y 1 -direction, as well as for the loop in y 1 y 2 -plane, is not linear for small asymptoting to the linear dependence only for large . It should be noted that the deviation from linearity strengthens with increasing ν. We also see that for the configuration in the transverse plane, the deviation of the thermalization time for the anisotropic cases from the isotropic one the thermalization time is much stronger than for the other orientations. In Fig.10(b) the comparison of thermalization times for different orientations in the case ν = 4 is plotted. This plot shows that the dependence on the orientation is crucial, varying the orientation we change the order of thermalization time. This means that characteristic scale depends on the orientation. The similar behavior of the thermalization time on was observed for the thermalization time of two-point correlators in [53].

Thermalization times of different observables
It is interesting to compare the thermalization times of different observables. In our work [53] we have studied two-point correlation functions and the holographic entanglement entropy in the Lifshitz-like backgrounds. By virtue of the spatial anisotropy of the metric we had two different configurations of the correlators and entropy with respect to the longitudinal and transversal directions. We have observed that the entanglement entropy for a subsystem delineated in the transversal direction thermalizes faster then the two-point correlator and Wilson loop in the longitudinal one. In [53] we also have calculated thermalization times for two-point correlators. In Appendix B the additional computations for the thermalization time of the holographic entanglement entropy are given. In Figure 11 we show the comparison of the thermalization times for two point correlation functions, holographic entanglement entropy and Wilson loops for different configurations. We see that the order of the thermalization process in a certain direction is similar to the case in the isotropic and ordinary Lifshitz backgrounds [18]. The two-point correlator is the observable that thermalizes first, then we observe the thermalization of the Wilson loops, and the entanglement entropy is the observable that thermalizes last. One should be noted that the thermalization process of the Wilson loop and the entanglement entropy in the anisotropic background is faster even in the longitudinal direction than thermalization of the same observables in the isotropic case. From Figure 11 (a) it is also interesting to see that the curves for the two-point correlator in the x-direction in the anisotropic background and the correlator in the isotropic spacetime match.

Conclusions
In this paper, we have explored the holographic scenario of the formation of the quark-gluon plasma using the bulk backgrounds (12), which possess spatially anisotropy. To probe the formation of the quark-gluon plasma we have used the rectangular spatial Wilson loops located on the boundary of our background (12). We have considered three possible configurations of Wilson loops on the boundary: the infinite rectangular strip located on the plane including one longitudinal and one transverse directions, the xy 1 -plane, with a short extent in the x-or y 1 -direction, and the infinite rectangular strip located on the transverse y 1 y 2 -plane.
We have analyzed Wilson loops both for static and time-dependent cases using the static black brane (9), and Vaidya solutions (12), respectively. The results obtained in this paper show how the expectations of Wilson loops are modified in the presence of anisotropy in the strong coupling limit.
We have found, that at small distances, the pseudopotential derived from the Wilson loop located in the xy 1 -plane has a nontrivial dependence on the parameter ν. Namely, for ν > 1 a breaking of the Coulomb phase has been observed. For Wilson loops lying in the transverse plane the Coulomb phase is unbroken, and all dependence for small is encoded in the ν dependent constant. At large all pseudopotentials are linear growing functions. Also we have found that the magnetic string tension is also affected by the anisotropic parameter ν. For the contour located on the transverse directions, the dependence of the string tension on the temperature is suppressed by anisotropy, so the magnetic string tension becomes close to a constant value. This effect is clearly seen for large ν. We have also observed interesting results for Wilson loops in the Vaidya backgrounds (12). The effect of anisotropy parameter on the thermalization time for different orientations also has been investigated. In the x, y 1(∞) -case of orientation the dependence of thermalization time on scale is linearly growing function for all ν. For other orientation the dependence is not linear. In the transverse orientation this can be seen very clear. Until some critical value the dependence is slowly growing function, and after this critical value it shows the linear growth.
A common feature in the behavior of the pseudopotential is the tendency of achievement of the saturation for large values of the boundary time. The dynamical exponent also influences to the thermalization of Wilson loops, so that the value of ν increases this behavior strengthens. We have seen that the thermalization is much faster than in the AdS case. It is worth noticing that the approach to the saturation also depends on the orientation of the Wilson loop. The configuration on the transversal directions the system saturates quicker than for the contours on the xy 1 -plane.
Comparing to results for the evolution of the holographic entanglement entropy in the backgrounds (12) in [53] we have found that the thermalization process both of the Wilson loops and the entanglement entropy is faster in the transverse direction. We have seen the similar behavior of the thermalization time as a function of for the two-point correlators in [53].
We also calculated the dependence of the jet quenching parameter on the orientation and on the anisotropic parameter ν. We also noticed that for special orientations the string action defined the jet quenching parameter becomes complex that leads to a suppression of the jet quenching parameter. This phenomena requires more elaborations as well as the study of the time evolution of the jet quenching parameter modeling by the light-like Wilson loop in the time-depending background (12).
It would be interesting to generalize our results to the case of the modified backgrounds that describe confinement/deconfinement phase transition [69].

Acknowledgments
We are grateful to Oleg Andreev and Oleg Teryaev for useful discussions.

A. Asymptotics for static pseudopotentials
Here we collect the expressions for asymptotics of static pseudopotentials. In Figure A.12 we present the comparison of the pseudopotentials given by exact formula (25), (26) and the approximated formula (A.4).
Inverting (A.1) for arbitrary ν we obtain The pseudopotential V x,y 1(∞) can be evaluated approximately For large ν one can expand the coefficients and obtain ). where .
For large ν we have the expansion for (A.7) where The pseudopotential for general ν can expressed as 4V y1, y 2(∞) (z * ) √ πνz B. Thermalization times of holographic two-point correlators and entanglement entropy

B.1. Thermalization time of two-point correlators
Under the holographic approach one can find the thermalization time t therm of the two-point correlator at the scale using the Vaidya background. For this, one should consider a geodesic of a bulk particle with equal time endpoints located at the distance and find the time when the geodesic covered by the shell. In the Vaidya-Lifshitz background (12)-(13) we should study the thermalization in both longitudinal and transversal directions.
For the thermalization in the longitudinal direction we have the following relation for the length where w = z/z * and the turning point is assumed to lie above the horizon, i.e. z h > z * .
The distance in the transversal direction is given by The thermalization time of the two-point correlator in both directions is defined by .

B.2. Thermalization time of entanglement entropy
To study the thermalization of the entanglement entropy we should also consider configurations in the longitudinal and transversal directions.
The longitudinal length scale is given by The length for a subsystem delineated in the transversal direction is Here we are also interested in the value of the boundary time when the surface is covered by the shell, i.e. the thermalization time of the entanglement entropy has the same expression as for the two-point correlator (B.1).