Thermal stability of phenolic resin: new insights based on bond dissociation energy and reactivity of functional groups

Abstract

Density functional theory (DFT) was applied to model molecules of phenolic resin (PR) to interpret the relationship between the atomistic structure and thermal properties of the cured PR. The bond dissociation energy (BDE) of C–C (in methylene bridges) and C–O bonds (in hydroxyls), as well as the Fukui function of hydroxyl, benzene and methylene groups, were calculated. The isomers of bisphenol-F and the methyl substituents have a slight effect on the BDEs of C–C and C–O bonds, while the oxidized structures, such as p-benzoquinone and aldehyde groups, lead to a drastic decrease in the C–C bond BDEs. The high reactivity of the carbon atoms in the benzene groups and the oxidized structures results in an increased possibility of being attacked by free radicals and protects the methylenes from being attacked, but it will also lead to ring-opening reactions and weight loss. These results provide a great opportunity to understand the relationship between the atomistic structure and the thermal stability of the cured PR, which plays a pivotal role in the design and optimization of thermal stable polymers.

---

Introduction

Phenolic resins (PRs) are widely used in thermal protection and ablative systems, owing to their low cost, excellent processability and good thermal stability.^1,2 Even so, the rapid development of the aerospace industry challenges the traditional PRs because of the requirement for higher thermal stability.³ Currently, improved thermal stability of PRs is achieved through modifications, for example, by heteroatoms (B, Si, P and Ti, etc.^4–7), nanomaterials (nano-clay, carbon nanotubes and graphene, etc.^8–10) and heat-resistant polymers (polyimide, bismaleimide and melamine, etc.^11–13). However, due to the insoluble and infusible properties of these materials and their structural complexity, the modification of PR is mostly costly and time consuming. A better understanding of the pyrolysis mechanism of PR, especially the initial stage of pyrolysis, will be helpful for the synthesis, process and modification of PR.

The pyrolysis of PR is a multi-scale and multi-stage process that contains several types of physical and chemical changes (see Fig. 1).^2,14 The thermal stability of PR largely depends on the weight loss during the initial stage of thermal decomposition () and volatilization ().^15–17 Thermal decomposition temperature (T_d,5% or T_d,10%), which is a widely used parameter on the evaluation the polymer thermal stability and depends on the difficulty of backbone bond cleavage and the reactivity of functional groups in the unpyrolyzed region ().¹⁸ With the exception of the C–C bonds in the benzene rings, all the other bonds (C–C, C–O, C–H and O–H) can break during the thermal decomposition (or homolytic cleavage) stage () and lead to the production of free radicals such as H˙, HO˙, R˙ and RO˙. It is these radicals which attack the unpyrolyzed structures ( and ) owing to their high mobility. The thermal stability of PR is further determined by the reactivity of the hydroxyl, methylene and benzene groups in the unpyrolyzed regions.^17,19 Therefore, two aspects we should be focused on in studying the structure–thermal stability relationship of PR, (i) the bond strength of C–O, C–H, O–H and C–H; (ii) the reactivity of hydroxyl, methylene and benzene groups with radicals. Bond dissociation energy (BDE) is often used to evaluate the bond strength. However, limited with the insoluble and infusible properties of the cured PR, the BDE of the bonds in the backbone and the reactivity of the main groups can hardly be measured or evaluated through the experimental approaches.^16,17,20

Fig. 1 Multi-stage pyrolysis model of PR². Thermal stability of PR mainly depends on the homolytic reactions in region and and radical reactions in region and .

Molecular simulation techniques are effective and efficient in revealing the structural–thermal stability relationship of PR. For example, Jiang et al.²¹ and Qi et al.²² used the ReaxFF molecular dynamics method to track the initial pyrolysis paths, Bauschlicher et al.²³ and Bian et al.²⁴ calculated the BDE of the backbone bonds through density functional theory (DFT), Khavryuchenko et al.²⁵ adopted the semi-empirical methods (PM7) to model the carbonization of PR. Among the molecular simulation methods, DFT has gained immense popularity over the past few decades due to its high computational efficiency and accuracy. For example, Wu et al.²⁶ studied the thermal stability of boron containing PR composites based on graphene-like model molecules, the results showed that the BDEs of C–C bonds calculated on the B3LYP/6-31G(d) level were in excellent agreement with the results on the HF/6-31G(d) level. The electrostatic potential, Fukui function and frontier orbitals are effective in describing the reactivity of the functional groups during pyrolysis, as proposed by Hemelsoet et al.²⁷ Conner et al.²⁸ revealed that the reactivity of phenolic compounds was correlated with the electrostatic charges at the reaction sites in the phenolic ring. Our previous study²⁴ also indicated that the boron substitution had a significant effect on the reactivity of carbon atoms next to the boron atoms, which will further influence the thermal stability of boron-containing PRs.

