On the time evolution of holographic n-partite information

We study various scaling behaviors of n-partite information during a process of thermalization for n disjoint system consisting of n parallel strips whose widths are much larger than the separation between them. By making use of the holographic description for entanglement entropy we explore holographic description of the n-partite information by which we show that it has a definite sign: it is positive for even n and negative for odd n. This might thought of as an intrinsic property of a field theory which has gravity dual.


Introduction
Entanglement entropy for a spatial region A in a quantum field theory is defined by the von Neumann entropy of the corresponding reduced density matrix, S A = −Tr(ρ A log ρ A ).
Here ρ A is the reduced density matrix given by ρ A = TrĀ ρ withĀ being the complement of A and ρ is the total density matrix describing the state of the corresponding quantum field theory. In general, for a local field theory, the entanglement entropy is UV divergent and the coefficient of the most divergent term for spatial dimensions bigger than one is proportional to the area of the entangling region [1], while for the spatial dimension equal to one the divergent term is logarithmic (see for example [2,3] for two dimensional CFT). Entanglement entropy for two disjoint regions has been studied in [4][5][6][7]. We note, however, that for two disjoint regions A and B, it is more natural to compute the amount of correlations (both classical and quantum) between these two regions which is given by the mutual information. It is actually a quantity which measures the amount of information that A and B can share. In terms of the entanglement entropy it is given by property of the entanglement entropy, it is evident that the mutual information is always non-negative and it is zero for two uncorrelated systems. More generally one may want to compute entanglement entropy for a subsystem consisting of n disjoint regions A i , i = 1, · · · , n. Following the notion of mutual information for a system of two disjoint regions, it is natural to define a quantity, n-partite information, which could measure the amount of information or correlations (both classical and quantum) between them. Intuitively, one would expect that for n un-correlated systems the n-partite information must be zero. Moreover, for n separated systems it should be finite. Actually for a given n disjoint regions, there is no a unique way to define n-partite information and indeed, it can be defined in different ways; each of them has its own advantage. In particular in terms of the entanglement entropy one may define an n-partite information as follows [8] S(A i ∪A j ∪A k )−· · · · · ·−(−1) n S(A 1 ∪A 2 ∪· · ·∪A n ), (1.2) where S(A i ∪A j · · · ) is the entanglement entropy of the region A i ∪A j . . . with the rest of the system. Note that with this definition the 1-partite information and 2-partite information are, indeed, entanglement entropy and mutual information, respectively. It is clear that, for this definition, n-partite information for n ≥ 2 is finite. It is worth noting that the n-partite information (1.2) may be re-expressed in terms of (n − 1)-partite information as follows (1.3) Therefore the n-partite information I [n] may be thought of a quantity which measures the degree of extensivity of the (n − 1)-partite information. Moreover, in terms of the mutual information, the n-partite information (1.2) may be recast into the following form I [2] (A 1 , A i ∪ A j ) + n i=2<j<k I [2] (A 1 , A i ∪ A j ∪ A k ) − · · · +(−1) n I [2] (A 1 , A 2 ∪ A 2 · · · ∪ A n ). (1.4)

JHEP09(2015)165
Note that this quantity is finite for a system with n disjoint regions and is zero for n uncorrelated regions. In this paper we will mainly study the n-partite information defined in equation (1.2).
Entanglement entropy (or other measures defined above) may provide a useful quantity in studying the process of thermalization under which an excited state could thermalize to an equilibrated state. Indeed, being out of equilibrium the thermodynamical quantities such as temperature, entropy, pressure, etc. are not well defined during the process of thermalization. Evolution of a system after a global quantum quench [10] is an example of the thermalization. In a field theory, global quantum quench is a sudden change in the system which might be caused by turning on/off a parameter in the Hamiltonian of the system in an interval δt → 0. This change takes the system to an excited state for a new Hamiltonian with non-zero energy density that could eventually thermalize to an equilibrium state.
Although the entanglement entropy could be useful to probe the thermalization process, in general, for a generic quantum system it is difficult to compute it. We note, however, that for those strongly coupled systems which have gravitational duals [11], in order to compute the entanglement entropy one may employ its holographic description [12,13]. Of course when the system is time-dependent, one should use the covariant proposal of entanglement entropy [14]. For completeness and further use, we have reviewed the holographic computations of the entanglement entropy in an appendix.
In the context of AdS/CFT correspondence, the thermalization process after a global quantum quench may be described by the process of a black hole (brane) formation due to a gravitational collapse of a thin shell of matter. The corresponding metric of a collapsing shell of neutral matter in d + 1 dimensions is given by the AdS-Vaidya metric where ρ is the radial coordinate, x i 's (i = 1, · · · , d−1) are the spatial boundary coordinates and v is the null coordinate, note that the AdS radius is set to be one. Here, θ(v) is the step function, and hence, one is dealing with an AdS geometry for v < 0 while for v > 0 the geometry is an AdS-Schwarzschild black brane whose horizon is located at ρ H = m −1/d .
The above background having a theta function on it (m(v) = m θ(v)) could provide a gravitational description for a sudden change in a strongly coupled field theory which might be caused by turning on a source for an operator in an interval δt → 0. This change can excite the system to an excited state with non-zero energy density that could eventually thermalize to an equilibrium thermal state. Since we are considering a sudden change in the theory, it is then natural to think of the process as a thermalization after a global quantum quench.
Therefore in order to study entanglement entropy during the process of thermalization after a global quantum quench one needs to compute the holographic entanglement entropy in the above time-dependent background using its covariant proposal [14]. Indeed, using the above metric, time dependent behaviors of entanglement entropy has been studied in JHEP09(2015)165 several works including [15][16][17][18][19][20][21][22][23][24][25][26][27] 1 where it was shown that the entanglement entropy grows linearly with time then it saturates to its equilibrium value (see also [30]).
The above consideration may be compared with the results of [10] where the behavior of the entanglement entropy after a global quantum quench for a two dimension CFT was studied. Although the quench considered in this case is different from that studied holographically, there is an agreement between the results of these two different setups. Of course this agreement might be understood from the fact that in both cases one is considering the evolution of an excited state in a CFT. Keeping this distinction in our mind, in what follows we will refer to our setup as a gravitational description for a thermalization process after a global quantum quench.
Using the gauge/gravity correspondence the mutual information has also been studied in [31][32][33] where it was shown that the mutual information exhibits a phase transition from a positive value to zero as one increases the distance between two regions. Time dependent behavior of mutual information in a global quantum quench has been studied in [34][35][36][37] where it was numerically shown that the mutual information in a time dependent background is always non-negative if the solution satisfies null energy condition.
Tripartite information I [3] during a global quantum quench for a strongly coupled field theory has been also studied in [34,35] using its holographic description. It was shown, numerically, that the tripartite information for a strongly coupled field theory which has gravity description is always non-positive. Actually, it was observed in [8] that the holographic mutual information is monogamous. Therefore one may consider the monogamy condition of mutual information for a strongly coupled field theory as a necessary condition for a theory to have a gravity description.
The main aim of this article is to study different scaling behaviors of the n-partite information in the thermalization process of a strongly coupled field theory undergoing a global quantum quench using the holographic description. To do so, we will consider n disjoint parallel strips with the same width ℓ separated by distance h. Motivated by the study of mutual information, we will consider the case where ℓ ≫ h. In this case and under certain assumptions the expression for n-partite information will be simplified significantly so that we could study its scaling behaviors analytically. Indeed taking into account that n-partite information may be expressed in terms of the entanglement entropy of different entangling regions, one may utilize the procedure of [38,39] to compute the corresponding entanglement entropy and thereby to study the evolution of n-partite information during a global quantum quench. By making use of this procedure we will show that the holographic n-partite information has definite sign, which following [8] might be thought of as a necessary condition for a theory to have a gravity description.
The paper is organized as follows: in the next section, we will review the computations of holographic mutual information in static backgrounds. In order to explore the procedure we will first study the scaling behaviors of mutual information in the thermalization process after a global quantum quench in a strongly coupled field theory, in section 3. The time evolution of the n-partite information will be discussed in section 4 and 5. In section 6 we will present numerical results to examine our analytical results. Finally section 7 is devoted to conclusions and discussions. We have enclosed the paper with some further details of calculations in an appendix.

