Factors influencing thermal solidification of bent-core trimers

Bent-core trimers are a simple model for many small organic molecules whose bulk liquids exhibit both crystallization and glass-formation in experiments.^1–4 Their shape can be characterized using three dimensionless parameters (Fig. 1): the bond angle θ₀, the ratio $r̃$ of the end-monomer radius to the center-monomer radius, and the ratio R of the bond length to the center-monomer diameter. It is well known that the properties of systems composed of such molecules depend strongly on all three of these parameters. For example, the three terphenyl isomers, which can be modeled as trimers with the same $r̃$ and R but different θ₀, form very differently structured bulk solids under the same preparation protocol.¹ Analogous θ₀-dependent differences occur between the three xylene isomers.² Two more classes of small molecules, the diphenylcycloalkenes and cyclic stilbenes, which can each be modeled as having the same θ₀ but different $r̃$ and R, show similarly complex and poorly understood dependence of crystallizability on molecular shape.^3,4

One of the reasons why our understanding of such phenomena and hence our ability to engineer crystallizability/glass-formability at the molecular level remains very limited is that only a few theoretical studies have isolated the role played by molecular shape using simple models. Reference 5 reported the densest packings of 2D R = 1/2 trimers as a function of $r̃$ and θ₀. References 6–8 examined the thermal solidification of Lewis-Wahnstrom-like models (⁠$R=2−1/6,r̃=1,θ0=105○$⁹) and reported several nontrivial effects of trimeric structure, e.g., that its enhancement of the interfacial energy between crystalline and liquid phases promotes glass-formation. The tangent-sphere (⁠$r̃=R=1$⁠) case shown in Fig. 1(b) is of considerable interest because it allows straightforward connection to results obtained for monomers—and hence isolation of the role played by the bond and angular constraints—while remaining a reasonable minimal model for trimeric molecules. References 10 and 11 reported the densest packings of 2D and 3D trimers and contrasted them to the jammed packings they form under dynamic athermal compression. Reference 12 examined how the dynamics of supercooled trimer liquids vary with θ₀.

In this paper, we examine the thermal solidification of bent-core tangent-sphere Lennard-Jones trimers as a function of their bond angle θ₀ using molecular dynamics simulations. We obtain three key results: (i) the sharply increasing solidification temperature T_s(θ₀) for θ₀ ≲ 20° is driven by the same gradual loss of effective configurational freedom that drives athermal trimers’ decreasing ϕ_J(θ₀);¹¹ (ii) θ₀ that allow formation of both FCC and HCP order crystallize, while θ₀ that only allow formation of HCP order glass-form; (iii) measurements of local cluster-level structure via characteristic crystallographic element (CCE) and topological cluster classification (TCC) analyses^13,14 are as predictive of whether trimers will crystallize or glass-form as they are for monomeric systems.

Each simulated trimer contains three monomers of mass m. These trimers are rigid; their bond lengths and angles are held fixed by holonomic constraints. Monomers belonging to different trimers interact via the truncated and shifted Lennard-Jones potential,

where ϵ is the energy scale of the pair interactions, σ is the monomer diameter, and r_c = 2^7/6σ is the cutoff radius.

The Lewis-Wahnstrom model⁹ is a good glass-former largely because its equilibrium backbone bond length ℓ₀ is incommensurable with its equilibrium nearest neighbor distance for nonbonded neighbors r₀; specifically, it has ℓ₀ = σ and r₀ = 2^1/6σ. In contrast, the current model makes these lengths commensurable (ℓ₀ = r₀ = σ). Any frustration against crystallization is driven primarily by the bond angle θ₀.^15,16 Six θ₀ are commensurable with 3D close-packing: 0°, cos⁻¹(5/6) ≃ 33.5°, 60°, cos⁻¹(1/3) ≃ 70.5°, 90°, and 120°.¹¹ All of these θ₀ allow formation of HCP crystals, whereas only 0°, 60°, 90°, and 120° allow formation of FCC crystals. We will argue below that this distinction is critical.

