Generalization of the Hall-Petch and inverse Hall-Petch behaviors by tuning amorphous regions in 2D solids

: The strength σ y ( D ) of a polycrystal can decrease or increase with the grain diameter D , i


INTRODUCTION
Crystalline and amorphous materials are ubiquitous in nature and are widely used in the industry [1]. However, their crossover regime is much less explored and poorly understood [2]. Conventional crystals, such as silicon [3] or water ice [4] crystal can be catastrophically compressed into an amorphous solid (i.e., a glass) without forming intermediate structures of crystalline and amorphous mixtures. In recent years, glass-crystal composites (GCCs), also known as dual-phase materials, have been fabricated especially for amorphousnanocrystalline alloys [5]. They possess the advantages of both crystals and glasses and exhibit exceptional properties, such as ultrahigh strength [6], super plasticity [7,8], high reversible strain, world-record fatigue life [9], and less strength-ductility tradeoff [2,10], thus resulting in their broad industrial applications. These super alloys are typically fabricated by severe plastic deformation [11,12], whose microscopic kinetics are unclear. The size and shape of crystalline and amorphous regions can hardly be controlled; thus, their effect on material properties has been poorly explored and understood. Furthermore, while great efforts have been made to improve GCC's properties, their trial-and-error fabrication processes have limited theoretical rationale. Simulations show that thick amorphous regions better absorb dislocation [8] and suppress dislocation formation [13], but studies about the impact of l on material properties are very limited [8,[14][15][16].
The effect of mean grain diameter D has been well-studied in polycrystals, such as ceramics, metals and alloys with normal thin GB (i.e., l ≃ 1 particle) [17]. The strength, or yield stress σ y , increases as D decreases, i.e., the famous Hall-Petch (HP) behavior [18,19] discovered in 1950s. This trend reverses in 3 D 10-15 nm, i.e., the inverse Hall-Petch (IHP) behavior [20][21][22]. The HP and IHP behaviors generally hold in atomic and molecular polycrystals. However, when the GB thickness l is not as thin as about one particle in traditional polycrystals, how does D affect material properties? Moreover, how does σ y vary with the GB thickness l and does it exhibit HP-and IHP-like behaviors? These questions have rarely been investigated. Recent simulations found that σ y of Cu-based GCCs is higher at l = 1.0 nm than those at l = 0.5 and 1.5 nm [14,15], suggesting a nontrivial σ y (l). Similar behavior has also been reported by a simulation of quasi-two-dimensional high entropy GCC alloys, which focuses on the deformation mechanism of thick GB materials [16].
Motivated by the above basic questions, in this study, we measure the strength of solids by systemically changing D and l using two-dimensional (2D) molecular dynamics simulations. The traditional HP and IHP behaviors of σ y (D, l = 1) are generalized to σ y (D, l). The maximum σ y (D, l) exceeds the maximum σ y (D), thus providing a new approach to enhance the strength of solids.
The properties of a GCC can only be properly measured when it contains a sufficient number of crystalline grains, which is computationally expensive for three-dimensional (3D) samples, especially at large D and l. Therefore, we simulate 2D samples that are also poorly explored. 2D materials have attracted great interest in recent years because of their novel properties thanks to the new fabrication techniques. For example, 2D polycrystalline graphene [23,24], silicene [25], 2D amorphous cobalt-anadium hydroxide [26], and vanadium pentoxide [27] show exceptional chemical, physical, electronic, and surface properties [28]. Fabricating 2D solids with thick GBs is challenging. Monodispersed isotropic particles at high densities in 2D only form single crystals or polycrystals without amorphous domains, because the local favorable packing (i.e., equilateral triangle) is compatible with the global favorable packing (i.e., triangular lattice) without geometrical frustrations. Thus, we use binary-sized particles to form the amorphous GB regions and monodispersed particles in the crystalline regions. This should be one of the simplest model systems to achieve GCCs. Such composition distribution is common in polycrystals (e.g., carbon atoms segregating on the GBs of stainless steels) and in GCCs (e.g., binary Cu-Zr [13] or Ni-Co [29] amorphous regions with monodispersed Cu or Ni crystalline regions).