The computational study of polymer or amorphous organic system is usually based on simplified model molecules, hence we need to design simplified model molecules for the cured PR. The multistep procedures of synthesis and curing of PR lead to the rather complex structure of the cured PR. The crosslinked network structure formed by phenols and methylenes, and a small amount of carbonyls, ethers and benzoquinones is found in the cured PR. Based on the studies of Trick et al.,¹⁶ Jackson et al.¹⁹ and Qi et al.,²² here we classified the structures in cured PRs as three types (see Fig. 2), (i) bisphenol-F-like isomers with different connection patterns between phenolic rings; (ii) bisphenol-F-like structures with different methylene substituted sites and numbers on phenol rings; (iii) bisphenol-F-like structures containing oxidized structures, such as quinones, ketones, ethers, aldehydes and alcohols. Such model construction strategy is effective and reliable for BDE calculation in organic systems, as Younker et al.²⁹ suggested.

Fig. 2 Bisphenol-F-like model molecules for cured PR. The model molecules are classified as (i) isomers with different connection patterns between phenol pairs; (ii) bisphenol-F-like model molecules with different methyl substituent status; (iii) bisphenol-F-like model molecules containing oxidized structures.

The main goals of this study are to obtain the BDEs of the skeleton of PR and the reactivity of main functional groups (hydroxyl, benzene and methylene) in PR, and illustrate the structural–thermal stability of PR based on these results. Based on the simplified model molecules of PR as mentioned above, DFT methods were applied to calculate the BDEs and bond orders (BOs) of main chemical bonds, as well as the electrostatic potentials, Fukui functions and molecular orbitals.

Computational details

Bond order evaluation of the influence of model size

In order to study the influence of model size on the bond strength, we calculated the BOs of C–C bonds in methylene bridges and C–O bonds in hydroxyls based on model molecules with different size (Fig. S1†). The simplified model molecules of PR with different connection patterns contain 2, 4 and 6 phenol rings. Because the C–C bonds in the methylene bridges can freely rotate, there will be many conformations in model molecules containing 4 and 6 phenol rings. To traverse as many conformations as possible for each models molecule, 10 frames were screened out from 100 randomly generated fames, which are extracted from a 10 ns molecular dynamics simulation. The bond order of C–C bonds in the methylene bridges and C–O bonds in the hydroxyls were calculated with the aid of Gaussian 09 (ref. 30) and Multiwfn³¹ softwares. Details of model screening and bond order calculation can be found in the ESI (page 2†).

Bond dissociation energy calculation

Considering the computational accuracy, the primary issue is to find a suitable theoretical method for calculating the BDEs. The simplified model molecules for cured PR (Fig. 2) were optimized with the hybrid functionals of BPW91, B3LYP, B3PW91 and M06-2X,^32–34 combined with the basis sets of 6-311G(d,p) and def2-TZVP.³⁵ Bauschlicher et al.³⁶ found that the BDEs of O–H and C–O in bisphenol-F calculated with B3LYP were only slightly less than that with CCSD(T). Hence all these hybrid functionals we chose are credible for BDE calculation. Owing to its high percent of Hartree–Fock exchange (54%), M06-2X is better than the other functionals in describing the long range interaction, which will benefit BDE calculation.^37,38 The non-bond interactions, such as hydrogen bonds, π–π and p–π interactions, are also playing an important role of the bond strength of phenolic model molecules.³⁹ Considering the good description of the weak interaction and comparable precision with CCSD(T), the DFT-D3 dispersion correction method⁴⁰ was used to improve the accuracy of BDE calculation. With good quality and high efficiency, def2-TZVP and 6-311G(d,p) have already been widely used for BDE calculation.^35,41 The vibrational frequencies were calculated at the same level for characterizing the nature of structures and for zero-point energy corrections. All the above calculations were performed using the Gaussian 09.

The BDE of chemical bond A–B in molecule AB defined as enthalpy change when the bond is cleaved (reaction (1)) at certain temperature (such as 298 K, see Fig. 3),⁴²

A–B → A˙ + B˙ (1)

ΔH_T(A–B) = H_T(A) + H_T(B) − H_T(A–B) (2)

where H denotes as the enthalpy of molecules or free radicals at certain temperature. For this study, the BDEs were calculated at 298 K. The enthalpy of both the homolytic products and the reactants can be obtained based on the total electronic energy and vibrational frequencies, which were obtained from optimization and frequency calculation jobs. The homolytic products for each model were optimized with unrestricted open-shell method. The scale factor used for frequency correction of BPW91/6-311G(d,p), B3PW91/6-311G(d,p), B3LYP/6-311G(d,p), B3LYP-D3/6-311G(d,p), B3LYP/def2-TZVP, M06-2X/6-311G(d,p), M06-2X-D3/6-311G(d,p) and M06-2X/def2-TZVP were 0.986, 0.986, 0.991, 0.967, 0.999, 0.947, 0.979 and 0.983, respectively.^43–45 Further details of the BDE calculation can be found in ESI (page 3†).

Fig. 3 Procedures for the calculation of bond dissociation energy, bond order, Fukui function for radical attack, electrostatic potential and HOMO/LUMO energy levels.

Fukui function for free radical attack

