Accurate conductivity predictions of KOH(aq) are crucial for electrolysis applications. OH– is transferred in water by the Grotthuss transfer mechanism, thereby increasing its mobility compared to that of other ions. Classical and ab initio molecular dynamics struggle to capture this enhanced mobility due to limitations in computational costs or in capturing chemical reactions. Most studies to date have provided only qualitative descriptions of the structure during Grotthuss transfer, without quantitative results for the transfer rate and the resulting transport properties. Here, machine learning molecular dynamics is used to investigate 50,000 transfer events. Analysis confirmed earlier works that Grotthuss transfer requires a reduction in accepted and a slight increase in donated hydrogen bonds to the hydroxide, indicating that hydrogen-bond rearrangements are rate-limiting. The computed self-diffusion coefficients and electrical conductivities are consistent with experiments for a wide temperature range, outperforming classical interatomic force fields and earlier AIMD simulations.

Hydrate-based H2 storage is based on the mechanism of trapping H2 in water-based structures that are environmentally friendly and cost-efficient. Understanding the effects of common promoters on hydrate-based H2 storage at the molecular level is crucial for designing efficient storage systems, and for discovering novel promoters. Here, a series of μs-scale molecular dynamics simulations are performed to investigate the nucleation of binary H2 hydrates from gas-liquid two-phase solutions in the presence of various promoters, i.e., CH4, CO2, C2H6, C3H8, C5H10, and THF. The simulation results indicate that the H2 and promoter molecules first dissolve in water from the gas phase and then are absorbed on the cage-faces, promoting the nucleation and growth of binary H2 hydrates. THF is the most effective promoter for hydrate-based H2 storage, exhibiting high performance in converting H2 from the gas phase to hydrates. It is followed by CH4, C2H6, and CO2; C3H8 and C5H10 molecules are less effective H2 hydrate promoters. The presence of large promoter molecules enhances multi-occupied cage formation. The molecular insight into the nucleation of binary H2 hydrates with various promoters provided here not only contributes to a broader understanding of hydrate-based H2 storage but is expected to motivate further experimental and computational studies.

One of the most promising energy carriers for transport applications are hydrogen-based energy carriers. NaBH₄ is a hydrogen energy carrier and produces hydrogen bubbles when it is dissolved in water. The formation of hydrogen bubbles hinders experimental measurements of the thermodynamic and transport properties of aqueous NaBH₄ solutions at elevated temperatures. Accurate knowledge of these properties is essential for the NaBH₄ hydrolysis reactor modeling and design. Molecular dynamics (MD) simulations provide the option to study the thermodynamic and transport properties of NaBH₄ aqueous solutions without hindering hydrogen bubble formation. In this work, a new force field is developed for BH₄^-, namely, the Delft force field of BH₄^- (DFF/BH₄^-), which, combined with additional force fields, can accurately describe experimental densities and viscosities of 0 to 5 m (mol salt/kg water) NaBH₄, 0 to 3 m NaB(OH)₄, and 1 m NaOH aqueous solutions at 295 K within 1.8% and 10.8% maximum deviation, respectively. Empirical fitting correlations are created for densities, viscosities, and self-diffusivities obtained from the MD simulations of 0 to 5 m NaBH₄, 0 to 5 m NaB(OH)₄, and 0 to 1 m NaOH aqueous solutions at 295-363 K for NaBH₄ hydrolysis reactor modeling and design purposes.

The use of reactive working fluids in thermodynamic cycles is currently being considered as an alternative to inert working fluids, because of the preliminarily attested higher energy-efficiency potential. The current needs to simulate their use in thermodynamic cycles, which may operate in liquid, vapour or vapour-liquid state, are an accurate real-fluid equation of state and ideal gas thermochemical properties of each molecule constituting the mixture, to calculate the equilibrium constant. To this end, the appeal to a multi-scale theoretical methodology is paramount and its definition represents the objective of the present work. This methodology is applied and validated on the system N₂O₄ ⇌ 2NO₂. Firstly, the equations solved for simultaneous two-phase and reaction equilibrium are presented. Secondly, ideal gas thermochemical properties of N₂O₄ and NO₂ are computed at atomic scale by quantum mechanics simulations. Then, to apply the selected cubic equation of state, pure-component properties of the species forming the reactive mixture (critical point coordinates and acentric factor) are required as input. However, these properties are not measurable, since NO₂ and N₂O₄ do not exist in nature as pure components. To get around this difficulty, the methodology relies on molecular Monte Carlo simulations of the pure N₂O₄ and NO₂, as well as on the reactive N₂O₄ ⇌ 2NO₂, enabling the determination of those missing pure-component properties and thus the calculation, on a macroscopic scale, of the reactive mixture properties. Finally, the comparison of calculated mixture properties with available experimental data leads to validate the accuracy of the proposed methodology.

