Interfacial Characteristics of Boron Nitride Nanosheet / Epoxy Resin Nanocomposites: A Molecular Dynamics Simulation

: The interface between nanoﬁllers and matrix plays a key role in determining the properties of nanocomposites, but the interfacial characteristics of nanocomposites such as molecular structure and interaction strength are not fully understood yet. In this work, the interfacial features of a typical nanocomposite, namely epoxy resin (EP) ﬁlled with boron nitride nanosheet (BNNS) are investigated by utilizing molecular dynamics simulation, and the e ﬀ ect of surface functionalization is analyzed. The radial distribution density (RDD) and interfacial binding energy (IBE) are used to explore the structure and bonding strength of nanocomposites interface. Besides, the interface compatibility and molecular chain mobility (MCM) of BNNS / EP nanocomposites are analyzed by cohesive energy density (CED), free volume fraction (FFV), and radial mean square displacement (RMSD). The results indicate that the interface region of BNNS / EP is composed of three regions including compact region, bu ﬀ er region, and normal region. The structure at the interfacial region of nanocomposite is more compact, and the chain mobility is signiﬁcantly lower than that of the EP away from the interface. Moreover, the interfacial interaction strength and compatibility increase with the functional density of BNNS functionalized by CH 3 –(CH 2 ) 4 –O– radicals. These results adequately illustrate interfacial characteristics of nanocomposites from atomic level.


Introduction
Epoxy resin (EP) is a kind of thermosetting material with excellent performances including electrical insulation, mechanical properties, and chemical stability, which have drawn wide application in power and industrial equipment [1,2]. However, with the rapid development of power and industrial fields, epoxy materials equipment with superior performances is required. One common approach is to add nanomaterials such as carbon nanotubes (CNT), graphene, oxides, nitrides, and their mixture into various polymers, which not only retains the original excellent performance of pristine polymer, but also greatly enhances electrical, mechanical, and thermal properties [3][4][5][6][7][8][9][10]. Especially, boron nitride nanosheet (BNNS) as nanofiller for epoxy composites has higher bandgap and better insulating performance when compared with graphene oxide (GO)/EP or graphene (GNP)/EP composites [11][12][13][14][15][16][17]. Since the physicochemical properties of nanomaterials and polymers are distinctly different, a clear interface is formed between them which plays a significant role in the macroproperties of nanocomposites [18][19][20]. Consequently, the knowledge and understanding of the interface characteristics are critical for design and manufacturing of novel polymer nanocomposites. The BNNS model is with size of 18.67 × 16.94 Å and includes 60 N atoms and 60 B atoms. The -O-(CH2)4-CH3 group is used to functionalize the BNNS [53]. The oxygen atom prefers to covalently bond with B atoms on BNNS, leaving one of the adjacent N atoms (after breakage of the B-N bond) saturated with H [54]. Two kinds of BNNS with different functionalization density were adopted in this simulation. One included four functional groups (4FBNNS), and another one included eight functional groups (8FBNNS). The BNNS and functionalized BNNSs are shown in Figure 2.

The Process of Crosslinking
The principle of cross-linking between EP and curing agent follows the two main reaction formulas, as shown in Figure 3. Two epoxy functional groups can theoretically react with one secondary amine functional groups. Each EP molecule has two epoxy functional groups, and each DETDA molecule includes two secondary amine functional groups. The cross-linked degree of the actual system is generally in the range of 80-95% [55]. Therefore, the cross-linked degree was selected as 85% in this study.  The BNNS model is with size of 18.67 × 16.94 Å and includes 60 N atoms and 60 B atoms. The -O-(CH 2 ) 4 -CH 3 group is used to functionalize the BNNS [53]. The oxygen atom prefers to covalently bond with B atoms on BNNS, leaving one of the adjacent N atoms (after breakage of the B-N bond) saturated with H [54]. Two kinds of BNNS with different functionalization density were adopted in this simulation. One included four functional groups (4FBNNS), and another one included eight functional groups (8FBNNS). The BNNS and functionalized BNNSs are shown in Figure 2.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 3 of 15 Then the cross-linked process of EP and curing agent is introduced in detail. Finally, the method of building compound models is illustrated.