Mutual information for static backgrounds
In this section we review the computation of holographic mutual information for two parallel strips in static backgrounds. The backgrounds we will be considering is either an AdS geometry or an AdS black brane which could provide gravitational descriptions for a conformal field theory at the ground state or a thermal state, respectively. To fix our notation, let us consider two parallel infinite strips with the equal width ℓ separated by a distance h in a d-dimensional field theory as depicted in figure 1.
It was argued in [31] that the holographic mutual information undergoes a first order phase transition as one increases the distance between two strips. Indeed, there is a critical value of h ℓ above which the mutual information vanishes. This peculiar behavior has to do with the definition of entanglement entropy of the union A ∪ B. Holographically this phase transition may be understood from the fact that for a given two strips with the widths ℓ and distance h, there could be two minimal hypersurfaces associated with the entanglement entropy S(A ∪ B) and thus the corresponding entanglement entropy behaves differently. More precisely one gets where S(l) is the entanglement entropy of a strip with width l. From the above expression and the definition of the mutual information (1.1), it is then clear that in the case of h ≫ ℓ, the mutual information becomes zero, while for h ≪ ℓ, one finds

JHEP09(2015)165
In what follows we will consider h ≪ ℓ. Therefore to find the mutual information of two parallel strips depicted in figure 1, one essentially, needs to compute the entanglement entropy of three strips 2 with widths h, ℓ and 2ℓ + h. To make the paper self contained we have reviewed the computation of holographic entanglement entropy in appendix A.1. By making use of the holographic entanglement entropy of a strip in a d-dimensional CFT whose gravity dual is provided by the AdS geometry (see (A.4)) the mutual information for two parallel strips for the vacuum state is found . Note that h ≪ ℓ condition guarantees the positivity of the resultant mutual information.
Let us consider the mutual information of the same strips for a thermal state whose gravity dual is provided by an AdS black brane metric The corresponding mutual information may be analytically expressed in certain limits and it is illustrative to study it in these limits. For example, one may assume that ℓ ≪ ρ H in which all entanglement entropies 3 involving in the computation of mutual information, (2.2), may be expanded as (A.6) leading to where On the other hand using the holographic renormalization one finds that the dual excited state has non-zero energy which is given by Therefore combining equations (2.5) and (2.7) one arrives at In the light of the first law of entanglement thermodynamics [40][41][42][43] one may think of the 2 If the widths of two parallel strips in figure 1 are not the same, one should compute four entanglement entropies corresponding to ℓ1, ℓ2, h and ℓ1 + ℓ2 + h. 3 Note that since h ≪ ℓ one also has h ≪ ρH .

JHEP09(2015)165
above equation as the first law for mutual information. Note that due to the minus sign it is obvious that as one increases the energy by ∆E the mutual information decreases by ∆I such that equation (2.8) holds. In other words it indicates that the mutual information of two static regions is maximal when the system is in the vacuum state.
On the other hand for the case of h ≪ ρ H ≪ ℓ, the corresponding entanglement entropy for the region h should be approximated by equation (A.6), while for those of ℓ and 2ℓ + h one has to use the large entangling region expansion given by equation (A.12). Therefore one arrives at [23] Here c 2 is a positive number which can be found numerically using, e.g. Mathematica (see appendix A.1). For example for d = 3, 4 it is c 2 = 0.88, 0.33, respectively. It is worth noting that although the term containing c 2 is subleading in the expression of the entanglement entropy at large entangling region, it plays an important role in the expression of mutual information and indeed it might be as important as the other terms. Finally for ρ H ≪ h and ρ H ≪ ℓ the mutual information is identically zero [23].

Time evolution of mutual information
In this section we shall study the scaling behavior of the holographic mutual information during the process of thermalization. We will consider a case where the quench occurs in a time interval δt → 0 so that the corresponding gravitational description of the process may be provided by AdS-Vaidya metric given by (1.7). We will compute the holographic mutual information for the parallel strips depicted in figure 1. We must emphasis that in the following studies the resultant semi-analytic expansion just gives us a piece wise and not a smooth function for mutual information. Following our discussions in the previous section we will assume h ≪ ℓ so that to find the mutual information, essentially, one needs to compute holographic entanglement entropy of three strips with widths h, ℓ and 2ℓ + h in the AdS-Vaidya metric (1.7). To do so, one should use the covariant proposal of the entanglement entropy which has been reviewed in appendix A.2. Note that in the present case we will have to deal with three hypersurfaces. For each of these hypersurfaces we denote the crossing point and the turning point by (ρ i c , ρ i t ) with i = 1, 2, 3.
The system has several scales and therefore as time passes one should look for different behaviors of the corresponding entropies in different scales. In the present case where two strips have the same width there are four scales given by ρ H , h, ℓ and 2ℓ + h. As a matter of fact, having noted that h ≪ ℓ, there are four main possibilities for the order of scales as follows which we will study them separately.

JHEP09(2015)165
Note that in all cases for v < 0 the step function in (1.7) is zero and the dual geometry is an AdS solution which is a static background. Therefore the mutual information of the vacuum state before the quench is given by equation (2.3). Note also that as one increases the width of strips there is an upper limit for the mutual information given by which is the absolute value of the finite part of the entanglement entropy for a strip with the width h.
In this case the width of the all entangling regions involved in the computation of the mutual information are much greater than the horizon radius and, consequently, the corresponding co-dimension two hypersurfaces can penetrate the horizon. To study the time scaling behavior of the mutual information, one may distinguish between five time intervals as stated below.

Early time: t ≪ ρ H
In this time interval, the co-dimension two hypersurfaces in the bulk associated with the entanglement entropies appeared in equation (2.2) cross the null shell almost at same point which is very close to the boundary, Therefore the holographic entanglement entropy for all regions may be well approximated by equation (A.33) Here and throughout this section we use a notion in which l 1 = h, l 2 = ℓ and l 3 = 2ℓ + h. Plugging these expressions into equation (2.2), one finds (3.5) One observes that the mutual information starts from its value in the vacuum, I vac , and remains fixed up to order of O(t 2d ) at the early times.

Steady behavior
The system reaches a local equilibrium at t ∼ ρ H when it has ceased production of thermodynamic entropy, though the entanglement entropy still increases. In this time interval all entanglement entropies appearing in (2.2) grow linearly with time (see (A.36)). In other

JHEP09(2015)165
words, one has wheref (ρ) and ρ m are defined in appendix A.2. The mutual information is then obtained from (2.2) as follows Since we are dealing with the large entangling regions, the turning points of all hypersurfaces are large, and therefore from equation (A.39) one can deduce ρ i m = ρ m = (2(d − 1)/(d − 2)) 1/d ρ H . As a result, the second term in the above equation vanishes leading to a constant mutual information in this time interval too. Thus starting from a static solution one gets almost constant mutual information all the way from t = 0 to t ∼ h 2 .

Linear growth:
Using (A.12) and (A.36) one can show that the entanglement entropy associated with the entangling region h will be saturated to its equilibrium value at t ∝ h 2 − c 2 ρ H + c 0 . Therefore in this time interval the entanglement entropy S(h) is given by (A.12) (3.8) On the other hand entanglement entropies associated with the entangling regions ℓ and 2ℓ + h are still increasing linearly with time. Thus from equation (A.36) one has Plugging these results into equation (2.2) one finds (3.10) which can be recast into the following form

JHEP09(2015)165
where v E is given by (A.40). Here we have used the fact that the entangling regions are large so that ρ im = ρ m with ρ m is given by (A. 39). Note that the above mutual information is positive as long as ρ H ≪ h 2 and h 2 ≪ t. It is also clear that the mutual information in this time interval is always bigger than I vac and grows linearly with time. It is also worth noting that to get a positive mutual information it was crucial to keep the subleading term Assuming to have a linear growth all the way up to t ∼ ℓ 2 − c 2 ρ H + c 0 where the entanglement entropy associated with the entangling region ℓ saturates to its equilibrium value, the mutual information reaches its maximum value during the thermalization process. More precisely setting one finds Here t (1) max is the time when the mutual information reaches its maximum value I max .