Focused ultrasound has experimentally been found to enhance the diffusion of nanoparticles; our aim with this work is to study this effect closer using both experiments and non-equilibrium molecular dynamics. Measurements from single particle tracking of 40 nm polystyrene nanoparticles in an agarose hydrogel with and without focused ultrasound are presented and compared with a previous experimental study using 100 nm polystyrene nanoparticles. In both cases, we observed an increase in the mean square displacement during focused ultrasound treatment. We developed a coarse-grained non-equilibrium molecular dynamics model with an implicit solvent to investigate the increase in the mean square displacement and its frequency and amplitude dependencies. This model consists of polymer fibers and two sizes of nanoparticles, and the effect of the focused ultrasound was modeled as an external oscillating force field. A comparison between the simulation and experimental results shows similar mean square displacement trends, suggesting that the particle velocity is a significant contributor to the observed ultrasound-enhanced mean square displacement. The resulting diffusion coefficients from the model are compared to the diffusion equation for a two-time continuous time random walk. The model is found to have the same frequency dependency. At lower particle velocity amplitude values, the model has a quadratic relation with the particle velocity amplitude as described by the two-time continuous time random walk derived diffusion equation, but at higher amplitudes, the model deviates, and its diffusion coefficient reaches the non-hindered diffusion coefficient. This observation suggests that at higher ultrasound intensities in hydrogels, the non-hindered diffusion coefficient can be used.

H₂-CO₂ mixtures find wide-ranging applications, including their growing significance as synthetic fuels in the transportation industry, relevance in capture technologies for carbon capture and storage, occurrence in subsurface storage of hydrogen, and hydrogenation of carbon dioxide to form hydrocarbons and alcohols. Here, we focus on the thermodynamic properties of H₂-CO₂ mixtures pertinent to underground hydrogen storage in depleted gas reservoirs. Molecular dynamics simulations are used to compute mutual (Fick) diffusivities for a wide range of pressures (5 to 50 MPa), temperatures (323.15 to 423.15 K), and mixture compositions (hydrogen mole fraction from 0 to 1). At 5 MPa, the computed mutual diffusivities agree within 5% with the kinetic theory of Chapman and Enskog at 423.15 K, albeit exhibiting deviations of up to 25% between 323.15 and 373.15 K. Even at 50 MPa, kinetic theory predictions match computed diffusivities within 15% for mixtures comprising over 80% H₂ due to the ideal-gas-like behavior. In mixtures with higher concentrations of CO₂, the Moggridge correlation emerges as a dependable substitute for the kinetic theory. Specifically, when the CO₂ content reaches 50%, the Moggridge correlation achieves predictions within 10% of the computed Fick diffusivities. Phase equilibria of ternary mixtures involving CO₂-H₂-NaCl were explored using Gibbs Ensemble (GE) simulations with the Continuous Fractional Component Monte Carlo (CFCMC) technique. The computed solubilities of CO₂ and H₂ in NaCl brine increased with the fugacity of the respective component but decreased with NaCl concentration (salting out effect). While the solubility of CO₂ in NaCl brine decreased in the ternary system compared to the binary CO₂-NaCl brine system, the solubility of H₂ in NaCl brine increased less in the ternary system compared to the binary H₂-NaCl brine system. The cooperative effect of H₂-CO₂ enhances the H₂ solubility while suppressing the CO₂ solubility. The water content in the gas phase was found to be intermediate between H₂-NaCl brine and CO₂-NaCl brine systems. Our findings have implications for hydrogen storage and chemical technologies dealing with CO₂-H₂ mixtures, particularly where experimental data are lacking, emphasizing the need for reliable thermodynamic data on H₂-CO₂ mixtures.