Molecular Structure
In this work, diglycidyl ether of bisphenol A (DGEBA) and diethyl toluene diamine (DETDA) molecule were selected as EP and curing agent, as shown in Figure 1. The polymerization degree of EP was set as 0, which can reduce the calculating burden of the simulation and have no effect on the simulation results [52]. Atoms on the monomers were defined as reactive atoms by assigning a special name to the atoms (R1 and R2). The BNNS model is with size of 18.67 × 16.94 Å and includes 60 N atoms and 60 B atoms. The -O-(CH2)4-CH3 group is used to functionalize the BNNS [53]. The oxygen atom prefers to covalently bond with B atoms on BNNS, leaving one of the adjacent N atoms (after breakage of the B-N bond) saturated with H [54]. Two kinds of BNNS with different functionalization density were adopted in this simulation. One included four functional groups (4FBNNS), and another one included eight functional groups (8FBNNS). The BNNS and functionalized BNNSs are shown in Figure 2.

The Process of Crosslinking
The principle of cross-linking between EP and curing agent follows the two main reaction formulas, as shown in Figure 3. Two epoxy functional groups can theoretically react with one secondary amine functional groups. Each EP molecule has two epoxy functional groups, and each DETDA molecule includes two secondary amine functional groups. The cross-linked degree of the actual system is generally in the range of 80-95% [55]. Therefore, the cross-linked degree was selected as 85% in this study.

The Process of Crosslinking
The principle of cross-linking between EP and curing agent follows the two main reaction formulas, as shown in Figure 3. Two epoxy functional groups can theoretically react with one secondary amine functional groups. Each EP molecule has two epoxy functional groups, and each DETDA molecule includes two secondary amine functional groups. The cross-linked degree of the actual system is generally in the range of 80-95% [55]. Therefore, the cross-linked degree was selected as 85% in this study. Then the cross-linked process of EP and curing agent is introduced in detail. Finally, the method of building compound models is illustrated.

Molecular Structure
In this work, diglycidyl ether of bisphenol A (DGEBA) and diethyl toluene diamine (DETDA) molecule were selected as EP and curing agent, as shown in Figure 1. The polymerization degree of EP was set as 0, which can reduce the calculating burden of the simulation and have no effect on the simulation results [52]. Atoms on the monomers were defined as reactive atoms by assigning a special name to the atoms (R1 and R2). The BNNS model is with size of 18.67 × 16.94 Å and includes 60 N atoms and 60 B atoms. The -O-(CH2)4-CH3 group is used to functionalize the BNNS [53]. The oxygen atom prefers to covalently bond with B atoms on BNNS, leaving one of the adjacent N atoms (after breakage of the B-N bond) saturated with H [54]. Two kinds of BNNS with different functionalization density were adopted in this simulation. One included four functional groups (4FBNNS), and another one included eight functional groups (8FBNNS). The BNNS and functionalized BNNSs are shown in Figure 2.

