# Note on the Kaplan–Yorke dimension and linear transport coefficients

## 1. INTRODUCTION

The Kaplan–Yorke (KY) or Lyapunov dimension was introduced$`^{\text{[1]}}`$ as a conjecture relating the Hausdorff (H) dimension and the Lyapunov exponents of the invariant measure of a given dynamical system.<sup>5</sup><sup>5</sup>5In the original paper only the natural invariant measure, i.e. the weak limit of the Lebesgue measure or Sinai-Ruelle-Bowen measure, was considered. In  the result is formulated in a more general way such that it holds for every invariant measure. This allows a computation of the H–dimension of the attractor of a dissipative dynamical system in terms of its Lyapunov exponents. Its validity has been proven for two dimensional dynamical systems$`^{\text{[2]}}`$ and for a rather large class of stochastic systems.$`^{\text{[3]}}`$ In general one only knows that the KY–dimension is an upper bound for the H–dimension.$`^{\text{[4]}}`$ Although it is not hard to construct rather artificial counter examples, it is generally believed that the conjecture holds for “generic” dynamical systems.$`^{\text{[5]}}`$

In this paper we will discuss a relation (see Eq. (4.10) below) between the KY–dimension in large thermostatted systems and their physical properties such as the transport coefficients in a nonequilibrium stationary state. Such a relation allows us to estimate the KY–dimension from a measurement of the transport coefficient. In doing so we obtain a new relation between a dynamical quantity (the KY–dimension) and a physical quantity (the transport coefficient).

A difficulty in doing this is that while the dynamical quantities are usually defined for any finite number of particles $`N`$, the physical quantities usually refer to systems of very large $`N`$, so that one can meaningfully define intensive quantities, which only depend on intrinsic parameters like the number density $`n=N/V`$ rather than on $`N`$ and $`V`$ separately, where $`V`$ is the volume of the system. We would like to state this by saying that strictly speaking, a thermodynamic limit has to be taken, i.e. $`N\mathrm{}`$, $`V\mathrm{}`$ with the number density ($`N/Vn`$) and other intensive quantities, such as in particular the shear rate, considered in this paper, held constant. This is straightforward if the linear transport coefficients are required and the limit of the external field $`F_e0`$ is taken before the limit $`N\mathrm{}`$, as in linear response theory. However, for finite fields, changes in the behavior of the system can occur, when the thermodynamic limit is approached - as e.g. the onset of turbulence in a sheared system.$`^{\text{[6]}}`$ In this paper we are interested only in the behavior of systems before such a transition takes place, like the laminar flow of a sheared fluid considered in section 5. We think, nevertheless, that our results can be usefully formulated for large systems using expressions like “for sufficiently large $`N`$” without taking the mathematical limit (see comment after Eq. (2.4) for a more precise discussion). Although this expression is not mathematically well defined, we think that its meaning will be clear in any practical application (see note 7 on page 14 for an attempt to clarify this point). Moreover, because errors in intensive quantities are typically $`O(N^1)`$ where $`N10^{23}`$, the approach is physically reasonable.

The above mentioned connection between a dynamical and thermodynamic treatment requires the usually discrete Lyapunov spectrum to be effectively replaced by a continuous intensive spectrum and an intrinsic version of the KY–dimension to be introduced. In section 4 we show how this can be implemented, after having introduced the basic equations which connect the dynamical and physical quantities in section 2, and deriving a new exact relation for the linear transport coefficients in section 3.

## 2. BASIC RELATIONS AND DEFINITIONS

As has been shown before,$`^{\text{[7, 8]}}`$ there is a direct relationship between the sum of the non-zero Lyapunov exponents $`\lambda _i`$ with $`\lambda _i\lambda _{i+1}`$, $`1i2dNf`$ (where $`N`$ is the number of particles, $`d`$ the Cartesian dimension and $`f`$ the number of zero Lyapunov exponents) and the phase space contraction rate in a thermostatted system, subject to an external force $`F_e`$, in a non equilibrium stationary state, of the form:

$$\underset{i=1}{\overset{2dNf}{}}\lambda _{i,N}(F_e)=\mathrm{\Lambda }_N(F_e).$$
(2.1)   
Here the subscript $`N`$ indicates the $`N`$-dependence of the various quantities in Eq. (2.1), $`2dNf`$ is the effective number of degress of freedom in phase space of the system and $`\mathrm{\Lambda }_N(F_e)=\frac{}{𝚪}\dot{𝚪}`$ is the phase space contraction rate, where $`𝚪`$ stands for the collection of the coordinates and momenta of the $`N`$ particles and $`\dot{𝚪}`$ for its time derivative. For macroscopic systems, i.e. systems with very large $`N`$, one can use the equality of the dynamical phase space contraction and the physical entropy production$`^{\text{[8, 9]}}`$ and obtain from Eq. (2.1):

