CHARMM-GUI Nanomaterial Modeler for Modeling and Simulation of Nanomaterial Systems
Molecular modeling and simulation are invaluable tools for nanoscience that predict mechanical, physicochemical, and thermodynamic properties of nanomaterials and provide molecular-level insight into underlying mechanisms. However, building nanomaterial-containing systems remains challenging due to the lack of reliable and integrated cyberinfrastructures. Here, we present Nanomaterial Modeler in CHARMM-GUI, a web-based cyberinfrastructure that provides an automated process to generate various nanomaterial models, associated topology, and configuration files to perform state-of-the-art molecular dynamics simulations using most simulation packages. The nanomaterial models are based on the interface force field (IFF), one of the most reliable FFs. The transferability of nanomaterial models among the simulation programs was assessed by single-point energy calculations, which yielded 0.01% relative absolute energy differences for various surface models and equilibrium nanoparticle shapes. Three widely-used Lennard-Jones (LJ) cut-off methods are employed to evaluate the compatibility of nanomaterial models with respect to conventional biomolecular FFs: simple truncation at r = 12 \xc3\x85 (12 cut-off), force-based switching over 10 to 12 \xc3\x85 (10\xe2\x80\x9312 fsw), and LJ particle mesh Ewald with no cut-off (LJPME). The FF parameters with these LJ cut-off methods are extensively validated by reproducing structural, interfacial, and mechanical properties. We find that the computed density and surface energies are in good agreement with reported experimental results, although the simulation results increase in the following order: 10\xe2\x80\x9312 fsw < 12 cut-off < LJPME; nanomaterials in which LJ interactions are a major component show relatively higher deviations (up to 4% in density and 8% in surface energy differences) compared to the experiment. Nanomaterial Modeler\xe2\x80\x99s capability is also demonstrated by generating complex systems of nanomaterial-biomolecule and nanomaterial-polymer interfaces with a combination of existing CHARMM-GUI modules. We hope that Nanomaterial Modeler can be used to carry out innovative nanomaterial modeling and simulation to acquire insight into the structure, dynamics, and underlying mechanisms of complex nanomaterial-containing systems.
INTRODUCTION<!>Workflow of Nanomaterial Modeler<!>Step 1 - Building bulk crystal(s)<!>Step 2 - Treatment of unbalanced atoms and surface modification<!>Model building<!>Computational details<!>Metals: gold surfaces with various Miller indices and its equilibrium particle shapes<!>Clay minerals: Kaolinite, Pyrophyllite, Montmorillonite, and Muscovite<!>Silica: \xc9\x91-cristobalite silica slab with different surface chemistry<!>Phosphate minerals: hydroxyapatite with different surfaces and Wulff construction<!>Transition Metal Dichalcogenide: MoS2<!>Complex system<!>CONCLUSIONS
<p>Understanding the dynamical evolution of biological and materials systems at the atomic scale is essential for groundbreaking advances in health science, materials science, energy conversion, sustainability, and overall quality of life.1–5 However, progress is limited as current experimental techniques alone cannot provide complete information about structures and dynamical processes at the nanometer to micrometer scale. Classical molecular modeling and simulation using force fields (FFs) and complex configuration databases are playing an increasingly important role in explaining experimental data, elucidating design principles, and making transformative property predictions for unknown biomolecular and nanomaterial structures and dynamics, as well as their interfaces. These computational methods are suitable for harnessing big data and accelerating discovery from the quantum scale to the microscale.6–8</p><p>State-of-the-art FFs for biomolecular systems (e.g., proteins, lipids, nucleic acids, and carbohydrates) have matured to the point where they can explain experiments and accelerate experimental discovery via testable hypotheses.9, 10 The same is true for the growing number of surface models for nanomaterials and nanomaterial-biological (nano-bio) interfaces.11, 12 In particular, the surface model database and parameters for metals, clay minerals, silica, apatites, layered materials (e.g., MoS2 and graphite), cement minerals, and gas molecules in the INTERFACE force field11 (IFF) routinely achieve more accurate predictions of surface energies, binding energies, and molecular recognition than typical density functional theory (DFT) methods and is compatible with biomolecular FFs.13–15 These recent developments present a unique opportunity for the life and material sciences to harness the predictive power of computer simulation methods to explore a broad range of nano-bio interfaces and complex electrolytes.</p><p>A recent round-robin study showed that distinct user groups working with different simulation FFs and programs yielded inconsistent results even for calculating simple thermodynamic properties such as density and potential energy.16 The lack of a reliable and unified cyberinfrastructure to build complex nano-bio interfaces poses major challenges to the molecular modeling and simulation community in terms of steep learning curves, risks of choosing unsuitable FF and faulty interface models, as well as mistakes in file conversion and input scripts that render simulations less useful or invalid. Several programs have been developed to help users to build nanomaterial model systems, including web applications, such as NanoModeler17 and PubVINAS,18 and stand-alone software packages such as Atomic Simulation Environment (ASE),19
pysimm,20
Molecular Simulation Design Framework (MoSDeF),21 and NanoMaterialCAD.22
NanoModeler supports 16 gold nanocluster models with ligand grafting function and PubVINAS provides 11 material types with corresponding physicochemical properties and/or bioactivities. ASE, pysimm, and MoSDeF provide methods to prepare various nanomaterial systems and API to integrate different features of existing software packages using Python-based scripting. NanoMaterialsCAD offers a graphical user interface (GUI) to build and manipulate a nanomaterial system. However, all aforementioned software requires significant preprocessing to prepare structures, topologies, and parameters of nanomaterials for simulation and/or is limited to use specific prebuilt nanostructures. Commercial packages also exist, such as Material Studio,23
Schrödinger,24 and Amsterdam Modeling Suite,25 but, are not freely available to everyone. Furthermore, the models from these programs are not transferable to other molecular dynamics (MD) simulation packages. Therefore, building models of functional nano-bio materials such as nanoparticle therapeutics, imaging agents, biomineral structures, bioinspired composites, and biosensors has been challenging to accomplish. Moreover, the simulation input preparation for nanomaterials and nano-bio interfaces currently involves multiple operations by researchers. This process includes the choice of building tools or scripts and file conversion and reassignment of FF parameters depending on the chosen simulation platform. Unfortunately, no single, user-friendly cyberinfrastructure is available to accomplish these crucial tasks.</p><p>CHARMM-GUI (https://www.charmm-gui.org), a web-based GUI, provides a well-designed workflow to interactively construct various complex biomolecular systems and seamlessly handles complicated internal data structures and simulation input files for CHARMM,26 NAMD,27 GROMACS,28 AMBER,29 GENESIS,30 LAMMPS31, TINKER,32 Desmond,33 and OpenMM34 with several biomolecular FFs.7 The simulation protocols are optimized6 following the principles of the original FF development.35, 36 Therefore, CHARMM-GUI is widely utilized in studies on the role of biomolecular motions, conformational changes, and thermodynamic relationships in biological function.</p><p>This work presents Nanomaterials Modeler, an important extension of CHARMM-GUI to a broad range of nanomaterials based on the IFF, which can bridge the gap between biomolecular and material simulations through compatibility with multiple simulation platforms.8, 37
Nanomaterial Modeler enables researchers to build nanomaterials models with up to 5,000,000 atoms and addresses aforementioned needs by merging IFF and CHARMM-GUI in an easy-to-use and state-of-the-art platform. The following sections discuss the methods, workflow, available nanomaterials models, the user interface, supported simulation engines, validation of the models, and example applications. The paper ends with brief conclusions.</p><!><p>Currently, CHARMM-GUI Nanomaterial Modeler supports ten classes of nanomaterials, including fcc metals, clay minerals, calcium sulfates, cement minerals, calcium silicate hydrate, silica, phosphate minerals, transition metal dichalcogenides (TMDC), and carbonaceous materials (Table 1).</p><p>Figure 1A shows an overall nanomaterial system building process that has been generalized and automated in two subsequent steps. Each step is designed to incorporate a user's specific options through a web interface and run CHARMM input files. Individual input and output files, including the generated structure and an archive of all created files, are available at each step. Video demonstrations on how to use Nanomaterial Modeler are available on the CHARMM-GUI website (http://www.charmm-gui.org/demo/nanomaterial). Nanomaterial Modeler adopts a GUI that allows researchers to quickly check and design nanomaterials (i.e., size and shape of nanomaterials, chemical modification of the surface, and periodicity along each axis).</p><!><p>In step 1, users can set a material type, shape, Miller index, size, and periodicity along X, Y, and Z directions. A unit cell structure of a selected nanomaterial (Figure 1B) is used to generate a user-specified nanomaterial system through unit cell duplication and translation (Figure 1C).</p><!><p>For specific nanomaterials, bonds between the primary and neighboring cells (i.e., image bonds) are required to model an infinite surface or molecule along with the periodic images. To facilitate such image bonds, the "patch information" necessary to create such connections has been defined for the currently supported nanomaterials (Figure 1D). The patch information includes all bonds, angles, dihedrals, and partial charge and atomic type information that are created when the primary cell connects to the 26 image cells. After setting periodicity, Nanomaterial Modeler performs necessary surface ionization, defect generation, and/or surface chemical modification (Figure 1E). At the end, users can obtain a nanomaterial system with the desired structure, topology, FF parameters, and simulation configuration files for further simulation. Furthermore, the generated structure can be used in CHARMM-GUI Multicomponent Assembler to model nano-bio systems as described below.</p><p>Nanomaterial Modeler provides validated all-atom simulation inputs for various MD programs, including CHARMM, GROMACS, NAMD, LAMMPS, AMBER, GENESIS, and OpenMM, enabling researchers to employ the package of their choice (see Supporting Information, Table S1, and S2). In the original IFF, electrostatic interactions are calculated using the particle mesh Ewald (PME) method and two Lennard-Jones (LJ) potentials (i.e., the 12-6 and 9-6 forms) are adapted. Nanomaterial Modeler supports a 12-6 LJ potential since the conventional biomolecular FFs such as AMBER,46 CHARMM,47–49 GROMOS,50 and OPLS-AA51 have been developed with a 12-6 LJ potential:

VLJ(r)=4ε[(σ/r)12−(σ/r)6]

where r is the distance between two interacting particles, ε is the depth of the potential well, and σ represents the distance at which the particle-particle potential energy is zero. In practice, the VLJ interactions are negligible at large distances, and therefore an energy cut-off is introduced at a certain interparticle distance. A 12 Å cutoff was chosen for IFF due to the rapid loss of electron-electron correlation of London dispersion interactions with distance after few layers of nearest neighbors in contrast to long-range Coulomb interactions between permanently charged atoms.52, 53 Also, a 12 Å cut-off is generally accepted for various FFs. However, different cut-off methods are used in conventional biomolecular FFs for different FF forms. The most widely-used methods are: i) simple truncation at r = 12 Å (12 cut-off), ii) force-based switching over 10 to 12 Å (10–12 fsw54), and iii) LJ particle mesh Ewald (LJPME). To compare the implementations of all supported MD programs and LJ methods, single point energy calculations of hydroxyapatite and α-quartz were performed (Table S3). The results show that all simulation programs are in excellent agreement with a maximum difference of 0.01% in the total energy, indicating that the nanomaterials IFF is correctly implemented within CHARMM-GUI. Note that the small observed differences are inevitable due to the usage of slightly different conversion factors within the programs and 10–12 fsw function implemented in GROMACS are slightly different with other software.55 In this work, the mechanical, physicochemical, and thermodynamic properties were analyzed using these LJ cut-off methods to assess the transferability of IFF for different simulation packages and ability to model complex nano-bio systems in combination with various biomolecular FFs.</p><!><p>The unit cell structures of all nanomaterials were obtained from X-ray diffraction data56–58 and the interface force field (IFF) database8. Starting from the conventional unit cell, a series of lattice vector transformations are performed to create an 'oriented' unit cell (OUC) where the a and b lattice vectors are parallel to the plane with Miller indices (hkl).59 Note that the c lattice vector is not necessarily perpendicular to the plane, although an orthogonal vector obtained within a reasonable cell size is used. For example, 18 fcc metals with 13 Miller indices have α = β = γ = 90°. All lattice parameters of OUC of supported nanomaterials are summarized in Table S4 and Table S5. Model building methods for larger models from multiple unit cells are described in Figure 1. All model building procedures are controlled by CHARMM scripts and individual input and output files, including the CHARMM scripts and an archive of all the created files (i.e., structure, topology, and configuration files), are available on the Nanomaterial Modeler (http://charmm-gui.org/input/nanomaterial). Details of complex nano-bio and nano-polymer system building are summarized in Supporting Information.</p><!><p>All simulations used the interface force field (IFF) for nanomaterials, the CHARMM force field for proteins, peptides, and lipids, and the CHARMM generalized force field for polymer.9, 36, 60 The TIP3P water model was utilized for water-containing systems. To achieve better sampling and check convergence, three independent molecular dynamics (MD) simulations were performed for each system. Periodic boundary conditions (PBCs) are employed for all simulations, and the particle mesh Ewald (PME) method61 is used for long-range electrostatic interactions. A leap-frog algorithm was used to integrate Newton's equation of motion. The simulation time step is set to 1 fs for equilibration and 2 fs for the production run in conjunction with the SHAKE algorithm62 to constrain the covalent bonds involving hydrogen atoms for all programs except GROMACS in which the LINCS algorithm63 is used instead. All simulation trajectories were recorded every 10 ps except vibrational frequency calculation (see Supporting Information). For each nanomaterial model, all structure and corresponding parameter files are generated in the CHARMM format (i.e., rtf and prm files). For GROMACS, NAMD, and LAMMPS, FF-Converter in CHARMM-GUI64 is used for format conversion from CHARMM data format to corresponding program readable formats. For specific nanomaterials, bonds between the primary and neighboring cells (image-bond) are required to model infinite surfaces or molecules along with the periodic images. Nanomaterial Modeler supports CHARMM, NAMD, GROMACS, LAMMPS, OpenMM, Amber, and Genesis for non-image-bond systems such as fcc metals. Four simulation programs, OpenMM, Gromacs, LAMMPS, and NAMD, are supported for image-bond systems. Three types of LJ cutoff methods, including 10–12 fsw, 12 cut-off, and LJPME are employed for the LJ interactions to investigate the cut-off method effect on the structural, physical, and mechanical properties. For LJPME, grid spacing and interpolation order were set to 1.2 Å and 4, respectively. Details of computational methods such as thermostat, barostat, and coupling constants for each simulation program are described in Supporting Information S1. Computational details of density, surface energies, mechanical properties, and vibrational spectra calculation for various nanomaterials are provided in Supporting Information S2–7. Details of computational methods, including density, surface energies, mechanical properties, and vibrational spectra, are provided in the Supporting Information.</p><!><p>Nanomaterial Modeler supports 18 fcc metals (Table 1): Ac, Ag, Al, Au, Ca (α), Ce (γ), Cu, Es (β), Fe (γ), Ir, Ni, Pb, Pd, Pt, Rh, Sr (α), Th (α), and Yb (β). One can generate i) bulk crystal, ii) different cleavage surfaces with Miller indices ranging from (100) to (332) with any combination of periodicity, and iii) a variety of shapes such as a sphere, cylinder, rod, polygon, box, and Wulff construction. Figure 2A shows the user interface of Nanomaterial Modeler for building a (111) gold surface with a size of 51.9 × 59.9 × 35.3 Å3; n.b., the input dimension and the final model dimension can be different because of the unit cell size. The selected (111) Miller index plane is perpendicular to the Z-axis. With user-specified X, Y, and Z dimensions, Nanomaterial Modeler displays an estimated system size from the unit cell information (see Supporting Information for lattice parameters of all unit cells). To build an infinite surface parallel to the XY plane with vacuum layers (i.e., a slab), periodicity can be selected in the X and Y directions and vacuum selected under "System Type". This simple user interface enables easy system building for various fcc metals. With these options, a final all-atom model of (111) gold surface is illustrated in Figure 2B. In this work, experimental densities and surface energies are employed for model validation. Computed densities of 18 fcc metals with different LJ cut-off methods (i.e., 12 cut-off, 10–12 fsw, and LJPME) are in excellent agreement with the experimental data for all supported simulation programs (Table S6 and Figure S1). Solid-vapor interface tensions (γSV(111)) for 18 fcc metals with 13 Miller indices were computed (Figure S2). Figure 2C shows comparison of computed and experimental γSV(111) of 18 fcc metals with three LJ cut-off methods. The 10–12 fsw, 12 cutoff, and LJPME reproduce surface energy for all fcc metals with less than −6.5%, 0.9%, and 6.5% deviations, respectively (Table S4 and S7).</p><p>Nanomaterial Modeler provides versatile methods to generate nanoparticle structures and surface modeling. One of the classic models to describe a particle shape is the Wulff construction65 based on its orientation-dependent surface free energy. The computed surface energies of all fcc metals are presented in Table S4, which were used to generate the Wulff constructions (Figure 2D). One can select any combination of Miller indices and adjust surface energy values in the user interface. Figure 2D shows predicted equilibrium shapes using the Wulff construction and corresponding all-atom models of gold nanoparticles. Selecting (100) and (110) surfaces yield a cuboid (Shape 1). Changing from (110) to (111) produces a truncated octahedron (Shape 2). As the (100) surface energy increases, (100) surface area decreases and the shape changes to an octahedron (Shape 3). When a (110) surface is added to Shape 2, edges are covered with (110) surface and overall nanoparticle shape is closer to a spherical shape (Shape 4).</p><!><p>Clay minerals have been widely employed not only for industrial applications,66 but also as additive biomaterials for drug delivery67. Nanomaterial Modeler provides facile modeling methods for four different clay minerals: kaolinite, pyrophyllite, montmorillonite, and muscovite. Unlike fcc metals, an infinite surface model of clay minerals has image bonds across the periodic boundary conditions (PBCs), and 27 image bonds are constructed. In addition, some clay minerals may contain ions between layers, which can be measured by cation exchange capacity (CEC) states. For example, montmorillonite exists in various CEC states with different amounts of Mg replaced by Al. This substitution leaves net negative charges, which attracts cations for charge neutrality. Nanomaterial Modeler provides options to control the CEC states by randomly replacing Al with Mg and randomly distributing either Na+ or K+.</p><p>Figure 3A shows the user interface for building a montmorillonite surface model with 0.333 mmol/100g CEC, a size of 60 × 60 × 30 Å3, and XY periodicity. Users can select desired CEC states by setting the ratio of the defect (i.e., a ratio of Al and Mg) and ion types. In addition, X and Y periodic options are selected to build an infinite surface along the XY plane. Representative snapshots of three clay minerals (pyrophyllite, montmorillonite, and muscovite) with various CEC states are shown in Figure 3B. The first row displays top views, where layers below the top layer are hidden for clarity. The second row illustrates a side view of each material. While pyrophyllite does not have ions in between layers, montmorillonite and muscovite do, and thus, the number of ions increases as the CEC states become larger. As a result, the spacing between layers of pyrophyllite is ~2.9 Å that is smaller than the other clay minerals due to the lack of ions between the layers (e.g., ~3.9 Å for montmorillonite and muscovite). Supported Miller indices and corresponding unit cell parameters are summarized in Table S5.</p><p>The densities, bulk modulus, and Young's modulus of the clay minerals are calculated with three different LJ cut-off methods to validate the model. Computed densities of clay minerals with other LJ cut-off methods show good agreement with experimental data with less than 4% deviation for all supported simulation programs (Table S8 and Figure S3). Figure 3C shows the bulk modulus of muscovite as a function of the applied pressure, which is in good agreement with a previously reported bulk modulus.68 Also, all of the LJ cut-off methods yield consistent results; 10–12 fsw and LJPME have deviations of −2.6% and 3.7% relative to the 12 cut-off. In experiments, solid-vapor interface tensions (γSV) show a broad range from 0.050 to 0.200 J/m2 because preparation process of Montmorillonite surface with specific CEC value is difficult, and even if the surface is made, the surface is very sensitive to environmental conditions such as humidity. The computed γSV of Montmorillonite are ranging from 0.087 to 0.113 J/m2, which is consistent with the range observed in experiments (Table S9). For Muscovite, computed γSV reproduce experimental observation with less than 4 % deviation. Note that the energy differences of computed γSV among the LJ cut-off methods are negligible for all clay minerals.</p><p>Figure 3D represents Young's modulus of muscovite calculated along the Cartesian directions. The calculated moduli (132.0 GPa, 113.4 GPa, and 28.6 GPa for X, Y, and Z direction, respectively) show consistent trends and similar values to previous experimental and DFT results69, although the values are about 10% lower. Identical to the reported study, the modulus along the X direction is approximately 16.4% larger than Y direction, and the modulus along the Z direction is smaller by at least 50% than either X or Y direction.</p><!><p>Silicon dioxide and silica are widely available in nature and biologically enriched in various organisms.70, 71 These are also important materials in the semiconductor industry. The simulation of bulk and surface properties of silica has been of great interest since the emergence of computational modeling methods.43, 72–74
Nanomaterial Modeler provides diverse modeling capabilities of silica that cover various surface chemistry and pH (i.e., the surface density of silanol and siloxide groups, and degree of ionization) for two shapes (box and sphere).</p><p>Figure 4A shows the user interface for building a surface of α-cristobalite with a 4.7 nm−2 density of silanol groups. The surface chemistry of silica depends on the surface characteristics (e.g., cleavage plane, particle size, and porosity), heat treatment, and environmental pH.43 Various forms of silica at high pH contain Q2 surface environments, which correspond to two silanol groups per superficial silicon atom (=Si(OH)2), and mixed Q2/Q3 surface environments where Q3 surface represents one silanol group per silicon atom (≡Si(OH)). The area density of silanol groups is in the range of 9.4 to 4.7 nm−2 (Figure 4B). Most silica glasses and medium size nanoparticles (~100 nm) contain 70–90% Q3 environments on the surface (Figure 4C). The silica surfaces after thermal treatment consist of a high portion of Q4 environments in which siloxide bridges without silanol groups (Figure 4D). All silanol groups on the silica surface are subject to deprotonation or/and protonation upon environmental conditions such as the area density of silanol groups, pH, ionic strength of the solution, and the type of ions present in solutions. In the physiological conditions with an ionic strength of 0.1~0.3 M of sodium ions, the ionization degree of silanol groups ranges from 0% to 20%. Nanomaterial Modeler supports ionization degrees up to 50% (Figure 4A). Figures 4E–G shows a Q3 surface with a total Si-O(H,Na) density of 4.7 nm−2 with 5%, 15%, and 25% ionization, respectively.</p><p>The silica models show good agreement between computed and experimental properties of (1) densities, (2) water contact angles on the silica surfaces, (3) heat of immersion of silica surfaces in water, and (4) vibrational properties. Figure 5A shows the water contact angle, θc, on charge-neutral silica surfaces ranging from Q3 to Q4 environments in simulation and experiment.75 The details of contact angle calculations are summarized in the Supporting Information and in Figure S4. The Q3 surface exhibits a contact angle of θc = 0. This surface is strongly hydrophilic due to the formation of hydrogen bonds between surface Si-OH groups and water molecules, thus fully wetted with water. Heat treatment from 200 to 1000°C decreases the area density of surface silanol groups due to condensation of adjacent silanol groups in Q2 and Q3 environments, yielding Q4 environments. The transition in surface chemistry from Q3 to Q4 is modeled with 4.7, 2.35, 1.2, and 0 silanol groups per nm2. Note that the silanol groups are non-ionized to represent interfaces with deionized water. The agreement between simulation and experimental measurements falls within ±3°. The heat of immersion (ΔHimm) represents the enthalpy released upon immersion of clean particles or surfaces into the water and provides insight into silica-water interactions. Figure 5B shows a computational procedure to compute ΔHimm of Q3 silica surfaces in water. ΔHimm of Q3 silica obtained from calorimetric measurements is 160 ± 5 mJ/m2 at 300 K, and the calculated results show good agreement within the uncertainty for all LJ cut-off methods. The model also reproduces the vibrational frequency of silica with infrared and Raman measurements (Figure 5C and Figure S5). Strong bands at 960~1200 cm−1 and 600~800 cm−1 correspond to asymmetric and symmetric Si-O-Si stretching vibrations, respectively. The O-Si-O bending vibrations at 400~500 cm−1 and O-H stretching vibrations near 3700 cm−1 are also reproduced in the simulation. Although the computed value for one of the bands near 800 cm−1 is shifted approximately 100 cm−1 lower and the intensities of vibrational spectrum are not reproduced due to the lack of full electronic structure, most peak positions of the vibrational modes are clear and sufficient to monitor peak shifts in chemically different environments.</p><!><p>Hydroxyapatite (HAP), Ca10(PO4)6(OH)2, is of great importance for human health because it is the major component in human bone and teeth and plays a central role in maladies such as osteoporosis.44, 76 Atomistic models of HAP could facilitate a better understanding of the complex surface chemistry and provide insight into the interaction between biomolecules and HAP. Recent simulation results have shown that specific interactions between HAP surfaces and proteins are strongly dependent on pH, type of facet, surface defect, and specific details of the amino acid arrangement.44, 77 To facilitate the research of HAP surface chemistry, surface models of various facets and nanocrystals as a function of pH are introduced in Nanomaterial Modeler.</p><p>Figure 6A shows the user interface for building a (001) HAP surface with a size of 60 × 60 × 30 Å3 at pH ~10. Nanomaterial Modeler provides two shapes (box and Wulff construction). The number of possible (hkl) cleavage plane is in principle unlimited, but cleavage preferably occurs in (hkl) directions with weaker nonbonded interactions. Nanomaterial Modeler offers the common cleavage planes of HAP, (001), (010), (020), and (101), and surface models of HAP at various pH conditions, including pH ~5, ~10, and > 14 (see details in Supporting Information). MD simulations were performed to investigate the cleavage energies of the common low index planes of HAP (Figure 6B). The computed cleavage energies of PO43- terminated surfaces range from 897.0 to 1223.3 mJ/m2 and increase in the following order: (001) < (101) < (010) < (020). Note that the cleavage energy differences among the simulation packages are negligible (Table S9). For the LJ cut-off methods, 10–12 fsw and LJPME produce only −2.4% and 3.2% deviations relative to the 12 cut-off method (Table S10). Such a small difference compared to the deviations observed for fcc metals arises from the fact that electrostatic interactions dominate the HAP surface energy (see details in Supporting Information). In an experiment, the preparation of ideal cleaved surfaces of HAP is complex as they are hygroscopic and thus sensitive to environmental conditions such as humidity. To our knowledge, direct measurements of cleavage energies in vacuum have not been reported. Nevertheless, the cleavage energies of minerals with similar chemical composition and previously calculated cleavage energies for HAP are in good agreement with the current simulation results (see Table S10).</p><p>Nanomaterial Modeler enables modeling of HAP nanocrystals (i.e., Wulff construction) and various surface models. Figure 6C illustrates predicted shapes of HAP nanocrystal and final all-atom models according to the combination of Miller indices and corresponding surface energies. Like the Wulff construction of gold nanoparticles in Figure 2, one can freely change the surface energy values in the text box and add/delete surfaces by clicking the ± button. In addition, the environment pH option can be specified during Wulff construction. Since cleavage energies are highly affected by environment pH, Nanomaterial Modeler automatically updates energy values and predicts corresponding shapes based on the selected pH range.</p><!><p>Two-dimensional (2D) materials, including graphene and transition metal dichalcogenides (TMDC), such as molybdenum disulfide (MoS2), have received significant attention due to their unique structural and electronic properties.78, 79
Nanomaterial Modeler provides a facile modeling capability for graphene and 2H-MoS2.</p><p>Figure 7A shows the user interface for building four 2H-MoS2 layers with a size of 50 × 50 × 24 Å3. In the Unit Cell Info section, one can set the number of MoS2 layers to be generated according to the Z-length value in the Box Options. Note that MoS2 is limited to XY-PBC as there is no information available for surface-end modification. The cleavage energy of the basal plane of a layered material is a key property for its applications. To evaluate the cleavage energy of MoS2, a potential energy difference of two systems (i.e., a box of surface slab separated with a 60 Å vacuum layer (Esep) and an equivalent periodic bulk model of the same number of atoms without a vacuum layer (Ebulk)) is computed with 10–12 fsw, 12 cut-off, and LJPME using four simulation packages (Figure 7B and Supporting Information). The results agree with an experimental observation of 150 ± 10 mJ/m2 for 12 cut-off. The use of different cut-off methods yields 137 ± 2 mJ/m2 (10–12 fsw) and 156 ± 1 mJ/m2 (LJPME), respectively. These relatively high deviations are also observed in the fcc metal cases as LJ interactions are the major component for their cleavage energy. For MoS2, the contributions of LJ and electrostatic interactions to the cleavage energy are 93% and 7%, respectively. This indicates that, when LJ interaction is the main contributing factor to the cleavage energy, the LJ parameter adjustments for 10–12 fsw and LJPME methods are required to achieve the level of accuracy at experiment.</p><p>Structural and mechanical properties, as well as surface properties are essential factors to validate the MoS2 model.15
Figure 7C shows that the computed XRD pattern from experiment and simulation. The characteristic peak of MoS2 is observed at 14.4°, corresponding to the (002) plane, and agrees well with experimental observation. For mechanical properties, the compressibility is computed with 12 cut-off method by recording the volume change at different pressures (1, 2, 3, and 4 GPa) and compared to the experiment (Figure 7D). The simulation results are identical with experiment up to about 2 GPa. Even when the pressure exceeds 2 GPa, the difference between the computed and experimental data remain below a 2% deviation.</p><!><p>Several complex nano-bio systems are modeled and simulated to illustrate the potential use of Nanomaterial Modeler in combination with other CHARMM-GUI modules. Note that these illustrative systems are chosen to show the new capability of Nanomaterial Modeler but not for FF validation. As a first example, a platinum (Pt) nanoparticle system with T7 peptides (acyl-TLTTLTN-amide) is built using three modules in CHARMM-GUI: Nanomaterial Modeler for the Pt nanoparticle, PDB Reader & Manipulator7, 80 for T7 peptides, and Multicomponent Assembler for assembly of all models and solvation (Figure 8A). A cuboctahedron Pt nanoparticle is generated using Wulff construction (Shape 2 in Figure 2 and Figure S6) and located at the system center. T7 peptides are randomly distributed in the system with 3, 6, 12, 18, and 36 peptides corresponding to concentrations of 5, 10, 20, 50, and 100 μg/ml (see Methods and Supporting Information for details of the system building and simulation methods). After MD simulations of 200 ns duration, all T7 peptides are adsorbed onto the nanoparticle surface at all concentration levels (Figure 8B). At low concentration (< 20 μg/ml), the binding of the T7 peptides to the Pt nanoparticle occurred preferentially near edges as opposed to the inner portions of the (100) facets. This is due to the substantially reduced binding of water at the edges, which is consistent with previously reported data.81 These simulations provide insight into the mechanisms of nanocrystal growth and the spatial distribution of facet-specific ligands as a function of concentration.</p><p>As a second example, 10th type III modules of fibronectin (FN-III10, PDB ID: 1TTF82) on a HAP surface are generated (Figure 8C). Fibronectin (FN) is known to regulate cell adhesion, growth, differentiation, or survival of osteoblasts, and support osteogenic cell responses in vitro.83 This complex nano-bio interface system can be readily generated using a combination of the Nanomaterial Modeler, PDB Reader, and Multicomponent Assembler modules. This system can be used to investigate the effects of HAP surface properties on the adsorption of FN modules, which is important because the surface topology may modulate the biological activity of FN and the corresponding cell adhesion process.</p><p>The third example consists of a supported lipid bilayer (SLB) that is a popular model of cell membranes with potential biotechnological applications.84 Many experimental techniques such as atomic force microscopy (AFM), quartz crystal microbalance with dissipation monitoring (QCM-D), and ellipsometry have been employed to investigate the structure and physical properties of lipid bilayers reconstituted with membrane proteins.85, 86 However, most of these methods only provide superficial information and/or lack molecular-level insight into underlying mechanisms. Figure 8D shows an SLB system built using a combination of the Nanomaterial Modeler, Membrane Builder,87–89 and Multicomponent Assembler modules (see Methods and Supporting Information).</p><p>In addition to various nano-bio interfaces, Nanomaterial Modeler can also be used to build polymer-containing nanomaterial systems (i.e., nano-polymers). As a final example, a nano-polymer complex system that consists of a cement mineral (in this case, tobermorite) and poly(acrylic acid) (PAA) was investigated. It is known that PAA strongly binds to the cement surface via ionic or possibly chelate binding and induces the cement to harden. Three modules in CHARMM-GUI are employed: Nanomaterial Modeler for tobermorite, Polymer Builder60 for PAA, and Multicomponent Assembler for integration and solvation. A tobermorite (004) slab is located in the simulation box, and a 60% ionized PAA chain is placed 10 Å above the surface (Figure 8E). At an early simulation stage, carboxylates in the PAA side chains begin to interact with the surface via strong ionic bonds with Ca2+ ions, and finally, the entire polymer chain is adsorbed on the surface.</p><!><p>This work presents Nanomaterial Modeler in CHARMM-GUI, a web-based cyberinfrastructure for building all-atom models of various nanomaterials and providing all necessary FF and configurational files for MD simulations. The model building workflow is generalized and automated in two steps: (i) building a bulk crystal through duplication and translation of a unit cell structure, and (ii) applying image patches and capping unbalanced atoms based on periodicity and chemical environment. Nanomaterial Modeler's versatile and efficient modeling features are illustrated by building various nanomaterial surface models and equilibrium nanoparticle shapes. Moreover, the transferability of nanomaterial models among the simulation programs was assessed by single-point energy calculations, which yielded 0.01% relative absolute energy differences for various surface models and equilibrium nanoparticle shapes. The significance of this work is that Nanomaterial Modeler provides a convenient modeling capability for various nanomaterial systems. Generated nanomaterial models can be used to model complex systems with other CHARMM-GUI modules, as demonstrated here with selected representative test cases.</p><p>To assess the transferability of our models and the IFF, we have investigated the effect of LJ cutoff methods on the structural, mechanical, and thermodynamic properties of nanomaterials. We have determined that the different LJ cut-off methods exhibit overall consistent results for most nanomaterial cases. The exceptions are nanomaterials for which LJ interactions are a major component for their cleavage energy. Such systems exhibit relatively larger deviations (up to 8%) compared to the electrostatic-driven materials if an LJ cut-off other than 12 Å is used. However, in most cases deviations are minor and consistent cut-offs or adjustments in the ε and σ parameters can be explored. Future work will include testing transferability for other nanomaterial properties and comparison with reactive and machine learning FFs through the OpenKIM framework.90 We hope that Nanomaterial Modeler can be useful for carrying out innovative and novel nanomaterial modeling and simulation research to acquire insight into structures, dynamics, and underlying mechanisms of complex nano-bio and nano-polymer interface systems.</p>