Thermodynamic factors for diffusion connect the Fick and Maxwell-Stefan diffusion coefficients used to quantify mass transfer. Activity coefficient models or equations of state can be fitted to experimental or simulation data, from which thermodynamic factors can be obtained by differentiation. The accuracy of thermodynamic factors determined using indirect routes is dictated by the specific choice of an activity coefficient model or an equation of state. The Permuted Widom’s Test Particle Insertion (PWTPI) method developed by Balaji et al. enables direct determination of thermodynamic factors in binary and multicomponent systems. For highly dense systems, for example, typical liquids, it is well known that molecular test insertion methods fail. In this article, we use the Continuous Fractional Component Monte Carlo (CFCMC) method to directly calculate thermodynamic factors by adopting the PWTPI method. The CFCMC method uses fractional molecules whose interactions with their surrounding molecules are modulated by a coupling parameter. Even in highly dense systems, the CFCMC method efficiently handles molecule insertions and removals, overcoming the limitations of the PWTPI method. We show excellent agreement between the results of the PWTPI and CFCMC methods for the calculation of thermodynamic factors in binary systems of Lennard-Jones molecules and ternary systems of Weeks-Chandler-Andersen molecules. The CFCMC method applied to calculate the thermodynamic factors of realistic molecular systems consisting of binary mixtures of carbon dioxide and hydrogen agrees well with the NIST REFPROP database. Our study highlights the effectiveness of the CFCMC method in determining thermodynamic factors for diffusion, even in densely packed systems, using relatively small numbers of molecules.

A microscopic insight into hybrid CH₄ physisorption-hydrate formation in halloysite nanotubes (HNTs) is vital for understanding the solidification storage of natural gas in the HNTs and developing energy storage technology. Herein, large-scale microsecond classical molecular dynamics simulations are conducted to investigate CH₄ storage in the HNTs via the adsorption-hydration hybrid (AHH) method to reveal the effect of gas-water ratio. The simulation results indicate that the HNTs are excellent nanomaterials for CH₄ storage via the adsorption-hydration hybrid method. The CH₄ physisorption and hydrate formation inside and outside of the HNTs are profoundly influenced by the surface properties of the HNTs and the kinetic characteristics of CH₄/H₂O molecules. The outer surfaces of the HNTs exhibit relative hydrophobicity and adsorb CH₄ molecules to form nanobubbles. Moreover, CH₄ molecules adsorbed on the outer surface are tightly trapped between the hydrate solids and the outer surface, inhibiting their kinetic behavior and favoring CH₄ storage via physisorption. The inner surface of the HNTs exhibits extreme hydrophilicity and strongly adsorbs H₂O molecules; thus, CH₄ hydrate can form inside of the HNTs. It is more difficult for CH₄ and H₂O molecules inside of the HNTs to convert into hydrates than for those outside of the HNTs. A moderate gas-water ratio is advantageous for CH₄ physisorption and hydrate formation, whereas excessively high or low gas-water ratios are unfavorable for efficient CH₄ storage. These insights can help to develop an efficient CH₄ solidification storage technology.

Vapor-Liquid Equilibria (VLE) of hydrogen (H₂) and aqueous electrolyte (KOH and NaCl) solutions are central to numerous industrial applications such as alkaline electrolysis and underground hydrogen storage. Continuous fractional component Monte Carlo simulations are performed to compute the VLE of H₂ and aqueous electrolyte solutions at 298-423 K, 10-400 bar, 0-8 mol KOH/kg water, and 0-6 mol NaCl/kg water. The densities and activities of water in aqueous KOH and NaCl solutions are accurately modeled (within 2% deviation from experiments) using the non-polarizable Madrid-2019 Na⁺/Cl⁻ ion force fields for NaCl and the Madrid-Transport K⁺ and Delft Force Field of OH⁻ for KOH, combined with the TIP4P/2005 water force field. A free energy correction (independent of pressure, salt type, and salt molality) is applied to the computed infinite dilution excess chemical potentials of H₂ and water, resulting in accurate predictions (within 5% of experiments) for the solubilities of H₂ in water and the saturated vapor pressures of water for a temperature range of 298-363 K. The compositions of water and H₂ are computed using an iterative scheme from the liquid phase excess chemical potentials and densities, in which the gas phase fugacities are computed using the GERG-2008 equation of state. For the first time, the VLE of H₂ and aqueous KOH/NaCl systems are accurately captured with respect to experiments (i.e., for both the liquid and gas phase compositions) without compromising the liquid phase properties or performing any refitting of force fields.