The Process of Crosslinking
The principle of cross-linking between EP and curing agent follows the two main reaction formulas, as shown in Figure 3. Two epoxy functional groups can theoretically react with one secondary amine functional groups. Each EP molecule has two epoxy functional groups, and each DETDA molecule includes two secondary amine functional groups. The cross-linked degree of the actual system is generally in the range of 80-95% [55]. Therefore, the cross-linked degree was selected as 85% in this study.  The structure of the above molecular model was optimized, then the energy minimization calculation was performed with Forcite module until the geometry of them reached convergence. After that, a stoichiometric mixture of 64 DEGBA molecules and 32 DETDA molecules was established by employing the amorphous cell module. Next, the above model was balanced and relaxed until the structure and energy became quite steady. This model was balanced once by canonical ensemble (NVT), in which the temperature of ensemble maintained 300K, the step of time was 0.5 fs, and this process kept 100 ps. Then, the model was subsequently relaxed twice by isothermal-isobaric ensemble (NPT) with pressure 1.0 GPa and 0.1 MPa (simulation time 100 ps and time step 0.5 fs). Afterwards, the structure of this model was optimized and balanced again. Through the above process, the uncross-linked model was obtained. The force field of this molecular dynamics (MD) simulation is universal force field (UFF), and the method of thermostat and barostat is Andersen and Berendsen, respectively [56].
Based on the established model, the cross-linked process was performed by using Perl scripting language according to the procedures in [55]. The main steps are described with four steps.
Step 1: Define the values of simulation parameters. The initial and final cross-linked cutoff distances are 3.5 Å and 12.5 Å, respectively. The step length is selected as 0.5 Å, and the reaction temperature keeps 300K. The number of iterations is 3. The ultima cross-linked degree is set as 85%.
Step 2: The distance between the unreacted atoms (R1, R2) is identified. If the reaction atoms of epoxy groups are located within the reaction cutoff distance, the ring of epoxy groups is opened. Then the reaction atoms of epoxy groups react with unreacted secondary amine to form new C-N chemical bonds. The reactive principle is shown in Figure 3.
Step 3: If three times of iterations are finished or there are no unreacted atoms within the reaction cutoff distance, the cutoff distance increases with steps of 0.5 Å. The geometry is optimized and balanced until this model reaches the most stable status. The simulating condition of this process is the same as that of building the uncross-linked model.
Step 4: The program repeats Steps 2-3, and when the cross-linked degree or the maximum reaction cut-off distance is reached, this process stops.
After above cross-linked process, the structure of this model was optimized and balanced again. The built cross-linked EP model is with density of 1.12 g/cm 3  The structure of the above molecular model was optimized, then the energy minimization calculation was performed with Forcite module until the geometry of them reached convergence. After that, a stoichiometric mixture of 64 DEGBA molecules and 32 DETDA molecules was established by employing the amorphous cell module. Next, the above model was balanced and relaxed until the structure and energy became quite steady. This model was balanced once by canonical ensemble (NVT), in which the temperature of ensemble maintained 300K, the step of time was 0.5 fs, and this process kept 100 ps. Then, the model was subsequently relaxed twice by isothermal-isobaric ensemble (NPT) with pressure 1.0 GPa and 0. 1 MPa (simulation time 100 ps and time step 0.5 fs). Afterwards, the structure of this model was optimized and balanced again. Through the above process, the uncross-linked model was obtained. The force field of this molecular dynamics (MD) simulation is universal force field (UFF), and the method of thermostat and barostat is Andersen and Berendsen, respectively [56].
Based on the established model, the cross-linked process was performed by using Perl scripting language according to the procedures in [55]. The main steps are described with four steps.
Step 1: Define the values of simulation parameters. The initial and final cross-linked cutoff distances are 3.5 Å and 12.5 Å, respectively. The step length is selected as 0.5 Å, and the reaction temperature keeps 300K. The number of iterations is 3. The ultima cross-linked degree is set as 85%.
Step 2: The distance between the unreacted atoms (R1, R2) is identified. If the reaction atoms of epoxy groups are located within the reaction cutoff distance, the ring of epoxy groups is opened. Then the reaction atoms of epoxy groups react with unreacted secondary amine to form new C-N chemical bonds. The reactive principle is shown in Figure 3.
Step 3: If three times of iterations are finished or there are no unreacted atoms within the reaction cutoff distance, the cutoff distance increases with steps of 0.5 Å. The geometry is optimized and balanced until this model reaches the most stable status. The simulating condition of this process is the same as that of building the uncross-linked model.
Step 4: The program repeats Steps 2-3, and when the cross-linked degree or the maximum reaction cut-off distance is reached, this process stops.
After above cross-linked process, the structure of this model was optimized and balanced again. The built cross-linked EP model is with density of 1.12 g/cm 3

Compound Models
In order to explore the interfacial characteristics between BNNS and EP, BNNS/EP and FBNNS/EP compound models were built. Firstly, BNNS and FBNNS were put into empty cells. Then 64 DGEBA and 32 DETDA molecules were packed into the above cells. Next, the optimization, relaxation, and balance processes were performed with the same procedures as in the previous section. Considering the effect of periodic boundary on the interface, a vacuum layer of 10 Å was added to the upper and lower sides of the pure EP, BNNS/EP, and FBNNS/EP along the z direction, as shown in Figure 5.