Initial states are generated by placing n_tri = 1333 trimers randomly within a cubic cell of volume V₀ = πn_triσ³. Periodic boundary conditions are applied along all three directions and Newton’s equations of motion are integrated with a time step δt = 0.005τ, where the unit of time is $τ=mσ2/ϵ$⁠. Fixed covalent bond lengths and bond angles are maintained using a standard method.¹⁷ Systems are equilibrated at temperature k_BT = ϵ and pressure P = 0 until the intertrimer structure has converged, then cooled to T = 0 at a rate Ṫ = −10⁻⁶ϵ/(k_Bτ) while maintaining zero pressure. All MD simulations are performed using LAMMPS.¹⁸ 

During the cooling runs, we monitor the cell volume V(T) and, hence, the monomer number density ρ(T) = 3n_tri/V(T) and packing fraction ϕ(T) = πρ(T)σ³/6.¹⁹ We then identify the solidification temperature T_s of each system using one of the two methods. For systems exhibiting crystallization, we simply locate the center of the first-order-transitionlike jump in ϕ(T) as T decreases. For systems exhibiting glass-formation, we identify T_s as the intersection points of linear fits to the high- and low-T portions of their ϕ(T) curves. Consistent with our focus on trimers’ solidification under dynamic cooling, these T_s(θ₀) are below the equilibrium solid-liquid transition temperatures T_melt(θ₀):T_s(θ₀) = T_melt(θ₀) − ΔT(θ₀, Ṫ), where ΔT increases with |Ṫ| and with glass-formability. The cooling rate employed here approaches the lowest rate that is currently computationally feasible for the n_tri and number of different θ₀ we study. When mapped to real units using the Lewis-Wahnström parameter values for ortho-terphemyl (OTP) (m = 1.27 × 10⁻²⁵ kg, σ = 4.83 × 10⁻⁹ m, and ϵ = 8.28 × 10⁻²¹ J⁹), |Ṫ| = 10⁻⁶ϵ/(k_Bτ) = 5.27 × 10⁴ K/s. Although this rate is quite large, it is well within the range of |Ṫ| achievable in experiments²¹ and is low enough that ΔT(θ₀) is small for our best crystal-formers. A further discussion of cooling-rate effects is given in the  Appendix.

We characterize systems’ local structure using the characteristic crystallographic element (CCE)¹³ and topological cluster classification (TCC)¹⁴ methods. Recent studies employing CCE^22–24 or TCC^25–29 analyses have led to much progress in our understanding of thermal solidification. These studies have found that the presence of energetically stable amorphous clusters in a system’s liquid state (at temperatures slightly above its T_melt) strongly promotes glass-formation. Long-lived clusters that are fivefold-symmetric and/or are subsets of icosahedra are particularly effective glass-promoters.^22,23,27,28

CCE employs descriptors known as “norms” that quantify the orientational and radial similarities of a given monomer’s local environment to those of various reference structures such as HCP and FCC crystals.¹³ These norms are built around sets of point symmetry groups that uniquely characterize the reference structure. Criteria based on these norms are applied to determine whether the monomer can be associated with structure X. Here, we identify monomers as FCC-like, HCP-like, or fivefold-like if their respective norms n_X (where X = FCC, HCP, or 5f) are less than 0.21. Monomers that satisfy none of these conditions are classified as “other”; these typically possess locally amorphous order. The relevant mathematical formulas and the algorithmic implementation are described in Ref. 13.

TCC is similar in spirit to CCE, but identifies differently structured clusters by their differing bond topology. Here, we employ TCC to track how the populations of various microstructural motifs within our systems vary with θ₀ and T, using the same procedures detailed in Ref. 14. We focus on the four cluster types shown in Fig. 2. 6A is the octahedron, a common motif in close-packed crystals. 6Z, a polytetrahedral structure, has higher energy but also much higher entropy than 6A³⁰ and is a common motif in glassy and jammed systems.³¹ 9B is a partial icosahedron; recall that icosahedral order has long been known to promote glass-formation.³² Finally, 9K is a subset of the HCP lattice.