SIMULATION SYSTEM
The molecular dynamics simulations are conducted using large-scale atomic/molecular massively parallel simulator (LAMMPS) [30]. Two particles at separation r have the 12-6 Lennard-Jones potential:  Figure 1 The initial 2D GCC is composed of crystalline grains with diameter D = 100 in an amorphous matrix with thickness l = 40.
The zoom-in of the local amorphous and crystalline areas are shown in the left and right panels, respectively. The amorphous region is a 50%:50% mixture of large (blue) and small (grey) particles. The 16 crystalline grains have random lattice orientations and are composed of small particles.
U(r > r c ) = 0, where the cutoff r c = 2.5d. U(r) is shifted to 0 at r c . The diameters d of small and large particles are 1 and 1.3, respectively. Such size ratio has high glass-forming ability and has been widely used in 2D glass studies [31]. We set (d, U 0 ) = (1.0, 3), (1.15, 2), and (1.3, 1) for interactions between small-small, small-large and large-large particles, respectively, because smaller particles with short separations are often set to have a stronger interaction against phase separation [32].
Bidispersed particles are packed into amorphous GBs with thickness l, and monodispersed small particles are packed into hexagonal lattice domains with diameter D (Figure 1). The sample contains 16 crystalline grains with the periodic boundary conditions. The number of particles ranges from 5000 to 1.4 million in samples with small and large (l, D). To prepare a sample, we first produce a glass by rapidly quenching a liquid composed of 50%:50% binary spheres to temperature T = 0.2T * during 5τ time period and well relax it for 200τ so that the structure does not change over time. Then the 16 hexagonal regions are replaced by single crystalline grains composed of monodispersed particles. The lattice orientation difference between adjacent grains is randomly set in the range of 10 • to 15 • . Such sample is relaxed in the NPT ensemble at 0.5T * and hydrostatic pressure 40P * for 200τ, and then changed to 0.2T * at 40P * during 200τ and further equilibrates for 200τ before the measurement. The units All particles have the unit mass m = 1.