Compound Models
In order to explore the interfacial characteristics between BNNS and EP, BNNS/EP and FBNNS/EP compound models were built. Firstly, BNNS and FBNNS were put into empty cells. Then 64 DGEBA and 32 DETDA molecules were packed into the above cells. Next, the optimization, relaxation, and balance processes were performed with the same procedures as in the previous section. Considering the effect of periodic boundary on the interface, a vacuum layer of 10 Å was added to the upper and lower sides of the pure EP, BNNS/EP, and FBNNS/EP along the z direction, as shown in Figure 5.
Based on above models, they were imported to LAMMPS, and the ReaxFF was employed. The element parameters of EP in ReaxFF were set according to Ref. [57,58], and those of BN were set according to Ref. [59,60]. Besides, minimize command was used to optimize structure, and iterations were 10,000. Fix nvt command was used to equilibrate the system, and the simulation conditions were temperature at 300K, time step of 0.5 fs, and total time of 100 ps. The Gaussian distribution was employed for initial velocity of the atoms. What's more, the parameters of interfacial characteristic were also calculated in LAMMPS based on the above models.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 5 of 15 Based on above models, they were imported to LAMMPS, and the ReaxFF was employed. The element parameters of EP in ReaxFF were set according to Ref. [57,58], and those of BN were set according to Ref. [59,60]. Besides, minimize command was used to optimize structure, and iterations were 10,000. Fix nvt command was used to equilibrate the system, and the simulation conditions were temperature at 300K, time step of 0.5 fs, and total time of 100 ps. The Gaussian distribution was employed for initial velocity of the atoms. What's more, the parameters of interfacial characteristic were also calculated in LAMMPS based on the above models.

Results and Discussion
The interface is the transitional region from the nanomaterial to the substrate, and its physicochemical properties are different from those of the matrix. In this section, the interfacial characteristic of BNNS/EP composites is discussed based on the established models. Firstly, the interfacial thickness and structure are obtained by calculating the RDD. After that, the IBE and interfacial shear strength (ISS) of BNNS/EP composites are computed, which can characterize the interaction strength of the interface. Moreover, the difference of MCM between all composite models is gained by simulating their RMSD. In order to further analyze the reason for changes of chain mobility at the interface, the FFV and CED of different composites are calculated, and the interfacial compatibility is also discussed.

Interfacial Structure
The structural feature is described by the RDD of BNNS/EP composites that is defined as the average per-atom densities in slices with a specified thickness along a specified axis. The RDD is computed by using the fix/ave/chunk command in LAMMPS [61]. The interfacial position and thickness are characterized by analyzing RDD. The RDD along the z-axis of four kinds of composites was computed, in which the slice thickness was chosen as 0.6 Å, as shown in Figure 6.

Results and Discussion
The interface is the transitional region from the nanomaterial to the substrate, and its physicochemical properties are different from those of the matrix. In this section, the interfacial characteristic of BNNS/EP composites is discussed based on the established models. Firstly, the interfacial thickness and structure are obtained by calculating the RDD. After that, the IBE and interfacial shear strength (ISS) of BNNS/EP composites are computed, which can characterize the interaction strength of the interface. Moreover, the difference of MCM between all composite models is gained by simulating their RMSD. In order to further analyze the reason for changes of chain mobility at the interface, the FFV and CED of different composites are calculated, and the interfacial compatibility is also discussed.

Interfacial Structure
The structural feature is described by the RDD of BNNS/EP composites that is defined as the average per-atom densities in slices with a specified thickness along a specified axis. The RDD is computed by using the fix/ave/chunk command in LAMMPS [61]. The interfacial position and thickness are characterized by analyzing RDD. The RDD along the z-axis of four kinds of composites was computed, in which the slice thickness was chosen as 0.6 Å, as shown in Figure 6.
From the RDD curves of BNNS/EP composites, the high density of about 4.0 g/cm 3 is at 40 Å at the position of BNNS. According to the RDD of pure EP, the average density is 1.12 g/cm 3 , which is consistent with the characteristics of epoxy-based amorphous materials. The mass density near BNNS is significantly higher than that of pure EP, which can explain the existence of interface at nanocomposites, which is formed under the action of Van der Waals long-range attraction on the surface of the BNNS. Furthermore, the interfacial thicknesses of BNNS/EP, 4FBNNS/EP, and 8FBNNS/EP are about 8 Å, 9.5 Å, and 10.5 Å, respectively, as shown in Figure 6. Besides, the interface region of BNNS/EP is composed of three regions, namely compact region, buffer region, and normal region. The mass density of the compact region is 15~25% higher than that of the normal region. The mass density of the buffer region decreases with the distance from the BNNS increasing, but the mass density is still higher than that of the normal region. The mass density of the normal region is the same as that of pure EP. Considering the interphase density may be related to the nanomaterials-matrix interaction, interface interaction strength is analyzed in the next section to further explain this phenomenon. From the RDD curves of BNNS/EP composites, the high density of about 4.0 g/cm 3 is at 40 Å at the position of BNNS. According to the RDD of pure EP, the average density is 1.12 g/cm 3 , which is consistent with the characteristics of epoxy-based amorphous materials. The mass density near BNNS is significantly higher than that of pure EP, which can explain the existence of interface at nanocomposites, which is formed under the action of Van der Waals long-range attraction on the surface of the BNNS. Furthermore, the interfacial thicknesses of BNNS/EP, 4FBNNS/EP, and 8FBNNS/EP are about 8 Å, 9.5 Å, and 10.5 Å, respectively, as shown in Figure 6. Besides, the interface region of BNNS/EP is composed of three regions, namely compact region, buffer region, and normal region. The mass density of the compact region is 15~25% higher than that of the normal region. The mass density of the buffer region decreases with the distance from the BNNS increasing, but the mass density is still higher than that of the normal region. The mass density of the normal region is the same as that of pure EP. Considering the interphase density may be related to the nanomaterialsmatrix interaction, interface interaction strength is analyzed in the next section to further explain this phenomenon.