During the cooling runs, we monitor the fractions f_X(T) of particles identified as X-like (for CCE analyses) or belonging to at least one cluster of type X (for TCC analyses). We will show below that θ₀-dependent differences in how these quantities evolve with decreasing T predict whether a given system will crystallize or glass-form.

Figure 3 summarizes the most basic features of bent-core trimers’ response to slow thermal cooling. Panel (a) shows the packing fraction ϕ(T) for the six values of θ₀ that are commensurable with 3D close-packing and contrasts these to results for monomers. Three features are immediately apparent. First, some θ₀ produce sharp, first-order-transitionlike jumps in ϕ(T) that indicate rapid crystallization, while others produce typical glassy behavior where a smooth crossover regime centered at T ≃ T_g(θ₀) connects roughly linear behavior at high and low T. It is noteworthy that some θ₀ that are commensurable with the formation of close-packed crystals glass-form even at the low cooling rate (|Ṫ| = 10⁻⁶/τ) employed here. Second, even for those systems that do crystallize, the sharpness of the transition depends strongly on θ₀; systems with sharper transitions (e.g., θ₀ = 120°) have faster crystallization kinetics. Third, trimers always solidify at higher temperatures and lower densities than their monomeric counterparts.

As shown in Fig. 3(b), trimers’ solidification temperatures remain relatively constant as θ₀ decreases from 120° to ∼20°, then increase rapidly as θ₀ → 0°. This trend presumably results from the loss of effective configurational freedom as trimers approach linearity. Specifically, the middle monomer in a bent trimer can relax away from obstacles by rotating about the line connecting the end monomers even if the end monomers are held fixed, whereas the middle monomer in a straight trimer cannot. The same decrease in effective configurational freedom appears to drive the decreasing ϕ_J in athermal bent-core trimers¹¹ as well as the increasing fragility of bidisperse bent-core trimer liquids¹² as θ₀ → 0°. See the  Appendix for a discussion of how this decrease might influence trimers’ θ₀-dependent liquid-state dynamics.

Figure 3(c) shows the packing fractions ϕ(T) for all θ₀ at three representative temperatures: T = 0, T = T_s(θ₀), and T = T_s(θ₀) + 0.025ϵ/k_B. The final densities reached at T = 0 do not depend very strongly on θ₀; the densest system (θ₀ = 120°) is only ∼6% more tightly packed than the least-dense system (θ₀ = 102.5°). The broad minimum around θ₀ = 102.5° is consistent with—and to some extent³⁴ supports—Lewis-Wahnström-like models’ choice of θ₀ = 105° as a bond angle suitable for modeling glassforming trimerlike molecules such as OTP.^6–9 Similar trends with θ₀ are apparent at T = T_s; the θ₀-dependence of ϕ(0) − ϕ(T_s) is fairly weak. In contrast, the much stronger θ₀-dependence of ϕ(T_s) − ϕ(T_s + 0.025ϵ/k_B) reflects the trends shown in panel (b), i.e., small-θ₀ trimers crystallize from much-lower-density supercooled liquids than their large-θ₀ counterparts.

The wide range of thermal solidification behavior highlighted in Fig. 3 results entirely from changing the bond angle θ₀, i.e., from changing trimers’ shape. This suggests that it can be understood at a microscopic level by applying tools and concepts like those applied to model atomic/colloidal systems in Refs. 22–29.