Linear decreasing: ℓ
the entanglement entropy S(ℓ) saturates to its equilibrium value. Therefore in this time interval both entanglement entropies S(ℓ) and S(h) have to be approximated by their equilibrium values as follows while the one associated with the entangling region 2ℓ + h still grows linearly with time (see (A.36)) Therefore, in this case the mutual information is which may be simplified as follows max , also note that I declines linearly with time and is positive for t < ℓ + h 2 .

Saturation
If one waits enough the entanglement entropy associated with the entangling region 2ℓ + h will also saturate to its equilibrium value at t ∼ ℓ + h 2 − c 2 ρ H + c 0 So that the mutual information will also saturates to its equilibrium value studied in the previous section. Of course generally it is not correct to plug just the equilibrium values of the corresponding entanglement entropies into equation (2.2) to find the mutual information. Indeed if one naively do that, in the present case, the resultant mutual information would become negative. Therefore the mutual information must reach its equilibrium value at a saturation time t (1) To find the saturation time and the equilibrium value of mutual information we note that at the end of the thermalization process the resultant background will be an AdS black brane. On the other hand as we have already mentioned in the previous section when both the width of strips ℓ and distance between them h are large compared to the radius of the horizon namely, ρ H ≪ ℓ and ρ H ≪ h, the mutual information is zero. Therefore in the present case one would expect that the mutual information becomes zero at the end of the thermalization process. Using this fact, one may estimate the saturation time as follows.
Indeed assuming the mutual information decreases all the way till it becomes zero, from (3.16), one should set so that the saturation time reads 19) which shows that the mutual information saturates long before ℓ + h 2 − c 2 ρ H + c 0 which would be the saturation time of the entanglement entropy of a strip with width 2ℓ + h. Indeed this result is consistent with the numerical results of [34][35][36][37].
Let us summarize the results of this subsection. We have found that for the case where ρ H ≪ h 2 the mutual information starts from its value in the vacuum and remains almost constant up to t ∼ h 2 , then it starts growing with time linearly. It reaches its maximum value at t (1) max after that it decreases linearly with time till it becomes zero at the saturation time which takes place approximately at t figure 2).

Second case:
In this case, similar to the previous subsection, one can study the behavior of the mutual information in five time intervals. We note, however, that since we have h 2 ≪ ρ H condition, the co-dimension two hypersurface corresponding to the entangling region h cannot probe the v < 0 region.  In this time interval the behavior of the entanglement entropies appearing in equation (2.2) at the early times are the same as that in the previous case. Thus in this time interval, the mutual information is essentially given by equation (3.5), which means that it remains constant for t ≪ h 2 .
3.2.2 Quadratic growth: h 2 ≪ t ≪ ρ H In this time interval the co-dimension two hypersurface corresponding to the entangling region h remains all the time in the region of v > 0 which is, indeed, a static AdS black brane. Therefore the corresponding entanglement entropy S(h) reaches its equilibrium value which, in the present case, is given by equation (A.6). On the other hand since we are still in the range of t ≪ ℓ 2 , the entanglement entropies associated with the entangling regions ℓ and 2ℓ + h are still at the early times so that they should be approximated by (A.33). Therefore one gets and Plugging these expressions into equation (2.2) one arrives at showing that the mutual information has a quadratic growth up to t ∼ ρ H .

JHEP09(2015)165
In this case the entanglement entropy S(h) is still given by equation (A.6), while since the system has reached a local equilibrium and moreover ρ H ≪ ℓ 2 , equation (A.36) should be used to describe the entanglement entropies associated with the entangling regions ℓ and 2ℓ + h, This leads to the following expression for mutual information Here, again, we have used the fact that the entangling regions are large so that ρ i m = ρ m . Moreover, in this time interval, the conditions h ≪ ρ H and ρ H ≪ t guarantee that the resultant mutual information will be positive and bigger than I vac .
The linear growth lasts all the way up to when S(ℓ) saturates to its equilibrium value. By making use of equation (3.12) one may also estimate the maximum value of the mutual information as follows In this time interval both the entanglement entropies S(h) and S(ℓ) are saturated to their equilibrium values, though because of their size of entangling regions, the corresponding equilibrium values are given by different expressions. Indeed, although the equilibrium value of S(h) is given by equation (A.6), for that of S(ℓ) one should use (A.12). The entanglement entropy S(2ℓ + h) is still growing with time as (A.36). Therefore, the mutual information in this time interval linearly decreases as time goes on and it is given by Note that the mutual information is positive and also I < I (2) max .

Saturation
As we have already mentioned, the final state of our system after a global quench we are considering is a thermal state whose gravity dual is provided by an AdS black brane. On Figure 3.

JHEP09(2015)165
Schematic behavior of mutual information during the thermalization process for h max , I max , I sat and t the other hand for this static background the mutual information of two strips depicted in figure 1 with the condition h 2 ≪ ρ H ≪ ℓ 2 is given by (2.9). Therefore in the present case the equilibrium value of the mutual information is It is then possible to estimate the saturation time by assuming that the mutual information decreases linearly with time till it reaches its equilibrium value (3.28). Indeed equating equations (3.27) and (3.28) one finds Moreover, from (3.28) it is obvious that which shows that in this case with the condition h 2 ≪ ρ H ≪ ℓ 2 the expression in the parentheses is always positive leading to the fact that I (2) sat < I vac . Let us summarize the results of this subsection. In fact we have found that the mutual information starts from its value in the vacuum and remains almost constant up to t ∼ h 2 , then it grows with time quadratically till t ∼ ρ H . After that a linear behavior starts and it reaches its maximum value at t (2) max . Finally it decreases linearly with time till it saturates to a constant value at the saturation time which takes place approximately at t (2) In this case the entanglement entropies associated with the entangling regions h and ℓ saturate to their equilibrium values before the system reaches a local equilibrium. Therefore the corresponding co-dimension two hypersurfaces cannot probe the region near and behind the horizon. In other words, the entanglement entropies S(h) and S(ℓ) do not exhibit linear growth, though S(2ℓ + h) could still grow linearly with time before it reaches its equilibrium value. Actually in this case the behavior of the mutual information for early times is almost the same as that in the previous case. Namely it starts from its value in the vacuum and remains fixed up to t ∼ h 2 then it begins to grow quadratically with time Let us assume that the entanglement entropy associated with the entangling region ℓ grows quadratically with time till it reaches its equilibrium value. Then using the fact that in this case the corresponding equilibrium value is given by (A.6), one may estimate the time when the mutual information becomes maximum as follows by which the maximum value of the mutual information reads Let us now study the other time intervals in more details.
3.3.1 Quadratic decreasing ℓ 2 < t < ρ H In this time interval, the entanglement entropies associated with the entangling regions h and ℓ are saturated to their equilibrium values given by (see (A.6)) On the other hand since we are in the regime of t < ρ H < ℓ + h 2 the entanglement entropy S(2ℓ + h) is still at the early times and should be approximated by equation (A.33) (3.35) Plugging these results into equation (2.2), one finds Since t > t max it is clear that I < I max .

JHEP09(2015)165
In this time interval, the entanglement entropies S(h) and S(ℓ) are the same as that in the previous case. On the other hand since in this time interval the system is locally equilibrated the entanglement entropy S(2ℓ + h) exhibits a linear growth. Therefore one has Plugging these results into equation (2.2), one finds Here we have used the large entangling region approximation for ρ 3 m . Note that in this time interval because of t < ℓ + h 2 , one obtains a positive value for mutual information.