Interfacial Bonding Strength
The interfacial binding energy reflects the bonding strength and compatibility between the matrix and the nanomaterial, and the interfacial bonding strength of composite materials can be described using IBE and ISS. A greater binding energy means that the interface strength of the composites is higher and compatibility is better. The calculation method of IBE is described as follows. The total energy (E) of the molecular system is the sum of kinetic energy (Ekinetic) and potential energy (Epotential). Molecular potential energy includes intramolecular potential energy and intermolecular potential energy and interaction energy. The intermolecular potential energy (Enon-bond) includes Van der Waals energy (Evdw), electrostatic energy (Ecoul), and hydrogen bonding energy (EHbond).
kinetic potential The IBE of different BNNS/EP composites can be calculated by the Equation (3): where EBNNS, EEP, and ET are the potential energy of BNNS, EP, and BNNS/EP composite, respectively.

Interfacial Bonding Strength
The interfacial binding energy reflects the bonding strength and compatibility between the matrix and the nanomaterial, and the interfacial bonding strength of composite materials can be described using IBE and ISS. A greater binding energy means that the interface strength of the composites is higher and compatibility is better. The calculation method of IBE is described as follows. The total energy (E) of the molecular system is the sum of kinetic energy (E kinetic ) and potential energy (E potential ). Molecular potential energy includes intramolecular potential energy and intermolecular potential energy and interaction energy. The intermolecular potential energy (E non-bond ) includes Van der Waals energy (E vdw ), electrostatic energy (E coul ), and hydrogen bonding energy (E H-bond ).
The IBE of different BNNS/EP composites can be calculated by the Equation (3): where E BNNS , E EP , and E T are the potential energy of BNNS, EP, and BNNS/EP composite, respectively. The interfacial shearing strength is another way to characterize interfacial bonding strength. The ISS can be deduced by E pullout that is the energy required to completely extract the nanomaterials from the composites and equals to the value of E interface . The ISS is calculated by Equations (4) and (5).
where S, h, and L are the surface area, the width, and length of the BNNS, respectively. z is the coordinate of the drawing direction. τ i is the shear stress. Since BNNS of this paper is not a square, τ i in the different direction as X and Y are calculated and averaged, shown as Equation (6).
The potential energy curves of each model for 100 ps of the balance process were taken, as shown in Figure 7. It can be seen that the curves are basically stable, indicating that the model is already in equilibrium. Furthermore, the potential energy of the BNNS/EP composites is significantly higher than that of the pure EP. The IBE and ISS of the BNNS/EP composites with different functional densities were calculated, as shown in Figure 8. It can be found that the IBE of BNNS/EP, 4BNNS/EP, and 8BNNS/EP is −397 kcal/moL, −584 kcal/moL, and −692 kcal/moL, respectively, and corresponding ISS is 9.5 MPa, 13.9 MPa, and 17.7 MPa. In other words, the interfacial strength of BNNS/EP composites increases with the functional density rising. Furthermore, the higher IBE and ISS values of the interface imply that the bonding strength of the interface between BNNS and EP is larger, so that the region near BNNS is more compact within a certain range, and the mass density in the vicinity of BNNS is significantly higher than that of the matrix. Besides, existing researches [37,38,51,52] showed that interfacial compatibility is positively correlated with interfacial bonding strength, so the compatibility of BNNS/EP also becomes better with the increase of functional density according to Figure 8, which means that the functional group of CH 3 -(CH 2 ) 4 -O-is appropriate for BNNS. And the favorable interfacial compatibility will contribute to molecular entanglement near BNNS, which is beneficial for enhancing interfacial bonding strength.
where S, h, and L are the surface area, the width, and length of the BNNS, respectively. z is the coordinate of the drawing direction. τi is the shear stress. Since BNNS of this paper is not a square, τi in the different direction as X and Y are calculated and averaged, shown as Equation (6).
The potential energy curves of each model for 100 ps of the balance process were taken, as shown in Figure 7. It can be seen that the curves are basically stable, indicating that the model is already in equilibrium. Furthermore, the potential energy of the BNNS/EP composites is significantly higher than that of the pure EP. The IBE and ISS of the BNNS/EP composites with different functional densities were calculated, as shown in Figure 8. It can be found that the IBE of BNNS/EP, 4BNNS/EP, and 8BNNS/EP is −397 kcal/moL, −584 kcal/moL, and −692 kcal/moL, respectively, and corresponding ISS is 9.5 MPa, 13.9 MPa, and 17.7 MPa. In other words, the interfacial strength of BNNS/EP composites increases with the functional density rising. Furthermore, the higher IBE and ISS values of the interface imply that the bonding strength of the interface between BNNS and EP is larger, so that the region near BNNS is more compact within a certain range, and the mass density in the vicinity of BNNS is significantly higher than that of the matrix. Besides, existing researches [37,38,51,52] showed that interfacial compatibility is positively correlated with interfacial bonding strength, so the compatibility of BNNS/EP also becomes better with the increase of functional density according to Figure 8, which means that the functional group of CH3-(CH2)4-Ois appropriate for BNNS. And the favorable interfacial compatibility will contribute to molecular entanglement near BNNS, which is beneficial for enhancing interfacial bonding strength.

