Molecular Simulations of Thermomechanical Properties of Epoxy-Amine Resins.

All-atom molecular dynamics (MD) simulations were performed with the CHARMM force field to characterize various epoxy resins, such as aliphatic and bisphenol-based resins. A multistep cross-linking algorithm was established, and key properties such as density, glass temperature, and elastic modulus were calculated. A quantitative comparison was made and was proven to be in good agreement with experimental data, with average absolute deviations between experiments and molecular simulation comprised between 2% and 12%. Additional findings on structure-property relationships were highlighted such as the effect of the cross-linking rate and oligomerization of the resin.

Introduction

Epoxies are some of the most prominent thermosetting polymers valued nowadays due to their wide range of applications. This class of polymers is of interest for instance for structural components in aerospace applications, as well as electronics packaging or various types of coatings. This is made possible by their interesting versatility. The final material can be easily modulated by the chemical structure of the epoxy or the nature of the hardener, responsible for the striking thermal and mechanical properties of the resulting material.[1]

Therefore, it is of great importance to be able to characterize them and comprehend the influence of the structural network on their properties. Aside from the usual experimental testing, computational chemistry has proven to be a great asset in this task. Indeed, predictive simulations are playing an increasing role in understanding the atomistic origin of the properties of polymers as well as in speeding up the process of selecting materials, thus complementing experiments in the search for tailored materials. However, one of the greatest challenges is to be able to reproduce reliably the physical, mechanical, and thermodynamic properties of such materials from atomistic models. Molecular dynamics (MD) has proven to be a powerful method to generate equilibrated atomistic configurations on which macroscopic properties are averaged.[2]

One of the first steps to reproduce thermosetting polymers in molecular simulations is to be able to cross-link the monomers, making it possible to characterize the polymeric material. Two approaches are typically used for the cross-linking of the monomers: a “back-mapping” procedure[3−6] and a “reactive” one.[7−15] The former is based on the cross-linking of coarse-grained monomers and then the all-atom model is later back-mapped in lieu of the coarse grains, while the latter already considers an all-atom model and the cross-linking is done between the experimental reactive pairs. These cross-linking methods can also be single or multiple steps but exhibit the same topology overall.[16]

Regarding the epoxy systems investigated by molecular simulations, the most extensive research effort was found in resins resulting from the cross-linking between the diglycidyl ether of bisphenol F (DGEBF) and diethyl toluenediamine (DETDA),[9,11,12,17,18] triethylenetetramine (TETA)[10,19−21] or diethylenetriamine (DETA),[22] especially in the fields of modern aeronautics. The cross-linking of the diglycidyl ether of bisphenol A (DGEBA) resin was also widely studied with different hardeners, such as isophorone diamine (IPDA),[23−25] trimethylene glycol di-p-aminobenzoate (TMAB),[26] diethyltoluenediamine (DETDA),[4,27−29] triethylenetetramine (TETA),[27] ethylenediamine (EDA),[7] diaminodiphenyl sulfone (DDS),[30] methylenedianiline (MDA),[31] poly(oxopropylene) diamines (POP),[16,32,33] 4,4′-methylenebis(cyclohexylamine) (MCA),[33] diethylenetriamine (DETA),[13,34−36] and polyetheramine JEFFAMINE D-230.[37−40]

To the best of our knowledge, novel epoxy resins based on the diglycidyl ether of vanillyl alcohol (DGEVA), diglycidyl ether of 1,4-butanediol (DGEBU), triglycidyl ether of phloroglucinol (PHTE), 4,4′-methylenebis(2,6-diethylaniline) (MDEA), m-xylylenediamine (MXDA), 4,4′-methylenebis(2-isopropyl-6-methylaniline) (MMIPA), tetraethylenepentamine (TEPA), and 4,7,10-Trioxa-1,13-tridecane diamine (TTDA) have not yet been studied by molecular simulations. Due to the broad variety of chemical structures, a comparison between experimental and simulation results on several key properties such as density, glass transition temperature, and elastic moduli is made possible.

Through this, an evaluation of the performance and the reliability of general atomistic force fields and the capacity of simulation methods to reproduce the thermomechanical properties of these resins is carried out. Indeed, predicting the physical, mechanical, and thermodynamic properties of resins by using full atomistic models remains challenging because of the complexity of the interactions involved and the algorithm needed to post-treat the data. Such molecular simulations can open the way to interpreting the thermomechanical properties of resins on the basis of their structure and the molecular interactions involved.