STRENGTH σ y AS A FUNCTION OF D AND l
A uniaxial compression is applied along the x direction ( Figure 1), deforming the sample with a compression strain ε. The yield stress and flow stress are very close and usually regarded as the same as the virial stress in simulations [33,34]. The measured stress-strain curves at different l in Figures 2A and S1 (Supporting Information) all exhibit the standard shape with three regimes: a linear increase (i.e., elastic regime), a nonlinear increase (i.e., strain-hardening), and a plateau (i.e., steady flow). The sawtooth shape of the plateau reflects the short quasi-elastic events, such as dislocation nucleation and movement in small systems [35]. We measure the mean flow stress by the plateau height at 7.5% < ε < 10% in three trials of simulations as σ y , which is sufficient to average out most noises [36]. The measured σ y (D, l) are shown in Figures 2B, 3A, Figure 2B increases and reaches a plateau at l = 1 and monotonically increases at l > 1. Both the plateau and increases are considered as the IHP behavior [17] as they deviate from the HP behavior of the monotonic decrease of σ y (D). By contrast, 3D polycrystals exhibit both the HP and IHP behaviors [17]. Studies on σ y (D) in 2D are very limited and focus only on polycrystals with norm thin GBs (i.e., l = 1) [23,24]. These studies only observed the IHP behavior without HP behavior, which has been attributed to the lack of dislocations in 2D [37]. A perfect 2D crystal lacks a source to generate dislocations according to the dislocation pile-up model [38]. Dislocations in 2D crystals are point defects that cannot entangle and can thus be easily absorbed by GBs. By contrast, dislocations in 3D crystals are topological line defects that can entangle into complex structures and cannot be easily absorbed by GBs. Thus, the HP behavior caused by dislocation motions is often absent in 2D polycrystals with normal thin GBs [23,24]. Thick GBs further suppress dislocations and their motions, thus the HP behavior is absent in Figure 2B. As shown in Figure 3A, σ y (l) monotonically decreases at D < 50, and changes non-monotonically with a peak at D > 50. They can be fitted by σ y (D, l) = e a 1 −l/a 2 − a 3 l −1 + 44. The constant 44 is the yield stress directly measured from the completely amorphous sample. l −1 term has been predicted by ref. [39] as the effect of thick grain boundaries with strain gradient. The fittings in Figures 3C-3E further show σ y (D, l) = 0.17D 0.35 e − l 0.08D+1.34 − 0.204D/l + 44, and the competition between the exponential term and D/l term gives the peak of σ y . The strength peak σ * y corresponds to D * and l * , which can be fitted by a power law D * (l * ), as shown by the dashed curve in Figure 3F. Since D * /l * is not a constant, σ * y is not solely determined by A GB /A tot . In addition, σ * y should exhibit a maximum rather than the monotonic increase with D * and l * in Figure 3F, but it is beyond the system size of our simulation. For normal polycrystals with l = 1, the fitted D * = 50 in Figure 3F for our 2D crystals is analogous to the peak position at the HP-IHP boundary of 15 nm (i.e., about 100 atoms) for typical 3D polycrystals. Interestingly, when σ y at fixed D is plotted as a function of the area fraction of GB, A GB /A tot , all curves collapse at A GB /A tot > 20% ( Figure 3B). This indicates that σ y is dictated by the total fraction of amorphous or crystalline regions, but is not sensitive to their spatial distribution. When A GB /A tot > 50%, σ y equals to that of the completely disordered glass, i.e., the crystalline regions do not affect σ y . When A GB /A tot < 20%, larger grains yield at larger σ y . The simulation in ref. [40] reports that σ y increases with the fraction of crystalline regions, i.e., decreases with A GB /A tot as shown in most regimes in Figure 3B. Nevertheless, ref. [40] measures only one curve without fixing D or l, thus are not related to the generalization of the HP and IHP behaviors.