Molecular Chain Mobility
The MCM is an important parameter to characterize various properties of nanocomposites. The movement ability of molecular chain can be described by mean square displacement (MSD) in polymer. The MSD represents the change in the position of atoms relative to their initial space during

Molecular Chain Mobility
The MCM is an important parameter to characterize various properties of nanocomposites. The movement ability of molecular chain can be described by mean square displacement (MSD) in polymer. The MSD represents the change in the position of atoms relative to their initial space during the dynamics simulation, and it is calculated by Equation (7).
where r i (t) is the position of i atom at time t, and r i (0) the initial position of i atom.
A major distinction of the performance between the interface layer and the matrix is the difference in the atoms mobility. The MSD of the four models in NVT ensemble for 200 ps was calculated, as shown in Figure 9. It can be seen that the molecular mobility of the BNNS/EP composite is lower than that of the pure EP, indicating that the molecular chain mobility of nanocomposites is reduced by incorporating BNNS. In order to further characterize the mobility of molecular chains at different positions of models, the compute MSD/chunk command in LAMMPS was used to calculate the RMSD between different layers of nanocomposites. The space was divided into 50 layers along z direction, and each layer was set to be 1.2 Å in this paper. This result is shown in Figure 10. The MCM of the region with a distance to BNNS less than 10 Å is clearly lower, that is, the movement of atoms is restricted, which is a manifestation of the cage effect. The molecular mobility of the region with distance greater than 10 Å is almost same as that of EP. Through comparing the RMSD curves of four models, it can be found that the MSD value of the interface region is relatively volatility, and the value of RMSD near the interface region is significantly lower than other locations, which indicates that MCM is reduced in the interfacial region. The strong interaction between the molecules and reduction of chain movement space are the main causes of this phenomenon. To further explore the reasons for the decrease of MCM of BNNS/EP composites, the CED and valid free volume fraction (VFFV) of compound models were calculated, as shown in Figure 11. The intermolecular forces of polymers have a prominent influence on many properties of the nanocomposites and is usually characterized by CED. A higher CED represents a stronger interaction force between molecules. And the glass transition temperature and flexibility of the molecular chains are also related to the CED [23,24,33,36,62]. The CED value is the cohesive energy per unit volume. Since the cross-linked EP is a networked thermoset polymer, the cohesive energy and CED can be calculated by Equations (8) and (9). The intermolecular forces of polymers have a prominent influence on many properties of the nanocomposites and is usually characterized by CED. A higher CED represents a stronger interaction force between molecules. And the glass transition temperature and flexibility of the molecular chains are also related to the CED [23,24,33,36,62]. The CED value is the cohesive energy per unit volume.
Since the cross-linked EP is a networked thermoset polymer, the cohesive energy and CED can be calculated by Equations (8) and (9).
where E coh and V m represent cohesive energy and volume of composites, respectively. The intermolecular forces of polymers have a prominent influence on many properties of the nanocomposites and is usually characterized by CED. A higher CED represents a stronger interaction force between molecules. And the glass transition temperature and flexibility of the molecular chains are also related to the CED [23,24,33,36,62]. The CED value is the cohesive energy per unit volume. Since the cross-linked EP is a networked thermoset polymer, the cohesive energy and CED can be calculated by Equations (8) and (9). (9) where Ecoh and Vm represent cohesive energy and volume of composites, respectively. The CED of different kinds of nanocomposites were compared, as shown in Figure 11. The CED of the BNNS/EP composites is significantly higher than that of the pure EP, which indicates that the incorporation of BNNS enhances the interaction force between molecules of nanocomposites. That is, the movement of molecular chain is restricted with the increase of CED, which reduces the MCM. Moreover, it can be observed from the IBE and CED that the intermolecular interaction force is enhanced with the increase of functional density. As a result, the compatibility and flexibility of The CED of different kinds of nanocomposites were compared, as shown in Figure 11. The CED of the BNNS/EP composites is significantly higher than that of the pure EP, which indicates that the incorporation of BNNS enhances the interaction force between molecules of nanocomposites. That is, the movement of molecular chain is restricted with the increase of CED, which reduces the MCM. Moreover, it can be observed from the IBE and CED that the intermolecular interaction force is enhanced with the increase of functional density. As a result, the compatibility and flexibility of molecules become better, resulting in the molecular chain being easier to wrap around the BNNS [37,47]. These conclusions can also explain that the functional group of CH 3 -(CH 2 ) 4 -O-is beneficial for enhancing the interfacial compatibility between BNNS and EP.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 10 of 15 molecules become better, resulting in the molecular chain being easier to wrap around the BNNS [37,47]. These conclusions can also explain that the functional group of CH3-(CH2)4-O-is beneficial for enhancing the interfacial compatibility between BNNS and EP. Free volume refers to the space inside the material that is not occupied by molecules, which plays a crucial role in understanding the movement ability of atoms in the amorphous polymer. The free volume within the model was calculated using the Open Visualization Tool (OVITO) software. This method uses a probe with a specific radius to probe the polymer, and the probe radius was selected as 0.4 Å here [63]. Due to the presence of vacuum layer in the models, the calculation formula of VFFV is calculated from Equation (10). Free volume refers to the space inside the material that is not occupied by molecules, which plays a crucial role in understanding the movement ability of atoms in the amorphous polymer. The free volume within the model was calculated using the Open Visualization Tool (OVITO) software. This method uses a probe with a specific radius to probe the polymer, and the probe radius was selected as 0.4 Å here [63]. Due to the presence of vacuum layer in the models, the calculation formula of VFFV is calculated from Equation (10).
where V occupy is occupied volume of molecular or atomic, V free is free volume, and V vacuum is vacuum volume. Free volume is an important factor to affect the movement ability of molecular chains. The VFFV of four models is shown in Figure 11. It can be seen that the VFFV of the BNNS/EP nanocomposites is significantly lower than that of pure EP, and the VFFV values become smaller with the functional density rising. Consequently, the movement space of molecular chains in BNNS/EP composite is limited, leading to a decrease in chain mobility. Besides, this conclusion also illustrates that the incorporation of BNNS will reduce the free volume of composites and the functional group of CH 3 -(CH 2 ) 4 -O-is effective to modify the surface of BNNS.