Figure 4 shows snapshots of θ₀ = 70.5° and θ₀ = 120° systems at T = T_s + 0.025ϵ/k_B, T_s, and 0, with monomers color-coded by their CCE norms. The systems appear similar at T_s + 0.025ϵ/k_B, while a few fivefold-symmetric sites are present, most lack any distinguishable order. For the θ₀ = 70.5° systems, the only obvious changes as T decreases are that more fivefold-symmetric sites appear, as do a very small number of HCP-ordered sites. The behavior of θ₀ = 120° trimers is very different. A large degree of mixed FCC/HCP order is present by T = T_s, and a relatively well-ordered crystal develops as T decreases further; the T = 0 snapshot clearly shows crystalline grains separated by stacking faults and the amorphous interphase. Since these two θ₀, respectively, have the lowest and highest fraction of close-packed monomers at T = 0, we characterize them as the best glass-former and best crystal-former.

We now examine these differences more quantitatively. Figure 5 shows f_X(T) for these systems. Panel (a) illustrates how the crystal-former shows sharp (and typical^22,23) upward jumps in both f_FCC and f_HCP and sharp downward jumps in f_5f and f_oth ≡ 1 − (f_FCC + f_HCP + f_5f) at T = T_s, whereas the glass-former shows no such jumps. Intriguingly, the total fraction of locally close-packed monomers f_cp ≡ f_FCC + f_HCP is about an order of magnitude higher in the crystal-former at temperatures slightly above T_s than it is in the glass-former, presumably indicating the crystal-forming liquid’s greater population of subcritical nuclei. The results from TCC analyses [panel (b)] show analogous trends; the crystal-former shows sharp (and typical^26,29) upward jumps in both f_6A and f_9K and sharp downward jumps in f_6Z and f_9B at T = T_s, whereas the glass-former does not. Remarkably, the formation propensity of midsize crystalline clusters, e.g., 9K (midsize amorphous clusters, e.g., 9B), is significantly higher (lower) in the crystal-former even at T ≃ 1.5T_s.

These differences between θ₀ = 70.5° and θ₀ = 120° systems are consistent with the trends seen in numerous model colloidal systems.^22–29 Considered in isolation, they are not surprising. What is surprising is that they are so stark despite the fact that both θ₀ = 70.5° and θ₀ = 120° trimers can form 3D close-packed crystals.¹¹ We now attempt to understand them further by comparing CCE and TCC results for all θ₀.

Figure 6 shows the θ₀-dependence of f_X(T) for T = 0, T_s, and T_s + 0.025ϵ/k_B. Panel (a) summarizes the essential features of the bond angle’s effect on crystallizability. f_cp(0) exhibits four broad maxima centered at θ₀ ≃ 0°, 60°, 90°, and 120°. These four θ₀ all allow formation of both FCC and HCP crystals. θ₀ = 0°, 60°, and 120° trimers all have f_cp(0) > 0.5. θ₀ = 90° trimers crystallize far less well: their f_cp(0) ≃ 0.1. A potential explanation for this difference is suggested by the fact that θ₀ = 0°, 60°, and 120° trimers can tile the 2D triangular lattice, whereas θ₀ = 90° trimers cannot. Specifically, if crystal nuclei preferentially grow by trimers attaching to their triangular ({1 1 1}) lattice planes, their growth kinetics will necessarily be slower for θ₀ = 90° than for θ₀ = 0°, 60°, and 120°.

Figure 6(a) also illustrates a second key difference between crystal-forming and glass-forming systems. Most crystal-formers exhibit comparable f_FCC and f_HCP for T ≲ T_s. This is not surprising; formation of RHCP-crystalline order is typical for pair interactions ranging from Lennard-Jones to hard-spherical.³⁶ What is surprising is that the formation of FCC (but not HCP) order is strongly suppressed in all glass-forming systems. The reason for this suppression is clear for θ₀ = 33.5° and 70.5°—these trimers can tile the HCP but not the FCC lattice—and the full f_FCC(T_s) dataset suggests that this suppression of FCC relative to HCP ordering extends to θ₀ well away from these special values. This is a nontrivial result because trimers with any θ₀ ≤ 120° (θ₀ ≤ 60°) can form close-packed bilayers (trilayers).¹¹ 