Continuous Fractional Component Monte Carlo (CFCMC) and molecular dynamics (MD) simulations are performed to calculate the solubilities and self-diffusion coefficients of four light n-alkanes (methane, ethane, propane, and n-butane) in aqueous NaCl solutions as well as the thermodynamic properties of their corresponding hydrate crystals. Correction factors k_ij to the Lorentz-Berthelot combining rules for alkane groups (CH₃) and water are optimized (k_ij = 1.04) by fitting excess chemical potentials to experimental data at 1 bar and 298.15 K. Using these values of k_ij, we calculate the solubilities of the four alkanes in aqueous NaCl solutions with different molalities (0-6) mol/kg at different temperatures (278.15-308.15) K and pressures (1, 100, 200, 300) bar. The diffusion coefficients of the four alkanes in NaCl solutions (0-6) mol/kg are calculated at different temperatures (278.15-308.15) K and 1 bar and corrected for the finite-size effects. The lattice parameters of the corresponding hydrates with different guest molecules are computed using MD simulations at different temperatures (150-290) K and pressures (5-700) MPa. Isothermal compressibilities at 287.15 K and thermal expansion coefficients at 14.5 MPa for the corresponding hydrates are calculated. We present an extensive collection of thermodynamic data related to gas hydrates that contribute to a fundamental understanding of natural gas hydrate science.

We explore the impact of force field parameters and reaction equilibrium on the scaling behavior towards the critical point in reactive binary systems, focusing on NO ₂/N ₂O ₄. This system can be considered as a special single-component system since NO ₂ and N ₂O ₄ are in chemical equilibrium via the chemical reaction 2NO ₂⇌N ₂O ₄. We simplify the system by representing both components as single LJ particles, achieving excellent agreement with densities computed using molecular simulations in which all-atom force fields were used. We investigate the effect of force field parameters (ɛ and σ) on phase behavior and show that the critical exponent β remains constant, which means that intermolecular interactions do not affect the scaling to the critical point when the chemical reaction takes place. We also investigate the sensitivity of the reaction equilibrium constant and show that even small changes in isolated molecule partition functions lead to large differences in chemical equilibria. We show that the critical exponent β is different for systems with different reaction equilibrium constants, so a careful parameterization of β is needed for an accurate computation of critical temperatures of reactive mixtures. We perform a screening of reactive binary mixtures for a wide range of ideal gas reaction equilibrium constants, revealing key insights into the thermodynamic behavior and critical properties. Thereby we facilitate the efficient screening of reactive binary mixtures for various applications. Our results emphasize the importance of accurately parameterizing β and provide valuable insights into the critical scaling behavior of complex reactive systems.

Both CH₄ hydrate accumulation and hydrate-based CO₂ sequestration involve hydrate formation in mixed clay sediments. The development of realistic clay models and a nanoscale understanding of hydrate formation in mixed clay sediments are crucial for energy recovery and carbon sequestration. Here, we propose a novel molecular model of pseudo-hexagonal montmorillonite nanoparticles. The stress-strain curves of tension, compression, and shear of pseudo-hexagonal montmorillonite nanoparticles exhibit linear characteristics, with tension, compression, and shear moduli of ∼435, 410, and 137 GPa, respectively. We perform microsecond molecular dynamics simulations to study CH₄ and CH₄/CO₂ hydrate formation in montmorillonite-illite mixed clay sediments with surface defects. The results indicate that hydrate formation in mixed clay sediments is significantly influenced by the presence of clay defects. CH₄ and CH₄/CO₂ mixed hydrates are challenging to form at the junction between the inside and outside clay defects. CH₄ and CH₄/CO₂ mixed hydrates exhibit a preference for forming outside the clay defects rather than inside the clay defects. Some CH₄ and CO₂ molecules from the inside clay defect migrate to the outside clay defect, thereby promoting CH₄ and CH₄/CO₂ mixed hydrate formation outside the clay defects. This molecular insight advances the development of clay particle models and expands an understanding of natural gas hydrate accumulation and hydrate-based CO₂ sequestration.

Knowledge of the microscopic behavior of CO₂ hydrates in oceanic sediments is crucial to evaluate the efficiency and stability of hydrate-based CO₂ sequestration in oceans. Here, systematic molecular dynamics simulations are executed to investigate the growth and dissociation of CO₂ hydrates, and the mechanical instability of CO₂ hydrate-Illite interface in the brine-urea-Illite system. Simulation results show that the CO₂ hydrate growth is jointly affected by the confined space, Illite surface properties, and presence of urea. Specifically, the interfacial H₂O and the ion layer on the Illite surface hinder the growth of CO₂ hydrate crystals toward Illite surfaces. Urea molecules can bind salt ions and increase CO₂ concentrations in the water, thus kinetically promoting CO₂ hydrate growth. The dissociation of the CO₂ hydrate is affected by Illite surface properties and the CO₂ hydrate structure. CO₂ hydrate starts from the regions where hydrate particles are minimally in contact and extends on both sides. The mechanical tension and compression of the CO₂ hydrate-Illite interface exhibit nonlinear characteristics by changing the hydrogen bonds and the CO₂ hydrate structure. The molecular insight into the microscopic behavior of CO₂ hydrates in the brine-urea-Illite system contributes to a broader understanding of hydrate-based CO₂ sequestration.