Saturation
At this stage let us consider the situation that takes place after a long time when all the entanglement entropies appeared in (2.2) saturate to their equilibrium values. Since the entanglement entropies S(h) and S(ℓ) have already saturated (given by equation (A.6)), they do not change as the system evolves with time, though the one associated with entangling region 2ℓ + h will be saturated whose equilibrium value is given by (A.12). This would lead to the following mutual information which may be recast to the following form (3.41) From this expression it is clear that in the range we are interested in (i.e. h 2 ≪ ℓ 2 < ρ H < ℓ + h 2 ) the expression in the parentheses is always negative and therefore one gets I sat < I vac . 4 To estimate the saturation time, with the assumption that the mutual information JHEP09(2015)165 Figure 4. Schematic behavior of mutual information during the thermalization process for h sat and t To summarize the results of this subsection, one has observed that the mutual information starts from its value in the vacuum and remains almost constant up to t ∼ h 2 , then it starts growing with time quadratically till t  figure 4).

Fourth case:
In this case due to the fact that all entangling regions h, ℓ and 2ℓ + h are smaller than the radius of horizon, the corresponding entanglement entropies saturate to their equilibrium values before the system reaches a local equilibrium. Therefore neither the entanglement entropies nor the mutual information exhibit linear growth with time during the process of thermalization.
In fact to study the behavior of the entanglement entropies S(h), S(ℓ) and S(2ℓ+h) one should use either equation (A.6) or (A.33) depending on whether they have been saturated or not. Actually the situation is very similar to the third case studied above. Namely the mutual information starts from its value in vacuum and remains fixed up to t ∼ h 2 when it starts growing quadratically with time

JHEP09(2015)165
Assuming to have this quadratic growth up to its maximum value given by one may estimate a time when the mutual information becomes maximum Then it starts decreasing quadratically with time which is positive as long as t < ℓ+ h 2 . Finally the mutual information reaches its equilibrium value as system evolves with time. Indeed the saturation takes place, when the values of the entanglement entropies appearing in (2.2) become that of an AdS black brane given by (A.6). Therefore from equation (2.2), the saturated mutual information is obtained as Assuming to have the quadratic decreasing all the way to the saturation point, one may estimate the saturation time by equating equations (3.46) and (3.47) which leads to These behaviors are summarized in figure 5.

n-partite information for static backgrounds
In this section by making use of the AdS/CFT correspondence we will study n-partite information of a subsystem consists of n disjoint regions A i , i = 1, · · · , n in a d-dimensional CFT for the vacuum and thermal states whose gravity duals are provided by AdS and AdS black brane geometries, respectively. The n disjoint regions are given by n parallel infinite strips of equal width ℓ separated by n − 1 regions of width h, as depicted in figure 6. Following our discussions in the introduction we shall define the n-partite information as follows [8] max , I . n disjoint entangling regions A i , i = 1, · · · , n for computing n-partite information.
order to compute the holographic mutual information there are two possibilities to get minimal surface in the bulk associated to the entanglement entropy of the union S(A ∪ B).
In the present case where we are dealing with parallel strips with h ≪ ℓ, taking the minimal surface leads to Similarly for the union of three regions one uses

JHEP09(2015)165
and more generally for arbitrary integer numbers k, m and j > 1 one has By making use of these expressions, equation (4.1) evaluated for the system depicted in figure 6 may be simplified significantly as follows Interestingly enough, one observes that among various co-dimension two hypersurfaces only three of them corresponding to nℓ + (n − 1)h, (n − 1)ℓ + (n − 2)h and (n − 2)ℓ + (n − 3)h contribute to the n-partite information. 5 Therefore in order to compute the npartite information one should essentially redo the same computations we have done for the mutual information. In what follows using the AdS/CFT correspondence we will computẽ I [n] which as we will see it is always positive. In other words the holographic n-partite information has definite sign: for even n it is positive and for odd n it is negative. Let us consider the holographic n-partite information for the vacuum state of a CFT whose gravity dual is given by an AdS background. Indeed from equation (A.4) one finds (4.6) Using a numerical calculation one can show that for fixed h ℓ the above quantity for all values of d and n is positive and approaches zero in large ℓ limit.
For a thermal case whose gravity dual is provided by an AdS black brane geometry, and in the limit of ℓ ≪ ρ H , utilizing equation (A.6) one arrives at On the other hand, by making use of equation (2.7) one finds This relation shows that by increasing the temperature, the n-partite information is increased (decreased) for n even (odd).
On the other hand for the case of ρ H ≪ ℓ since all entangling regions appearing in the definition of n-partite information (4.5) contains a factor of ℓ, the corresponding JHEP09(2015)165 entanglement entropy should be approximated by equation (A.12). We note, however, that since the resultant n-partite information vanishes.

Time evolution of n-partite information
In this section we would like to study the scaling behavior of n-partite information during the process of thermalization. This is indeed a generalization of the mutual information studied in section 3. Again, the thermalization process we are considering is holographically modelled by the AdS-Vaidya metric (1.7). Therefore one, essentially, needs to study different scaling behaviors of three entanglement entropies appearing in the n-partite information (4.5) in the AdS-Vaidya metric. To do so, we will utilize the results reviewed in appendix A.2.
In general for the system we are considering there are four time scales given by the radius of the horizon ρ H and three entangling regions appearing in equation (4.5) which are (n − 2)ℓ + (n − 3)h, (n − 1)ℓ + (n − 2)h and nℓ + (n − 1)h. Since we are interested in h ≪ ℓ situation, one may recognize four possibilities for the order of these scales as follows which could be studied separately. We note, however, that since the situation is very similar to the mutual information, one would expect to get the same qualitative behaviors for the n-partite information. In what follows we will explore the first case listed above in more detail and will briefly present the results of other cases.

First case: 2ρ
In this case since all entangling regions involved in the computations of n-partite information are larger than the radius of horizon, the corresponding co-dimension two hypersurfaces in the bulk would have a chance to penetrate the horizon. Indeed in this case there are five time intervals in which the n-partite information behaves differently. We will study these intervals separately. It is worth noting that before the thermalization process, the system is in the vacuum state and therefore the corresponding n-partite information is given by equation (4.6).

Early time: t ≪ ρ H
At the early times all co-dimension two hypersurfaces cross the null shell almost at the same point that is very close to the boundary. Therefore all entanglement entropies appearing

JHEP09(2015)165
in the n-partite information (4.5) should be approximated by equation (A.33). Thus at the early times one findsĨ showing that the n-partite information starts from I [n] vac in the vacuum and remains fixed up to order of O(t 2d ).

Steady behavior
The system reaches a local equilibrium at t ∼ ρ H after which it does not produce thermal entropy, though the entanglement entropies associated with the entangling regions appearing in the n-partite information still increase with time. Actually since all the entangling regions are larger than the radius of horizon, the corresponding entanglement entropies grow linearly with time (see (A.36)). Therefore one finds where Note also that in order to find the above expression we have used the large entangling region limit by which ρ im = ρ m with ρ m is given by (A.39).
From the above equation it is clear that in this time intervalĨ [n] is bigger than its value in the vacuumĨ [n] vac , and grows linearly with time. Actually the linear growth continues until the entanglement entropy S((n − 1)ℓ + (n − 2)h) saturates to its equilibrium value at whichĨ [n] reaches its maximum at

JHEP09(2015)165
with the maximum value given bỹ (5.6) It is important to note that sinceĨ [n] is positive the actual sign of the n-partite information is given by the prefactor (−1) n in equation (4.5). Therefore for odd n the n-partite information has, actually, a minimum, though for even n it has a maximum. (5.7)