$$\frac{1}{N}\underset{i=1}{\overset{2dNf}{}}\lambda _{i,N}(F_e)=\frac{\sigma _N(F_e)}{nk_B}=\frac{J_N(F_e)F_e}{nk_BT}=\frac{L_N(F_e)F_e^2}{nk_BT}.$$
(2.2)   
Here $`\sigma _N(F_e)`$, $`J_N(F_e)`$ and $`L_N(F_e)`$ are the entropy production rate per unit volume, the dissipative flux and the transport coefficient respectively, induced in the system in the stationary state by the external force $`F_e`$, where a nonlinear constitutive relation $`J_N(F_e)=L_N(F_e)F_e`$ has been used. The subscript $`N`$ indicates the $`N`$-dependence for finite systems. The kinetic temperature $`T`$ is determined by the relation:

$$\frac{1}{dNd1}\underset{i=1}{\overset{N}{}}\frac{𝐩_i^2}{m}k_BT,$$
(2.3)   
where $`m`$ is the particle mass, $`\{𝐩_i,i=1,N\}`$ are the peculiar momenta and $`T`$ is the kinetic temperature.

For systems at equilibrium the thermodynamic temperature appearing in Eq. (2.3) should strictly only be calculated in the thermodynamic limit. However, for nonequilibrium systems, as mentioned above the application of the limit $`N\mathrm{}`$ is not straightforward.

Still, for sufficiently large values of $`N`$, Eq. (2.2) can be interpreted as:

$$\frac{1}{N}\underset{i=1}{\overset{2dNf}{}}\lambda _{i,N}(F_e)=\frac{\sigma (F_e)}{nk_B}+O(N^1)=\frac{L(F_e)F_e^2}{nk_BT}+O(N^1),$$
(2.4)   
where by $`O(N^1)`$ we mean that - at least at equilibrium - the finite size corrections can be bounded by a function of the form $`CN^1`$ with $`C`$ of order 1. Although this condition on the constant $`C`$ is not mathematically precisely defined, we discuss in section 5 numerical experiments, which will give an indication of the magnitude of the $`O(N^1)`$ corrections in the relationship between the KY–dimension and the viscosity. In what follows when we write $`O(N^1)`$, we will always intend it to carry the particular meaning given by the above discussion.

Since we are interested here in the linear ($`F_e`$) regime one can use that $`L_N(F_e)=L_N+O(F_e^2)`$,<sup>6</sup><sup>6</sup>6Here and in what follows we assume that the transport coefficients are even in $`F_e`$. where $`L_N=L_N(0)0`$ is the linear transport coefficient. We note that the left hand side of Eq. (2.2) is not restricted to small fields, due to its dynamical origin, and that the right hand side is related for small fields to the usual entropy production rate per unit volume of Irreversible Thermodynamics.$`^{\text{[10]}}`$

We now introduce the Kaplan–Yorke (KY) dimension. If $`N_{KY}`$ is the largest integer for which $`\underset{i=1}{\overset{N_{KY}}{}}\lambda _{i,N}(F_e)>0`$, the KY–dimension, $`D_{KY,N}`$,$`^{\text{[1]}}`$ for a finite system with a discrete Lyapunov spectrum, is defined by$`^{\text{[1, 7]}}`$:

$$D_{KY,N}=N_{KY}+\frac{\underset{i=1}{\overset{N_{KY}}{}}\lambda _{i,N}(F_e)}{|\lambda _{N_{KY}+1,N}(F_e)|}$$
(2.5)   
where we have not indicated the $`N`$-dependence of $`N_{KY}`$.

## 3. SMALL PHASE SPACE REDUCTION

In case the phase space dimension reduction is smaller than one, Eq. (2.5) can be reduced to the exact equation:

$$\underset{i=1}{\overset{2dNf}{}}\lambda _{i,N}(F_e)=\lambda _{min,N}(F_e)(2dNfD_{KY,N}(F_e))=\frac{\sigma _N(F_e)V}{k_B},$$
(3.1)   
where the minimum Lyapunov exponent, $`\lambda _{min,N}(F_e)=\lambda _{2dNf,N}(F_e)`$. From Eqs. (3.1) and (2.2) we can then trivially calculate the linear (i.e. the limiting zero field) transport coefficient in linear response theory as,

