While performing NEMD simulations of thermostatted shear flow for hard-sphere fluids, Erpenbeck (1984) observed that at very high shear rates, fluid particles organised themselves into strings. This was an early observation of a nonequilibrium phase transition. This organisation of particles into strings reduces the rate at which entropy is produced in the system by the external field. This effect is in competition with the kink instability of the strings themselves. If the strings move too slowly across the simulation cell, thermal fluctuations in the curvature of the strings lead to their destruction. A snapshot of a string phase is shown in Figure 6.11. The velocity gradient is vertical and the streaming velocity is horizontal. The system is 896 soft discs at a state point close to freezing and a reduced shear rate of 17.
The string phase is in fact, stabilised by the use of a thermostat which assumes that a linear velocity profile, (implicit in equation (6.44)), is stable. Thermostats which make some assumption about the form of the streaming velocity profile are called Profile Biased Thermostats (PBT). All the thermostats we have met so far are Profile Biased. At equilibrium there can be little cause for worry, the streaming velocity must be zero. Away from equilibrium we must be more careful.
Any kink instability that might develop in Erpenbeck’s strings, leading to their breakup, would necessarily lead to the formation of large scale eddies in the streaming velocity of the fluid. The Profile Biased Thermostat would interpret any incipient eddy motion as heat, and then thermostat would try to cool the system by suppressing the eddy formation. This in effect stabilises the string phase (Evans and Morriss, 1986).
Profile Biased Thermostats for shear flow assume that the kinetic temperature TB, for a system undergoing planar Couette flow can be defined from the equation,
(6.50)
In this equation d is the number of dimensions and N is the number of particles. The term iγyi is the presumed streaming velocity at the location of particle i. Once the form of the streaming velocity profile is established it is a simple matter to use peculiar velocity scaling, Gaussian isokinetic or Nosé methods to thermostat the shearing system.
At small shear rates and low Reynolds number, the Lees-Edwards shearing periodic boundary conditions do indeed lead to a planar velocity profile of the form assumed in (6.50). In Erpenbeck’s (1984) simulations the Reynolds numbers, (Re=ρmγL2∕η), were very large (103-105). The assumption of a linear streaming velocity profile under these conditions is extremely dubious. Suppose that at high Reynolds number the linear velocity profile assumed in (6.50) is not stable. In a freely shearing system with Lees-Edwards geometry, this might manifest itself in an S-shaped kink developing in the velocity profile. If (6.44) is used to maintain the temperature, the thermostat will interpret the development of this secondary flow as a component of the temperature. This increase in temperature will be continuously removed by the thermostat, leading to a damping of the secondary flow.
If we rewrite the SLLOD equations in terms of laboratory momenta,
(6.51)
then the momentum current, J,
(6.52)
satisfies the following continuity equation,
(6.53)
The derivation of this equation is carried out by a simple supplementation of the Irving-Kirkwood procedure (§3.7,8). We have to add the contribution of the thermostat to equations (3.115) and (3.123). Comparing equation 6.53) with the momentum conservation equation (2.12) we see that the thermostat could exert a stress on the system. The expected divergence terms (ρuu+P), are present on the right hand side of (6.53). However the term involving α, the thermostatting term, is new and represents the force exerted on the fluid by the thermostat. It will only vanish if a linear velocity profile is stable and,
(6.54)
At high Reynolds number this condition might not be true. For simulations at high Reynolds numbers one needs a thermostat which makes no assumptions whatever about the form of the streaming velocity profile. The thermostat should not even assume that a stable profile exists. These ideas led to development (Evans and Morriss, 1986), of Profile Unbiased Thermostats (PUT).
The PUT thermostat begins by letting the simulation itself define the local streaming velocity u(r,t). This is easily done by replacing the delta functions in (6.52) by microscopically small cells in the simulation program. The temperature of a particular cell at r, T(r,t), can be determined from the equation,
(6.55)
where n(r,t) is the number density at r,t (the delta function has unit volume). The number of degrees of freedom in the cell is dn(r,t)-d, because d degrees of freedom are used to determine the streaming velocity of the cell.
The PUT thermostatted SLLOD equations of motion can be written as,
(6.56)
The streaming velocity, u(r,t), is not known in advanced but is computed as time progresses from its definition, (6.52). The thermostat multiplier α, could be a Gaussian multiplier chosen to fix the peculiar kinetic energy (6.55). Equally well the multiplier could be a Nosé-Hoover multiplier. The momentum equation for the PUT thermostatted system reads,
(6.57)
From the definition the streaming velocity of a cell we know that, n(r,t)u(r,t)=∑(pi∕m)δ(ri(t)-r).
We also know that, n(r,t)u(r,t)=∑u(r,t)δ(ri(t)-r)=u(r,t)∑δ(ri(t)-r). Thus the thermostatting term in (6.57), vanishes for all values of r.
In terms of practical implementation in computer programs, PUT thermostats can only be used in simulations involving large numbers of particles. Thus far their use has been restricted to simulations of two dimensional systems. At low Reynolds numbers where no strings are observed in Profile Biased simulations, it is found that Profile Unbiased simulations yield results for all properties which are indistinguishable from those computed using PBT methods. However at high strain rates the results obtained using the two different thermostatting methods are quite different. No one has observed a string phase while using a PUT thermostat.