Determination of interaction parameters in a bottom-up approach employed in reactive dissipative particle dynamics simulations for thermosetting polymers.

The limitations in previous dissipative particle dynamics (DPD) studies confined simulations to a narrow resin range. This study refines DPD parameter calculation methodology, extending its application to diverse polymer materials. Using a bottom-up approach with molecular dynamics (MD) simulations, we evaluated solubility parameters and bead number density governing nonbonded interactions via the Flory-Huggins parameter and covalent-bonded interactions. Two solubility parameter methods, Hildebrand and Krevelen-Hoftyzer, were compared for DPD simulations. The Hildebrand method, utilizing MD simulations, demonstrates higher consistency and broader applicability in determining solubility parameters for all DPD particles. The DPD/MD curing reaction process was examined in three epoxy systems: DGEBA/4,4'-DDS, DGEBA/MPDA and DGEBA/DETA. Calculations for the curing profile, gelation point, radial distribution function and branch ratio were performed. Compared to MD data for DGEBA/4,4'-DDS, the maximum deviation in secondary reactions between epoxy and amine groups according to DPD simulations with Krevelen-Hoftyzer was 14.8%, while with the Hildebrand method, it was 1.7%. The accuracy of the DPD curing reaction in reproducing the structural properties verifies its expanded application to general polymeric material simulations. The proposed curing DPD simulations, with a short run time and minimal computational resources, contributes to high-throughput screening for optimal resins and investigates mesoscopic inhomogeneous structures in large resin systems.

Molecular dynamics simulations for interfacial structure and affinity between carboxylic acid-modified Al2O3 and polymer melts.

Controlling the dispersion state of nanoparticles in a polymer matrix is necessary to produce polymer nanocomposites. The surface modification of nanoparticles is used to enable their dispersion in polymers. Moreover, molecular dynamics (MD) simulations are useful for revealing the interfacial properties between nanoparticles and polymers to aid in the design of materials. In this study, the effect of surface coverage, modifier length, and polymer species on the interfacial structure and affinity between surface-modified Al2O3 and polymer melts were investigated using all-atom MD simulations. Hexanoic, decanoic, and tetradecanoic acids were used as surface modifiers, and polypropylene (PP), polystyrene (PS), and poly (methyl methacrylate) (PMMA) were used as polymers. The work of adhesion Wadh and the work of immersion Wimm were selected as quantitative measures of affinity. Wadh was calculated using the phantom-wall approach, and Wimm was calculated by simply subtracting the surface tension of polymers γL from Wadh. The results showed that Wadh and Wimm were improved by surface modification with low coverage, owing to a good penetration of the polymer. The effect of modifier length on Wadh and Wimm was small. Whereas Wadh increased in the following order: PP < PS < PMMA, Wimm increased as follows: PMMA < PS < PP. Finally, the trend of Wadh and Wimm was organized using the Flory-Huggins interaction parameter χ between the modifier and the polymer. This study demonstrates that the interfacial affinity can be improved by tuning the surface coverage and modifier species depending on the polymer matrix.

Thermoset resin curing simulation using quantum-chemical reaction path calculation and dissipative particle dynamics.

Thermoset resin, which is commonly used as a matrix in carbon-fiber-reinforced plastic, requires curing procedures. We propose a curing simulation technique involving a dissipative particle dynamics (DPD) simulation, which can simulate a larger system and longer time period than those of conventional all-atom molecular dynamics (AA-MD) simulations. The proposed curing DPD simulation can represent the thermoset resin exothermic reaction process precisely by considering each reactivity according to the reaction types calculated via quantum-chemical reaction path calculations. The cure reaction process given by the curing DPD simulation agrees well with that given by a conventional curing AA-MD simulation, but with run-time and computational-resource reductions of 1/480 and 1/10 times, respectively. We also conduct reverse mapping, through which the AA-MD system can be reconstructed from the DPD system, to evaluate the structural and thermomechanical properties. The X-ray diffraction pattern and thermomechanical properties of the reconstructed system agree well with those of the systems derived from the curing AA-MD simulation and experimental setup. Therefore, a cured-resin AA-MD system can be obtained from a curing DPD simulation at an extremely low computational cost, and the thermomechanical properties can be evaluated precisely using this system. The proposed curing simulation technique can be applied in high-throughput screening for better materials properties and in large system calculations.