$$L_N=\underset{F_e0}{lim}\frac{(2dNfD_{KY,N}(F_e))\lambda _{max,N}(F_e)nk_BT}{NF_e^2}.$$
(3.2)   
A similar relation has been obtained for the periodic Lorentz gas on the basis of periodic orbit theory.$`^{\text{[11]}}`$ We remark that from the point of view of linear response theory, i.e. Eq. (3.2), a phase space dimension reduction smaller than one occurs for any $`N`$, including $`N\mathrm{}`$. In that case one can let $`N\mathrm{}`$ in Eq. (3.2) and obtain a new exact relation for the linear transport coefficients, $`L=lim_N\mathrm{}L_N`$, equivalent to the Green-Kubo formulae. The corresponding expression for the KY–dimension of the steady state attractor for sufficiently small fields, is given by:

$$D_{KY,N}(F_e)=2dNf\frac{LF_e^2N}{\lambda _{max,N}nk_BT}+O(F_e^4).$$
(3.3)   
To obtain the Eqs. (3.2) and (3.3), we have used that for systems which are symplectic at equilibrium, (i.e. all Hamiltonian equilibrium systems), one can write for small $`F_e`$: $`\lambda _{max,N}(F_e)=\lambda _{max,N}+O(F_e^2)=\lambda _{min,N}+O(F_e^2)`$, where $`\lambda _{max,N}\lambda _{max,N}(0)`$ and $`\lambda _{min,N}\lambda _{min,N}(0)`$.

For systems which satisfy the Conjugate Pairing Rule$`^{\text{[12]}}`$ (CPR) the sum of each conjugate pair, $`i,i^{}=2dNfi+1`$ of Lyapunov exponents is

$$\lambda _{i,N}(Fe)+\lambda _{i^{},N}(Fe)=\frac{2\sigma _N(Fe)V}{k_B(2dNf)},i.$$
(3.4)   
We note that in nonequilibrium systems, the Conjugate Pairing Rule is expected to hold only in systems that are thermostatted homogeneously. In some systems there is numerical evidence that the maximal exponents satisfy the Conjugate Pairing Rule, while the other pairs do not (that is Eq. (3.4) is true for i=1).$`^{\text{[13]}}`$ These systems are said to satisfy the weak Conjugate Pairing Rule (WCPR). By combining the Eqs. (3.1) and (3.4), one obtains for sufficiently small fields:

$$\frac{D_{KY,N}(F_e)}{(dNf/2)}=1\frac{\lambda _{max,N}(F_e)}{\lambda _{min,N}(F_e)}.$$
(3.5)   
Substituting Eq. (3.4) into Eq. (2.2) and using Eq. (3.3), one obtains another expression for the limiting KY–dimension, for sufficiently small fields:

$$\frac{D_{KY,N}(F_e)}{(dNf/2)}=3+\frac{\lambda _{min,N}(F_e)}{\lambda _{max,N}(F_e)}+O(F_e^4)$$
(3.6)   
In Eqs. (3.5) and (3.6), we have chosen to use the maximal Lyapunov exponents as the conjugate pair in Eq. (3.4), and therefore these equations are valid provided WCPR is obeyed. A similar looking formula, as Eq. (3.6) with the first two terms on the right hand side only, has been quoted in .

As mentioned before, all the results in this section hold under the hypothesis that the phase space dimension reduction is smaller than unity. For this to be true for any given small $`F_e`$, $`N`$ is constricted to be of $`O(F_e^2)`$, or equivalently, for any given large $`N,F_e`$ to be of $`O(N^{1/2})`$, as can be seen from Eq. (3.3), so that $`F_e`$ and $`N`$ are coupled. <sup>7</sup><sup>7</sup>7One might think that a phase space dimension reduction of one could hardly have any practical consequence in a macroscopic system whose phase space dimension is of the order of $`10^{23}`$. As we will discuss in Section 6, such a very small phase space dimension reduction is expected to occur under physically realizable conditions.

However, from a general physical point of view, we would like to have a theory for $`D_{KY,N}`$, which holds uniformly in $`N`$, i.e. with Eq. (3.3) valid for every $`N`$ with a small but fixed $`F_e`$, so that $`N`$ and $`F_e`$ are independent variables. Now, it is trivial to generalize Eq. (3.1) to the case in which the dimensional reduction is greater than one. In fact one then obtains that:

$`{\displaystyle \underset{i=1}{\overset{2dNf}{}}}\lambda _{i,N}(F_e)`$ $`=`$ $`\lambda _{N_{KY}+1,N}(F_e)(2dNfD_{KY,N}(F_e))`$ (3.7)   
$`+{\displaystyle \underset{i=N_{KY}+2}{\overset{2dNf}{}}}(\lambda _{i,N}(F_e)\lambda _{N_{KY}+1,N}(F_e)).`$   
If we assume, as is usually done, that for sufficiently large $`N`$, $`\lambda _{2dNfj}=\lambda _{2dNf}+O(N^1)`$ for fixed $`j`$ not varying with N; then for any fixed dimensional reduction, Eq. (3.7) simply becomes Eq. (3.1) with a correction of $`O(1)`$ in $`N`$. Keeping the phase space dimension reduction fixed as $`N`$ increases still imposes a condition on $`F_e`$ of the form discussed above. To better control the errors when $`2dNfN_{KY}`$ becomes large ($`O(N)`$), a more careful treatment of the Eq. (3.1) is needed in order to obtain an expression uniform in $`N`$ for the phase space dimension reduction.

## 4. LYAPUNOV SPECTRUM FOR VERY LARGE $`N`$ AND LARGE PHASE SPACE REDUCTION

For large $`N`$ one can indeed reformulate the definition of $`D_{KY,N}`$ in a more analytical way. Consider thereto the stepwise continuous function of a continuous variable $`0<x2dNf`$: $`\stackrel{~}{\lambda }_N(x,F_e)=\lambda _{i,N}(F_e)`$ for $`i1<xi`$. If one introduces the integral:

$$\stackrel{~}{\mu }_N(x,F_e)=_0^x\stackrel{~}{\lambda }_N(y,F_e)𝑑y,$$
(4.1)   
Eq. (2.5) can be rewritten in the form:

$$\stackrel{~}{\mu }_N(D_{KY,N}(F_e),F_e)=0.$$
(4.2)   
Since $`D_{KY,N}(F_e)`$ as well as $`\stackrel{~}{\mu }_N(x,F_e)`$ are expected to be extensive quantities, i.e. they are proportional to $`N`$ for large $`N`$, it is natural to define the quantity

$$\delta _N(F_e)=\frac{D_{KY,N}(F_e)}{2dNf}$$
(4.3)   
where $`\delta _N(F_e)`$ is the dimension per effective degree of freedom ($`2dNf`$) in phase space. We now use again that, as $`N`$ grows, the difference $`\lambda _{i,N}\lambda _{i+1,N}`$ is expected to go to zero as $`N^1`$. This suggests a possible rescaling of the variable $`x`$ in $`\stackrel{~}{\lambda }_N(x,F_e)`$ to define the function $`\lambda _N(x,F_e)=\stackrel{~}{\lambda }_N(xN,F_e)`$, i.e. $`\lambda _N(x,F_e)=\lambda _{i,N}(F_e)`$ for $`\frac{i1}{2dNf}<x\frac{i}{2dNf}`$. The function $`\lambda _{i,N}(F_e)`$ can then be expected to be well approximated by a continuous function of the variable $`x`$ when $`N`$ is sufficiently large. Thus rewriting Eq. (4.1) as

$$\mu _N(x,F_e)=_0^x\lambda _N(y,F_e)𝑑y$$
(4.4)   
where $`0x1`$, Eq. (4.2) is equivalent to:

$$\mu _N(\delta _N(F_e),F_e)=0$$
(4.5)   
where $`\delta _N(F_e)`$ and $`\mu _N(x,F_e)`$ can be assumed to be intensive quantities. We can make this more precise by the following crucial assumption on the nature of the function $`\lambda _N(y,F_e)`$, i.e. of the Lyapunov spectrum:

Smoothness Hypothesis: If $`N`$ is sufficiently large, one can write:

$$\lambda _N(x,F_e)=l(x,F_e)+O(N^1)$$
(4.6)   
with $`l(x,F_e)`$ a smooth function<sup>8</sup><sup>8</sup>8It is enough to assume that $`l(x,F_e)L^1`$ in $`x`$ and $`C^2`$ near $`x=1`$. Moreover we will require that $`l(x,F_e)`$ is $`C^2`$ in $`F_e`$ for every $`x`$. This may not be the case when a phase transition occurs.$`^{\text{[15]}}`$ of $`x`$ and $`F_e`$.