In order to establish this comparison, a cross-linking algorithm was first made based on previous studies using the reactive method.[7−15] This was followed by several procedures to ensure the reliability of the properties obtained, thus making it possible to enlighten key properties of materials and their relationship to their structure. To some extent, this study also aims at evaluating the transferability of the molecular models on resins not yet investigated by molecular simulations.

The paper is organized as follows. Section describes the experimental and computational procedures. Section presents the main results of this work in terms of density, glass transition temperature, and Young’s modulus of various epoxy-amine resins. We conclude in Section .

Experimental and Simulation Methods

Experimental Details

Hardeners and Resins Studied

Various epoxy resins are studied, exhibiting a large panel of thermomechanical properties. These resins are widely used in industry. The formulations presented in this study include the compounds DGEBA, DGEBF, DGEBU, PHTE, and DGEVA as well as the hardeners IPDA, DETA, EDA, TEPA, TTDA, MDEA, MMIPA, and MXDA. DGEVA and PHTE can be derived from bio-based raw materials and were synthesized by Specific Polymers. These components are presented in Figure . The average number of repeating units (n) of DGEBA (Sika, Biresin CR170) was determined by 1H NMR (Figure ).

Figure 1

Molecular representation of the epoxy resins (a–e) and of the hardeners (f–m). Carbon, oxygen, nitrogen, and hydrogen atoms are shown in gray, red, blue, and white, respectively. (n, o) Snapshots of cross-linked systems, where orange and gray phases represent respectively the epoxy resins and hardeners.

Figure 2

Depiction of the meaning of n in terms of the oligomerization of DGEBA.

Table gives for the epoxies studied the degree (n) of oligomerization of the resin, the resin:hardener ratio, the molecular weight of the constitutive repeat unit, and the empirical density (detailed in the following subsection) as well as the number of atoms present and the size of the box for the simulations.

Table 1

Description of the Characteristics of Resins and Hardeners Studied, Breaking down the Chemical Structure of the Repeat Unit: Composition, Mass, and Expected Densitya

resin	hardener	n	R:H ratio	mass of CRU (g/mol)	density (g/cm3) (eq2)	no. of atoms	side of the box (Å)
DGEBU	IPDA	0.32	2:1	574	1.040	7938	41.9
DGEBA	DETA	0.32	5:2	2368	1.141	11267	48.7
DGEBA	EDA	0.32	2:1	924	1.145	10919	47.8
DGEBA	IPDA	0.32	2:1	1034	1.117	12701	50.5
DGEBA	MMIPA	0.32	2:1	1174	1.128	14240	52.6
DGEBA	MXDA	0.32	2:1	1000	1.152	11729	49.1
DGEBA	TEPA	0.32	7:2	3404	1.137	11639	49.2
DGEBA	TTDA	0.32	2:1	1084	1.129	13106	50.8
DGEBF	IPDA	0	2:1	1008	1.140	12150	49.8
DGEVA	DETA	0	5:2	1616	1.174	11500	48.5
DGEVA	IPDA	0	2:1	734	1.134	8910	44.1
DGEVA	MDEA	0	2:1	874	1.146	12900	50.3
DGEVA	MXDA	0	2:1	700	1.187	9800	45.7
PHTE	IPDA	0	3:4	1686	1.119	12900	50.3
PHTE	MDEA	0	3:4	2106	1.136	15750	54.1
PHTE	MXDA	0	3:4	1584	1.186	11100	48.0

Calculated from eq as the number of atoms present and the size of the box for the simulations.

Determination of the Experimental Density and Glass Transition Temperature

The measurement of density has been performed on materials fully cured using a Ohaus density determination kit for a scale. Density determination is performed by Archimedes’ principle. The density of the material is determined with the aid of a liquid whose density ρ0 is known (water was used as an auxiliary liquid here). The solid is weighed in air (A) and then in the auxiliary liquid (B). The density ρ can be calculated from the two weights as follows ρ = [A/(A – B)](ρ0 – ρ) + ρL with ρ0 being the density of the auxiliary liquid and ρL the air density (0.0012 g/cm3). The temperature of the auxiliary liquid was measured and taken into account in all density determinations.