Evaluation of the work of adhesion at the interface between a surface-modified metal oxide and an organic solvent using molecular dynamics simulations.

Advancing the practical applications of surface-modified nanoparticles requires that their dispersion in solvents can be controlled. The degree of dispersion depends on the affinity between surface-modified nanoparticles and solvents, which can be quantified using the work of adhesion at the interface. Herein, the affinity between a surface-modified inorganic solid and an organic solvent was evaluated by calculating the work of adhesion at the interface. The phantom-wall method, which is a thermodynamic route for evaluating the work of adhesion at an interface using molecular dynamics simulations, was applied to the decanoic acid-modified Al2O3/hexane interface. Molecular dynamics simulations were performed for flat interface systems to focus on the interactions between substances that affect the affinity on the surface. As a result, the surface coverage of decanoic acid was found to affect the work of adhesion, with a maximum value of 45.66 ± 0.75 mJ/m2 at a surface coverage of 75%. An analysis of the mass density profiles of Al2O3, decanoic acid, and hexane in the vicinity of the interface showed that the increase in the work of adhesion with the surface coverage was due to the penetration of hexane molecules into the decanoic acid layer on the Al2O3 surface. At a surface coverage of 75%, some hexane molecules were trapped in the layer of oriented decanoic acid molecules. These results suggested that the interfacial affinity can be enhanced by controlling the surface modification so that the solvent can penetrate the layer of the modifier.

Cross-Plane and In-Plane Heat Conductions in Layer-by-Layer Membrane: Molecular Dynamics Study.

A material with anisotropic heat conduction characteristics, which is determined by molecular scale structure, provides a way of controlling heat flow in nanoscale spaces. As such, here, we consider layer-by-layer (LbL) membranes, which are an electrostatic assembly of polyelectrolyte multilayers and are expected to have different heat conduction characteristics between cross-plane and in-plane directions. We constructed models of a poly(acrylic acid)/polyethylenimine (PAA/PEI) LbL membrane sandwiched by charged solid walls and investigated their anisotropic heat conduction using molecular dynamics simulations. In the cross-plane direction, the thermal boundary resistance between the solid wall and the LbL membrane and that between the constituent PAA and PEI layers decrease with increasing degree of ionization (solid surface charge density and the number of electric charges per PAA/PEI molecule). When the degree of ionization is low, the cross-plane thermal conductivity of a constituent layer is higher than that of the bulk state. As the degree of ionization increases, however, the cross-plane thermal conductivity of PAA, a linear polymer, decreases because of the increase in the number of in-plane oriented polymer chains. In the in-plane direction, we investigated the heat conduction of each layer and found the enhancement of effective in-plane thermal conductivity again due to the in-plane oriented chain alignment. The heat conduction in the LbL membrane is three-dimensionally enhanced compared to those in the bulk states of the constituent polymers because of the electrostatic interactions in the cross-plane direction and the molecular alignment in the in-plane direction.

Interpretation of Young's equation for a liquid droplet on a flat and smooth solid surface: Mechanical and thermodynamic routes with a simple Lennard-Jones liquid.