BOs, Fukui functions and frontier orbitals were calculated at B3LYP/6-311G(d,p) and M06-2X-D3/6-311G(d,p) levels, with the aid of Multiwfn. Because it is preferable for the reactivity analysis to be regional,⁴⁶ the introduction of diffusion function may lead to unreliable results, thus the wavefunction analyses were not performed on the triple zeta basis set of def2-TZVP. To calculate f⁰(r) (Fukui function for free radical attack),⁴⁷ the wavefunction of these model molecules with neutral (Q = 0), cationic (Q + 1) and anionic (Q − 1) charges were also obtained.where f⁰(r) denotes the Fukui function for free radical attack of a molecule with N electrons, f⁰_k(r) denotes the Fukui indice for free radical attack of atom k; ρ_N+1(r) and ρ_N−1(r) denote electron density of a molecule with N + 1 and N − 1 electrons, respectively; ρ_k(N + 1) and ρ_k(N − 1) denote atomic charge of atom k in molecule with N + 1 and N − 1 electrons, respectively.

Results and discussion

The influence of model size on the bond strength

The complexity of cured PR originates from the diversity of connection patterns between phenolic groups and the rotational freedom of C–C bonds in methylene bridges. The simplest structure in the cured PR is the bisphenol-F isomers (2-1, 2-3 and 2-5, see Fig. S1†).² There are more than 20 minimum points on the potential energy surface of the methylene bridge rotation in o–o′ type bisphenol-F (2-1),⁴⁸ with the rotational barriers of approximately 5–10 kcal mol⁻¹ that can be easily overcome at the room temperature.⁴⁹ Except for the o–o′ type phenol pairs, in which hydrogen bonds form between hydroxyls and lead to restrained motion of the phenolic pairs, the other types of phenol pairs are almost statistically distributed in the cured PR. In order to study the influence of isomers and conformation of bisphenol-F on the bond strength of C–C and C–O bonds in PR, the BO of these bonds were calculated with different model size.

The BO of C–C and C–O bonds (see Fig. 4) reflected that the model size has little effect on the bond strength. Taking 2-1, PR-1, 4-1 and 6-1 as examples, the BO of C–C and C–O does not change obviously with the increasing model size. This means that bond strength evaluation based on bisphenol-F like model molecules will be reasonable. Note that owing to the limitations imposed by the BO definition, the influence of benzyl substituents on the bond strength of C–C bond can hardly be inferred from the obtained BO values.

Fig. 4 Bond orders of C–C (a) and C–O (b) in model molecules of PR with different model size calculated at the B3LYP/6-311G(d,p) level. The C–C bonds refer to the bonds in the methylene bridge, while the C–O bonds refer to the bonds in the hydroxyls. The dark yellow and blue bars represent for the Mayer bond order (MBO) and Laplacian bond order (LBO), respectively.

The influence of the isomers and conformation on the bond strength of C–C and C–O was analyzed based on the Mayer BO (MBO) and Laplacian BO (LBO), which show good correlation with bond strength.⁵⁰ By comparing the BOs of C–C bonds in 2-1, 2-3, 2-5; 2-2, 2-4, 2-6; 4-1, 4-3, 4-5; 6-1, 6-3, 6-5; it can be found that the C–C bond strength in the methylene bridges is in the order of o–o′ < o–p′ < p–p′. A possible reason for this maybe that the hydrogen bonds formed in the o–o′ type phenol pairs lead to high strain and low bond strength of methylene bridge. The bond strength of the C–C bonds at the two sides of the methylene in o–p′ type model molecules is almost equivalent, as can be found by the comparison of the 2-3, 2-4, 4-3, 4-4, 6-3, 6-4 BO values. The similar rules can be found in the p–p′ type model molecules. The BOs of C–O bonds in most of the model molecules do not change obviously with the chemical environment, except for those of the o–o′ type structures (2-1, 2-2, 4-1, 4-2, 6-1 and 6-2) in which the C–O bond strength are relatively weaker than the others. The hydrogen bonds formed between hydroxyls in o–o′ type structures are the main reason for this effect.^48,49

The error bars in Fig. 4 can be regarded as the ranges of calculated BOs caused by conformation change, rather than the errors of BO calculation. Take the model molecules formed by two phenols (2-1 to 2-6) as examples, the error bars of both C–C and C–O of o–o′ type model molecules (2-1 and 2-2) are higher than those of the o–p′ and p–p′ type model molecules (2-3 to 2-4), the main reason for this phenomenon is that the formation of hydrogen bonds will significantly change the BO of C–C and C–O in 2-1 and 2-2. There are a large number of conformations in the model molecules containing four (4-1 to 4-6) and six phenols (6-1 to 6-6), which consequently lead to significant higher BO error bars (about ±0.02, for both C–C and C–O), compared with those obtained from the model molecules containing two phenols.

Bond dissociation energy

While most of the chemical bonds in cured PR could rupture at high temperature, the stability of C–C and C–O bonds are different in different chemical environment. Based on the BO analysis results and previous studies,^17,51,52 we know that the rupture of C–C and C–O bonds will lead to significant weight loss. The C–C and C–O BDEs in the simplified model molecules of PR were calculated with different functional/basis set combination (Fig. 5, 6, Table S1 and S2†), the BOs were also calculated to further support the BDE analysis results (Fig. 5, 6 and S2†). The DFT-D3 dispersion correction method was adopted to improve the reliability of BDE calculations.