DEFORMATION MECHANISMS
For large-grain (D 50) samples with l 10, GBs emit dislocations that glide through one grain and are absorbed by nearby GBs, see the yellow line segments in Figure 4A. Note that dislocations in small 2D crystalline grains can be easily absorbed by GBs, thus they only temporarily exist when they glide through the grains. Such dislocation motion is the major mode of plastic deformation. We observe that dislocations are more difficult to be emitted and easily absorbed by thicker GBs, in accordance with previous studies on Cu-Zr GCCs [8,13]. As l increases, the dislocation motions give way to another mode of deformation, namely, GB deformations ( Figures 4B and 4C). Thick GBs suppress dislocations but promote GB sliding. These two opposite effects on plastic deformation result in a peak of σ y (l) at l = 3 for the samples with D = 100 ( Figure 4D). When l > 10 and D = 100 (i.e., A GB /A tot > 20%), all the applied plastic deformation is borne by the GB sliding in amorphous regions. Therefore, σ y is solely determined by A GB /A tot and insensitive to D, which explains the collapsed σ y at A GB /A tot > 20% shown in Figure 3B. Similarly, dislocations barely exist in small crystalline grains, and thus all the applied plastic deformation is borne within the amorphous regions for samples with D < 50 ( Figure S2). Therefore, σ y (l) monotonically decreases at D < 50.
Deformations are dominated by dislocation motions in crystalline regions [41] and by GB sliding [42,43], GB diffusion creep [44,45] and shear banding in amorphous regions. These basic modes of deformation widely exist in polycrystals, glasses and their hybrid structures [41]. They have been used to explain the HP and IHP behaviors of σ y (D) in normal polycrystals with l = 1. Here, we similarly use these modes to explain the generalized HP-and IHP-like behaviors in σ y (D, l). As D varies, the competition between the deformations via dislocations and GBs gives rise to the classical HP and IHP behaviors with a peak of σ y at D ≃ 15 nm for most 3D polycrystals. Such a competition of the two modes of deformations shown in Figure  4 similarly gives rise to a peak of σ y (l) when D 50 in Figure 3A because reducing D has a similar effect to increasing l, i.e., increasing the fraction of amorphous regions. The increasing and decreasing part of σ y (l) at a fixed D are analogous to the classical HP and IHP behaviors of σ y (D), respectively, at a fixed l = 1.
The fraction of large-strain (η Mises > 16%) particles in GBs changes monotonically with ε ( Figures 5A and  S4) as expected, but non-monotonically with l, as shown in Figure 5B. For example, the fractions at different ε all reach maxima at l = 10 in samples with D = 100, as shown in Figure 5B. The peaks in Figures 5B and 5C divide σ y (l) into three regimes, see Figure 5C as an example at D = 100.
Besides the shear strains of individual particles in Figure 4, the relative motions between particles are characterized using (1) a reference line ( Figures 6A-6F, S5, S6); and (2) a displacement field (Figures 6G-6I). A row of particles is chosen as a reference line to highlight the key parts of the deformation. It has been used in simulation [46] and colloid experiment [47], but rarely in atomic polycrystals because atomic motions can hardly be tracked under electron microscopy. In Figures 6G-6I affine deformation. Such vector field shows both the magnitude and direction of particles' displacements relative to the background; thus, it can better illustrate the collective motions. The reference line has been used to illustrate the shear-coupled GB migration [46,47] and the displacement field has often been used in fluid flows. Figures 6D and 6E show that the compression along the x direction deforms the initially straight yellow reference line into small steps. Each step is along the [01] lattice directions inside the grains and random directions in GBs, see ellipses in Figure 6E. A newly formed step inside a crystalline region is caused by a dislocation gliding through the reference line. A gliding dislocation produces a straight line in the displacement field, as shown in Figure 6D. Deformations are mainly due to the dislocation sliding at l < 10 and GB sliding at l 10 ( Figure 6). At l 10, the bending of the reference line is localized to one large step near the middle of two triple junctions (ellipse in Figure 6F), that is, the center of a swirl ( Figure  6I). Particles with large displacements are near the triple junctions ( Figures 6E and 6F). The compression along the x direction expands l to l 2 for GBs along the x direction and reduces l to l 1 for GBs along ±120 • directions ( Figure 6E) because the amorphous regions have a positive Poisson ratio. The total amorphous area remains the same.
Interestingly, the displacement field in each crystalline grain exhibits a flower shape with a bright center, that is, a region with zero displacements ( Figure 7B). The two-in-two-out flow field (blue arrows in Figure  7C) around each fixed point (i.e., bright spot) demonstrates that the bright spot is a topological defect with charge -1 commonly observed in liquid crystals [48]. Such flow pattern reflects the compression and elongation along the x and y directions, respectively. A "flower" has more than four bright "petals" because there are more directions for vectors to be parallel, and thus they strongly overlap and produce brighter regions (Figure  7B). Usually, a flower has six or twelve bright petals, which reflect the six-fold symmetry of the hexagonal lattice. As a topological defect, a dislocation glides through a flower and cuts it into two topological defects with a domain wall ( Figures 7D and 7E).

