Noura Dawass

4 70 S URFACE E FFECTS temperature T = 2, dimensionless densities ρ ranging from 0.2 to 0.8 and using number of particles N equals to 100, 500, 1000, 5000, 10000, 30000, and 50000. For each system size, the length of the simulation box L was set according to the required density. The desired temperature was imposed using the Andersen thermostat [22] . All MD simulations started from a randomly-generated configuration for which an energy minimization was used to eliminate particle overlaps. A suffi- cient number of time steps was used to initialize the system. After initialization, RDFs were sampled every 100 time steps. For both, initialization and produc- tion, a dimensionless time step equal to 0.001 was used. The simulation length was chosen depending on the size of the system and the available computational resources. For instance, for systems with N = 100, 1x10 9 production time steps were carried out, while for the maximum size N = 50000, 7x10 5 steps were used. Multiple independent simulations were performed for each point ( ρ , N ). The re- sulting RDFs were then averaged and used to compute G ∞ and F ∞ . At high den- sities ( ρ > 0.4), RDFs from at least 10 runs are used. At lower densities, at least 20 runs are performed to enhance statistics. 4.3. R ESULTS 4.3.1. E STIMATION OF K IRKWOOD –B UFF I NTEGRALS KBIs in the thermodynamic limit G ∞ are obtained using the three different ap- proaches discussed earlier. To compare the estimation methods, WCA systems were studied while fixing temperature and density. These parameters define the thermodynamic state of the system. Values of KBIs, computed using different methods, for other densities for LJ and WCA fluids are provided in the Support- ing Information of Ref. [128] . After comparing estimation methods of KBIs, the relation between density of the system and KBIs for LJ and WCA system is dis- cussed. Figure 4.1 shows RDFs for systems of different sizes of a WCA fluid at T = 2 and ρ = 0.6 (dimensionless units). Figure 4.1 (b) shows that using small system sizes, specifically N = 100 and N = 500, results in RDFs with higher oscillations than large systems, where N equals to or is larger than 1000. As will be shown later, this causes implications in the computation of G ∞ . In Figure 4.2, the scal- ing of KBIs of finite subvolumes G V with 1/ L is presented. For large systems, where N > 500, a linear range is identified which can be extrapolated to the limit 1/ L → 0. Instead of computing G V , G ∞ can be directly estimated from RDFs us- ing Eqs. (4.1) and (4.2) . Figure 4.3 shows the estimation of G 2 for systems with varying sizes. When plotted as a function of L , the values of the integrals G 2 show a plateau at a constant value which corresponds to G ∞ . However, Fig-