Eq. (4.6), as Eq. (2.4), is to be interpreted that for sufficiently large $`N`$, $`|\lambda _N(x,F_e)l(x,F_e)|<C^{}N^1`$ with $`C^{}`$ of order 1. Similarly, as mentioned already above, all the $`O(N^1)`$ terms appearing in forthcoming equations (see e.g. Eq. (4.12)) have to be interpreted in this way and the constants obtained then can be directly expressed in term of $`C^{}`$.

We now define, in analogy with Eqs. (4.4) and (4.5):

$$m(x,F_e)=_0^xl(y,F_e)𝑑y$$
(4.7)   
and $`d(F_e)`$ through the equation:

$$m(d(F_e),F_e)=0$$
(4.8)   
respectively. Clearly our Hypothesis implies that:

$$\delta _N(F_e)=d(F_e)+O(N^1)$$
(4.9)   
where $`d(F_e)`$ and $`m(F_e)`$ are now intensive quantities. The function $`m(x,F_e)`$ is sketched in Fig. 1 1, especially near $`x=1`$.

From this figure one easily deduces that:

$$d(F_e)=1\frac{m(1,F_e)}{m^{}(1,F_e)}+O(m(1,F_e)^2)=1+\frac{1}{l(1,F_e)}\frac{\sigma (F_e)}{2dnk_B}+O(F_e^4),$$
(4.10)   
where $`m^{}(F_e)`$ is the derivative of $`m(x,F_e)`$ with respect to x at $`x=1`$. Here we used that $`m^{}(1,F_e)=l(1,F_e)=\lambda _{min}(F_e)`$ and that $`m(1,F_e)`$, is the sum over all Lyapunov exponents, divided by $`(2dNf)`$, i.e. the phase space contraction (or entropy production) rate per (effective) degree of freedom in phase space.<sup>9</sup><sup>9</sup>9In a more analytic way Eq. (4.10) follows from the inverse function theorem together with our Smoothness Hypothesis and the fact that $`m(0,0)=0`$. Using Eq. (2.2), Eq. (4.10) can be rewritten as:

$$d(F_e)=1\frac{1}{\lambda _{max}}\frac{LF_e^2}{2dnk_BT}+O(F_e^4)$$
(4.11)   
where the maximum and minimum Lyapunov exponents at equilibrium, (i.e. when $`F_e=0`$) are $`\lambda _{max}=l(0,0)`$ and $`\lambda _{min}=l(1,0)`$, respectively. Moreover we have used that $`\lambda _{max}=\lambda _{min}`$ and that $`l(1,F_e)=l(1,0)+O(F_e^2)`$. In terms of the extensive quantity $`D_{KY,N}(F_e)`$, Eq. (4.11) can be rewritten as:

$$\frac{D_{KY,N}(F_e)}{2dNf}=1\frac{1}{\lambda _{max}}\frac{LF_e^2}{2dnk_BT}+O(F_e^4)+O(N^1).$$
(4.12)   
where we kept the term $`f/(2dN)`$ of $`O(N^1)`$ on the left hand side of Eq. (4.12) to facilitate comparison in figures 4 and 5 of section 5 for $`N=32`$ particles. We observe that Eq. (4.12) is formally very similar to Eq. (3.3), except for the correction term $`O(N^1)`$ and the substitution of the asymptotic $`\lambda _{max}`$ for the finite $`N`$ value $`\lambda _{max,N}`$. It is clear that the same generalization can be performed for the Eqs. (3.5) and (3.6), if one notes that in the present context, the CPR, Eq. (3.4), becomes:

$$l(x,F_e)+l(1x,F_e)=\frac{\sigma (F_e)}{dnk_B}.$$
(4.13)   
From Eq. (4.11) it follows that in the linear regime the reduction in phase space dimension in large thermostatted dissipative systems is extensive. For those systems for which the Smoothness Hypothesis holds, this result is exact. The extensive nature of the reduction has been noted before.$`^{\text{[8, 16, 17, 18]}}`$

## 5. NUMERICAL TEST