MultiRegional Interface Model
In order to describe the interfacial characteristics of nanocomposites in detail, a multiregional interface model is presented based on the work in [23,33,47], as shown in Figure 12. The interfacial region is composed of compact region, buffer region, and normal region. The interfacial interaction strength and compatibility of nanocomposite are enhanced by filling BNNS, so that the molecular chains are easy to wrap around BNNS. Therefore, the mass density of the compact region is 15~25% higher than that of the normal region. Interaction strength reduces with the increase of distance from BNNS, which leads to the decrease of mass density in the buffer region compared with that of compact region. Furthermore, the decrease of free volume and the enhancement of interfacial interaction strength cause the molecular chain mobility to reduce. The molecular chain mobility in the compact region is significantly lower than that of other regions, and molecular chain mobility of the buffer region is higher than that of the compact region, but still lower than that of the normal region. Moreover, the molecular structure and chain mobility of the normal region is the same as the EP matrix. As a result, there are visible differences in characteristics between different regions of the nanocomposite interface. region. Moreover, the molecular structure and chain mobility of the normal region is the same as the EP matrix. As a result, there are visible differences in characteristics between different regions of the nanocomposite interface. In this work, BNNS is functionalized by CH3-(CH2)4-O-radicals, the interfacial interaction strength between EP and functionalized BNNS is enhanced compared with that of untreated BNNS/EP, which improves compatibility and reduces free volume of nanocomposite, so that the molecular chains are easier to twine around treated BNNS. Consequently, the mass density of the compact region and the buffer region is higher than that of untreated BNNS/EP. At the same time, the molecular chain mobility of functionalized BNNS/EP is lower than that of untreated BNNS/EP. These results show that the interfacial thickness, interaction strength, and compatibility of nanocomposite are enhanced by filling functionalized BNNS, and the interfacial characteristics such molecular structure and chain mobility are improved with the functional density increasing.