Fig. 5 BDEs (298 K, kJ mol⁻¹) and BOs of C–C bonds in methylene bridges of model molecules of PR. Red and olive bars represent the BDEs calculated with B3LYP-D3/6-311G(d,p) and M06-2X-D3/6-311G(d,p), respectively. The dark yellow solid squares and hollow circulars represent for Mayer bond orders (MBOs) and Laplacian bond orders (LBOs), respectively. The referenced model molecule is diphenylmethane.

Fig. 6 BDEs (298 K, kJ mol⁻¹) and BOs of C–O bonds in hydroxyls of model molecules of PR. The referenced model molecule is phenol.

C–C bonds in methylene bridges. Diphenylmethane, which is similar to the bisphenol-F-like PR model molecules in this study, is used as control model in the C–C bond strength evaluation (Fig. 5, Table S1†). The BDEs of C–C bond in diphenylmethane calculated at M06-2X/6-31(d,p) and M06-2X/def2-TZVP levels are closer to the reported value (378.2 ± 8.4 kJ mol⁻¹),⁵³ compared with the results at the other levels. The DFT-D3 corrected results, i.e., the BDEs at B3LYP-D3/6-311(d,p) and M06-2X-D3/6-311(d,p) levels, are more credible than the uncorrected results, but still underestimate (for B3LYP-D3/6-311(d,p)) or overestimate (for M06-2X-D3/6-311(d,p)) the BDE, compared to the reported value. This demonstrates the importance of dispersion interaction in BDE calculation.⁵⁴ The discussions of C–C BDEs below are mainly based on the results obtained at the B3LYP-D3/6-311(d,p) and M06-2X/def2-TZVP levels. Generally, for the same bond in similar chemical environments, a higher bond order reflects a higher bond strength.^55,56 Hence the BO analysis at the two levels will provide some in situ information about the bond strength of C–C bonds.

The comparison between PR-1, PR-2 and PR-3 shows that, the isomers of bisphenol-F has a negligible effect on the BDE of C–C bonds, in seeming contradiction to the BO analysis results described in Section 3.1 (Fig. 4 and 5). This can be explained by the hydrogen bonds formed in the o–o′ type structure (PR-1), which lead to a lower BO of C–C bond,^23,36,48 compared with o–p′ (PR-2) and p–p′ type (PR-3) structures. The BDEs are based on hypothetical homolytic reactions of C–C bond and are the enthalpy differences between the products and the reactant. The intra hydrogen bond will lead to the reduced enthalpy of the reactant and increased BDE of the C–C bond in o–o′ type structure (PR-1). Although the bond strength of C–C bond in o–o′ type structure is slightly weaker than that in the o–p′ and p–p′ type structures, the former will still be stronger than the latter two for the intra hydrogen bond formation. This result is helpful for explaining the experimental observation that the cured PR based on the high-ortho novolac resin is more thermally stable than the random linear type novolac resin.^6,57

The comparison among PR-4 to PR-10 shows that, the methyl substituents' site and number have a slightly negative effect on the BDE of C–C bonds. The BDE of C–C bonds adjacent to the two or one methyl substituted phenol are less than that next to the three methyl substituted phenol (PR-4 to PR-9). It is found that the ortho-substituted methyl will decrease the BDE of C–C bonds, while the effect of para-substituted methyl is negligible (PR-6 and PR-10). That is largely a result of the increased conjugative effect on the phenol ring and the decreased electron cloud on the methylene bridge, derived from the inductive effect of the ortho-substituted methyls.^23,58 Therefore, structures such as PR-5 and PR-6 will be beneficial for the thermal stability of PR, while structures like PR-4, PR-11 and PR-12 have a negative effect.

The influence of oxidized structures on the thermal stability of PR can be analyzed using PR-11 to PR-14. A small proportion of oxidized structures, such as ketone, aldehyde and quinone can be formed during the synthesis, storage, transportation and curing.^6,59,60 The HO˙ and H₂O generated during the pyrolysis of PR will also lead to the generation of such oxidized structures.^19,61 It can be found that the BDE of C–C bond next to the p-benzoquinone, which can be easily formed during the pyrolysis of PR,⁶² is significantly less than that next to the ketone and aldehyde. The p-benzoquinone structure is presumably either non or weakly antiaromatic,⁶³ hence it will be easily ruptured at high temperature, leading to the CO and CO₂ formation.^21,22 The BDE of C–C bonds adjacent to the aldehydes (PR-13 and PR-14), the hydroxymethyl (PR-15 and PR-16) and the ether linkage (PR-17) are higher than that in the methylene bridge in PR-1, hence it will be difficult for these groups to be dissociated from the resin skeleton during the homolytic stage. This is in agreement with our recent experimental studies.^51,61

C–O bonds in hydroxyls. The previous study showed that the HO˙ peeled from the hydroxyls during the initial stage of pyrolysis had destructive effect on the skeleton of PR.^21,22,36,61 The calculation of the C–O BDE will be helpful for determining whether it can be influenced by the chemical environment of the hydroxyls. The calculated BDEs of C–O in phenol at M06-2X-D3/6-311G(d,p) (449.68 kJ mol⁻¹) and M06-2X/def2-TZVP (449.15 kJ mol⁻¹) level are closer to the reported value (445.67 ± 4.2 kJ mol⁻¹)⁵³ (Fig. 6) than the values at the other levels (Table S2†). The discussions of C–O BDE below will mainly base on the data at M06-2X-D3/6-311G(d,p) level.

