Bulk viscosity of molecular ﬂuids

The bulk viscosity of molecular models of gases and liquids is determined by molecular simulations as a combination of a dilute gas contribution, arising due to the relaxation of internal degrees of freedom, and a configurational contribution, due to the presence of intermolecular interactions. The dilute gas contribution is evaluated using experimental data for the relaxation times of vibrational and rotational degrees of freedom. The configurational part is calculated using Green-Kubo relations for the fluctuations of the pressure tensor obtained from equilibrium microcanonical molecular dynamics simulations. As a benchmark, the Lennard-Jones fluid is studied. Both atomistic and coarse-grained force fields for water, CO2, and n-decane are considered and tested for their accuracy, and where possible, compared to experimental data. The dilute gas contribution to the bulk viscosity is seen to be significant only in the cases when intramolecular relaxation times are in the μs range, and for low vibrational wave numbers (<1000 cm-1); This explains the abnormally high values of bulk viscosity reported for CO2. In all other cases studied, the dilute gas contribution is negligible and the configurational contribution dominates the overall behavior. In particular, the configurational term is responsible for the enhancement of the bulk viscosity near the critical point.


Full pressure tensor
The pressure tensor elements are given by [1,2] where the sum is over N a atoms, and p i α and p i β are the particle momenta in the respective direction. The second term is a sum over all forces between atoms i and j.
P αβ is expressed as: Thus it ia 0 for α = β. In the microcanonical ensemble (constant number, volume, energy) we have N = N and E = E , such that P αβ = δ αβ P.