CONCLUSION AND DISCUSSION
We systemically measured the properties of the 2D GCCs with different GB thicknesses l and mean grain diameters D. The classical HP and IHP behaviors of strength σ y (D) at l = 1 are generalized for the first time in two directions: (1) σ y (D) at l > 1 ( Figure 2B) and (2) σ y (l) at a fixed D ( Figure 3A). For direction (1), only IHP behavior is observed, in accordance with the absent HP behavior reported in 2D polycrystalline graphene [23,24] and silicene [25]. The absent HP behavior is attributed to the lack of dislocations because they cannot entangle in 2D and thus can be easily absorbed into GBs. For direction (2), σ y (l) exhibits a peak at l * (D) when D 50 due to the competition between the deformations via dislocations and GBs (Figure 4), and monotonically decreases when D < 50 because GB deformations dominate ( Figure S2). The increases and decreases of σ y (l) in Figure 3A are analogous to the classical HP and IHP behaviors of σ y (D), respectively, because both increasing l and decreasing D increase the amorphous fraction A GB /A tot .
As A GB /A tot increases, dislocation gliding gradually gives way to GB sliding or deformation inside GBs, that is, the mechanism of HP behavior giving way to the mechanism of IHP behavior as D decreases in normal 3D polycrystals. Therefore, the general HP and IHP behaviors in σ y (D, l) can be unified as σ y (A GB /A tot ), as shown in the collapse in Figure 3B. Although simulations for large (D, l) in 3D are challenging, our generalization of HP and IHP behaviors of σ y (D, l) in 2D could similarly exist in 3D because the competition between dislocation motion and GB deformation similarly exists in both 2D and 3D. Furthermore, the HP regimes due to dislocation motions should be wider in 3D compared with those in 2D because 3D solids have much more dislocations. Another prediction is that the HP behavior of σ y (D) should be absent in 3D at large l because the thick GBs, also known as the amorphous intergranular films in 3D [13], suppress the dislocation motions. It is worth noting that the deformation mechanism in real materials can be more complicated, such as caused by residual stress from the fabrication process [49] or via cracking in brittle materials [50], which may affect the HP and IHP behaviors. Moreover, the grain sizes are not uniform in real materials which may produce less sharp peak of σ y (D, l) in Figure 3 due to the mixed deformation mechanisms of thin and thick GBs.
The traditional σ y (D) classifies normal polycrystals with l = 1 into three regimes: the HP regime (D 15 nm), the IHP regime (3 nm D 15 nm), and unstable regime (D 3 nm) up to the amorphous limit D = 0. As an analogy, we classify solids into three regimes using σ y (l) at a fixed D, see the example at D = 100 in Figure 5C: (1) the HP-like regime at l 3; (2) the IHP-like regime which can be further divided by the peak in Figure 5B as (2a) the competition between dislocation motion and GB sliding (3 < l < 10) and (2b) only GB sliding (10 < l < 40, i.e. 20% < σ y (A GB /A tot ) < 50%); and (3) the glass regime. In regime (2b), σ y only depends on the total fraction of amorphous region, not depends on grain size distribution, see the collapse in Figure 3B. In regime (3), σ y (l) is independent to the crystalline regions at l > 40, i.e., σ y (A GB /A tot ) > 50% ( Figures 3A and 3B). The dislocation gliding mechanism in HP-like regimes and the GB deformation mechanism in IHP-like regimes are observed by strain fields (Figure 4), deformation steps on a reference line ( Figures 6A-6C), and displacement vector fields ( Figures 6D-6F). Steps in the reference line are along the lattice direction due to dislocation gliding in the crystalline regions, but localized at the center of the swirl in each thick GB. Particles in triple junctions have the maximum displacements. Particles' displacements in each crystalline grain exhibit a two-in-two-out flow around a fixed point, that is, a topological charge of -1 commonly observed in liquid crystals. A dislocation gliding through the field can cut the topological charge into two.
To the best of our knowledge, HP and IHP behaviors have only been generalized by tuning the twinboundary spacing inside grains of polycrystal in simulation [36]. In the current study, our generalizations along the other two directions bridge the crystal and amorphous solids beyond polycrystals. As a by-product of the generalization, both ref. [36] and our results provide ways to exceed the conventional maximum strength of σ y (D) by tuning new parameters. Which microstructure gives the maximum strength is an important question. The strength of normal polycrystals reaches the maximum at l ≃ 15 nm, that is, the boundary between the HP and IHP regimes. Here, we show that expanding GBs can further enhance the strength and reach the maximum at l * (D), as shown in Figure 3C. Besides strength, the ductility, fatigue life and other material properties can be similarly studied as a function of (D, l) in this system, which can guide the fabrication of high-strength materials. This system can also be used to study the polycrystal-glass transformation via a new approach by changing l, in contrast to changing D in ref. [51]. These results are of basic importance in materials science and have practical implications for GCC and super alloy fabrications.

Data availability
Data that support the plots within this paper are available from the corresponding authors upon reasonable request. Source data and simulation code are available for this paper.