Saturation
In the present case where ρ H ≪ ℓ, as we have seen in the previous section all entanglement entropies saturate to their equilibrium values and the n-partite information becomes zero. Then assuming to have linear decreasing all the way to the saturation, one may estimate the saturation time by setting equation (5.7) to zero, It is worth noting that the n-partite information saturates before n 2 ℓ + n−1 which is essentially the time when entanglement entropy S(nℓ + (h − 1)h) saturates to its equilibrium value.
As a result, we found that in the case of ρ H ≪ l i , the quantityĨ [n] starts from its value in the vacuum and remains almost constant up to t ∼ n−2 2 ℓ + n−3 2 h, then it grows linearly with time till it reaches its maximum value at t One observes that the quantityĨ [n] has the same behavior as the mutual information, though scaling behaviors occur at different time scales. Note that the to find the actual value of the n-partite information, the factor of (−1) n should also be taken into account. Therefore although the behavior should be the same, the n-partite information is either negative (for odd n) or positive (for even n). To illustrate the situation we have summarized the JHEP09(2015)165 I [3] t I [3] vac I [3](1) min I [3](1) sat = 0 ℓ 2 t [3](1) min t [3](1) s Figure 7. Schematic behavior of tripartite information during the thermalization process for the case ρ H ≪ ℓ 2 . Here I [3] vac , I [3](1) min and I [3](1) sat , up to a minus sign, are given by equations (4.6), (5.6) and (5.8), respectively. results in figure 7 for tripartite information (Note that in this case because n is odd we have t min , I min instead of t max , I max ).

Other cases
Since for the model we are considering the n-partite information (or more precisely the quantityĨ [n] ) has the same structure as the mutual information (three entanglement entropies have to be computed), the behavior ofĨ [n] should be the same as that of the mutual information. Indeed we have explicitly shown this in the previous subsection for the case where all entangling regions are bigger than radius of the horizon. Having reached to this conclusion in what follows we just briefly present the results of other cases.

Second case:
In this caseĨ [n] starts from its value and remains constant at the early times till the first entanglement entropy associated with the entangling region (n − 2)ℓ + (n − 3)h saturates. Then one gets a quadratic growth as follows in the time interval (n−2) 2 ℓ + (n−3) 2 h < t < ρ H . When system reaches a local equilibrium one gets linear growth. Indeed for time interval ρ H < t < (n−1) 2 ℓ + (n−2) 2 h one has linear growth, while for (n−1) 2 ℓ + (n−2) 2 h < t < n 2 ℓ + (n−1) 2 h one gets linear decreasing with time.

JHEP09(2015)165
I [3] t I [3] vac I [3](2) min I [3](2) sat ℓ 2 ρ H t [3](2) min t [3](2) s Figure 8. Schematic behavior of tripartite information during the thermalization process for vac , I [3](2) min and I [3] (2) sat , up to a minus sign, are given by equations (4.6), (5.11) and (5.13), respectively. Therefore it has a maximum value given bỹ to its equilibrium value given bỹ showing thatĨ  2 ℓ + (n−2) 2 h < ρ H < n 2 ℓ + (n−1) 2 h In this case the situation is exactly the same as the third case of the mutual information. Namely the quantityĨ [n] starts from its value in vacuum and remains fixed at the early times. Then it grows quadratically with time till reaches a maximum after that it decreases quadratically and then linearly with time up to the saturation point. The maximum and

Fourth case:
In this case which all entangling regions involving in the computation of the n-partite information (4.5) are smaller than the radius of the horizon, the corresponding entanglement entropies saturate to their equilibrium values before the system reaches a local equilibrium. Therefore during the process of thermalization the n-partite information does not exhibit linear growth. IndeedĨ The situation is illustrated in figure 10 for tripartite information.

Numerical results
So far following [38,39] we have analytically studied the behavior of n-partite information in a process of thermalization with certain assumptions. In order to examine our results in this section we will study the behavior of n-partite information numerically. In particular we will mainly focus on the mutual information and 3-partite information in more details and then briefly study 4-partite and 5-partite information.
It is worth mentioning that although such a numerical analysis has been already done in e.g. [34,35], in what follows our main interest is to explore various scaling regimes we have obtained in the previous sections. This could be used to examine the validity of our assumptions, approximations and results. To be concrete we will consider a 2+1 dimensional boundary theory, however, the result could be extended to higher dimensions. In this case the area of the extremal surface using (A.14) and setting L = 1 reads