Knowledge on the kinetics of gas hydrate dissociation in clay pores at static and dynamic fluid conditions is a fundamental scientific issue for improving gas production efficiency from hydrate deposits using thermal stimulation and depressurization respectively. Here, molecular dynamics simulations were used to investigate poly- and mono-crystalline methane hydrates in Na-montmorillonite clay nanopores. Simulation results show that hydrate dissociation is highly sensitive to temperature and pressure gradients, but their effects differ. Temperature changes increase thermal instability of water and gas molecules, leading to layer-by-layer dissociation from the outer surface. Under flow conditions, laminar flow predominates in nano-pores, and non-Darcy flow occurs due to clay-fluid interactions. Viscous flow disrupts hydrogen bonding at the hydrate surface, enhancing kinetic instability of water. Grain boundaries of polycrystalline hydrates are less stable compared to bulk phases and preferentially decompose, forming new dissociation fronts. This accelerates dissociation compared to monocrystalline hydrates. Fracture occurs at the grain boundaries of polycrystalline hydrate in the fluid, resulting in separate hydrate crystal grains. This fracture process further accelerates hydrate dissociation. In flow systems, methane nanobubbles form in fluid and readily transport with fluid flow. Unlike surface nanobubbles at static conditions, these liquid nanobubbles exhibit mobility. The findings of this study can contribute to a better understanding of the complex phase transition behavior of hydrate in confined environment, and provide theoretical support for improving production control technology.

Hydrogen is a clean-burning fuel that can be converted to other forms. of energy without generating any greenhouse gases. Currently, hydrogen is stored either by compression to high pressure (>700 bar) or cryogenic cooling to liquid form (<23 K). Therefore, it is essential to develop safe, reliable, and energy-efficient storage technology that can store hydrogen at lower pressures and temperatures. In this work, we systematically designed 2902 Mg-alkoxide-functionalized covalent-organic frameworks (COFs) and performed high-throughput (HT) computational screening for hydrogen storage applications at 111, 231, and 296 K. To accurately model the interaction between Mg-alkoxide sites and molecular hydrogen, we performed MP2 calculations to compute the hydrogen binding energy for different types of functionalized models, and the data were subsequently used to fit modified-Morse force field (FF) parameters. Using the developed FF models, we conducted HT grand canonical Monte Carlo (GCMC) simulations to compute hydrogen uptakes for both original and functionalized COFs. The generated data were subsequently used to evaluate the materials’ gravimetric and volumetric storage performance at various temperatures (111, 231, and 296 K). Finally, we developed machine learning (ML) models to predict the hydrogen storage performance of functionalized structures based on the features of the original structures. The developed model showed excellent performance with a mean absolute error (MAE) of 0.061 wt % and 0.456 g/L for predicting the gravimetric and volumetric deliverable capacities, enabling a quick evaluation of structures in a hypothetical COF database. The screening results demonstrated that the Mg-alkoxide functionalization yields greater improvements in volumetric H2 storage capacities for COFs with smaller pores compared to those with larger (mesoporous) pores.

An in-depth review of the available experimental and molecular simulation studies of CO₂ diffusion in H₂O, which is a central property in important industrial and environmental processes, such as carbon capture and storage, enhanced oil recovery, and in the food industry is presented. The cases of both bulk and confined systems are covered. The experimental and molecular simulation data gathered are analyzed, and simple and computationally efficient correlations are devised. These correlations are applicable to conditions from 273 K and 0.1 MPa up to 473 K and 45 MPa. The available experimental data for diffusion coefficients of CO₂ in brines are also collected, and their dependency on temperature, pressure, and salinity is examined in detail. Other engineering models and correlations reported in literature are also presented. The review of the simulation studies focuses on the force field combinations, the data for diffusivities at low and high pressures, finite-size effects, and the correlations developed based on the Molecular Dynamics data. Regarding the confined systems, we review the main methods to measure and compute the diffusivity of confined CO₂ and discuss the main natural and artificial confining media (i.e., smectites, calcites, silica, MOFs, and carbon materials). Detailed discussion is provided regarding the driving force for diffusion of CO₂ and H₂O under confinement, and on the role of effects such as H₂O adsorption on hydrophilic confining media on the diffusivity of CO₂. Finally, an outlook of future research paths for advancing the field of CO₂ diffusivity in H₂O at the bulk phase and in confinement is laid out.