DSC analyses were performed by using a Mettler Toledo DSC instrument. Samples having masses of approximately 10 mg were placed in 100 μL aluminum pans with pierced lids under a nitrogen atmosphere. The applied heating rate was 10 °C min–1 under a nitrogen atmosphere. The thermal behavior of the samples was investigated using two repeated heating–cooling cycles. The following thermal procedure was used: first ramp from −60 to 250 °C and then from 250 to 25 °C. The glass transition temperature (Tg) was determined from the second heating curve.

Computational Details

The CHARMM general force field (CgenFF)[41−43] is used to describe the intra- and intermolecular interactions of the thermosetting polymers. This class I force field has the following functional form:For completeness, it is worth noting that force field parameters are usually denoted using σ rather than Rmin, the relationship between the two being . Additionally, the CHARMM force field presents an additional component in the form of the Urey–Bradley term, which consists of a harmonic potential as a function of the distance between the nonbonded atoms 1 and 3 of a 1–2–3 angle.

The Lennard–Jones parameters between pairs of different atoms are obtained from the Lorentz–Berthelot combination rules, in which ϵ values are calculated by using the geometric mean and Rmin values by using the arithmetic mean. The cutoff distance for the van der Waals interactions is set to 8.0 Å, while the particle–particle particle mesh (PPPM) method with an accuracy of 1.0 × 10–4 and a cutoff of 8 Å is considered for the long-range electrostatic interactions. For both the Lennard–Jones and Coulombic interactions, a force-switching function was applied between 8 and 10 Å.

Establishing the Cross-Linking Algorithm

In the scope of this work, a multistep reactive cross-linking algorithm was used and is described in Figure , where the statistical NVT and NPT ensembles refer to the canonical and isothermal–isobaric ensembles, respectively. It includes the following steps, where the sequence of 3 + 4 + 5 steps was repeated until a predefined conversion rate was acquired (here 95% in accordance with the experimental conversion obtained after curing).This cross-linking algorithm was designed without regard to the actual polymerization reaction kinetics. Experimentally, the curing temperature of such systems is below 200 °C. However, here the cross-linking was done partially at 627 °C (above the thermal decomposition temperature). The need for a higher temperature was driven by the appeal of a quick rearrangement of the atoms within the span of the MD simulations. Additionally, a buck pair_style was used for the unreacted pairs. Only the attractive term of the Lennard–Jones potential (1/r6) was used. Reactive pairs are thus attracted to one another in order to facilitate bond creations. Without this potential between reactive pairs, the acquisition of a highly cross-linked resin as depicted in this work would be very time-consuming. Additionally, the cutoff range was chosen in agreement with Figure . Radial distribution functions, shown in Figure , were calculated on an uncured mixture of DGEBA+IPDA and DGEVA+MXDA between the reactive carbon of the epoxy resin and the reactive nitrogen of the amine hardener. The first maxima are close to 3.7 and 4.1 Å, respectively. This is why a cutoff of 2–3 Å is deemed relevant when it comes to the compromise between computational efficiency and the resulting cross-linked resins. Too high of a cutoff guarantees faster cross-linking but can also lead to the formation of metastable resins where bonds cannot come back to their equilibrium value (∼1.45 Å).

Figure 3

Multistep reactive cross-linking algorithm.

Figure 4

Radial distribution functions of uncured resin DGEBA+IPDA (red) and DGEVA+MXDA (green) between the reactive carbon of the epoxy resin and the reactive nitrogen of the amine hardener.

Stoichiometric packing of the unreacted monomers into a box of predefined density using Packmol.[44] The parameters (type, bond, angle, dihedral, charges) are those of the reacted mixture.

Relaxation of the initial mixture with an NPT ensemble during 500 ps at 300 K.

NPT simulation at 300 K to relax the system.

NVT simulation at 900 K to allow the quick movement of the atoms while searching for unreacted pairs within a defined cutoff distance. At each time step, the interatomic distances between unreacted pairs are calculated. The simulation was stopped once a pair was found within 2 Å.

Creation of the bonds for two pairs that are in agreement with the distance conditions. Additionally, a second cutoff of 3 Å was applied to allow the formation of multiple bonds per loop and to quicken the process.

Updating of the parameters for the pairs that did not react.

Quick NPT relaxation at 300 K on the fully cross-linked system during 1 ns.