JHEP09(2015)165
The minimization condition leads to the following equations of motion which should be solved with the following boundary conditions In order to study the equations numerically one should approximate the theta function appearing in f with a smooth analytic function. Actually in what follows we will consider the following function (see for example [15]) where m 0 is a measure of the horizon radius for the final static black-brane geometry, i.e. ρ H = m −1/2 0 and a is the parameter that controls the thickness of the null shell. In the limit of a → 0 this profile coincides with the step function, as illustrated in figure 11.
To find the profile of the extremal surface numerically one should solve equations (6.  Having found the extremal surface numerically one could plug the profile of the extremal surface into (6.1) to read the area of the extremal surface as a function of boundary time. It is, however, important to note that due to the large volume limit, the area (6.1) is divergent and needs to be regularized by introducing a UV cut-off at ρ = ǫ. In this case the finite part of the area is given by Here we have used the conservation law, ρ 2 L = ρ 2 t , to simplify the final expression. Evolution of the area of extremal surface for a strip in a thin shell limit for the large (ℓ > ρ H ) and small (ℓ < ρ H ) entangling regions is depicted in figure 13. 6 Actually in this figure we have plotted ∆A defined by Note that the actual value of the entanglement entropy has an extra factor of (4G N ) −1 in front of the area, though in this note we will neglect this factor.

Mutual information
In this section by making use of the numerical results of the holographic entanglement entropy we will numerically explore different scaling regimes of the holographic mutual information. To do so a non-trivial task is how to compute the entanglement entropy of a union of two subsystems. As we have already mentioned there are two minimal surfaces JHEP09(2015)165 S dis. : S con. : Figure 14. Two different configurations for computing S(A 1 ∪ A 2 ). associated with the entanglement entropy S(A 1 ∪ A 2 ) (see figure 14) indicating that there is a transition between connected and disconnected configurations as one increases h ℓ . Note also that the value of the disconnected configuration is independent of h.
Indeed it is easy to find the transition point between these two configurations which can be done by solving the S con. (h, ℓ) = S dis. (ℓ) for h. In particular for a four dimensional AdS background in which the corresponding expression for the entanglement entropy is given by equation (A.4) for d = 3 one finds that the transition occurs at h = 1 2 √ 5 − 1 ℓ ∼ 0.618 ℓ. This means that for the vacuum state and for h < 0.618 ℓ one must consider the connected configuration where the resulting mutual information will be a finite positive number, though for h > 0.618 ℓ, the disconnected configuration is favored and the resulting mutual information is zero. Since the results we have presented in the previous sections depend crucially on the assumption of whether the connected or disconnected configurations are favored in what follows for all cases we will compute the evolution of S(A ∪ B) too.
To proceed with the numerical computations we will set ρ H = 1 and a = 0.001. Having collected all information and the procedure of our numerical method, let us present our numerical results for the mutual information for all scaling cases we have considered in section 3.

First case
For this case we will fix the width of strips to be ℓ = 4.5 which is larger than the radius of horizon ρ H = 1. Then we will consider different values for h = 2.1, 2.2, 2.4, 2, 6. Note that for these cases we have the condition h < 0.618 ℓ and therefore the mutual information for the vacuum state (AdS geometry) is non-zero. The numerical results are given in figure 15. As one observes the numerical results are in good agreement with the analytical results we have obtained in section 3. In particular the left plot in this figure shows that saturation of the mutual information (which takes place at the crossing point of the dashed curve with others) happens long before the saturation of the HEEs. Also according to (3.19) as we increase the separation between the strips the saturation time decreases.
On the other hand for h = 3 one has h > 0.618 ℓ so that the mutual information in the vacuum state is zero. The corresponding behavior is shown in figure 16.

Second case
To study the second case we set ℓ = 3 and consider h = 0.2, 0.3, 0.4, 0.5 to make sure that the condition h 2 < ρ H < ℓ 2 is satisfied. The numerical results for this case are presented in figure 17. It is worth to mention that in this case, the numerical results indicate that the saturation value is independent of ℓ which is in agreement with that corresponding analytical results (see (3.28)).

Third case
To examine the third case we set ℓ = 1.18 and let h = 0.4, 0.42, · · · , 0.48. For these values of ℓ and h we have plotted the numerical results for the connected and disconnected configurations as well as the evolution of the mutual information in figure 18. Again the results are in a good agreement with that discussed in section 3.3.

Fourth case
For this case we consider ℓ = 0.45 and h = 0.23, 0.232, . . . , 0.24 and the corresponding numerical results are depicted in figure 19. These plots should be compared with figure 5 in subsection 3.4.

3-partite information
To further examine our analytical results, in this subsection, we will numerically study the behavior of 3-partite information during the process of thermalization. The corresponding system consists of three parallel strips with width ℓ separated by distances h as drawing in figure 20.
The corresponding 3-partite information is given by  As we have already mentioned in section 5 the main subtlety in evaluating the 3-partite is the way we compute the entanglement entropy of union of subsystems. In order to compute these quantities let us review the assumptions which led to a simple expression given in equation (4.5) for h ≪ ℓ.
Actually for the union of two subsystem one may consider different configurations for the extremal surfaces as depicted in figure 21. 7 S (2) dis. :

JHEP09(2015)165
dis. : More precisely one has h ≫ ℓ, i = 1 or 3 (6.9) and (6.10) Similarly one may also study the union of three subsystems where we could have different configurations for the extremal surface as given in figure 22.
Mathematically these configurations can be translated into the following expressions for the union of three strips A 1 , A 2 and A 3 h ≫ ℓ, S(3ℓ + 2h) + S(ℓ + 2h) + S(ℓ) ≡ S dis. : Putting these results into (6.8) and in the limit of h ≪ ℓ one arrives at which is the same as (4.5) for n = 3. Indeed to get the above expression for the 3-partite we have assumed that dis. , and S con. < min S (4) dis. , S dis. , S dis. , h ≪ ℓ. (6.13) Therefore to proceed with evaluating the 3-partite information, it is crucial to see in what extend our assumptions is reliable. Of course, in general, it is not an easy task to prove the above inequalities during the thermalization process even numerically. Nevertheless we have provided some numerical examples in figures 23 and 24 showing that in the desired range of parameters these inequalities are indeed hold.
Having explored the subtlety we are encountering when we are going to compute 3partite information, in the rest of this subsection we numerically study different scaling behaviors of 3-partite information during the process of thermalization. Actually for all cases which we would like study, one should first, numerically, check whether the conditions (6.13) are satisfied. If the conditions were satisfied, then one can evaluate 3-partite information using the expression (6.12). Indeed we have done these considerations and the results are as follows.

First case
In this case to meet the condition 2ρ H < ℓ < 2ℓ + h, we will set h = 0.2 and will consider different values for the width of strips as ℓ = 2, 2.2, · · · , 3.6. For these values the behavior of the 3-partite information is shown in figure 25. It is worth mentioning that for each case we have numerically checked that the conditions (6.13) are, indeed, satisfied.

JHEP09(2015)165
dis. and S (3) dis. from bottom to top. Right plot: S con. , S (4) dis. , S (5) dis. and S (6) dis. from bottom to top. These plots show that in this range of the parameters the conditions (6.13) satisfied.
dis. and S (3) dis. from bottom to top. Right plot: S con. , S (4) dis. , S (5) dis. and S (6) dis. from bottom to top. These plots show that in this range of the parameters the conditions (6.13) satisfied. From this figure it is clear that the saturation value of the 3-partite information tends to zero as one increases the size of the entangling regions, which this is in agreement with our analytic results in section 5. Also one can check that the saturation time increases when we increase the strip width as we expect from equation (5.9).

Second case
To study this case we will fix h = 0.1 and consider different values for the width of strips as ℓ = 0.96, 1, · · · , 1.16. It is clear that for these values one has ℓ < 2ρ H < 2ℓ + h. For these values the behavior of 3-partite is depicted in figure 26. Note for all results appearing in this figure we have numerically checked the conditions (6.13) too.

Third case
In this case to maintain the condition 2ℓ + h < 2ρ H we will set h = 0.1 and will study the behavior of the 3-partite information for different values of the width of strips given by ℓ = 0.4, 0.42, · · · , 0.52. The corresponding results are given in figure 27. Again, the conditions (6.13) are satisfied for our results.

Fourth case
The last case corresponds to the situation where all entangling regions appearing in the expression of 3-partite information are smaller than the radius of horizon. More precisely one has 3ℓ+2h < 2ρ H . In order to numerically study 3-partite information in this region we will set h = 0.1 and will evaluate the 3-partite for different widths ℓ = 0.24, 0.25, · · · , 0.32. The results are given in figure 28. For this case the conditions (6.13) are also numerically checked.

4-partite information
To extend the presented numerical computations beyond that which has been already considered in the literature in this subsection we will consider 4-partite information. Following our previous discussions the main point is to explore different configurations we may find for computing the entanglement entropy of union of different subsystems.
Indeed in the present case the different independent and nonintersecting configurations one should study are those appearing in S(   present case we just need to study 19 independent configurations one of which is connected and the others are disconnected, as shown in figure 29. Now in order to prove equation (4.5) for n = 4 and in the limit of h ≪ ℓ the following conditions must be satisfied dis. , and S (5) dis. < S (6) dis. < min S (7) dis. , S (8) dis. , S (9) dis. , S (10) dis. , and S con. < min S (11) dis. , · · · , S (18) dis. . (6.14) Actually we have numerically checked that these conditions are satisfied within the range of parameters we are interested in. For example figure 30 shows numerical results for certain values of h and ℓ.
Having check the validity of our assumptions on the minimal configurations one may proceed to compute the 4-partite information using equation (4.5). Of course as we have already mentioned, due to the relative values of the horizon radius and the width of the entangling regions, one recognizes four different cases. Since we have explored these possibilities in the cases of mutual information and 3-partite information (see also general argument in section 5) here we just present final numerical results in figure 31.

5-partite information
Similarly one could also work out the case of 5-partite information. Of course in this case one must compare different configurations corresponding to S( It is easy to see that in the present case one finds 46 independent nonintersecting configurations with different areas. One of them is connected and the others are disconnected. Indeed the situation is very similar to what we have done in the previous cases. Therefore we just shown the behavior of 5-partite information during the process of thermalization for the first and second cases (according to the notation in section 5) in figure 32.

Numerics vs. analytic expansions
In the previous subsections we have numerically studied n-partite information (for n = 2, · · · , 5) where we have found that the corresponding behaviors qualitatively are in agreement with the analytical results of sections three and five. To make the comparison clearer, in this subsection, we shall compare the actual numerical values which have been found for different quantities at distinguished points with those predicted analytically in sections 3 and 5. This might help us to better understanding of limitation and range of the validity of our study. Table 1 presents the numerical results for the special points of the mutual information and the corresponding saturation times for different cases following our notation in section 3. These numerical results are actually obtained by fixing the width ℓ while varying h within its allowed range (see different cases in section 3). The corresponding results for 3-partite information and those of 4 and 5-partite are, rather briefly, presented in the table 2 and 3, respectively. Note that for these cases we have set ρ H = 1 and a = 0.001.
Although the results given in tables 1-3 are enough to explore the level of agreement and range of validity of two different approaches which have been considered in this paper (numerical and semi-analytical), it is useful to visually present the results in different plots. To proceed it helps if one first studies the behavior of holographic entanglement entropy. Actually figure 33 shows the evolution of HEE for different values of ℓ, e.i. ℓ   11.2. In this figure the black circles show the numerical results while the dashed curves represent the results of our semi-analytic study. As we have already mentioned the semianalytic method gives just a piece wise plot which should be compared with certain regions of the numerical computations. From this figure one also observes that by increasing the width of the entangling region, the analytic expansions become more precise and smooth, as expected.
According to our semi-analytical expansions we would expect that the transition between quadratic growth and linear growth occurs at t trans. ∼ O(ρ H ) = 1, though our numerical results give actual value for the saturation times. For example for ℓ = 4.5 one has t trans. ∼ 1.15 while for ℓ = 11.2 one gets t trans. ∼ 1.5. It is worth noting that in both cases the saturation time obtained by the semi-analytic expansion is generally larger than the numerical ones, though by increasing ℓ they would converge at the same time. Since the n-partite information may be given in terms of entanglement entropy, such a difference would also affect the saturation times of the n-partite information.

Ana
Let us now consider mutual and n-partite information. To be specific in figure 34 we have depicted the results for mutual information in different cases for particular values of parameters. In these figures the numerical results are compared with the semi-analytical results. Note that for the latter approach the corresponding curves are plotted in different colors indicating different scaling behaviors we have found for the mutual information. In other words these curves have been plotted patch wise and the resulting function may not be a smooth function for whole time during the process of thermalization. Similarly one could graphically compare other cases. For example in the case of higher n, we have depicted the results for I [3](2) , I [3](3) , I [4](1) and I [5](1) in figure 35.  This may be understood as follows.
As far as the values of n-partite information at different distinguished points are concerned the small mismatch is related to our assumptions on the relative size of the entan-JHEP09(2015)165 Figure 35. Comparing analytical and numerical results for I [3](2) for h = 0.1 and ℓ = 1.12 (left up), I [3](3) for h = 0.1 and ℓ = 0.42 (right up), I [4](1) for h = 0.2 and ℓ = 2.8 (left down) and I [5](1) for h = 0.2 and ℓ = 2 (right down). In these plots the black circles show the numerical results, the dashed colored curves corresponds to different scaling behaviors we have found in our analytical studies.
gling widths, their separations and the horizon radius. Indeed to work out our analytical results we have assumed a stricken inequality such as ρ H ≪ h 2 ≪ ℓ 2 , though in our numerical computations our parameters satisfy ρ H < h 2 < ℓ 2 . As a result although in the analytical considerations we could drop all higher order corrections, in the numerical one, their effects are also taken into account. As an explicit example, in figure 36 we plot two parameters that appear in the analytic calculation for linear growth regime, i.e. ρ m and t s and compare them with the numerical data. In this figure the solid circles show the numerical results. The dashed blue line in the left plot shows the value of ρ m approximated by (A. 39). Note that in the analytical approximations that we have used, we always assume that in the large entangling region, ρ m is constant and according to (A.39) does not depend on ℓ. This plot shows that this assumption is more concrete when one considers ℓ > 6 (Note that we always consider ρ H = 1). The dashed red line in the right plot show the value of t s approximated by (A.41), where v E is given by (A. 40). Note that in this plot the slope of the curve is given by the entanglement velocity. In both plots in the large entangling region limit analytic expressions converge to numerical results. Actually as one increases the width of entangling regions our code becomes more unstable. Therefore we have some restrictions when we are considering the large entangling region limit.
On the other hand in order to compute the time scale of the distinguished points we have assumed a particular behavior for holographic entanglement entropy as it approaches JHEP09(2015)165  its saturation point. To be precise let us recall that in order to compute the n-partite information, with the assumption we made, one has to compute three entanglement entropies associated with the entangling regions ℓ 1 = nℓ + (n − 1)h, ℓ 2 = (n − 1)ℓ + (n − 2)h and ℓ 3 = (n − 2)ℓ + (n − 3)h. We note that although the distinguished points occur when one of these entanglement entropies saturates to its equilibrium value, there is a subtlety to compute the corresponding saturation time when the entangling region is in a shape of strip [39].
In order to further compare our numerical results with that of semi-analytic, in what follow we will consider another method for the comparison to explore the regime of validity of our analytic expansions. Indeed using "FindFit" command in Mathematica and assuming certain fit functions, we compare two approaches in specific examples. In figure 37 we present the evolution of HEE for ℓ = 2.2, 4.5 and 11.2 numerically together with certain piecewise fit functions given by  will consider the following fit functions (6.16) Table 5 shows the parameters {a 1 , b 1 , a 2 , b 2 } for both approaches (note that for analytic expansions we use eqs. (3.11) and (3.17)).
A similar analysis also works for I (2) (see figure 39). In this case using the following fit functions I (2) quad. = c 1 + d 1 t 2 , I lin.grow. = a 3 + b 3 t, I  Table 5. Comparing specific values of the parameters {a 1 , b 1 , a 2 , b 2 } for ℓ = 4.5.   Table 6. Comparing specific values of the parameters {c 1 , d 1 , a 3 , b 3 , a 4 , b 4 } for ℓ = 3. and utilizing eqs. (3.22), (3.24) and (3.27), one arrives at the results presented in table 6. One may go head to study I [n](i) , though the conclusions would be the same.
To conclude this section we observe that there is rather a good agreement between numerical and semi-analytical results. We note however that due to the limitation of the numerical computation as well as the semi-analytical approximations, the actual values of distinctive points may not be precisely the same.

Discussions
In this paper using the covariant prescription for computing the holographic entanglement entropy we have studied mutual information and n-partite information (defined by equation (1.2)) for a strongly coupled field theory whose gravitational description is provided by an AdS-Vaidya metric. We have computed the n-partite information for a system consisting of n parallel strips (two for mutual information) with the same width ℓ separated by distances h with the condition h ≪ ℓ. With this assumption the expression of n-partite information is simplified so that in order to study its behavior, one essentially needs to study entanglement entropy of three strips with different widths. Therefore it is possible to ex-

JHEP09(2015)165
plore the evolution of the n-partite information during the process of thermalization after a global quantum quench, by making use of the results for the entanglement entropy [38,39].
Of course time evolution of the n-partite information is sensitive to the size of three entangling regions appearing in the computation of n-partite information. Moreover the model has a distinctive time scale given by the horizon ρ H in which the theory reaches a local equilibrium. Then the behavior depends on relative size of the corresponding entangling regions and the radius of horizon: they could be larger or smaller than ρ H . Therefore in the intermediate region the n-pratite information could increase (decrease) linearly or quadratically with time.
An interesting observation we have made is that the holographic n-partite information has definite sign: it is positive for even n and negative for odd n, though for a generic field theory it could have either signs. Therefore following [8] one may suspect that having definite sign for the n-partite information is, indeed, an intrinsic property of a field theory which has gravity dual.
We also examined our analytical study by numerical computations. Actually for mutual information and 3-partite information some numerical computations have been performed in, e.g. [34,35]. Here, the corresponding numerical computations were studied in more details in order to explore different scaling behaviors of mutual information and 3-partite information. We have also considered 4 and 5-partite information, the generalization to higher n is straightforward. It should be mentioned that in the range of the parameters which we are interested in, the numerical computations confirm our results. We have also numerically checked our assumptions under which the expression of n-partite information simplified drastically. We also compare these two approaches for computing time evolution of n-partite information. The numerical results, patch wise, are best fitted with, quadratic, liner and constant curves, confirming the overall picture of the time dependent behavior of n-partite information.
Moreover we have studied the n-partite information (1.2) for a system consisting of n parallel strips with the same width separated by distances h. It is however, instructive to explore the results for the case where the system is not symmetric. In other words, one may consider n strips A i with width ℓ i , i = 1, · · · , n separated by h j , j = 1, · · · , n−1. Therefore the strips could have any size and are separated by arbitrary distances. Nevertheless, inspired by holographic mutual information, if for arbitrary numbers i, m, k and j > 1 one assumes then the n-partite information (1.2) may be recast into the following form 8

JHEP09(2015)165
On the other hand by making use of the holographic description of entanglement entropy and setting L = ℓ 2 + · · · + ℓ n−1 + h 2 + · · · + h n−2 one finds I [n] = (−1) n−1 S(L)−S(L+ℓ 1 +h 1 )−S(L+ℓ n +h n−1 )+S(L+ℓ 1 +ℓ 2 +h 1 +h n−1 ) , (7.3) where, as before, S(l) is the holographic entanglement entropy of a strip with the width l. It is then easy to follow our discussions to find the behavior of n-partite information during a process of thermalization. Note that in this case we have five different cases depending on whether the radius of horizon ρ H is larger or smaller than the width of strips appearing in (7.3). Indeed the general rule is as follows. If the width of the entangling region appearing in the expression of n-partite information is smaller than the radius of horizon, the corresponding entanglement entropy grows quadratically with time and saturates before the system reaches a local equilibrium where E is the energy density. Note that since the system has not reached a local equilibrium, the energy density is a proper quantity one may define. On the other hand if the width of the entangling region is larger than the radius of horizon, the corresponding entanglement entropy grows quadratically with time before the system reaches a local equilibrium, while it has linear growth after local equilibrium and then it saturates After local equ.S ∼ S vac + V d−1 S th t, Note that when the system is locally equilibrated the entanglement entropy may be given in terms of the thermal entropy. 9 To conclude we note that for a system of n parallel strips with different widths and distances between them, one would still get the same behavior though the corresponding behavior is less symmetric around the maximum or minimum points. As an explicit example we have numerically computed 3-partite information for h 1 = 0.1, h 2 = 0.5, ℓ 1 = 1, ℓ 2 = 2 and several values for ℓ 3 . The results are depicted in figure 40.
As we have already mentioned the holographic mutual information undergoes a first order phase transition as one increases the distance between two [31]. It is then natural to see whether such a transition would also occur for the n-partite information (1.2). Actually in this case it can be seen that the situation is very similar to that of mutual information. In other words the n-partite information vanishes as on increases the distance between the strips. More precisely, if one changes the distance between given two consecutive strips A i and A i+1 in the system such that S(A k , A k+1 ) = S(A k ) + S(A k+1 ), the 9 Note that in the above schematic expressions we have dropped the numerical factors. n-partite information (1.2) vanishes. More precisely using this identity and for 1 ≤ k < n, equation (7.2) reads
To conclude we have seen that if the mutual information of two consecutive strips of a system consisting of n parallel strip vanishes the n-partite information defined by (1.2) vanishes too. Therefore there could be a phase transition in n-partite information if one increases the distance h i . Since the n-partite information vanishes when the mutual information vanishes, the critical distance should be the same for both of them. More precisely for a given subsystem of n parallel strips specified by ℓ i , h i and ℓ i+1 , there is a critical h c i over which both mutual and n-partite information vanish. It is important to emphasis that this behavior is due to the facts that we are working with a field theory which has a holographic description and moreover the n-partite information is defined by equation (1.2). In general field theory it might not be true and moreover, as we will see, this is also not true if one uses another definition for n-partite such as that defined in equation (1.5). Actually our numerical results concerning the contributions of different hypersurfaces in evaluating 3-partite information, indeed, supports the existence of the above phase transition. Definitely the phase transition of the n-partite information deserved more investigations. We hope to further study this phase transition in near future.
In this paper we have only considered n-partite information based on the definition (1.2), though one could also study the behavior of multi-partite entanglement defined by equation (1.5). Indeed in this case for the system we have been considering (n strips with width ℓ separated by h with the condition ℓ ≫ h), equation ( It is then easy to show that vac . It is worth nothing that the above expression is the same as that of mutual information (2.8) up to a factor of n(n − 1). Actually the behavior of the above quantity in a process of thermalization is very similar to that of mutual information. In particular one can show that it is always positive during the process of thermalization.
It is also easy to see that unlike the previous case, in the present case when the mutual information of two consecutive strips becomes zero, the multi-partite information does not vanish and instead it breaks into two multi-partite informations. More precisely suppose I [2] (A k , A k+1 ) = 0, then from the definition of multi-partite (1.5) one gets (7.8)

JHEP09(2015)165
To proceed let us consider a strip with the width ℓ in a d-dimensional space time as follows where (t, x) are the space time coordinates.

A.1 Static background
Let us first compute the entanglement entropy for a d-dimensional static system whose gravitational dual is provided by the metric of (2.4). Following the holographic description of the entanglement entropy [12,13] one needs to minimize the area of a co-dimension two hypersurface whose boundary coincides with the boundary of the above strip. The profile of corresponding hypersurface in the bulk may be parametrized by x 1 = x(ρ), and therefore the area functional is where "prime" represents derivative with respect to ρ. It is then straightforward to minimize the above area to arrive at . For an excited state whose gravitational dual is provided by the black brane solution (2.4) the corresponding entanglement entropy may be found by minimizing the area when f = 1. In this case, in general, it is not possible to find an explicit expression for the entanglement entropy, though in certain limits one may extract the general behavior of the entanglement entropy. In particular in the limit of ml d ≪ 1, one finds which leads to the following expression for the entanglement entropy

JHEP09(2015)165
where S vac is the entanglement entropy of the vacuum solution given in (A.4), and .

(A.7)
On the other hand for mℓ d ≫ 1 the main contributions to the entanglement entropy comes from the limit where the minimal surface is extended all the way to the horizon so that ρ t ∼ ρ H . In this limit equation (A.3) for d > 2 reads . (A.8) Note that apart from the UV divergent term in S BH , due to the double zero in the square roots, the main contributions in the above integrals come from ξ = 1 point. Indeed around ξ = 1 the entanglement entropy S BH may be recast to the following form It is now clear that the first term in the above equation is divergent at ξ = 1 while the second one is finite, though the second term is UV divergent. Indeed the first term is exactly the one appears for ℓ. Therefore one has Now the aim is to compute the the second integral. Of course it can not been performed analytically, though one may solve it numerically to find its finite part. Indeed using "NIntegrate" command in the Mathematica one finds where c 2 is a positive number. For example for d = 3, 4 one gets c 2 = 0.88, 0.33, respectively. Therefore one arrives at [22] (A.12) Note that the first finite term in the above expression is proportional to the volume which is indeed the thermal entropy, while the second finite term is proportional to the area of the entangling region. Indeed this term plays a crucial role in our study.

A.2 Time dependent background
In this subsection we will review computations of the holographic entanglement entropy in the AdS-Vaidya background (1.7) for the case where the size of the entangling region is JHEP09(2015)165 larger than the radius of horizon [38,39]. More precisely the entangling region is given by a strip given in equation (A.1) with ℓ ≫ ρ H . As mentioned before for this time dependent background the covariant proposal for the holographic entanglement entropy is needed and the v(x) and ρ(x) may be used to parametrize the corresponding co-dimension two hypersurface in the bulk. Then the induced metric on the hypersurface is where "prime" represents derivative with respect to x. Therefore, the hypersurface's area can be obtained as (A.14) the corresponding entanglement entropy can then be found after evaluating (A.14) at the extremal surface as follows Note that (A.14) may be thought of as a one dimensional action for a dynamical system for the fields v(x) and ρ(x). Since the action is independent of x its corresponding Hamiltonian is a constant of motion This conservation law helps one to write equations of motion for v and ρ which read as where P 's are the momenta conjugate to v and ρ up to a factor of H −1 and are defined by These equations have to be solved by the following boundary conditions Note that with this boundary condition one obtains H = ρ d−1 t , where (ρ t , v t ) stands for the turning point coordinate of the extremal hypersurface in the bulk.
One should solve equations to find the extremal surface and the numerical method is actually needed, however, analytic solutions can still be found for some particular forms of m(v). It is known that in a quench there is a rapid change in the theory so that, one may assume that f (ρ, v) = 1 − θ(v)g(ρ), where θ(v) is the step function. This implies that f does not depend on v in most of time and hence ∂f (ρ,v) In the present case one has g(ρ) = ( ρ ρ H ) d where the horizon locates at ρ H with m = 1 ρ d H .
For v < 0 one has f = 1 and therefore the geometry is actually an AdS geometry which corresponds to the vacuum. In this case one gets which can also be used to find where (ρ t , v t ) is the extremal hypersurface turning point in the bulk and the crossing point where the hypersurface intersects the null shell is given by (ρ c , v c ). Since one is injecting the matter in v direction, one would expect that its corresponding momentum conjugate jumps once one moves from the initial phase to the final phase. While the momentum conjugate of ρ must be continuous. Therefore one gets v ′ (f ) = v ′ (i) . On the other hand by integrating equations of motion across the null shell one arrives at

JHEP09(2015)165
It is, then, easy to read the momentum conjugate of v in the final phase Now we have all ingredients to find the area of the corresponding extremal hypersurface in the bulk. In general the hypersurface could extend in both v < 0 and v > 0 regions of space-time. Therefore the width ℓ and the boundary time are found as follows where E = P (f )v which in the large entangling region limit becomes Finally, the area reads