In this study, we carried out molecular dynamics simulations of a cylindrical Lennard-Jones droplet on a flat and smooth solid surface and showed that Young's equation as the relation among solid-liquid, solid-vapor, and liquid-vapor interfacial tensions γSL, γSV, and γLV, respectively, was applicable only under a very restricted condition. Using the fluid stress-tensor distribution, we examined the force balance in the surface-lateral direction exerted on a rectangular control volume set around the contact line. As the mechanical route, the fluid stress integrals along the two control surfaces normal to the solid-fluid interface were theoretically connected with γSL and γSV relative to the solid-vacuum interfacial tension γS0 by Bakker's equation extended to solid-related interfaces via a thought experiment, for which the position of the solid-fluid interface plane was defined at the limit that the fluid molecules could reach. On the other hand, the fluid stress integral along the control surface lateral to the solid-fluid interface was connected with γLV by the Young-Laplace equation. Through this connection, we showed that Young's equation was valid for a system in which the net lateral force exerted on the fluid molecules from the solid surface was zero around the contact line. Furthermore, we compared γSL - γS0 and γSV - γS0 obtained by the mechanical route with the solid-liquid and solid-vapor works of adhesion obtained by the dry-surface method as one of the thermodynamic routes and showed that both routes resulted in a good agreement. In addition, the contact angle predicted by Young's equation with these interfacial tensions corresponded well to the apparent droplet contact angle determined by using the previously defined position of the solid-fluid interface plane; however, our theoretical derivation indicated that this correspondence was achieved because the zero-lateral force condition was satisfied in the present system with a flat and smooth solid surface. These results indicated that the contact angle should be predicted not only by the interfacial tensions but also by the pinning force exerted around the contact line.

Equivalence of the EMD- and NEMD-based decomposition of thermal conductivity into microscopic building blocks.

Thermal conductivity of a material can be comprehended as being composed of microscopic building blocks relevant to the energy transfer due to a specific microscopic process or structure. The building block is called the partial thermal conductivity (PTC). The concept of PTC is essential to evaluate the contributions of various molecular mechanisms to heat conduction and has been providing detailed knowledge of the contribution. The PTC can be evaluated by equilibrium molecular dynamics (EMD) and non-equilibrium molecular dynamics (NEMD) in different manners: the EMD evaluation utilizes the autocorrelation of spontaneous heat fluxes in an equilibrium state whereas the NEMD one is based on stationary heat fluxes in a non-equilibrium state. However, it has not been fully discussed whether the two methods give the same PTC or not. In the present study, we formulate a Green-Kubo relation, which is necessary for EMD to calculate the PTCs equivalent to those by NEMD. Unlike the existing theories, our formulation is based on the local equilibrium hypothesis to describe a clear connection between EMD and NEMD simulations. The equivalence of the two derivations of PTCs is confirmed by the numerical results for liquid methane and butane. The present establishment of the EMD-NEMD correspondence makes the MD analysis of PTCs a robust way to clarify the microscopic origins of thermal conductivity.

Hydrodynamic consideration of the finite size effect on the self-diffusion coefficient in a periodic rectangular parallelepiped system.

In the present study, we use molecular dynamics (MD) simulations to provide an insight into the system size effect on the self-diffusion coefficient of liquids in the periodic rectangular parallelepiped system, from the hydrodynamic perspective. We have previously shown that in the rectangular box system, the diffusivity exhibits anomalous behaviors, i.e., the diffusion tensor appears to be anisotropic despite the bulk liquid simulation and the diffusion component in the direction along the short side of rectangular box with a high aspect ratio exceeding the diffusivity in the infinite system [Kikugawa et al., J.Chem. Phys. 142, 024503 (2015)]. So far, the size effect on the diffusivity has been intensively studied in the cubic system and has been interpreted quite well by the theoretical considerations employing the hydrodynamic interaction. Here, we have extended the hydrodynamic theory to be applied to periodic rectangular box systems and compared the theoretical predictions with MD simulation results. As a result, the diffusivity predicted by the hydrodynamic theory shows good agreement with the MD results. In addition, the system size effect was examined in a rod-shaped rectangular box in which the two shorter side lengths were equivalent and a film-type rectangular box in which the two longer side lengths were equivalent. It is of interest that we found that the aspect ratio, at which the diffusivity coincides with that in the infinite system, is a universal constant independent of the cross-sectional area for the rod system or the thickness for the film system. By extracting the universal structure in the hydrodynamic description, we also suggested a simplified approximate model to accurately predict the size effect on the diffusivity over a practical range of aspect ratios.