NPT simulation during 2 ns at 600 K to further relax the system, following a ramp from 300 to 600 K for 1 ns.

Last NPT simulation during 2 ns at 300 K to obtain a fully relaxed system, following a ramp from 600 to 300 K for 1 ns.

Experimentally, primary amines go through an addition reaction with the epoxide groups to form by ring-opening a hydroxyl group and a secondary amine. The secondary amine can further react with an epoxide to form a tertiary amine and an additional hydroxyl group. A hydroxyl-epoxy addition can also occur but less frequently. This last step was ignored during the simulated cross-linking, as well as the difference in reactivity of the primary and secondary amines.

All potential reactive sites were activated equally by removing the hydrogen atoms from the amine of the curing agent and hydrating the oxygen atoms from the epoxide and thus opening up the oxirane. The force field parameters (bonds, angles, dihedrals, impropers, partial charges) were thus adjusted to take into account this reacted structure. This change was revoked for the molecules that did not react at the end of the cross-linking procedure (see step 6 in the list).

The sequence 7 + 8 + 9 in the list was the postcure annealing process performed to fully relax the system. Using a high-temperature annealing period after cross-linking increases the density due to the fact that this addition of thermal energy can disrupt metastable configurations and allow finding lower energy configurations closer to global minima. It was successfully demonstrated that the thermal agitation produced higher density configurations, even though the original densities had already plateaued with the first NPT equilibration.[22]

Eight simulation boxes were hence created for each system containing the reacted form of the liquid mixture. They differ in their structure and the initial positioning of the molecules as to produce a statistical sample. The multistep cross-linking procedure depicted in Figure was used to polymerize the resins at a conversion rate of 95% of the reactive end groups. The use of a high temperature of 900 K and a cutoff value of 2–3 Å grants a high rate of conversion while staying true to the interatomic distances. The full procedure can last between 6 and 20 ns.

To conclude on this methodological aspect, this method was employed to achieve the same degree of conversion as found experimentally within the time-limit constraints of the MD simulations. The accuracy of the structural, volumetric, and thermal properties of the resulting cross-linked polymer will be evaluated in the course of this paper.

Description of the MD Simulations

The velocity-Verlet integrator was used to integrate the equations of motion using a time step of 1 fs. The SHAKE algorithm was applied to constrain the length of covalent bonds to hydrogen atoms to their equilibrium values. In order to control the temperature and pressure for NVT and NPT simulations, a Nosé–Hoover thermostat and barostat were applied.[45] Unless stated otherwise, the MD simulations were performed at 300 K and 1 atm. The system was equilibrated during the first nanosecond, and data were acquired during the following nanosecond. Eight simulation boxes were used for each system. The molecular dynamics calculations were performed using the LAMMPS package.[46]

Results and Discussion

Measuring the Density

Before examining the density of cured resins, the uncured resins and hardeners were examined to appreciate the relevance of the CHARMM force field. Figure shows that the majority of the monomers are included in the 5% domain, with a warning concerning the hardener EDA.

Figure 5

Comparison of experimental and simulated densities of the uncured resins and hardeners. The dotted red lines correspond to a deviation of 5% from the experimental value.

Figure reports the values of density of different epoxy-amine resins obtained from experimental measurements, empirical calculations, and atomistic simulations. Results taken from other works are given inTable . An empirical power law was taken into account to provide a qualitative comparison for resins where the experimental density was not present at room temperature. This relation can be used to describe the density of linear and network polymers according to their average atomic mass (Ma) and constitutive repeat unit (CRU) as stated in eq (53)where k ≈ 320, ρ is in kg/m3, and Ma is in g/mol. Ma is the average atomic mass of the polymer. This empirical relationship remains valid for most organic polymers (1.1 < ρ < 1.4 g cm–3). Unfortunately, the maximum error of up to 10% endangers a quantitative comparison, thus requiring the knowledge of the experimental density for a more refined comparison.

Figure 6

Comparison of the simulated, empirical, and experimental densities. Results taken from other works are given inTable .

Table 2

Compiled Simulated and Experimental Data Sets for Cross-Linked Resins at 95% Confidencea

	density (g/cm3)		Tg (°C)		E (GPa)