Conclusions
A thorough investigation on the features of composite interface is of significance for nanocomposites modification. Molecular dynamic simulation for cross-linked EP and BNNS/EP with different CH3-(CH2)4-O-functionalized density was performed to explore the interfacial characteristics, in which several parameters including RDD, IBE, ISS, RMSD, CED, and VFFV were analyzed. The results indicated that the interface between BNNS and EP is composed of three regions including compact region, buffer region, and normal region. The interfacial interaction strength and compatibility of nanocomposite are enhanced by filling BNNS, so the mass density of the compact region is 15~25% higher than that of the normal region. The interaction strength reduces with the increase of distance from BNNS, which causes the decrease of mass density in the buffer region In this work, BNNS is functionalized by CH 3 -(CH 2 ) 4 -O-radicals, the interfacial interaction strength between EP and functionalized BNNS is enhanced compared with that of untreated BNNS/EP, which improves compatibility and reduces free volume of nanocomposite, so that the molecular chains are easier to twine around treated BNNS. Consequently, the mass density of the compact region and the buffer region is higher than that of untreated BNNS/EP. At the same time, the molecular chain mobility of functionalized BNNS/EP is lower than that of untreated BNNS/EP. These results show that the interfacial thickness, interaction strength, and compatibility of nanocomposite are enhanced by filling functionalized BNNS, and the interfacial characteristics such molecular structure and chain mobility are improved with the functional density increasing.

Conclusions
A thorough investigation on the features of composite interface is of significance for nanocomposites modification. Molecular dynamic simulation for cross-linked EP and BNNS/EP with different CH 3 -(CH 2 ) 4 -O-functionalized density was performed to explore the interfacial characteristics, in which several parameters including RDD, IBE, ISS, RMSD, CED, and VFFV were analyzed. The results indicated that the interface between BNNS and EP is composed of three regions including compact region, buffer region, and normal region. The interfacial interaction strength and compatibility of nanocomposite are enhanced by filling BNNS, so the mass density of the compact region is 15~25% higher than that of the normal region. The interaction strength reduces with the increase of distance from BNNS, which causes the decrease of mass density in the buffer region compared with that of the compact region. Furthermore, the decrease of free volume and the enhancement of interfacial interaction strength of nanocomposite result in the decrease of molecular chain mobility in the compact region, and molecular chain mobility of the buffer region is higher than that of the compact region. Moreover, these features of the normal region are the same as for the EP matrix. Besides, the interfacial interaction strength of functionalized BNNS/EP is stronger compared with that of untreated BNNS/EP. As a result, the interfacial thickness, interaction strength, and compatibility between BNNS and EP are enhanced with the increase of functional density of BNNS, which also implies that the functional group of CH 3 -(CH 2 ) 4 -O-radicals can improve the surface properties of BNNS. This work provides a novel insight to investigate the interfacial characteristics of nanocomposites from atomic level.