The isomers of bisphenol-F (PR-1 to PR-3) and the methyl substituents' number and site (PR-4 to PR-10) have a negligible effect on the BDE of C–O bonds in hydroxyls, while the oxidized structures have a negative effect. The ether linkage (in PR-17) is unstable due to the low BDE of C–O bond. Ether linkage exist extensively in cured PR, especially in the fast cure product. Because of the low BDE of C–O bond, the ether linkage is vulnerable and will rupture at low temperature. This can be one important reason for that the thermal stability of cured PR can be further improved by post-curing treatment.⁶⁴

The C–O bond adjacent to the carbonyl (oxidized from methylene, PR-11),⁵⁹ benzoquinone (PR-12)⁶⁵ and aldehyde (PR-13)¹⁹ are significantly weaker than that in the other chemical environment. This can be caused by the reduced electron cloud in hydroxyls, derived from the formation of hydrogen bonds between the hydroxyls and the neighboring carbonyls.^6,48 The BO analysis (Fig. 6) further verified this speculation. Combined with the analysis of C–C BDE in PR-11 and PR-12, it can be inferred that the p-benzoquinone structure has more negative effect on the thermal stability of PR, compared with the carbonyl.

Reactivity analysis

Relationship between reactivity and thermal stability. In addition to the BDEs of main bonds in the resin, the reactivity of phenolic hydroxyl, methylene and phenol ring are also vital to the thermal stability of PR. The free radical (H˙, HO˙, CH₃˙ and Ph˙) related reactions in the unpyrolyzed region dominate after the homolytic reactions (C–C, C–O and O–H bonds).^21,22 Considering the complexity of reaction paths resulting from the isomers, conformation, methyl substituents' number/site and oxidized structure, a detailed study of the reaction mechanism for the radical reactions is quite challenging.^36,66

Because this work mainly focused on the structural-thermal stability relationship of PR, here we adopted the electrostatic potential, Fukui function and frontier orbitals to analyse the reactivity of hydroxyl, methylene and benzene groups. During the initial stage of the free radical substitution reaction, the free radical are attracted by the electrostatic potential, hence the reactivity of the target groups can be evaluated by the electrostatic potential the model molecules. When the distance between the radical and target group is sufficiently small, the electron clouds begin to overlap and the electron population began to change significantly. Fukui function will be a more reliable tool for presenting the reactivity of functional groups during this stage, as proposed by Paul et al.⁴⁷ Further information of reactivity and thermal stability of these model molecules can also be obtained from frontier orbitals.^67,68

Electrostatic potential. For most model molecules, the most negative and most positive electrostatic potential distributed around the hydroxyls, while the electrostatic potential distributed around the methylene and benzene groups are close to zero (Fig. 7). This means that the hydroxyls are the most active sites in each model molecules. The reactivity of the methylenes can be influenced by the conformation of the model molecules. For example, the hydrogen bond formation in PR-1 leads to the methylene wrapped inside the range of negative potential of hydroxyls, consequently the methylene can hardly be attacked by radicals. However, owing to the low rotational barriers of the methylene bridges in PR-2 and PR-3, the possibility of the methylenes being attacked by radicals can be significantly increased greater than that in PR-1. This will be an important reason for the experimental observations that the PR cured from the high-ortho novolac resin has a better thermal stability than that from the random linear type novolac resin.^69,70

Fig. 7 Electrostatic potential maps of model molecules of PR calculated at B3LYP/6-311G(d,p) level.

For the model molecules with mono-, bis- and tri-methyl substituents (PR-6 to PR-10), the electrostatic potential around the hydroxyls decreased and led to increased possibility of methylenes being attacked, compared with PR-1. For example, the positive electrostatic potential around the methylene in PR-8 and PR-9 are obviously higher than that in PR-1. In other words, the increased crosslink density will lead to a decreased possibility of methylenes being attacked by radicals. Combined with the BDE analysis result, it can be inferred that the increased BDE and decreased reactivity of methylene bridges are the main reasons that can be used to explain the experimental observations that the increased crosslink density leads to improved thermal stability of the PR.^57,71,72

Frontier orbitals. The HOMO/LUMO analysis will be helpful to determine the distribution of reactive sites (Fig. 8, 9 and 10). We use the radical reactions between HO˙ and the model molecules as examples. The SOMO of HO˙ is characteristic of π orbital⁷³ (both α and β orbitals), which is symmetrical to the HOMO over the hydroxyl and benzene groups, but is asymmetrical to the HOMO over methylenes. The gaps between the SOMO of HO˙ and HOMO of the model molecules of PR are about 3 eV at both B3LYP/6-311G(d,p) and M06-2X/6-311G(d,p) levels (see Table S3†). Combined with the electrostatic potential analysis, it can be speculated that the HO˙ will mainly attack the hydroxyl and benzene groups. Meanwhile, the high concentration and high mobility of HO˙ will lead to high probability of HO˙ attack reactions.^21,61 The oxidized model molecules will be more reactive than the other model molecules due to their relatively low gaps (Table S1†). In short, based on the frontier orbital analyses, both the hydroxyl, benzene and methylene groups are the potential sites that can be attacked by the radicals.