system	exptl	sim	exptl	sim	exptl	sim
DGEBU/IPDA		1.0503	60[47]	558		2.2311
DGEBA/DETA	1.178	1.1113	119,[48] 132,[49] 142*	123,[34] 13310	2.5511[48]	2.5714
DGEBA/EDA		1.1253	70,[50] 118,[7] 115[51]	623, 117[7]		2.7915
DGEBA/IPDA	1.135	1.0845	150*, 155,[48] 159[49]	14921, 158[25]	2.655[48]	2.5613
DGEBA/MMIPA	1.118	1.0753		1367		2.1920
DGEBA/MXDA	1.185	1.1243		10110		3.4124
DGEBA/TEPA		1.1022	126[48]	12524	2.3413[48]	2.4323
DGEBA/TTDA		1.1162	69[47]	848		2.146
DGEBF/IPDA		1.1024		12910		1.9216
DGEVA/DETA		1.1792	79*	3210		2.5713
DGEVA/IPDA		1.1546	97,[52] 107*	898		3.9926
DGEVA/MDEA		1.1444	104*	9813		3.2318
DGEVA/MXDA		1.2193	74*	759		3.4124
PHTE/IPDA		1.1085	225*, #	22615		3.9933
PHTE/MDEA		1.1046	117*	1328		3.7129
PHTE/MXDA		1.2025	176*	17010		4.8818

Values with asterisks were obtained experimentally, while # indicates Tα. The subscripts give the accuracy of the last decimal(s): i.e., 1.0503 means 1.050 ± 0.003.

First, we observe that the agreement between empirical calculations (eq ) and molecular simulation is good, with an average absolute deviation of 2.2%. The predicted values by molecular simulations also reproduce the variations in densities among the various resins. Other molecular models have also been shown to reproduce accurately the density for DGEBA/IPDA[24,25] and DGEBA/DETA.[34] We can conclude that the CHARMM force field is able to yield simulated densities of epoxy-amine resins that match reasonably with experiments, even if this model is not exclusively designed for polymer and additives, unlike the COMPASS model. The good agreement between the different density values obtained by simulation and empirical calculations enables a validation of the choice of the force field and the methodology used to polymerize the cross-linked epoxy resins. The CHARMM force field can contribute for epoxy resins to the determination of a qualitative study but offers an error rate of up to 5.7% with the few experimental points. The simulated values reproduce well the variations between resins but give an underestimation of the density. On the other hand, the experimental density calculated with the Archimedes method can also present an underestimation of the real density of the polymer. Taking into account this bias, this first comparison of the density shows that the model and the methodology developed for the cross-linked epoxy resins perform well for this thermodynamic property. To go further in the validation of the methodology, we propose to extend the comparison to the glass transition temperature.

Assessing the Glass Transition Temperature

Modeling of glass transition temperature is based on the dilatometry principle; the density decreases with an increase in the temperature. The glass transition temperature is then characterized by a break in this slope, showing a change in the relaxation speed of the material. It separates the vitreous from the rubbery states. The modeling procedure consists of heating the polymer at Tg+200 K for 2 ns and cooling it at a rate of 20 K ns–1 down to Tg–200 K. This cooling is done stepwise to retrieve the density of the equilibrated system at each plateau of temperature. The value of the ramp is a good compromise between the rapidity of the simulation and the accuracy of the results. As previously observed by Godey,[54] a transient region persists between the rubbery and vitreous states. This is relevant to the enlargement of the glass transition domain in a confined environment.[55,56] Region II is delimited by two temperatures, T1 and T2, which correspond to a departure from the linear behavior of the free volume with respect to the temperature. Figure shows the three regions obtained, and the glass transition temperature is conventionally taken as the intercept of the linear fits of region I and III (denoted as the green dashed line).

Figure 7

Free volume of DGEBA+IPDA with respect to the temperature. The glass transition temperature is considered as the middle vertical line (green dashed line).

Usually, the coefficient of volumetric thermal expansion (CVTE) is adopted to determine the linear fit of the domains. Here, this coefficient is too noisy for the thermosets studied. A piecewise linear function is thus applied to fit the three domains accordingly; the script is provided by the Python library pwlf.[57] This brings some automation to the fit of the linear domains, often very complicated and subjective to do. For each system, eight independent simulations are performed and the change in free volume is measured as a function of temperature.