We tested our Smoothness Hypothesis as well as Eq. (4.12) and equations derived from it assuming that the WCPR is valid, for a system of 32 WCA particles$`^{\text{[19]}}`$ undergoing shear flow in two dimensions. Although this system does not satisfy the CPR, it does appear to satisfy WCPR to within $`0.7\%`$ when $`N=32`$$`^{\text{[13]}}`$ (which was within the limits of numerical error achieved). The equations of motion for particles in such a system are the so-called SLLOD equations,$`^{\text{[20]}}`$

$$\dot{𝐪}_i=𝐩_i/m+𝐢\gamma y_i,\dot{𝐩}_i=𝐅_i𝐢\gamma p_{yi}\alpha 𝐩_i.$$
(5.1)   
Here, at not too large Reynolds numbers, the momenta, $`𝐩_i`$ are peculiar momenta, $`𝐢`$ is the unit vector in the $`x`$ direction and $`𝐅_i`$ is the force exerted on particle $`i`$ by all the other particles, due to a Weeks-Chandler-Andersen pair interaction potential$`^{\text{[19]}}`$ between the particles. The value of $`\alpha `$ is determined using Gauss’ Principle of Least Constraint to keep the kinetic temperature fixed.$`^{\text{[20]}}`$ The SLLOD equations of motion given by Eq.(5.1) then model Couette flow when the Reynolds number is sufficiently small so that laminar flow is stable.$`^{\text{[20]}}`$ <sup>10</sup><sup>10</sup>10The $`\alpha `$ appearing in the Eq. (5.1) is related to the phase space contraction rate $`\mathrm{\Lambda }`$ in Eq. (2.1). It is given by the relation $`_{i=1}^{2dNf}\lambda _i=\mathrm{\Lambda }=dN\alpha +O(1)`$.

For this system the dissipative flux $`J`$ is just the xy element of the pressure tensor $`P_{xy}`$; the transport coefficient $`L`$ is the Newtonian shear viscosity $`\eta `$; and the external field $`F_e`$ is the shear rate $`\gamma =u_x/y`$,$`^{\text{[20]}}`$ where $`u_x`$ is the local flow velocity in the $`x`$-direction in the system. The calculations were carried out at a reduced kinetic temperature of unity and a reduced density n$`=N/V=0.8`$. All physical quantities occurring in the Eq. (5.1) as well as the temperature and density are made dimensionless by reducing them with appropriate combinations of molecular quantities. In particular the reduction factor for $`\gamma `$ amounts to about 1 ps<sup>-1</sup> for Argon. (c.f. ) We note that for the SLLOD equations for shear flow in 2 dimensions using non autonomous Lees-Edwards periodic boundary conditions,$`^{\text{[20, 21]}}`$ $`f=5`$ (due to the conservation of kinetic energy, momentum and position of the center of mass).

Figures 2 and 3 represent a direct test of our Smoothness Hypothesis. Figure 2 gives the Lyapunov spectrum for a system at a shear rate $`\gamma =0.05`$ and for a number of particles $`N=8,18,32`$. As can be easily observed the Lyapunov exponents for $`N`$=18 and 32 just fill the “open spaces” left by the exponents for $`N`$=8 and 18, respectively.

Figure 3 shows the behavior of $`\lambda _{max}(\gamma )`$ and $`\lambda _{min}(\gamma )`$ as functions of the externally applied field $`\gamma `$. Although no numerical experiment can in general confirm a mathematical hypothesis, the numerical results seem to agree very well with our Smoothness Hypothesis.

We now use the Lyapunov spectrum to compute $`D_{KY,N}(\gamma )`$. This is shown in figure 4a where the $`D_{KY,N}(\gamma )`$ is plotted for $`N=32`$, $`d=2`$ and $`0<\gamma <0.5`$. As can be easily seen, although from the definition Eq. (2.5) one would generally expect to see discontinuities in the first derivatives of this function for those values of $`\gamma `$ for which $`N_{KY}`$ change (c.f. ), the function appears very smooth for a rather large range of values of $`\gamma `$. This indicates that $`N=32`$ is already “big enough” to consider the Lyapunov spectrum as “effectively continuous”.<sup>11</sup><sup>11</sup>11More precisely one can say that for $`N=32`$ the numerical errors involved in computing the Lyapunov exponents are already larger than the errors introduced by neglecting the $`O(N^1)`$ correction. This observation provides a more precise meaning of expressions like “sufficiently large $`N`$” or “a constant $`C`$ of order 1”. This also permits us to check the validity of Eq. (4.12), and the large $`N`$ versions of Eqs. (3.5) and (3.6), which in this case take the form:

$`{\displaystyle \frac{D_{KY,N}(\gamma )}{2dNf}}=1{\displaystyle \frac{\eta _N\gamma ^2}{\lambda _{max,N}2dnk_BT}}+O(\gamma ^4)+O(N^1)`$ (5.2)   
$`{\displaystyle \frac{D_{KY,N}(\gamma )}{dNf/2}}=1{\displaystyle \frac{\lambda _{max,N}(\gamma )}{\lambda _{min,N}(\gamma )}}+O(\gamma ^4)+O(N^1)`$ (5.3)   
$`{\displaystyle \frac{D_{KY,N}(\gamma )}{dNf/2}}=3+{\displaystyle \frac{\lambda _{min,N}(\gamma )}{\lambda _{max,N}(\gamma )}}+O(\gamma ^4)+O(N^1),`$ (5.4)   
respectively. Figure 4a shows that these expressions can describe the system studied over the range of fields considered and that the correction terms are small. The deviations in the values obtained using the various expressions is within the numerical error in the data.

We observe here that $`D_{KY,N}`$ of Eq. (2.5) is well approximated by a fourth order even polynomial in $`\gamma `$ for fields up to $`\gamma 0.5`$ and that at $`\gamma 0.3`$ the terms of $`O(\gamma ^4)`$ are just $`5\%`$ of the terms of order $`O(\gamma ^2)`$, while the phase space dimension reduction is clearly greater than unity. We note that the quadratic dependence of $`D_{KY,N}`$ on $`\gamma `$ implies a linear dependence of the shear stress on $`\gamma `$,$`^{\text{[23]}}`$ so the linear regime for $`N=32`$ extends beyond strain rates where the phase space dimension reduction is smaller than one, which runs (for $`N=32`$) till approximately $`\gamma =0.175`$. At $`\gamma =0.3`$, $`\eta (\gamma )`$ is just $`8\%`$ smaller than $`\eta (0)`$.

For completeness, in Figure 4b, we expand Figure 4a in the regime where the phase space reduction is smaller than unity. In this regime, Eqs. (3.3), (3.5) and (3.6) are expected to be valid for the calculation of $`D_{KY,N}`$ and any deviations of the data from the value calculated using (2.5) are due to the limited numerical precision of the independent calculations of the viscosity and the Lyapunov exponents, $`O(F_e^4)`$ corrections, or to the assumption that the WCPR is obeyed. The numerical results indicate that the numerical error gives the most significant contribution to the deviations observed. The WCPR is found to be valid at least to within numerical error for the state points considered. This is consistent with previous work.$`^{\text{[13]}}`$

Using the fourth order fit to $`D_{KY,N}`$, mentioned above, i.e.  $`D_{KY,N}(\gamma )=2dNf+\frac{1}{2}D_{KY,N}^{^{\prime \prime }}\gamma ^2+\frac{1}{4}D_{KY,N}^{^{\prime \prime \prime \prime }}\gamma ^4`$ for the shear rate range $`[0.2,+0.2]`$, (c.f. figure 4b), one can calculate the zero strain rate or Newtonian shear viscosity from $`\eta =D_{KY,N}^{^{\prime \prime }}`$$`\lambda _{max}nk_BT/2N`$. Here $`D_{KY,N}^{^{\prime \prime }}`$ and $`D_{KY,N}^{^{\prime \prime \prime \prime }}`$ are the second and fourth derivatives of $`D_{KY,N}(\gamma )`$, respectively, with respect to $`\gamma `$, taken at $`\gamma =0`$. This leads to a value of $`\eta =2.54\pm 0.07`$, which agrees very well within statistical uncertainties with the Newtonian viscosity directly measured in the simulation and calculated from the defining constitutive relation $`\eta (\gamma )P_{xy}(\gamma )/\gamma `$, viz. $`\eta =2.52\pm 0.05`$ (see the + points in figure 4b). Note that although the precision of $`D_{KY,N}`$ is high (less than $`0.03\%`$ statistical error), the viscosity is related to the phase space contraction $`(2dNfD_{KY,N})`$ and thus calculation of the viscosity from the phase space dimension reduction involves, at small fields, a very small difference between two large numbers, resulting in a larger relative error in the viscosity. Furthermore since the dimensional contraction is a quadratic function of the shear rate, it results in a more difficult calculation of the viscosity from the phase space dimension reduction than the constitutive equation which is linear in the shear rate at low strain rates.