Fig. 8 SOMOs of H˙, HO˙, CH₃˙ and Ph˙ calculated at B3LYP/6-311G(d,p) and M06-2X/6-311G(d,p) levels. Both α and β SOMOs for each radical are considered.

Fig. 9 HOMO surfaces and energy levels of the model molecules of PR at B3LYP/6-311G(d,p) level.

Fig. 10 LUMO surfaces and energy levels of the model molecules of PR at B3LYP/6-311G(d,p) level.

Fukui function. The Fukui function for free radical attack (f⁰(r)) of the above model molecules is more reliable than the electrostatic potential and frontier orbitals in evaluating the reactivity of the radical reactive sites.⁷⁴ For clarity, the Fukui function for radical attack was mapped onto the total electron density of each model (Fig. 11, isovalue = 0.0028 a.u.).⁷⁵

Fig. 11 Fukui function for radical attack mapped onto the isosurface of the total electron density (isovalue = 0.0028 a.u., B3LYP/6-311G(d,p)) of model molecules of PR.

It can be found that except for the hydroxyl groups, the carbon atoms in the benzene groups are also the reactive sites (the dark blue regions), which is in accord with the HOMO/LUMO analysis results. The previous study showed that the C₂, C₄ and C₆ in phenol rings in PR are the active sites, and C₄ is more active that C₂ and C₆ (ref. 60 and 76) (the number rules for the mentioned atoms can be seen in Fig. S3†). However, based on the Fukui function analysis (Table S4 and S5†), it can be observed that the reactivity of the carbon atoms in the benzene groups are affected by the status of methyl substituents. For the phenol rings with two methyl substituents, the C₁, C₂, C₃ and C₄ are the active sites. For the phenol rings with one methyl substituents, C₂, C₃, C₅ and C₆ are the active sites when it is ortho-methyl substituents, while the C₁, C₂ and C₄ are the active sites when it is para-methyl substituted. The condensed Fukui indices for free radical attack (f⁰_k(r), see Table S1†) for each atom in the above molecules were calculated based on the Hirshfeld charges, quantitatively verified the isosurface analysis results.²⁷ The probable reactions on these carbon atoms are ring opening and radical substitution reactions, as reported Kovacevic et al.^77,78 All of these reactions have a negative effect on the thermal stability of PR.

The reactivity of the carbonyls and aldehydes in the oxidized structures (PR-11 to PR-14) are significantly higher than that of the hydroxyls, consequently leading to a decreased reactivity of hydroxyls and carbons in phenol rings. The carbon atoms in the carbonyl groups (PR-11 and PR-12) and aldehyde groups (PR-13 and PR-14) are particularly more active than the carbon atom in methylene (PR-1), further verifying the experimental speculation that the carbonyl and aldehyde structures in cured PR are vulnerable.^19,79,80 This result can hardly be revealed through the experimental approaches⁸¹ or BDE analysis, demonstrating the advantage of reactivity analyses. Considering the low activity and low BDE of the ether linkage, it can be inferred that the main reactions occurred in PR-17 are the homolytic reaction of C–O bond and the radical reaction on the carbon atoms in the phenol rings.

Conclusions

Based on simplified model molecules of cured PR, the BDEs of main bonds and the reactivity of main functional groups were obtained with DFT method. Our computational results provided the following clues for understanding the structure–thermal stability relationship of the cured PR.

(i) The bond strength of C–C bonds in oxidized model molecules is significantly decreased when the phenolic rings are connected to benzenequinones and ketones (PR-11 and PR-12), while the effects of model size, isomers of bisphenol-F and methyl substituent number/site on phenol rings are negligible.

(ii) Both hydroxyl and benzene groups are the potential reactive sites for radical reactions. In particular, hydroxyls are the most reactive sites in cured PR. The reactivity of methylene groups is effected by the conformation, while the reactivity of the carbons in the benzene groups is related to the methyl substituents on the benzene groups.

(iii) The formation of hydrogen bond in the o–o′ type structure has a beneficial effect on the thermal stability of the C–C bonds in methylene bridge, which can reduce the reactivity of methylene. Base on this conclusion, a reason for the experimental observations that the cured PR obtained from high-ortho type novolac resin is more thermally stable than that cured from random-linear type novolac resin is clear.

Acknowledgements

The authors would like to acknowledge the financial support to this work provided by the National Natural Science Foundation of China through grant no. 51473134 and no. 51273160.

References