Non-polarizable force fields fail to accurately predict free energies of aqueous electrolytes without compromising the predictive ability for densities and transport properties. A new approach is presented in which (1) TIP4P/2005 water and scaled charge force fields are used to describe the interactions in the liquid phase and (2) an additional Effective Charge Surface (ECS) is used to compute free energies at zero additional computational expense. The ECS is obtained using a single temperature-independent charge scaling parameter per species. Thereby, the chemical potential of water and the free energies of hydration of various aqueous salts (e.g., NaCl and LiCl) are accurately described (deviations less than 5% from experiments), in sharp contrast to calculations where the ECS is omitted (deviations larger than 20%). This approach enables accurate predictions of free energies of aqueous electrolyte solutions using non-polarizable force fields, without compromising liquid-phase properties.

Hydrate dissociation is often accompanied by the formation of nanobubbles. Knowledge of the effects of nanobubbles on hydrate dissociation is essential for understanding the dynamic behavior of the hydrate phase change and improving the gas production efficiency. Here, molecular dynamics simulations were performed to study the methane hydrate dissociation kinetics with and without a pre-existing methane nanobubble. The results show that the hydrate cluster in the liquid phase dissociates layer-by-layer. This process is shown to be independent of the temperature and nanobubble presence at the simulation conditions. Hydrate dissociation does not always lead to nanobubble formation because the supersaturated methane solution can be stable for a long time. A steep methane concentration gradient was observed between the hydrate cluster surface and the methane nanobubble, which can enhance the directional migration of methane and effectively minimize the methane concentration in the liquid phase, thereby increasing the driving force for the hydrate dissociation. Our findings indicate that the presence of a nanobubble near the hydrate surface does not decrease the activation energy of hydrate dissociation, but it can increase the intrinsic decomposition rate. The average hydrate dissociation rate is linearly correlated with the mass flow rate towards the nanobubble. The mass flow rate is determined by the nanobubble size and hydrate-nanobubble distance. Our findings contribute to the fundamental understanding of the dissociation mechanism of gas hydrates in the liquid phase, which is crucial for the design and optimization of efficient gas hydrate production techniques.

Available data from experiments and molecular simulations for the intra-diffusivities of H (Formula presented.) and O (Formula presented.) in H (Formula presented.) O, and for the self-diffusivity of pure H (Formula presented.) O (at pressure and temperature conditions in which the solvent is in the vapour phase) are compared against calculations based on the Chapman-Enskog theory or other semi-empirical/semi-theoretical methods. A novel methodology is proposed to extrapolate the intra-/self-diffusivities data computed from molecular dynamics simulations at low pressures. The extrapolated values are used to further refine the recently-proposed [Tsimpanogiannis et al., J. Chem. Eng. Data, 66, 3226-3244, (2021)], molecular simulation based correlation of intra-/self-diffusivities as a function of pressure and temperature with the solvent being in the vapour phase.

Water transport is one of the major factors responsible for moisture damage in asphalt pavements. To study the thermodynamics and kinetics of water transport in bitumen and to uncover microscale mechanisms of moisture-induced damage, molecular dynamics simulations were performed for up to 600 ns for water–bitumen systems with realistic water contents that varied from 0 to 1.76 wt%. Hydrogen bonding interactions and clustering of water molecules at various combinations of temperature and water content were investigated, and their effects on the self-diffusion coefficient of water and bitumen properties are computed and discussed. It is shown that the saturated water concentration in bitumen is small, especially at low temperatures, and additional water molecules tend to form large water clusters via hydrogen bonding, indicating micro-phase separation of the water and bitumen phases inside the simulation box. Hydrogen bonding and water clustering play a crucial role on the magnitude of the self-diffusion coefficient of water. Physical properties of bitumen that include viscosity and cohesive energy are affected by water. The presence of large water clusters is indicative of how degradation in cohesion is observed on the microscale.