The θ₀-dependence of the TCC cluster populations at these three characteristic temperatures reinforces the above observations. Figure 6(b) shows three clear maxima [minima] of f_6A(0) [f_6Z(0)] at θ₀ ≃ 0°, 60°, and 120°. These match the trends in f_cp shown in panel (a)—sensibly so, because small crystalline clusters like 6A support close-packing while small amorphous clusters like 6Z suppress it.^25–29 The higher-T results for f_6A and f_6Z are less conclusive; crystal formers show greater increases in f_6A as T approaches T_s from above, but the trends in f_6Z are not notably different for crystal-formers vs glass-formers.

Clearer predictive distinctions emerge at the nine-particle level [Fig. 6(c)]. For 60° ≲θ₀ ≲ 100°, glass-formers have notably higher f_9B at T = T_s + 0.025ϵ/k_B, and the largest f_9B for T ≥ T_s occurs in θ₀ = 70.5° systems. Although partial-icosahedral structures like 9B are well-known to suppress crystallization in atomic/colloidal systems,^25–29,32 few previous studies of molecular liquids have explored how these structures’ formation propensity varies with the constituent molecules’ shape, and this result would have been difficult to predict at the single-trimer level, especially since θ₀ = 70.5° trimers can close-pack. Trends in crystalline-cluster populations (e.g., 9K) are less clear, apart from their notably lower f_X(T ≳ T_s) for θ₀ ≲ 20°, which is consistent with these systems’ higher T_s and lower ϕ(T_s).³⁵ 

Thus far, our cluster-level analyses have treated all monomers equally and ignored dynamics. This approach is incomplete because we are considering molecular (i.e., trimeric) liquids and because cluster-level dynamics are critical in determining glass-formability.^27,28 Figure 7 illustrates the additional insights that can be gained at the six-particle level. In panel (a), $fXtri(T)$ is the fraction of trimers that lie entirely within a cluster of type X. Differences between the best crystal-former and the best glass-former are far more dramatic than those illustrated in Fig. 5(b). Although many monomers in the θ₀ = 70.5° liquid lie within 6A clusters, very few trimers do. In the θ₀ = 120° liquid, however, many trimers lie within 6A clusters even at T = 3T_s/2. Thermal diffusive motion of these trimers should make the 6A cluster lifetime much shorter in the θ₀ = 70.5° liquid than it is in the θ₀ = 120° liquid. Panel (b) confirms this hypothesis: $PXs(T)$ is the probability that a given cluster of type X that is present at T_s + 0.025ϵ/k_B is still present when the temperature has dropped to T. The data show that 6A (6Z) clusters are far more stable in θ₀ = 120° (θ₀ = 70.5°) liquids in this temperature range. Note that these dramatic differences in cluster stability are present despite the fact that the two systems’ T_s are almost identical. While one might argue that they should have been expected because 70.5° is very close to the ideal fivefold-symmetry-promoting bond angle (72°), we emphasize that they occur despite the fact that θ₀ = 70.5° trimers can close-pack.¹¹ 

We examined how thermal solidification of bent-core trimers is influenced by their bond angle θ₀. Their solidification temperature T_s(θ₀) is relatively constant for θ ≳ 20°, but increases rapidly with decreasing θ₀ ≲ 20° owing to the reduction in trimers’ effective configurational freedom as they straighten. This decreasing freedom also produces a sharp decrease in athermal bent-core trimers’ ϕ_J(θ₀) as θ₀ → 0°.¹¹ We therefore conclude that these systems’ thermal-solidification (not just glass-formation) and jamming transitions are intimately connected, in a manner consistent with the ideas of Liu and Nagel.³⁷ Similar connections have recently been demonstrated for polymers.^38,39

