Molecular Dynamics Simulation of Calcium-Silicate-Hydrate for Nano-Engineered Cement Composites—A Review

With the continuous research efforts, sophisticated predictive molecular dynamics (MD) models for C-S-H have been developed, and the application of MD simulation has been expanded from fundamental understanding of C-S-H to nano-engineered cement composites. This paper comprehensively reviewed the current state of MD simulation on calcium-silicate-hydrate (C-S-H) and its diverse applications to nano-engineered cement composites, including carbon-based nanomaterials (i.e., carbon nanotube, graphene, graphene oxide), reinforced cement, cement–polymer nanocomposites (with an application on 3D printing concrete), and chemical additives for improving environmental resistance. In conclusion, the MD method could not only compute but also visualize the nanoscale behaviors of cement hydrates and other ingredients in the cement matrix; thus, fundamental properties of C-S-H structure and its interaction with nanoparticles can be well understood. As a result, the MD enabled us to identify and evaluate the performance of new advanced nano-engineered cement composites.


Introduction
Historically, cementitious materials have been used as one of the most common and popular construction materials. The hydrated cement is composed of nanostructured multiple composite phases that include an amorphous phase, nano-/micro-size crystals, and bound water. Calcium-silicate-hydrate (C-S-H) plays a pivotal role in controlling mechanical, chemical, and transport properties of cement composites and their engineering performance as well [1]. Therefore, in order to improve and customize the macroscopic properties of cementitious materials, understanding of the structure of C-S-H at nano-/micro-level is prerequisite [2] because chemical process at molecular level eventually affects engineering performance of the composites in the bulk scale [1,[3][4][5].
Nano-engineering in cementitious materials encompasses alteration and modification of cement hydrates for enhancing and modifying properties and performance in macro scale [1,6]. In addition, it also deals with characterization and prediction techniques via atomic-or molecular-level modeling for better understanding of how chemical interactions correlate with the macro-level behaviors [7,8]. Via nano-engineering processes such as nano-particles, nano-reinforcements, and chemical admixtures, mechanical and durability properties and degradation processes of cement composites could be effectively controlled and enhanced; thus, novel and smart multi-functions can be incorporated into the cement composites [9]. Fundamental knowledge in the field of the nano-engineered cement composites has been advanced with high-resolution characterization tools [9]; however, there are still many challenges in understanding material behaviors at the atomic/molecular level, solving test repeatability, and realizing cost-benefit issues, etc. [9,10].
Molecular dynamics (MD) simulation is a powerful computational method that models physical movements of atoms and molecules accounting for potential energy in a given position and then numerically computes the atomic/molecular forces based on Newton's classical mechanics [1]. The MD simulation quantitatively determines statistical properties of multi-body systems because it models a large size of the molecular system with longer simulation durations for the required level of accuracy [11][12][13][14]. By this reason, the computational MD has been applied to cement composites to study mechanical behaviors and properties of C-S-H and other applications such as carbon-based nanomaterials, polymer-cement nanocomposites, and chemical treatments for the modification of material's properties, etc. [15][16][17]. There are still needs and opportunities for the MD in cement composites because of diverse molecular structures of cement hydrates and introduction of numerous nano-admixtures into the cement matrix.
This paper provides a comprehensive review on MD analyses on cement composites. Fundamentals of C-S-H including formation, classifications, and molecular modeling were also first reviewed. Thereafter, MD studies on C-S-H structures are reviewed, focusing on C-S-H MD modeling, water dynamics characteristics in pore structures, and mechanical properties of the C-S-H in nanoscale. Lastly, the recent applications of MD to nano-engineered cement composites were reviewed.

Formation and Classification of C-S-H
Tricalcium silicate (alite: C 3 S) and dicalcium silicate (belite: C 2 S) because constitute about 50~70% and 15~30% of Portland Cement (PC) by mass, respectively [18,19]. C 3 S dominates early hydration process of PC, which controls early-age properties of cement and concrete [20]. During that time, a large amount of calcium-silicate-hydrate (C-S-H) gel is formed by the alite hydration, whereas C 2 S hydration governs later composition of C-S-H [21]. A study on the early hydration process of alite shows the reaction begins immediately upon contacting the water. The initial reaction is a complex heterogeneous process that are congruent dissolution and incongruent dissolution with a formation of a silica-rich layer on the surface of C 3 S [22]. This alite-water reaction forms C-S-H and calcium hydroxide (CH) following the relationship by the schematic Equation (1). The hydration of C 2 S occurs similarly to the case of C 3 S and can be expressed by Equation (2). The rate of hydration is generally lower than alite. In addition, the slow process of C 2 S due to different structures of C 3 S results in low concentration of calcium ions; thus, CH crystals are generally large but the C-S-H form is similar [23]. In this reaction, roughly three Ca-O and two Si-O bonds are broken, including six O-H bonds and eventually forms C-S-H gel and portlandite (Ca(OH) 2 ). The energy required for the bond formation of C-S-H and portlandite is 138 kJ/mol less than the energy released due to the bond breaking of the reactants mentioned in the above reaction. As a note, it is known that some of the added C 3 S in PC triggers an accelerated hydration process of C 2 S, as the nucleation of CH could be accelerated and thus the C 2 S dissolution process is also accelerated [22].
It is important to note that the formation, composition, and structure of C-S-H can be influenced even by a slight variation of w/c, mix proportion, humidity, the degree of hydration, and curing temperature [21]. Those factors contribute the ratio of calcium to silicon (Ca/Si) [24][25][26] and the variability of the C-S-H structure comes from the varied Ca/Si ratio ranged from 0.7-2.3 [10,27].
The Ca/Si ratio may increase under higher hydration and better curing conditions, but the use of supplementary cementitious materials (SCMs) may decrease the ratio. Lower Ca/Si ratios may increase the length of silicate chains and the interlayer distance, resulting in changes of C-S-H structures [26]. Generally, C-S-H was firstly classified into two categories as tobermorite-like C-S-H(I) for Ca/Si < 1.5 and jennite-like C-S-H (II) for Ca/Si > 1.5 [28]. The C-S-H(I) was additionally divided into C-S-H(α) for Ca/Si < 1.0 and C-S-H (β) for 1 < Ca/Si < 1.5 [29]. Nowadays, depending on this ratio, C-S-H can be depicted as tobermorite, jennite, metajennite, and many other formations shown in Table 1 in the next section.
It was observed that the type of C-S-H largely depends on hydration times; for instance, low-density C-S-H is rapidly formed at the surface of C 3 S in the early stage of hydration [21,30]. As the hydration process continues, a group of C-S-H having higher density is subsequently formed, which mainly contributes to the development of strength. Powers and Brownyard firstly denoted the physical characterization of C-S-H at multiscale in the 1940s and 1950s [31,32]. They suggested inside and outside hydrated products at the cement clinker particles (e.g., C 3 S and C 2 S), which were denoted as inner (Ip) and outer (Op) of C-S-H by Taplin [33]. Diamond (1976) observed four different types of C-S-H morphologies including three types of Op-C-S-H and one type of Ip-C-S-H [34]. They are classified as Type I: fibrous particles look like partly rolled sheets, Type II: reticular network/interlocking structure but not in C3S or C2S pastes, Type III: "relatively nondescript" or "equant grain morphology", and Type IV: inner product in older pastes. Richardson and Groves investigated the formation and classification of C-S-H with transmission electron microscopy (TEM) and then adopted the terminologies of OP and IP [35,36]. Later, Jennings and Tennis have classified the C-S-H in density perspective by modeling colloidal particles that pack into two separate arrangements: high density (HD) and low density (LD) products [37,38]. They showed that the LD products are formed during the first or second day of hydration but those are more deformed during the drying stage due to its open structure than HD products [39]. Nicoleau described the formation of C-S-H being the outcome of cohesive forces during the hydration process, and its physical transformation could trigger coagulation and eventually form structures [40]. Bonding characteristics of the C-S-H formation were also investigated. The high-density C-S-H is a layered structure composed of calcium silicate sheets randomly connected by strong iono-covalent bonds [41].

Structural Models of C-S-H
Modeling the C-S-H structures needs to address the variation of chemical composition and its distribution in a cement matrix because engineering performance of cement hydrates is highly affected by structures and chemical compositions of C-S-H determined by w/c ratio and curing conditions. Therefore, specifically, the model should be able to account for: (1) Ca/Si ratio of C-S-H ranged from 0.7 to 2.3 [42], (2) compositional heterogeneity of C-S-H, and (3) age of hydration process [43,44]. In an early stage of C-S-H study, Taylor (1986) announced that C-S-H appears to be a disordered layer structure composed of structurally imperfect jennite (Ca 9 Si 6 O 32 H 22 ) and others similarly related to 14 Å tobermorite (Ca 5 Si 6 O 26 H 18 ) [45]. Tobermorite and jennite are layered structures that are calcium sheets flanked on both side by linear silicate chains called "dreierketten chain" [46,47]. Figure 1 shows a schematic diagram of C-S-H crystalline structure. The chains are composed in a repeating manner at three SiO 4 tetrahedra intervals. Two adjacent tetrahedrons are coordinated to the Ca 2+ s of the layer, and the third tetrahedron bridges two consecutive dimers, which is called bridging tetrahedra. Two oxygens from non-bridging tetrahedra are coordinated with Ca 2+ s in tobermorite [48,49], whereas one oxygen from the non-bridging tetrahedra is coordinated to the Ca 2+ on the sheet and the other oxygen is provided by hydroxide ions [50,51]. Taylor [45]. In the early research period of Taylor, he proposed that both tobermorite-like and jennite-like structures may simultaneously exist at the early age of hydration, and the distribution could finally merge into a jennite-type structure having an intermediate composition by interacting the regions of low and high Ca/Si ratio as time passes. Thereafter, however, he suggested that the tobermorite-like and jennite-like regions could be hardly determined and could merge into each other within individual layers [52]. Unlike Taylor's model, Cong and Kirkpatrick reported that the jennite-like C-S-H is rarely formed based on experimental studies, and they proposed tobermorite-like C-S-H model when Ca/Si < 1.5 [52,53]. There is a difference between the two models in that Tayler focused on the case of Ca/Si > 1.5 whereas Cong and Kirkpatrick focused on Ca/Si < 1.5. Nonetheless, both models involve quite disordered structures. Nonat and Lecoq later elaborated the model associated with the evolution of C-S-H in accordance with the Ca/Si ratio of 0.66~2.0 [29]. In their model, the main differences from the Taylor's and Cong and Kirkpatrick's models are that C-S-H structures are not necessarily to be disordered structure and the layers do not contain jennite-like regions. The X-ray diffraction (XRD) analysis of C-S-H showed similar patterns with the tobermorite over the wide range of Ca/Si ratios.    Selected studies of C-S-H structures and the corresponding models proposed by various researchers are summarized in Table 1. It is noticed that the layer structure is the most widely adopted model in which 11 Å and 14 Å tobermorite and jennite structures are used to model the structure of C-S-H. Accordingly, these C-S-H structures have been mainly adapted for the MD simulation of C-S-H to be presented in following sections.

MD Models of C-S-H
Since the 2000s, MD simulation of C-S-H has been explored and helped to investigate the structures, behaviors, and diverse properties of hydrated cement matrix at molecular level [2]. MD theoretically computes a set of molecular orbital phase spaces where each molecule follows Newton's laws of motion [68]. In other words, the MD simulation method calculate a typical trajectory of the molecular system based on the fact that the total energy is constant. As this is a deterministic algorithm, initial conditions use random velocities, and it is usually known as an initial configuration from Monte Carlo simulations. Accordingly, understanding properties and behaviors of hydrated cement in small scale will help understand the macroscopic behaviors. Pellenq et al. (2009) first suggested the C-S-H molecular model by MD method [69]. The model was developed based on the mean value of 1.7 for Ca/Si ratio repeating 14 Å tobermorite crystalline structure, which was suggested by Richardson [36]. Figure 3 shows the TEM image of C-S-H clusters [70] and the corresponding MD model by Pellenq et al. [69] The blue and white spheres are oxygen and hydrogen atoms of water molecules, respectively. The green and gray spheres are inter-and intra-layer calcium ions, respectively. Yellow and red sticks represent silicon and oxygen atoms in silica tetrahedral. The overall chemical composition of the computational model of hydrated C-S-H is (CaO) 1.65 (SiO 2 )(H 2 O) 1. 75 . The C-S-H model is based on a bottom-up atomistic simulation approach (i.e., atomic to composite scale) which considers only the chemical specificity of the system as a major constraint. This C-S-H pattern that involves the interaction of CaO, SiO 2 and H 2 O molecules allows the distributions of short silica chains such as monomers, dimers and pentamers, and thus, it could provide more realistic Ca/Si ratio and density values calculated by grand canonical Monte Carlo simulation of moisture adsorption at 300 K. . The blue and white spheres are oxygen and hydrogen atoms of water molecules, respectively; the green and gray spheres are inter-and intra-layer calcium ions, respectively; yellow and red sticks represent silicon and oxygen atoms in silica tetrahedral [69]. Table 2 summarizes the selected studies of MD simulation on C-S-H structures and experimental methods for validation. The table also summarizes the applications of MD analyses to C-S-H, looking into the structures of C-S-H, diffusion coefficients of water, density, elastic modulus, tensile/compressive/shear strength, chloride diffusion, and the most stable structures corresponding to minimum energy. One of the most important parameters in the simulation of C-S-H is Ca/Si ratio. Researchers used different Ca/Si ratio to simulate the closest structure to the C-S-H gel, which ranged from 0.66 to 2.0, however, the typical ratio of Ca/Si in hydrated cement composite is considered as 1.7. 9 Å, 11 Å, and 14 Å tobermorite layered C-S-H structures have been mainly modeled. Tobermorite structures could have several modifications under different levels of hydration and degree of cross-linking of silicate chains between two adjacent layers, and this could form different levels of separation of the interlayer. According to the spacing due to the separation, the structures are named by the 9 Å, 11 Å, and 14 Å tobermorite. On the other hand, only a few studies attempted to model amorphous C-S-H clusters [68]. MD simulation software (i.e., GULP, LAMMPS, and TREMOLO) involves NVT (Number of atoms, Volume and Temperature are constant) or NPT (Number of atoms, Pressure and Temperature constant). Appropriate potentials are used to simulate the real-life interactions (i.e., 2-body or 3-body) among the molecules. Reactive force field ReaxFF is useful for capturing bond breaking and reforming by considering fracture [71], electrostatic equilibrium, and water dissociation [72], and thus, it is widely adapted for simulating systems having dynamically changing bond topologies [73]. Whereas the ReaxFF reproduce the energy evolution due to the bonding and dissociating atoms, empirical force fields such as CLAYFF and CSHFF predefine the connections between atoms at a fixed state [74]. CLAYFF force field combines coulombic, van der Waals, bonded stretch and angle bend potentials altogether. On the other hand, CSHFF is more popular for investigating water dynamics and cohesion between C-S-H grains [75]. Accordingly, the empirical force fields have been employed for studying the nano-engineered cementitious materials and its applications such as nanomaterial-reinforced cement and polymer modified cement, which to be presented in the following sections.

Water Dynamics in C-S-H
The pore system having a wide range of size distribution from nano-to micro-meters significantly affect physical and chemical properties of hydrated cement composites such as strength, shrinkage, creep, and chemical reactivity [77]. Since 1970s, extensive studies about water interactions which include water bound into cement paste, "glassy water" interacted in gel pores, and unbound water in larger capillary pores in the near surface of cement hydrates have been made; for instance, dynamics of unbound water, physically bound and chemically bound water under variable conditions [88][89][90][91][92][93][94][95][96].  simulated C-S-H nanostructure under varied moisture conditions (from dry to saturated conditions) based on twelve (12) C-S-H gel structures with water/Ca ratios from 0 to 0.95. Specifically, they investigated the chemical bonds in the C-S-H structure using "CSHFF force field" [97][98][99][100]. Figure 4 illustrates the morphology change of silicate chains in C-S-H models under different water contents. In the dry state, the bridging silicate tetrahedron is connected to the surrounding monomers and adjacent tetrahedrons. When the water/Ca ratio is 0.3, two dimers and one monomer can be linked together and form a longer silicate chain. In the saturated state, the dimer structure develops orderly, and the water molecules block the silica chains from connecting to each other. The behavior of confined water in nanometer scale greatly differs than that of bulk-scale water. The viscosity of the nanoconfined water (≤ 1 nm separation) is seven orders of magnitude higher than the viscosity of bulk water at room temperature [102]. Solid surfaces can perturb the confined water up to several molecular diameters [77]. Kalinichev et al. (2007) and Youssef et al. (2011) simulated water dynamics in the pore systems of cement matrix. Kalinichev's MD analyses illustrate that (1) water molecules develop a three-dimensional H-bond network and (2) the structure of the H-bond network can be influenced by substrate structures (e.g., C-S-H structures and/or charge distribution) [77]. In addition, Youssef's research demonstrates that nano-pore space has hydrophilic characteristic due to non-bridging oxygen atoms on the disordered silicate chains acting as H-bond acceptor sites [103]. In their studies, the models were intentionally designed to simulate the solid-fluid interaction that is a strong attractive forces between water molecules and C-S-H chains, and used a "CLAYFF forced field" method to generate energetics between the solid-fluid interfaces that do not require the prior definition of chemical bonds. Since the disordered systems such as hydrous minerals, which is composed of a large number of atoms, could be kept relatively small and simple by modeling the number of interactions in the absence of defined chemical bonds, the O-H of H 2 O, the OH-group on the solid surface, and the bond at the aqueous oxyanion were only defined in the CLAYFF [77]. Moreover, they found that the diffusion coefficients of water molecules are much higher near the surface of solid than in the confined channels. However, the both values are much lower than the diffusivities of bulk water, which were validated using NMR technique. Consequently, it was concluded that the water at the surface of tobermorite and jennite is highly structured with reflecting the structure, composition, and charge distribution of the underlying C-S-H substrate. In these ways, MD methods offer significant potential for modeling fluid-solid interface and their interaction in cement systems.

Nanoscale Mechanical Properties and Performance of C-S-H
Understanding properties of C-S-H can help predict mechanical and chemical properties of hydrated cement materials, thus findings can be reflected on the design of cement/concrete composites. Manzano et al. (2007) [104] simulated different crystalline C-S-H models that have diverse Ca/Si ratios, and subsequently computed diverse elastic properties such as bulk (K), shear (G), and Young's Modulus (E). The results showed that the values of modulus slightly decreased when Ca/Si ratio of C-S-H increased and also when more water content was added. Moreover, the study also showed that mechanical properties of C-S-H with dimer or pentamer silicate chains were lower than those with infinite silicate chains. MD studies computing elastic properties of C-S-H have been widely performed under diverse force field types. It is shown that CLAYFF generally estimates higher modulus values than other types such as COMPASSFF, COMPASSIIFF, or Universal, but estimates Poisson's ratio in the middle range [105]. Table 3 summarizes the computed elastic properties of diverse C-S-H structures; however, nanoindentation tests show much lower values. The test results report bulk modulus (K) of 15-18 GPa, shear modulus (G) of 9.7 GPa, and Young's modulus (E) of 18-30 GPa [65,[106][107][108]. In other words, the computed moduli of C-S-H gel overestimates about 3-5 times larger than those of experimental test. Table 3. Elastic properties of C-S-H gels computed by MD simulations.  [78] studied tensile and compressive strength of the C-S-H structure from the uniaxial stress-strain data computed by MD. The results showed that the strength values of C-S-H in nanoscale is about three times of maximum compressive and tensile strength of hydrated cement at bulk scale. In their MD simulation, the tensile strength appeared 23% of the compressive strength. This study also addressed that attractive electrostatic forces and silicate bond (O-Si-O), which resist tensile stresses greater than 600 MPa, play a crucial role in determining the strength of hydrated cement, and the breakage of silicate chains is the primary cause of reduced tensile strength of C-S-H structure due to less O-Si-O bonding. The application area of MD simulation methods was also expanded as such from the conventional analysis of static performance and mechanical properties of C-S-H phases to investigate dynamic behavior characteristics such as propagation of stress waves under shock loading. For instance, the pattern of crack propagation may change as the strain rate increases due to higher shock wave; a higher stress level is required for the failure of specimen, and the energy absorption capacity increases [115][116][117]. Lin et al. (2018) investigated the propagation of stress waves in the C-S-H under the compressive shock loading by selecting the particle velocity of 0.2-3.0 km/s for the wave generation using MD simulations [118]. The elastic limit was calculated as 7.5 GPa, and the results also presented that 0.5 km/s is the critical particle velocity to transit the C-S-H structure from an elastic state to a plastic state. If the particle velocity is larger than 2.0 km/s, the C-S-H structure is collapsed due to the shock wave. Figure 5 illustrates the dynamic-deformation process of C-S-H structures under different particle velocities simulated by MD method. The black-dotted box in the figure indicates C-S-H layer, and the plastic deformation is related to the densification of the C-S-H layers. Figure 5b shows the C-S-H layer begins to compact, indicating that the system transits from elastic to plastic at the particle velocity of 0.6 km/s. Figure 5d demonstrates that the system lost the layered C-S-H morphology at the particle velocity of 3.0 km/s, and it is considered as structure collapse. Crack propagation mechanism in C-S-H gel is also studied via the MD method.  reported that the layered C-S-H gel shows dual types of crack growth natures under loading conditions [100]. In the x-y plane, since the ionic-covalent bonds of Si-O and Ca-O are stable, cracks are hardly coalesced in this direction, which may slow down the crack propagation and thus lead to more ductile characteristic. On the other hand, in the z-direction, cracks could be propagated in the interlayer region more frequently because of higher potential to break the H-bonds network under loading, and this could be exhibited as a brittle behavior in macroscale. They also simulated the development of crack under the tensile loading process in which a void is present at the center of C-S-H. The void could be a natural pore having diameter of 0-50 Å generated within hydration process rather than defects. Figure 6 demonstrates the crack evolution of C-S-H gel under the tensile load along the x-direction simulated through MD method. Figure 6a-f represent the strain at 0, 0.08, 0.16, 0.2, 0.3, and 0.4 Å/Å, respectively. As shown in Figure 6c, the crack is triggered from the boundary of central void by dissociating Si-O and Ca-O bonds as strain increases from 0.08 to 0.16. Consequently, it could be identified in Figure 6e that the cracks coalesce with small cracks and propagates to the direction of around ±45 • from the x-axis. This MD simulation provides considerable insight regarding the reduced plastic deformation during the cracking process of C-S-H gel. First of all, water molecules hinder the interflow of Ca and Si between adjacent calcium silicate sheets, thus weakens the plastic deformation of the C-S-H gel. Next, the water molecules instantly spread out to bond to the non-bridging oxygen site in the damaged region, and thus obstruct bridging by the interlayer calcium, resulting in the plastic deformation due to the fact that Ca, Si, and O diffusion may be reduced [100].

Carbon-Based Nanomaterials-Reinforced Cement Composite
It has been reported that carbon-based nanomaterials (e.g., carbon nanotubes/nanofibers (CNTs/CNFs) and graphene/graphene oxide (GO)) could enhance mechanical properties such as Young's modulus (E) and tensile strength (f t ); for instance, E = 10-950 GPa and f t = 11-150 GPa for CNT, E = 25-600 GPa and f t = 7-30 GPa for CNF, E = 1.5-2.8 TPa and f t = 130-195 GPa for pristine graphene, and E = 180-230 GPa and f t = 0.11-24 GPa for GO [119][120][121][122][123][124][125][126][127]. Cwirzen et al. indicated that adding CNTs to the cement paste improved the workability and increased the compressive strength up to 50% [128]. In addition, pristine CNTs could improve the flexural and compressive strength of the cement composites about 10-20% [129]. One of the reasons for the improvement of mechanical strengths is a nucleation effect of CNT allowing the growth of C-S-H gel in the early hydration process reported by Makar and Chan [130] (Figure 7a). They confirmed the CNT's nucleation effect by measuring the decrement of Ca(OH) 2 during the hydration process for the CNT-containing cement paste. A "bridging effect" of CNTs in the cement matrix is considered as another mechanism for the strength enhancement, which could increase the tensile strength of C-S-H as they can bridge the micro and nanocracks in the C-S-H (Figure 7b) [131][132][133], and hence can be termed as "bridging material" hereafter. Numerical investigation of Eftekhari et al. [134,135] presented that the addition of CNTs to the cement mixtures could significantly increase the fracture energy and tensile strength of the CNT-reinforced cement and delay the crack propagation. Manzur [136] addressed that the strength improvement could be adjusted by the CNT's specific surface area based on experimental results that the compressive strength of CNT-cement composite increases at the larger specific surface area of CNT, and this could be supported by that COOH or C-OH groups on the surface of CNT allow a possible interaction with C-S-H gel [137]. Graphene oxide (GO) has been considered as a viable alternative for reinforcing cementitious matrices with better economic advantages than CNT [138][139][140]. The interfacial strength of pristine graphene with -OH, -COOH and -NH 2 functionalized groups in the cement composites have been studied by AlKhateb et al., and the results presented that the electrostatic forces of the functional groups play a roles in determining the strength of GO-reinforced cement composites [141]. In addition, it is reported that the functionalized surface of GO may act as a seeding material accelerating the hydration process, thus contributing to increase the amount of C-S-H gel in the matrix [142][143][144][145][146]. Consequently, previous studies have shown the positive effects of GO on mechanical properties of the cementitious composites. In summary, for the cement paste and mortar, 0.02-0.08% of GO by weight of cement increased the compressive strength by 29-46% and flexural strength by 26-61% [147][148][149][150][151][152][153][154][155]. Moreover, Lv et al. performed advanced mechanical testing and reported the GO nanosheets significantly affect to reduce brittleness and enhance toughness of the cement based materials [156]. Effectivereinforcement potential was also investigated in nanoscale showing mitigation of crack propagation and lead to less sudden failure [157][158][159], and enhancement of Young's modulus of the GO-reinforced cement composites was evaluated [160,161]. The accelerated C-S-H formation due to GOs may also reduce porosity and hinder the ingress of water and chloride ions, which subsequently improve its durability [122,162]. For the unhardened fresh state of the GO-mixed cement mixtures, it is generally known that the GOs noticeably reduce the fluidity and viscosity of cement paste and mortar mixtures, and this is understood by the fast early cement hydration process due to the functionalized surface of GO [121,163].
In recent, a low-cost GO, edge oxidized graphene oxide (EOGO), has been introduced [164] and relevant investigations on both rheological and mechanical aspects (including long-term fatigue behavior) of GO-cement and concrete composites have been performed; the result demonstrated that the EOGO could improve not only static compressive/flexural strength but also flexural fatigue characteristics (i.e., fatigue life and ductility) [165][166][167]. The results also indicated the plastic strain of the GO-reinforced cement mortar increased under flexural cyclic loading condition, and this was interpreted as the EOGOs helped prevent or delay internal damage propagation in the material.
In addition, EOGO was also mixed with steel fibers in cement concrete mixtures, and the results showed the 0.1% GO by weight of cement improve the energy absorption capability of about 30% in the range of elastic behavior under the static flexural loading condition. It is assumed that the GOs extend the range of elasticity and "interact" with steel fibers in cement paste [168,169].

MD Simulation on CNT/Graphene/GO-Reinforced C-S-H
Eftekhari and Mohammadi (2016) studied the mechanical properties of CNT-reinforced C-S-H by MD simulation (see Figure 8) [170]. The MD model is firstly processed by that a hollow cylinder was placed in the C-S-H structure and expanded up to the CNT's diameter using "fix indent cylinder" function of LAMMPS. This pushes C-S-H atoms back without the breakage of their bonds. CNT was embedded in the hollow hole, and the C-S-H and CNT was then linked by Lennard-Jones (L-J) potential. After the system equilibrium, the interlayer space between C-S-H and carbon atoms by the L-J interaction was set as about 2.7 A. The simulation results confirm that the CNT embedment enhances diverse mechanical properties of C-S-H. The tensile strength of the CNT-reinforced C-S-H is substantially improved along the direction of CNTs. Especially the CNT increased the tensile strength along the Z direction (perpendicular to the silicate layer) up to 6 GPa. The results also demonstrated that the CNTs play a role in efficiently bridging ("bridge effect") cracks of C-S-H (see Figure 8b). In addition, the pullout simulation of CNT estimated the bilinear force-displacement response model, which can provide a viable mechanism to be used for the understanding of crack propagation and bridging effects at macro-scale. The study of MD simulation on graphene-cement composite was pioneered by Sanchez and Zhang [171], investigating the molecular-scale energetic, structural, and dynamic properties of the interface between the surface-functionalized graphitic structures and C-S-H. In this study, six types of functional groups on the surface of carbon nanosheet were considered: hydroxyl (OH), carboxyl (COOH), carboxylate (COO−), carbonyl (CO), and amine (NH 2 ) groups. As a C-S-H model, 9 Å tobermorite structure was used. The simulation results demonstrate that electrostatic forces play a pivotal role in interfacial interactions between the functionalized carbon sheet and C-S-H (i.e., 9 Å tobermorite). It was also identified that the polarity of functionalized groups could be used as an indicator of affinity for C-S-H. Later on, Alkhateb et al. investigated the interfacial strength between C-S-H and graphene nanoplates functionalized by different chemical groups, and the results presented that the functional groups such as OH, NH 2 , and COOH on the graphene nanoparticles enhance the interfacial strength of 13.5 GPa, 6.1 GPa, and 11.8 GPa, respectively, depending on the electrostatic forces of the functional groups, compared to that of the unfunctionalized graphene showing 1.2 GPa [141]. Another study investigating tensile and shear properties of 14 Å tobermorite reinforced with a pure graphene sheet (with no functional groups) [172] indicates that the graphene sheet contributes to the improvement of XY-plane tensile strength of 180-360% and shear strength of 90-225% and also increases stiffness and toughness of the graphene-reinforced tobermorite structures. In the meantime, water interface enhances intermolecular friction forces between graphene and C-S-H, resulting in improved tensile strength and toughness but lower shear toughness than that of the dried surface condition.
MDs has been used to investigate mechanical behavior of GO-reinforced cement composites, particularly looking into the interaction between the GOs and C-S-H phases. Hou et al. [173] investigated the effects of GO on the microstructure of hydrated cement and the corresponding mechanism. The test results showed that the GO increased flexural strength about 11.62%. The simulation studies indicated GOs enable greater hydration and cause nano-filling and crack-bringing effects (see Figure 9a). The functional hydroxyl groups in GOs provide non-bridging oxygen (NBO) sites that accept three types of hydrogen-bonds (i.e., C-O-OH-Si, C-OH-C-OH, and C-OH-OH-Ca) of interlayer water molecules in the C-S-H; thus, this may change polarity of GO surface and enhances the bonding with neighboring phases (see Figure 9b). Furthermore, Ca 2+ and Al 3+ ions near the surface of C-S-H help bridge the GO-silicate chains, thus a length of the silicate chain is increased (see Figure 9c,d). As a result, the defective GO structure is healed. It was also addressed that the aluminate-silicate chains, calcium ions, and functional hydroxyl groups in the cement matrix could establish stabilized connections between C-S-H and GOs. Kai et al., (2019) investigated chemical reactions, mechanical behaviors, and interfacial sliding of C-S-H and GO via MD simulation [174]. The simulation showed that chemical reactions such as turning epoxides into carbonyls and hydroxyls and deprotonating or detaching hydroxyls occur at the interface between GO and C-S-H matrices, resulting in high interfacial interaction energy and mechanical interlocking. Subsequently, mechanical test simulations showed that Young's modulus and strength of C-S-H are enhanced by 31.6-52.6% and 17.5-23.3%, respectively. Hou et al., (2018) explained more in-depth strengthening mechanisms of C-S-H and GO in molecular level [175]. This study investigated the structure, reactivity and interfacial bonding of C-S-H and GO, and the results showed that most COOH groups in GO nanosheet were de-protonated to COOgroups, which forms stable COO-Ca bonds with neighboring Ca ions. The de-protonated COO-groups could also generate H bonds from Si-OH in the C-S-H gel, which improves the interfacial connection. Moreover, uniaxial tensile testing simulation demonstrated that C-S-H reinforced with GO-COOH and GO-OH had increased interfacial cohesive strength and ductility.
In addition, Yang et al. reported the decrement of failure strength of C-S-H and GO composite under moisture condition [176]. The study investigated the effect of water molecules on the bonding in the interlayer region between C-S-H and GO nanosheet functionalized by deprotonated carboxyl (COO), hydroxyl (COH), epoxy (COC) groups. The results indicate that the interlayer bonds are weakened by water molecules as the moisture content increases, resulting in decreased tensile strength of the GO/C-S-H composite.

Cement-Polymer Nanocomposite
Polymers are widely used in cementitious materials to alter and/or improve the structures of cement matrix; thus improve mechanical and durability performance of the composites under severe environment [177][178][179][180]; thus, MDs have been used to understand chemical interactions between diverse types of polymers and cement hydrates. Hou et al., (2019) investigated the interfacial interaction between C-S-H and polymers, looking into structure, dynamics, energetics, and mechanical properties of the cement-polymer nanocomposite [181]. Polyethylene glycol (PEG), polyvinyl alcohol (PVA), and polyacrylic acid (PAA) were selected and simulated, and mechanical properties were computed for those polymer-C-S-H nanocomposites. In the MD model (see Figure 10a,b), H op , H cp , and H sp denote hydrogen atoms in hydroxyl groups, connected to carbon atoms, and in -COOH carboxyl groups, respectively. O hp , O fp , O dp , and O sp mean oxygen atoms in hydroxyl groups, in -CH 2 -, double-bonded in -COOH carboxyl groups, and single-bonded in the carboxyl groups, respectively. The simulation results show the Ca 2+ near the surface of C-S-H connects the functional groups in the polymers to oxygen in the silicate chains by forming O s -Ca-O p bond (where, O s and O p denote oxygen in silicate and polymer, respectively). In addition to the ionic bonding, functional groups in the polymers (i.e., bridging oxygen (C-O-C) in the PEG, hydroxyl (C-OH) in the PVA and carboxyl groups (-COOH) in the PAA) supply oxygen sites for the formation of H-bonds, which also enhances the connectivity of the interlayer region in C-S-H. Moreover, the computation results show that introducing polymers to C-S-H improves the H-bonds in the interface and heals defective silicate chains, resulting in delaying the propagation of crack under loading and ultimately enhancing cohesive strength and ductility of the C-S-H gels. Table 4 summarizes the results of MD simulations for the scenarios of PEG, PVA, and PAA. The PAA the most largely enhances the Young's modulus, tensile strength, and failure strain of C-S-H to 22.3%, 19.2%, and 66.7%, respectively, and the PVA and PEG followed. This result corresponds to the rank of interfacial binding energy (E), which follows the order of E(PAA) > E(PVA) > E(PEG).   Table 4. Mechanical properties of polymer-C-S-H nanocomposites. Reproduced from [181], with permission from Elsevier, 2019. Polymer-modified cement composite significantly improves adhesion in the interface between layers, which is a good alternative for 3D printing concrete industry. The benefits of 3D cement/concrete printing include energy-saving and material-efficient construction methods. Wang et al., (2020) used the MD and examined the interfacial mechanical properties of polymer-modified mortar layers that is used for 3D printing [182]. In this study, epoxy resin and chloroprene latex were adopted for enhancing the interlayer bonding. Direct tensile and shear strengths of the polymer-modified cement paste were evaluated by both experimental and MD numerical methods. The results show that weak interlaminar bonding derived from water molecules is offset and can be overcome by the electrostatic interaction (Coulomb force) between the epoxy resin and Ca 2+ from C-S-H. Consequently, enhancement of 219.55% and 201.41% of tensile and shear strengths were reported, respectively. Hosseini et al., (2019) introduced a new polymer-sand composites for enhancing the interface bonding of 3D printed cement layers [183]. They synthesized sulfur-black carbon (SBC) polymer as a binder and then mixed sulfur-black carbon-sand (SBCS) mortar for the interlayer bonding. This polymer-cement composite was evaluated by both experimental and MD methods. The results demonstrated that the CSH-CSH is mainly bonded by Van der Waals forces in the interlayer region but the CSH-SBC (mortar-SBCS-mortar) is bonded by electrostatic forces between Ca 2+ from C-S-H and the negatively charged SBC polymer. This different chemical bonding mechanism explains that the interlayer bonding strength of CSH-SBC is enhanced by more than 100% than that of CSH-CSH. A graphical summary of experimental and MD simulation studies performed by Hosseini et al. is presented in Figure 11. It is important to note that the MD method may help find out a practical feasible material solution to advance 3D concrete printing technology.

Chloride Ion Binding on Cement Hydrates
Chloride intrusion in concrete material is one of the primary deterioration mechanisms that result in corrosion of reinforcement steels in concrete infrastructure. If not well attended and maintained, it may lead to catastrophic failure. C-S-H and CH possess a diverse microstructural array which allows to compose a favorable pore networks for chloride ions to intrude concrete structures [18,19], and surface physico-chemical characteristics of these cement hydrates could potentially interact with chloride ions [79]. Pan et al., (2010) [79] simulated the chloride transport phenomenon through portlandite, tobermorite, and jennite channels filled with water molecules. The simulation found that the adsorption forces from solid surface could significantly reduce the movement rate of chloride ions in the near-surface domains. It was concluded that although portlandite strongly adsorbed the chloride ions through H-bonding, the C-S-H phases did not show noteworthy binding phenomena of the chloride ions due to lower oxygen atoms. Meanwhile, Zhou et al., (2018) investigated the chloride-ion adsorption capability of the C-S-H phases at variable Ca/Si ratios [184]. They reported that higher Ca/Si condition can provide more favorable environment for binding more chloride ions. Two potential mechanisms are addressed; (1) due to the H-bond network with -OH groups, water molecules may accumulate near the surface of C-S-H, which is beneficial for the Cl − adsorption, and (2) Ca 2+ ions are stably present in the interface layer, which can also promote the Cl − adsorption. Figure 12a shows experimental result of Cl − adsorption under different Ca/Si ratio of C-S-H, and Figure 12b illustrates the MD simulation of tobermorite model with aqueous solutions constituted of water molecules, chloride and calcium ions developed in the study [184]. The number of water molecules was determined based on a standard density of water (1 g/cm 3 ), and the concentration of chloride was selected as 0.44 mol/L that is the upper limit for the pore solution. This MD simulation results may be applied to provide a foundation to develop chemical additives and processing technology to slower the transport rate and adsorption of chlorides in reinforced concrete structures exposed to the marine environment by increasing Ca/Si ratio of the concrete mixture for the improvement of design.

Conclusions and Recommendations
This paper presents a comprehensive review on MD modeling on C-S-H and its applications in the nano-engineered cement composites. The MD is a promising tool for exploring new nano-additives because it can investigate and visualize fundamentals in a molecular level, which is difficult to be characterized by experimental methods. For instance, the MD can visually show the "bridging" effect of carbon-based nanomaterial (CNT, GO) and compute the bonding with C-S-H while experimental methods provide indirectly supporting results. The negative effect of water on the interlayer bonding between GO and C-S-H could be also investigated by the MD. For the cement-polymer nanocomposite, the MD enabled to identify the types of linking mechanisms between C-S-H and polymers; thus, computing the interfacial binding energy of polymers (i.e., PAA > PVA > PEG) was possible, which is matched well with the experimental result. The MD technique was even applied to the area of 3D-printed concrete research, and the simulation results well represented the experimental results. Lastly, the MD simulation could be effectively used in designing and assessing behaviors of chemical additives and treatment agents for the durability enhancement of cement composites under environmentally harsh conditions such as coastal areas. As nanomaterials and nanotechnology is becoming more popular and common in cement/concrete industry, the emergence of more efficient computational methods and tools capable of simulating at nanoscale are inevitable; therefore, the MD method will be an effective tool.