1. K. Hirano and M. Asami, React. Funct. Polym., 2013, 73, 256–269.
2. A. Gardziella, L. A. Pilato and A. Knop, Phenolic resins: chemistry, applications, standardization, safety and ecology, Springer Science & Business Media, 2013.
3. E. Fitzer and L. M. Manocha, Carbon Reinforcements and Carbon/Carbon Composites, Springer Science & Business Media, 2012.
4. S. Wang, X. Jing, Y. Wang and J. Si, Polym. Degrad. Stab., 2014, 99, 1–11.
5. C. Li, H. Fan, D.-Y. Wang, J. Hu, J. Wan and B. Li, Compos. Sci. Technol., 2013, 87, 189–195.
6. A. Knop and L. A. Pilato, Phenolic resins: chemistry, applications and performance, Springer Science & Business Media, 2013.
7. Y. Zhang, S. Shen and Y. Liu, Polym. Degrad. Stab., 2013, 98, 514–518.
8. J. Si, J. Li, S. Wang, Y. Li and X. Jing, Composites, Part A, 2013, 54, 166–172.
9. F. L. Nabil, A. Zaidon, M. Jawaid, U. M. K. Anwar, E. S. Bakar, M. T. Paridah, S. M. A. Ridzuan and G. M. Aizat, Int. J. Adhes. Adhes., 2015, 62, 124–129.
10. J.-M. Park, D.-J. Kwon, Z.-J. Wang, J.-U. Roh, W.-I. Lee, J.-K. Park and K. Lawrence DeVries, Composites, Part B, 2014, 67, 22–29.
11. T. Takeichi, Y. Guo and S. Rimdusit, Polymer, 2005, 46, 4909–4916.
12. C. Nair, Prog. Polym. Sci., 2004, 29, 401–498.
13. T. Takeichi, Y. Saito, T. Agag, H. Muto and T. Kawauchi, Polymer, 2008, 49, 1173–1179.
14. A. R. Bahramian, M. Kokabi, M. H. N. Famili and M. H. Beheshty, Polymer, 2006, 47, 3661–3673.
15. S.-S. Tzeng and Y.-G. Chr, Mater. Chem. Phys., 2002, 73, 162–169.
16. K. A. Trick and T. E. Saliba, Carbon, 1995, 33, 1509–1515.
17. Y. Chen, Z. Chen, S. Xiao and H. Liu, Thermochim. Acta, 2008, 476, 39–43.
18. J. Critchley, G. Knight and W. W. Wright, Heat-Resistant Polymers: Technologically Useful Materials, Springer Science & Business Media, 2013.
19. W. M. Jackson and R. T. Conley, J. Appl. Polym. Sci., 1964, 8, 2163–2193.
20. B. Yu and V. Kodur, Composites, Part B, 2014, 58, 510–517.
21. D.-e. Jiang, A. C. Van Duin, W. A. Goddard III and S. Dai, J. Phys. Chem. A, 2009, 113, 6891–6894.
22. T. Qi, C. W. Bauschlicher Jr, J. W. Lawson, T. G. Desai and E. J. Reed, J. Phys. Chem. A, 2013, 117, 11115–11125.
23. C. W. Bauschlicher Jr, J. D. Monk and J. W. Lawson, Polymer, 2015, 80, 265–274.
24. C. Bian, Y. Wang, S. Wang, Y. Zhong, Y. Liu and X. Jing, Polym. Degrad. Stab., 2015, 119, 190–197.
25. V. D. Khavryuchenko and O. V. Khavryuchenko, J. Phys. Chem. C, 2013, 117, 7628–7635.
26. X. Wu and L. R. Radovic, J. Phys. Chem. A, 2004, 108, 9180–9187.
27. K. Hemelsoet, V. Van Speybroeck, G. B. Marin, F. De Proft, P. Geerlings and M. Waroquier, J. Phys. Chem. A, 2004, 108, 7281–7290.
28. A. H. Conner, J. Appl. Polym. Sci., 2000, 78, 355–363.
29. J. M. Younker, A. Beste and A. C. Buchanan, Chem. Phys. Lett., 2012, 545, 100–106.
30. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, R. Scalmani, G. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. J. Montgomery, J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brother, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, N. J. Millam, M. Klene, J. E. Knox, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, O. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski and D. J. Fox, Gaussian 09, revision A.02, Gaussian, Inc., Wallingford, CT, 2009.
31. T. Lu and F. Chen, J. Comput. Chem., 2012, 33, 580–592.
32. A. D. Becke, Phys. Rev. A, 1988, 38, 3098.
33. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785–789.
34. E. G. Hohenstein, S. T. Chill and C. D. Sherrill, J. Chem. Theory Comput., 2008, 4, 1996–2000.
35. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
36. C. W. Bauschlicher Jr, T. Qi, E. J. Reed, A. Lenfant, J. W. Lawson and T. G. Desai, J. Phys. Chem. A, 2013, 117, 11126–11135.
37. J. L. Bao, R. Meana-Paneda and D. G. Truhlar, Chem. Sci., 2015, 6, 5866–5881.
38. Z. Kelemen, R. Streubel and L. Nyulaszi, RSC Adv., 2015, 5, 41795–41802.
39. R. M. Borges dos Santos and J. A. Martinho Simões, J. Phys. Chem. Ref. Data, 1998, 27, 707–739.
40. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
41. A. Baishya, V. Mundlapati, S. Nembenna and H. Biswal, J. Chem. Sci., 2014, 126, 1781–1788.
42. S. J. Blanksby and G. B. Ellison, Acc. Chem. Res., 2003, 36, 255–263.
43. I. Alecu, J. Zheng, Y. Zhao and D. G. Truhlar, J. Chem. Theory Comput., 2010, 6, 2872–2887.
44. Precomputed Vibrational Scaling Factors,http://cccbdb.nist.gov/vibscalejust.asp.
45. Database of Frequency Scale Factors for Electronic Model Chemistries,http://comp.chem.umn.edu/freqscale/version3b2.htm.
46. R. R. Contreras, P. Fuentealba, M. Galván and P. Pérez, Chem. Phys. Lett., 1999, 304, 405–413.
47. P. W. Ayers and R. G. Parr, J. Am. Chem. Soc., 2000, 122, 2010–2018.
48. S. Wang, C. Hao, Z. Gao, J. Chen and J. Qiu, J. Mol. Model., 2013, 19, 1913–1918.
49. N. R. Pillsbury, C. W. Müller and T. S. Zwier, J. Phys. Chem. A, 2009, 113, 5013–5021.
50. T. Lu and F. Chen, J. Phys. Chem. A, 2013, 117, 3100–3108.
51. S. Wang, Y. Wang, C. Bian, Y. Zhong and X. Jing, Appl. Surf. Sci., 2015, 331, 519–529.
52. C.-L. Liu, W.-S. Dong, J.-R. Song and L. Liu, Mater. Sci. Eng., A, 2007, 459, 347–354.
53. Y.-R. Luo, Handbook of Bond Dissociation Energies in Organic Compounds, CRC press, 2002.
54. S. Grimme, R. Huenerbein and S. Ehrlich, ChemPhysChem, 2011, 12, 1258–1261.
55. A. J. Bridgeman, G. Cavigliasso, L. R. Ireland and J. Rothery, J. Chem. Soc., Dalton Trans., 2001, 2095–2108.
56. A. J. Bridgeman and G. Cavigliasso, Faraday Discuss., 2003, 124, 239–258.
57. Y. Wang, S. Wang, C. Bian, Y. Zhong and X. Jing, Polym. Degrad. Stab., 2015, 111, 239–246.
58. A. L. R. Silva, Á. Cimas and M. D. M. C. Ribeiro da Silva, J. Chem. Thermodyn., 2014, 73, 3–11.
59. G. E. Maciel, I. S. Chuang and L. Gollob, Macromolecules, 1984, 17, 1081–1087.
60. X. Zhang and D. H. Solomon, Chem. Mater., 1999, 11, 384–391.
61. Y. Zhong, X. Jing, S. Wang and Q.-X. Jia, Polym. Degrad. Stab., 2016, 125, 97–104.
62. D. Rosu, C.-D. Varganici, L. Rosu and O. M. Mocanu, in Thermal Degradation of Polymer Blends, Composites and Nanocomposites, ed. M. P. Visakh and Y. Arao, Springer International Publishing, Cham, 2015, pp. 17–49.
63. A. Kayembe, E. Lewars and J. F. Liebman, J. Chem. Thermodyn., 2012, 52, 43–56.
64. L. Pilato, Phenolic resins: A century of progress, Springer, 2010.
65. R. Lee, Phenolic resin chemistry and proposed mechanism for thermal decomposition, NASA Marshall Space Flight Center, Washington, 2007, pp. 7–13.
66. M. A. V. Ribeiro da Silva, A. I. M. C. Lobo Ferreira and Á. Cimas, J. Chem. Thermodyn., 2011, 43, 1857–1864.
67. H. Zhang, Y. Shang, M. Li, H. Zhao, X. Wang and B. Han, RSC Adv., 2015, 5, 90343–90353.
68. B. Banerjee, D. Kundu, G. Pugazhenthi and T. Banerjee, RSC Adv., 2015, 5, 85280–85290.
69. M. Nomoto, Y. Fujikawa, T. Komoto and T. Yamanobe, J. Mol. Struct., 2010, 976, 419–426.
70. J. Gao, Y. Liu and F. Wang, Eur. Polym. J., 2001, 37, 207–210.
71. S. Koshibe and M. Nanjo, US Pat., 4297473 A, 1981.
72. A. Bandyopadhyay, P. K. Valavala, T. C. Clancy, K. E. Wise and G. M. Odegard, Polymer, 2011, 52, 2445–2452.
73. D. Minakata, W. Song, S. P. Mezyk and W. J. Cooper, Phys. Chem. Chem. Phys., 2015, 17, 11796–11812.
74. A. K. Chandra and M. T. Nguyen, Faraday Discuss., 2007, 135, 191–201.
75. J. Du, X. Sun, J. Chen and G. Jiang, RSC Adv., 2014, 4, 44786–44794.
76. C. Li, Z. Bu, J. Sun, H. Fan, J. Wan and B. Li, Thermochim. Acta, 2013, 557, 77–86.
77. R. Cartas-Rosado and M. Castro, J. Phys. Chem. A, 2007, 111, 13088–13098.
78. G. Kovacevic and A. Sabljic, J. Comput. Chem., 2013, 34, 646–655.
79. Z. Lu, K. Lin and J. Gan, Environ. Pollut., 2011, 159, 2546–2551.
80. R. Kaur and Vikas, RSC Adv., 2016, 6, 29080–29098.
81. C. A. Fyfe, M. S. McKinnon, A. Rudin and W. J. Tchir, Macromolecules, 1983, 16, 1216–1219.

---

Footnote

† Electronic supplementary information (ESI) available: Details for bond order and bond dissociation energy calculation, HOMO/LUMO energy and condensed Fukui indices of model molecules of phenolic resin. See DOI: 10.1039/c6ra07597e

---