On the other hand, our results also illustrate a fascinating contrast between thermal and athermal solidification. In athermal 3D systems, the competition between FCC and HCP ordering suppresses crystallization because it produces random-close-packed amorphous order.^40,41 In thermal trimer solidification, however, it seems that the absence of this competition suppresses crystallization. It is reasonable to suppose that HCP crystallites have a higher nucleation barrier (ΔG) than those of mixed FCC/HCP order owing to their lower entropy. A higher ΔG for trimers that can form only HCP crystals would make them less likely to crystallize (under cooling at a fixed rate Ṫ) than the trimers that can form both FCC and HCP order. This is indeed what we observe and is consistent both with Pedersen et al.’s finding that Lewis-Wahnstrom-like trimers glass-form more readily than monomeric Lennard-Jones systems primarily because their crystal nucleation barrier is higher^6–8 and with Russo and Tanaka’s conclusion that βΔG is the key quantity controlling glass-formability.⁴² 

We found that bent-core trimers exhibit a wide range of solidification behaviors; some are very good crystal-formers, while others are very good glass-formers. By examining all six of the θ₀ that are commensurable with 3D close-packing, we showed that these differences do not arise from geometric frustration. Instead, our CCE and TCC analyses showed that they arise from θ₀-dependent differences in trimer liquids’ tendencies to form local cluster-level structure that is stable and inhibits crystallization. Factors promoting greater stability of 6Z and other amorphous clusters (6A and other ordered clusters) unambiguously promote glass-formation (crystallization). The difference between the present study and previous studies with comparable findings^22–29 is that in trimer liquids the relative strength of these factors is determined by the quenched 2- and 3-body constraints (i.e., by the values of R and θ₀) rather than by details of the pairwise interactions or preparation protocol.

One of the principal goals of soft materials science is developing materials with tunable solid morphology. This can be done either by trial-and-error or by developing theories that predict the optimal molecular structure for obtaining a desired morphology and then synthesizing molecules possessing this structure. The latter approach, now commonly known as “molecular engineering,” has attracted great interest in recent years and has been applied to systems ranging from photonic crystals⁴³ to hydrogels.⁴⁴ The results presented in this paper should be a useful contribution to this effort because controlling bond angles in isomeric and near-isomeric small molecules is a molecular-engineering strategy that is often employed by experimentalists.^3,4,45,46

This material was based on work supported by the National Science Foundation under Grant No. DMR-1555242. We thank Austin D. Griffith for helpful discussions.

Precisely determining T_melt(θ₀) would require extensive simulations employing advanced techniques^6–8 and is beyond our scope. Here, to support our claim (Sec. II A) that ΔT(θ₀) is small for our best crystal-formers, we illustrate how T_s(θ₀) depends on Ṫ. Figure 8 shows ϕ(T) results for |Ṫ| spanning a factor of 100—from 10 times slower to 10 times higher than that employed above. The |Ṫ|-dependence is extremely weak for the glass-forming systems (θ₀ = 33.5° and 70.5°), indicating that their critical crystal nucleation rates are well below 10⁻⁷/τ. As expected, the crystal-forming systems exhibit larger (smaller) jumps in ϕ(T) at higher (lower) T, and hence higher (lower) T_s, for lower (higher) |Ṫ|. However, these differences in T_s are small, especially for monomers and θ₀ = 120° trimers.

We asserted that the dramatic increase in T_s(θ₀) [and decrease in ϕ_J(θ₀) in athermal systems¹¹] with decreasing θ₀ ≲ 20° is caused by the reduction of effective configurational freedom as trimers straighten. This idea could be tested by comparing the decays of the correlation functions

and

where $R→ee$ is a trimers’ end-end vector, $c→$ is the cross product of its bond vectors, and the averages are taken over all trimers and all t′, for various θ₀ at equal values of T − T_melt(θ₀) and/or T/T_melt(θ₀). An increasing decoupling of the decay times of $R(t)$ and C(t) with decreasing θ₀ ≲ 20° would directly support our assertion. Since performing this measurement in a rigorous fashion will require a high-precision calculation of T_melt(θ₀),^6–8 we leave it as a challenge for future work.