Effects of molecular structure on microscopic heat transport in chain polymer liquids.

In this paper, we discuss the molecular mechanism of the heat conduction in a liquid, based on nonequilibrium molecular dynamics simulations of a systematic series of linear- and branched alkane liquids, as a continuation of our previous study on linear alkane [T. Ohara et al., J. Chem. Phys. 135, 034507 (2011)]. The thermal conductivities for these alkanes in a saturated liquid state at the same reduced temperature (0.7Tc) obtained from the simulations are compared in relation to the structural difference of the liquids. In order to connect the thermal energy transport characteristics with molecular structures, we introduce the new concept of the interatomic path of heat transfer (atomistic heat path, AHP), which is defined for each type of inter- and intramolecular interaction. It is found that the efficiency of intermolecular AHP is sensitive to the structure of the first neighbor shell, whereas that of intramolecular AHP is similar for different alkane species. The dependence of thermal conductivity on different lengths of the main and side chain can be understood from the natures of these inter- and intramolecular AHPs.

Effect of the computational domain size and shape on the self-diffusion coefficient in a Lennard-Jones liquid.

In the present study, molecular dynamics (MD) simulations on the monatomic Lennard-Jones liquid in a periodic boundary system were performed in order to elucidate the effect of the computational domain size and shape on the self-diffusion coefficient measured by the system. So far, the system size dependence in cubic computational domains has been intensively investigated and these studies showed that the diffusion coefficient depends linearly on the inverse of the system size, which is theoretically predicted based on the hydrodynamic interaction. We examined the system size effect not only in the cubic cell systems but also in rectangular cell systems which were created by changing one side length of the cubic cell with the system density kept constant. As a result, the diffusion coefficient in the direction perpendicular to the long side of the rectangular cell significantly increases more or less linearly with the side length. On the other hand, the diffusion coefficient in the direction along the long side is almost constant or slightly decreases. Consequently, anisotropy of the diffusion coefficient emerges in a rectangular cell with periodic boundary conditions even in a bulk liquid simulation. This unexpected result is of critical importance because rectangular fluid systems confined in nanospace, which are present in realistic nanoscale technologies, have been widely studied in recent MD simulations. In order to elucidate the underlying mechanism for this serious system shape effect on the diffusion property, the correlation structures of particle velocities were examined.

Dissipative Particle Dynamics Simulation for Reaction-Induced Phase Separation of Thermoset/Thermoplastic Blends.

Reaction-induced phase separation occurs during the curing reaction when a thermoplastic resin is dissolved in a thermoset resin, which enables toughening of the thermoset resin. As resin properties vary significantly depending on the morphology of the phase-separated structure, controlling the morphology formation is of critical importance. Reaction-induced phase separation is a phenomenon that ranges from the chemical reaction scale to the mesoscale dynamics of polymer molecules. In this study, we performed curing simulations using dissipative particle dynamics (DPD) coupled with a reaction model to reproduce reaction-induced phase separation. The curing reaction properties of the thermoset resin were determined by ab initio quantum chemical calculations, and the DPD parameters were determined by all-atom molecular dynamics simulations. This enabled mesoscopic simulations, including reactions that reflect the intrinsic material properties. The effects of the thermoplastic resin concentration, molecular weight, and curing conditions on the phase-separation morphology were evaluated, and the cure shrinkage and stiffness of each cured resin were confirmed to be consistent with the experimental trends. Furthermore, the local strain field under tensile deformation was visualized, and the inhomogeneous strain field caused by the phase-separated structures of two resins with different stiffnesses was revealed. These results can aid in understanding the toughening properties of thermoplastic additives at the molecular level.

Uncovering the Mechanism of Size Effect on the Thermomechanical Properties of Highly Cross-Linked Epoxy Resins.