Influence of the Lennard-Jones cut-off on transport properties
The influence of the potential cut-off, r cut , is examined briefly in order to choose a cut-off value which does not affect results, as done in previous studies [3,4]. A nearcritical isotherm (T * = 1.35) was chosen to investigate the effects on the enhanced bulk viscosity near the critical point. As can be seen in Fig. S1a, the shear viscosity, η * , remains unaffected by the chosen cut-off. The bulk viscosity, however, is larger for larger cut-off values. The increase of the ratio, κ * /η * , for larger cut-offs is significant (Fig. S2). Therefore we can conclude that the long-range interactions contribute to the bulk viscosity, κ * , but not the shear viscosity. For larger r cut the effect of this becomes smaller.
S  Figure S1: Shear and bulk viscosity for an LJ fluid at T * = 1.35 for several cut-off parameters, r cut . The filled markers present data from this work with the crosses being obtained from Meier et al. (2004Meier et al. ( ,2005 [3,5] for r cut = 5.5. κ * /η * r cut = 2.5 r cut = 3.5 r cut = 4.5 r cut = 5.5 Figure S2: Viscosity ratio, κ * /η * , for an LJ fluid at T * = 1.35 for several cut-off parameters, r cut .

Simulation parameters for water, CO 2 and n-decane
Tables S1, S2, S3, S4 and S5 show the simulation parameters used in this work.    Table S3: Simulation parameters for different CO 2 models.

Additional data: Water at different temperatures
In addition to the data presented in the paper we also examine water using both atomistic (SPC/E, TIP4P/2005) and coarse-grained models (two versions of the SAFT-γ Mie model) at 393K and the near critical 613K.
The SAFT-γ Mie force field for water [6] not only has temperature dependent parameters σ and but there are also two different parametrizations performed over different temperature ranges. SAFT-ift is applicable in the lower temperature range (293 − 493 K) and was parametrized to fit to interfacial tension data whilst SAFT-vle was parametrized for 343 − 613 K using vapour-liquid equilibrium data.
We examine the pressure, shear viscosity and bulk viscosity data for water at two different temperatures; 393 K (SPC/E, SAFT-ift, SAFT-vle) and 613 K (SPC/E, TIP4P/2005 and SAFT-vle). The former is to compare the two coarse-grained parametrizations whilst the latter is to obtain data near the critical point. As SPC/E performs best in the liquid regime [7,8] this is further a crude test for the applicability of SPC/E in the vapour phase. As no rotational or vibrational relaxation time data is available for these temperatures, only the configurational contribution of κ is presented. Once again both SAFT models perform poorly in predicting bulk viscosity values of the order of the shear viscosity, with κ conf /η well below 1. SPC/E, however, obtains values between 3-4, similar to those at room temperature. These values are higher than that of 2.32 observed at 300 K. The increase of κ conf /η is contrary to an observed decrease with increase temperature in experiment [9,10].
In conclusion we can assert that, despite its shortcomings in thermodynamic quantities, transport properties are still well predicted by the SPC/E model. The SAFT-ift model is superior to SAFT-vle both in terms of the thermodynamic predictions and the shear viscosity. The performance of SAFT-ift is much improved over room temperature simulations. Bulk viscosities, however, are still not well predicted.

Conclusions
The SAFT-ift model outperforms the SAFT-vle model both in terms of thermodynamic and transport properties. Attention should be drawn to the relatively good prediction of the shear viscosity at higher temperatures (SAFT-ift at 393 K and SAFT-vle at 613 K) showing that coarse-grained models for water may be successfully used in transport applications in some physical situations. An increase in the viscosity ratio near the critical point is observed similar to other fluids. However, as this is not present in generally superior atomistic models, further investigations near the critical point should be performed.

Additional data: Supercritical CO 2
We study two atomistic (EPM2, TraPPE) and two coarse-grained models (SAFT dimer, SAFT monomer) for the supercritical isotherm at 323 K. All models studied show qualitatively the same behaviour for supercritical CO 2 as at 300 K ( Figure S6). The SAFT monomer model once again overestimates the shear viscosity at high densities. EPM2 is the most accurate in predicting η whilst TraPPE and the SAFT dimer model have comparable results. The viscosity ratio (where only the configurational term is considered for κ) is overestimated in the lower density range, in particular in the vicinity of the critical point. With the exception of the high density region, where there is a S-11 larger amount of scatter in the data together with bigger error bars, EPM2, TraPPE and SAFT dimer agree well with each other for the viscosity ratio. In the high density region the SAFT dimer model gives the lowest values compared to the atomistic and the monomer models. The dilute gas contribution, κ dilute , is not considered here as it is the same for all models.
6 Additional data: CO 2 in the two-phase region In order to investigate the simulation box size dependence of transport properties in the two-phase regions, we have simulated CO 2 at ρ = 0.6 g/cm 3 using different numbers of molecules: 500, 1372, 2048 (as presented in the paper), 4000, and 6912. The results are presented in Fig. S7. A clear dependence of the quantities on the simulation box size/number of molecules is observed. Notably, for all but the smallest box, κ conf /η agree relatively well (almost within their error). Therefore we conclude that the results presented in the paper may still serve as a guide to the general behaviour in a two-phase region.

Additional data: Influence of constraints on n-decane at 300K
The constraints placed on the bonds and angles in alkanes such as n-decane are not usually treated with much care. In the work presented in the paper we contrained both the bonds involving hydrogens and the H-C-H angle for the OPLS force field, whilst this force field is often employed without any constraints whatsoever. However, it has been stated in the literature that for certain properties these constraints can have a significant impact [11].
We present here data for n-decane, where either both the hydrogen bonds and the H-C-H angle is constrained (OPLS, as presented in the paper) or both are governed by a harmonic potential (OPLS (flexible)). Our work shows that for the shear viscosities, while there are differences, they are small (see Fig. S8a). It should be noted that the rigid model achieves slightly better accuracy with regard to NIST data than its flexible counter part. The difference is more notable in the viscosity ratio (Fig. S8b). The ratio with the experimental value [12]. For L-OPLS a similar investigation was conducted with no notable difference observed. As only one experimental observation was found to date, this may in practice lie within the error of the measurement. However, in any case it highlights an issue with how force fields for alkanes are used and presented.
Further analysis needs to be done to confirm which constraints produce the results most in line with the true physical molecule.

Compressibilities for different water models
In order to investigate the differences in predicting the bulk viscosity between the SAFT monomer and atomistic models (SPC/E and TIP4P/2005), we calculate the compressibilities, β c , for all three models at 300 K. As a relation between the bulk viscosity and the bulk modulus, K = 1/β c , is proposed in [13], a difference in the models could give an insight into the reasons for the varying performance. The isothermal bulk viscosity is defined through where V is the volume, P the pressure and V 0 the volume at the desired pressure. In

Benchmark testing for κ dilute
Benchmark testing of κ dilute was performed by the example of CO 2 (i.e. using values of c v for CO 2 ) in order to assess the importance of the relaxation times, τ rot and τ vib , and the wave number, k, on the rotational and vibrational contributions to κ, κ rot and κ vib .
The main conclusions are presented in the paper. The dependence on the parameters is illustrated in Figs. S9 and S10.

Complete data for the viscosities at different state points
Below is a collection of the data for water (T = (300, 393, 613) K, models = (SPC/E,       Table S15: Shear and bulk viscosity data for the SAFT1-vle water model at T = 613 K.