Figure 5 compares the values of $`D_{KY,N}`$ determined from Eqs. (5.2), (5.3) and (5.4) with those calculated from its definition in Eq. (2.5), as a function of $`N`$ for fields of $`\gamma =0.15`$ and $`\gamma =0.5`$. In Figure 5a the results for the field of $`\gamma =0.15`$ are shown, and at this strain rate the dimensional contraction will be less than unity for systems size up to approximately $`N=50`$. Therefore, for $`N\stackrel{<}{}50`$, Eq. (3.3), which contains no $`O(N^1)`$ corrections, will be valid, as will Eqs. (3.5) and (3.6) if the WCPR is obeyed. All the numerical results are consistent with the theory and with the assumption that the WCPR is obeyed.$`^{\text{[13]}}`$ At $`\gamma =0.15`$, the values determined using the various methods have a maximum difference of $`0.2\%`$ and are within the numerical errors at each particle number. In Figure 5b, the results for a strain rate of $`\gamma =0.5`$ are shown. In this case, a dimensional contraction of less than unity is only obtained for $`N\stackrel{<}{}5`$. Again, the deviations between the results calculated using Eqs. (5.2)- (5.4) and the definition of $`D_{KY,N}`$ given by Eq. (2.5), are small (at most $`1\%`$), and within the limits of error for all particle numbers considered. This confirms that the coefficients of the $`O(N^1)`$ and $`O(F_e^4)`$ terms are small for this system.

## 6. CONCLUSIONS

We mention here a few implications of the results presented in this paper.

* The extensivity of the phase space reduction for large $`N`$ and small fields is here, to the best of our knowledge, demonstrated for the first time, on the basis of the Smoothness Hypothesis of the Lyapunov spectrum and the extensivity of the total entropy production.
* The relationships given by Eqs. (3.2), (3.3) and (5.2) also apply to systems where not all particles are thermostatted. That is, they can be applied to systems where the Gaussian thermostat operates on selected particles, say those in the boundaries, while the remaining particles evolve under Newtonian dynamics, supplemented by a dissipative field.$`^{\text{[16, 24]}}`$ We note that the Eqs. (3.5), (3.6), (5.3) and (5.4) can only be assumed to apply to homogeneously thermostatted systems in general, since only for such systems can the WCPR be expected to hold.
* A simple calculation shows that for a typical case as that of one mole of Argon at its triple point, sheared at the rate of 1 Hz, the difference of the Kaplan-Yorke dimension and the phase space dimension $`(O(10^{23}))`$ is tiny, namely $`3`$. This follows from Eq.(3.3), which shows that the dimension loss, when measured in moles, is equal to the product of the total entropy production rate of the system and the reciprocal of the largest Lyapunov exponent. Since the largest Lyapunov exponent is controlled by the most unstable atomic processes, it is always very small $``$ 1 ps<sup>-1</sup> whether for atomic, molecular and even polymeric systems.

We note that this smallness of the phase space dimension reduction in irreversible processes near equilibrium could well be the reason that linear Irreversible Thermodynamics provides such a good description of nonequilibrium systems close to equilibrium. This is because the thermodynamic properties are insensitive to the high order distribution functions - including the full N-particle distribution function of the entire system - since they are determined by a few low order distribution functions, which “do not know” that the dimension of the steady state attractor is only a few dimensions smaller than the $`10^{23}`$ of the phase space of the system.
* For the system studied, the $`O(N^1)`$ corrections to $`D_{KY,N}/(2dNf)`$ that appear in Eqs. (5.2), (5.3) and (5.4), due to the smoothness hypothesis are small (see Figure 5) and less than $`1.0\%`$ even for system sizes of $`N=6`$.
* In conclusion, although the new relations involving $`D_{KY,N}`$ are simple consequences of the Eq. (4.10), they could nevertheless be useful for applications. In particular the equations give a simple expression for the Kaplan–Yorke dimension of the attractor of a class of many particle systems close to equilibrium i.e. in the regime of linear dissipation near equilibrium.

## 7. ACKNOWLEDGEMENTS

EGDC gratefully acknowledges the hospitality of the Research School of Chemistry of the Australian National University as well as financial support from the Australian Research Council, the Australian Defence Forces Academy and the Engineering Research Program of the Office of Basic Energy Sciences of the US Department of Energy under Grant No. DE-FG02-88-ER13847. DJS and DJE thank the Australian Research Council for support of this project.