Epoxy resins are widely used as matrix resins, especially for carbon-fiber-reinforced plastic, due to their outstanding physical and mechanical properties. To date, most research into cross-linking processes using simulation has considered only a distance-based criterion to judge the probability of reaction. In this work, a new algorithm was developed for use with the large-scale atomic/molecular massively parallel simulator (LAMMPS) simulation package to study the cross-linking process; this new approach combines both a distance-based criterion and several kinetic criteria to identify whether the reaction has occurred. Using this simulation framework, we investigated the effect of model size on predicted thermomechanical properties of three different structural systems: diglycidyl ether of bisphenol A (DGEBA)/4,4'-diaminodiphenyl sulfone (4,4'-DDS), DGEBA/diethylenetriamine (DETA), and tetraglycidyl diaminodiphenylmethane (TGDDM)/4,4'-DDS. Derived values of gel point, volume shrinkage, and cross-linked resin density were found to be insensitive to model size in these three systems. Other thermomechanical properties, i.e., glass-transition temperature, Young's modulus, and yield stress, were found to reach stable values for systems larger than ∼40 000 atoms for both DGEBA/4,4'-DDS and DGEBA/DETA. However, these same properties modeled for TGDDM/4,4'-DDS did not stabilize until the system size reached 50 000 atoms. Our results provide general guidelines for simulation system size and procedures to more accurately predict the thermomechanical properties of epoxy resins.

Heat conduction in chain polymer liquids: molecular dynamics study on the contributions of inter- and intramolecular energy transfer.

In this paper, the molecular mechanisms which determine the thermal conductivity of long chain polymer liquids are discussed, based on the results observed in molecular dynamics simulations. Linear n-alkanes, which are typical polymer molecules, were chosen as the target of our studies. Non-equilibrium molecular dynamics simulations of bulk liquid n-alkanes under a constant temperature gradient were performed. Saturated liquids of n-alkanes with six different chain lengths were examined at the same reduced temperature (0.7T(c)), and the contributions of inter- and intramolecular energy transfer to heat conduction flux, which were identified as components of heat flux by the authors' previous study [J. Chem. Phys. 128, 044504 (2008)], were observed. The present study compared n-alkane liquids with various molecular lengths at the same reduced temperature and corresponding saturated densities, and found that the contribution of intramolecular energy transfer to the total heat flux, relative to that of intermolecular energy transfer, increased with the molecular length. The study revealed that in long chain polymer liquids, thermal energy is mainly transferred in the space along the stiff intramolecular bonds. This finding implies a connection between anisotropic thermal conductivity and the orientation of molecules in various organized structures with long polymer molecules aligned in a certain direction, which includes confined polymer liquids and self-organized structures such as membranes of amphiphilic molecules in water.

Energy dissipation in non-isothermal molecular dynamics simulations of confined liquids under shear.

Energy is commonly dissipated in molecular dynamics simulations by using a thermostat. In non-isothermal shear simulations of confined liquids, the choice of the thermostat is very delicate. We show in this paper that under certain conditions, the use of classical thermostats can lead to an erroneous description of the dynamics in the confined system. This occurs when a critical shear rate is surpassed as the thermo-viscous effects become prominent. In this high-shear-high-dissipation regime, advanced dissipation methods including a novel one are introduced and compared. The MD results show that the physical modeling of both the accommodation of the surface temperature to liquid heating and the heat conduction through the confining solids is essential. The novel method offers several advantages on existing ones including computational efficiency and easiness of application for complex systems.

Heat Conduction Performance over a Poly(ethylene glycol) Self-Assembled Monolayer/Water Interface: A Molecular Dynamics Study.

This study examined the interfacial heat condition between a poly(ethylene glycol) (PEG) self-assembled monolayer (SAM) and water using molecular dynamics simulation. It was found that the PEG SAM has higher thermal boundary conductance (TBC) than the traditionally used alkane-based SAM. The TBC conditionally varied with the length of the PEG molecules, where interfacial thermal resistance was a key factor. Our results reveal that the TBC of the PEG SAM/water interface is greatly influenced by its structural properties rather than the matching of vibrational properties between the SAM terminal and water. The structural analysis shows that the water structure around the terminal oxygen atom of the SAM plays a vital role in controlling the TBC. In this study, the concept of free volume has also been exploited, and the result suggests that the reduction of the free volume fraction accommodates a higher TBC. The model was precisely validated against experimental data by calculating the tilt angle and dihedral angle of the PEG SAM, the persistence length of the PEG chain in the water medium, and the sulfur position of the PEG SAM headgroup on the gold surface using quantitative scanning transmission electron microscopy image simulation.