However, in order to compare the simulated to the experimental values, the relation of Williams–Landel–Ferry (eq ), which relies on the principle of time–temperature equivalence, is applied by using the empirical constants C1 = 17.44 and C2 = 51.6 K.It was deemed acceptable to use the universal constants because of the structure of the studied resins and the time constraint of using fitted constraints for each.[58] In the scope of this work, several systems have been done both experimentally and by molecular simulations; this is mainly the case for the bio-based resins (DGEVA and PHTE). First, the correlation between these experimental and simulated glass transition temperatures is acceptable, with an average absolute deviation of 12%. The difference can be as minimal as 0.4% for PHTE/IPDA and as important as 59.5% for DGEVA/DETA. Second, there is a disparity between experimental results if a comparison is made between values found in this work and those from references. As a matter of fact, a mean difference of ±11 K is found, with a disparity of up to ±27 K for the case of DGEBA/EDA. Nevertheless, the results in Figure show a general great correlation between simulated and experimental data. The average absolute deviation is about 7.0% without taking into account the three identified outliers and 12.2% with. Indeed, the simulated values of the glass transition temperature are well within the error range of the experiments, which confirms the interest in using molecular simulations to predict this property. The tendency is surely well-preserved even though generic WLF coefficients are applied. In cases where the prediction is less efficient, a tuning of those coefficients might be of interest. However, it will require additional time-consuming simulations of various cooling speeds for each system. The gap between experimental and simulated results regarding the DGEVA/DETA, DGEBA/EDA, and DGEBA/MMIPA systems has to be further investigated.

Figure 8

Comparison of simulated and experimental glass transition temperature Tg. Results taken from other works are given in Table .

Finding Young’s Modulus

Experimentally, Young’s modulus is typically obtained by traction testing. The results are widely dependent on the stress–strain applied during the experiment. Nevertheless, these types of tests are normed, and the classic speed used is on the order of mm/min. However, it is practically impossible to measure the modulus at experimental strain rates in an MD system. One way to bypass this constraint is by using an innovative method based on the fluctuation of local density[59] inside the solid epoxy materials at equilibrium. Consequently, an NVT simulation is performed during 2 ns at 300 K for each system. The structural factor S is determined from the atomistic positions and is a function of the wave vector q. The following extrapolation in eq is done to obtain the bulk and shear moduli, respectively B and G. The factor 2 is added here to take into consideration the two transversal directions compared to the expression established in ref (59) for a two-dimensional system (see the Supporting Information for a demonstration of the following expressions).

The first member of the preceding equations is obtained graphically by linearly fitting the structure factor S and acquiring the intercept, which corresponds to the limit when the wave vector q tends to zero. A visual representation of this graphical method is found in Figure . However, the property of interest is Young’s modulus. We suppose a homogeneous isotropic linear elastic material so that its elastic properties are uniquely determined by any two moduli. Thus, having two moduli grants access to any other of the elastic moduli. The following equation is used to determine Young’s modulus from the bulk and shear moduli computed previously:

Figure 9

Graphical extrapolation of the structure factors.

Figure depicts the simulated Young’s modulus determined here and a few experimental points taken from a previous experiment.[48] As can be seen, the simulation fits roughly the few experimental data points, which is promising. Knowing that the simulated values exhibit a significant correlation with an average absolute deviation of 2.7%, this seems like a fair assumption. This thus raises interest in the hope that this trend will be transferable to all modeled systems. However, a wider sample is needed in terms of experimental data sets and range of Young’s moduli to positively conclude.

Figure 10

Comparison of the simulated with the experimental Young’s modulus. Results taken from other works are given in Table .

Influence of the Molecular Mass Mn

One of the advantages of using MD simulations is to be able to straightforwardly customize the molecules used in the epoxy resins. As well as easily switching from one hardener to another, one can play with the possible oligomerization of the epoxies. The study done here enlightens the role that the molecular weight plays in the resulting thermomechanical properties of the resins. This investigation corresponds to a DGEBA+IPDA resin and includes six different oligomerizations of the DGEBA epoxy with values of n equal to 0.08, 0.32, 1, 3, 5, and 10 (see Figure ). Values of n equal to 0.08 and 0.32 were taken from commercial references Sicomin, SR GreenPoxy 28 and Sika, Biresin CR170, respectively. Their population of n = 0, n = 1, etc. were determined by 1H NMR. Otherwise, integer values of n were taken as a unique population. An additional experimental point at n = 5.3 is proposed, corresponding to the Epikote 1004 commercial reference. This sample’s density was experimentally measured, but a large amount of air was imprisoned inside the material, leading to an underestimation of the density.

A focus is made on the influence of the molecular mass on density, glass transition temperature, and Young’s modulus (see Figure ). Previous studies on DGEBA and Jeffamines show the link between the length of the cross-linker and the resulting properties.[39,40] Similarly, this relationship is also displayed here: an increase in molecular mass causes an increase of the density and a decrease of the glass transition temperature and Young’s modulus. This outcome can be explained first by the flexibility that the increase of the distance between cross-linking points brings, having a higher degree of oligomerization of the epoxies. By having this oligomerization, it is then easier to find a more relaxed structure, thus bringing the molecules closer together and increasing the density. The length of a cross-link is often one of the main focuses when it comes to tuning thermomechanical properties of materials. One could want to limit the brittleness of a material by lengthening a cross-link to bring some flexibility to a material, but it needs to be taken into account that this change is typically detrimental to the thermal stability. This study enlightens the dual effect of the cross-link length and how a middle ground needs to be reached when designing a new material, and how molecular simulations can help to that effect.

Figure 11

Influence of the molecular weight on the thermomechanical properties of a DGEBA/IPDA system.

Influence of the Conversion Rate

In the same spirit, a study on the influence of the conversion rate was done to enlighten another facet of the complexity of the structure–property relationship (Figure ). In the scope of this work, a DGEBA/IPDA resin was studied with various curing rates of 30%, 50%, 70%, and 90%. When it comes to conversion rate, this link has been shown previously on Jeffamine systems.[39] First, a slight decrease in density was observed here while it was not the case for the Jeffamine systems. This could be explained by the fact that an increase in cross-links could further constrain the system. The limited opportunity for monomers to be optimally stacked could lead to an increase in free volume and thus result in a decrease in density. Second, the trend regarding the glass transition temperature and Young’s modulus is coherent with that of the Jeffamine systems. The link between the structure and the resulting properties is more readily accessible than for the density observations. Since a low conversion rate accounts for fewer cross-links, it is coherent to have a higher elastic response for this material. Similarly, less thermal energy is needed to go from a rigid amorphous thermoset to a more flexible one when fewer covalent bonds are present.

Figure 12

Influence of the conversion rate on the thermomechanical properties of a DGEBA/IPDA system.

Conclusion

This work enlightens a broader study of the modeling of high-performance materials such as epoxy thermosets in terms of the chemical structure and thermomechanical properties. Novel resins were investigated, including bio-based ones, both experimentally and by molecular dynamics. The implementation of a “reactive” cross-linking algorithm allows for the creation of a highly cross-linked resin, of up to 99%. The use of the CHARMM force field, usually reserved for biomolecules, was proven to be suitable for these thermosetting resins. Indeed, this work exhibits a good correlation between experimental and simulated results for a variety of resins with backbones fluctuating from aliphatic to aromatic. Thermomechanical properties in the name of density, Young’s modulus, and glass transition temperature were measured. In order to acquire these properties by simulation, innovative methodologies and post-treatment of the data were performed and compared with the experimental data to fully implement a robust and generalized approach for each property studied. The most compelling comparison is found in the glass transition temperature, as it is the property with the most data points both experimentally and by molecular simulation. Also, the novelty of this work resides in the large range of glass transition temperatures covered by the resins studied. Overall, because of the consistency with the experimental data, with average absolute deviations between experiments and molecular simulation comprised between 2% and 12%, it is possible to conclude that all of the methods presented are transferable and successful. Thus, it is possible to efficiently determine the density, the glass transition temperature, and Young’s modulus for epoxy thermosets by atomistic methods. From an industrial point of view, the development of these methods is the starting point of a more ambitious goal. The definition of these first properties allows a quick characterization of a system and an evaluation if it is interesting for an experimental study according to the targeted specifications. As well as being a predicting tool, this investigation demonstrates the interest in using molecular simulations to extract important structural analyses, often difficult to obtain experimentally. Atomistic simulations will thus have a crucial role in the structure–property understanding of materials by visualizing the system and acquiring key properties. Building from this, other types of resins such as polyurethanes are currently being considered and studied with this approach.