A molecular dynamics study on heat conduction characteristics in DPPC lipid bilayer.

In this paper, nonequilibrium molecular dynamics simulations were performed on a single component 1,2-dipalmitoyl-sn-glycero-3-phosphatidylcholine lipid bilayer in order to investigate the thermal conductivity and its anisotropy. To evaluate the thermal conductivity, we applied a constant heat flux to the lipid bilayer along and across the membrane with ambient water. The contribution of molecular interaction to the heat conduction was also evaluated. Along the bilayer plane, there is little transfer of thermal energy by the interaction between lipid molecules as compared with the interaction between water molecules. Across the bilayer plane, the local thermal conductivity depends on the constituents (i.e., water, head group, and tail group of lipid molecule) that occupy the domain. Although the intramolecular transfer of thermal energy in the tail groups of lipid molecules works efficiently to promote high local thermal conductivity in this region, the highest thermal resistance appears at the center of lipid bilayer where acyl chains of lipid molecules face each other due to a loss of covalent-bond and low number density. The overall thermal conductivities of the lipid bilayer in the directions parallel and perpendicular to the lipid membrane have been compared, and it was found that the thermal conductivity normal to the membrane is higher than that along the membrane, but it is still smaller than that of bulk water.

Membrane attachment facilitates ligand access to the active site in monoamine oxidase A.

Monoamine oxidase membrane enzymes are responsible for the catalytic breakdown of extra- and intracellular neurotransmitters and are targets for the development of central nervous system drugs. We analyzed the dynamics of rat MAOA by performing multiple independent molecular dynamics simulations of membrane-bound and membrane-free forms to clarify the relationship between the mechanics of the enzyme and its function, with particular emphasis on the significance of membrane attachment. Principal component analysis of the simulation trajectories as well as correlations in the fluctuations of the residues pointed to the existence of three domains that define the global dynamics of the protein. Interdomain anticorrelated movements in the membrane-bound system facilitated the relaxation of interactions between residues surrounding the substrate cavity and induced conformational changes which expanded the active site cavity and opened putative pathways for substrate uptake and product release. Such events were less pronounced in the membrane-free system due to differences in the nature of the dominant modes of motion. The presence of the lipid environment is suggested to assist in decoupling the interdomain motions, consistent with the observed reduction in enzyme activity under membrane-free conditions. Our results are also in accordance with mutational analysis which shows that modifications of interdomain hinge residues decrease the activity of rat MAOA in solution.

Application of MDGRAPE-3, a special purpose board for molecular dynamics simulations, to periodic biomolecular systems.

We describe the application of a special purpose board for molecular dynamics simulations, named MDGRAPE-3, to the problem of simulating periodic bio-molecular systems. MDGRAPE-3 is the latest board in a series of hardware accelerators designed to calculate the nonbonding long-range interactions much more rapidly than normal processors. So far, MDGRAPEs were mainly applied to isolated systems, where very many nonbonded interactions were calculated without any distance cutoff. However, in order to regulate the density and pressure during simulations of membrane embedded protein systems, one has to evaluate interactions under periodic boundary conditions. For this purpose, we implemented the Particle-Mesh Ewald (PME) method, and its approximation with distance cutoffs and charge neutrality as proposed by Wolf et al., using MDGRAPE-3. When the two methods were applied to simulations of two periodic biomolecular systems, a single MDGRAPE-3 achieved 30-40 times faster computation times than a single conventional processor did in the both cases. Both methods are shown to have the same molecular structures and dynamics of the systems.