Probing the Effective Young's Modulus of ‘Magic Angle’ Inspired Multi‐Functional Twisted Nano‐Heterostructures

Two‐dimensional (2D) materials are crucially important nanomaterials because of their exciting multi‐functional properties. However, a single layer of 2D materials may not have a certain property adequately, or multiple application‐specific properties simultaneously to the desired and optimal level. For mitigating this lacuna, a new trend has emerged recently to develop nano‐scale engineered heterostructures by stacking multiple layers of different 2D materials, wherein each of the layers could also be twisted. The vast advantage of combining single layers of different 2D materials with different twisting angles has dramatically expanded this research field well beyond the scope of considering a 2D material mono‐layer, leading to a set of multifunctional physical properties corresponding to each possible combination of number of layers, different 2D materials therein, stacking sequence and the twisting angle of each layer. Effective mechanical properties such as Young's moduli are generally of utmost importance for analyzing the viability of such engineered nano‐heterostructures in various nanoelectromechanical applications. Efficient closed‐form generic formulae are proposed for the effective Young's moduli of twisted multi‐layer heterostructures. Based on this physics‐based analytical approach, a wide range of insightful new results are presented for twisted heterostructures, covering mono‐planar and multi‐planar configurations with homogeneous and heterogeneous atomic distributions.


Introduction
Mono-layer 2D materials and their multilayer Van der Waals heterostructures have been a topic of immense interest in the recent time since the discovery of graphene and its extraordinary mechanical, electronic, optical, thermal and chemical properties. The notion of superconductivity was predicted in the early stages of graphene research. After the recent exotic discoveries of correlated insulators and unconventional superconductivity in twisted bilayer graphene (TBG) by Cao et. al., [1] it has been getting a significantly increased attention from the scientific community. These phenomena can be observed in a narrow range of twisting angle (i.e., in-plane relative rotation between the two layers) near 1.05 • , which is referred as the first magic angle with isolated and relatively flat band near the chemical potential. Although in the last 2 years the origin of magic angles has not become clear, this can certainly be said that the magic angle superconductivity would have a range of realistic and exciting applications across different domains. As there is a twist involved to generate superconductivity of the two layers of graphene, the multi-physical properties of the heterostructure would become strikingly different from a graphene monolayer depending on the angle of twist. The vast possibilities of having exotic electronic structures in such engineered systems reflect the fact that TBG can be a playground for various kinds of rich physics, leading to tunable electronic properties depending on the application-specific requirements. It can be noted here that such twisting can be realized in multiple layers of a graphene based nanostructure, resulting in a different set of multi-physical properties (such as mechanical, optical, chemical, and thermal, besides electronic) corresponding to each possible combination while constructing the multi-layer structure with same material. Not only that, multiple 2D materials (such as graphene, hBN, MoS 2 , and stanene) could be stacked together to form a heterostructure where each layer could have different twisting angle. This opens up an entirely new domain of research for creating multi-functional metamaterials at nano-scale.
Superconductivity has been studied theoretically in graphene and other bi/multilayer heterostructures. [2][3][4][5] Guinea [2] showed that graphene heterostructures can achieve superconductivity by lowering of the average long wavelength interactions between layers in combination with a significant local repulsion at experimentally accessible temperatures, carrier densities above 5 × 10 13 cm −2 , and interlayer distances in the range of a few nanometers. Xu et. al. [6] demonstrated a strongly coupled graphene/2D superconducting Mo 2 C vertical heterostructure by chemical vapour deposition showing the possibility of superconductivity. Along with the superconductivity, various other novel photonic and optoelectronic properties are found to emerge in transition metal dichalcogenides (TMDCs) depending on the interlayer distances and alignments. [7] In general, the research track of 2D materials since the inception through the experimental extraction of a single layer of graphene is quite fascinating. As logically expected, research concerning pure graphene has reached to a matured stage after multiple decades of intensive investigation. Thus, research on other 2D materials [8,9] started receiving the due attention few years ago. While exploration of different monolayers of the wide range of 2D materials is still at the center of interest to the scientific community, researchers have recently started investigating the combination of these monolayers where a multi-layer nanoheterostructure [10][11][12][13][14][15][16] can be formed to enhance and modulate different physical properties like mechanical, electronic, optical, and thermal (refer to Figure 1). For instance, TMDCs such as MoS 2 exhibit rich electronic and piezoelectric properties, but their weak in-plane mechanical properties is a constraint for practical applications. In contrast, graphene is well-known for strong in-plane mechanical properties, while it is extremely soft in the out-ofplane direction with a very low bending stiffness. On the other hand, bending stiffness of MoS 2 is comparatively much higher due to significantly more single-layer thickness. With such complementary physical properties of graphene and MoS 2 , it is quite rational to combine these two materials for forming a graphene-MoS 2 heterostructure, that would possess the desired level of electronic properties and mechanical properties in the in-plane and out-of-plane directions.
The possibility of modulating multi-functional properties simultaneously has attracted significant amount of attention from the scientific community towards developing different nanoheterostructures. The advantage of combining single layers of different 2D materials has dramatically expanded the field of this research well beyond the scope of a single layer of 2D material including graphene or its derivatives. Similar to what was observed in case of graphene few years ago, the interest in such heterostructures is further growing rapidly with the advancement in scalable synthesizing capability of these materials. In the coming years, the attentiveness is expected to expand exponentially with the possibility to consider different tunable nanoelectromechanical properties of the prospective combination of so many two-dimensional materials. The impetus in this direction can be increased by many folds if the aspect of twisting each layer of the heterostructure is considered, leading to a set of multifunctional physical properties corresponding to each possible combination of number of layers, different 2D materials therein, stacking sequence and the twisting angle of each layer (refer to Figure 2 where the twisting of a two-layered nanostructure is explained).
Although there is a significant amount of scientific literature on the electronic, optical, thermal, and chemical properties, understanding of the mechanical properties of heterostructures are still fairly limited. Effective mechanical properties such as Young's moduli are of utmost importance for analysing the viability of such heterostructures in various nanoelectromechanical systems like nano-scale devices. The mechanical properties of single and multi-layer 2D materials along with their heterostructures have been studied without consideration of twist in any layer. [12,20,[30][31][32][33][34][35] Since the interest in twisted multi-layer 2D material heterostructures is picking up rapidly, there is a strong rationale to develop simple closed-form formulae that would be able to readily obtain the Young's moduli of the heterostructure corresponding to any combination of number of layers, different 2D materials therein, stacking sequence, and layer-wise twisting angles.
Our motivation of the current work initiated with the report of superconductivity in twisted "magic angle" bi-layer graphene. In such 2D material superlattices, it is also extremely important to characterize the mechanical properties so that the nanostructure could be utilized in a practical application. However, we soon realized that twisting in the layers of a generic multilayer heterostructure to modulate electronic, optical, mechanical, and thermal properties is just the beginning of an exciting research domain (such as the recent works [44][45][46][47][48] ), where various other individual and multiple of magic-like properties could potentially be identified corresponding to different combinations of twisting angle and stacking sequences of the 2D materials. In all such twisted heterostructures, Young's modulus is one of the most important mechanical features to access the viability of their adoption in many nanoelectromechanical devices and systems. In this article, we endeavor to systematically unfold the Young's moduli of twisted multi-layered heterostructures based on an efficient mechanics-based analytical approach. We would start with characterizing the Young's moduli of twisted bi-layer graphene at the 'magic angle,' followed by the Young's moduli of bi-layer graphene at any other possible twisting angles. Thereafter, the developed analytical formulae will be extended to a more general case of finding the Young's moduli for any multi-layer heterostructure corresponding to any possible combination of twisting angle of the constituting layers. The broad aim of this article is to show how the proposed efficient analytical models along with insightful new results will accelerate the domain of novel 2D multi-functional material innovation for achieving a single feature or multiple physical properties simultaneously to a desirable extent. In the following section, we would first provide the mechanics-based theoretical background for the effective Young's moduli of twisted heterostructures leading to their closed-form expressions, followed by numerical investigations considering various category of the nanostructures, stacking sequence, and twisting angle.

Closed-Form Analytical Expressions for Young's Moduli of Twisted Multi-Layered Heterostructures
Based on the placement of atoms in a single-layer and the nature of atomic distribution, 2D material lattices can be categorized for a better understanding of their mechanical behavior. In case of mono-planar structural configurations the atoms are placed in a single plane to form the lattice-like forms (such as graphene and hBN), while the atoms are placed in multiple planes in Figure 1. General overview of the nano-heterostructures and their structural classification. A) Top and side views of typical mono-layer nanostructures where a single atom (homogeneous atomic distribution) is used to form the structure in a single plane (such as the atomic structure of graphene). B) Top and side views of typical mono-layer nanostructures where more than one type of atom (heterogeneous atomic distribution) is used to form the structure in a single plane (such as the atomic structure of hBN and BCN). C) Top and side views of typical mono-layer nanostructures where a single atom (homogeneous atomic distribution) is used to form the structure in multiple planes (such as the atomic structure of silicene, germanene, and stanene). D) Top and side views of typical mono-layer nanostructures where more than one atoms (heterogeneous atomic distribution) are used to form the structure in multiple planes (such as the atomic structure of MoS 2 , WSS 2 , MoSeS 2 , WSeS 2 , and MoTeS 2 ). E) Typical 3D view and side views of a three-layer heterostructure made of two different 2D materials. It may be noted that different combination of 2D materials could be stacked to form such heterostructures following any stacking sequence with different number of layers. F,G) Most generalized view of a 2D material mono-layer with multi-planar structural form. The in-plane angles ( and ′ ) and out-of-plane angle ( ) are indicated in the figure, where = 90 • − ′ ∖2. It may be noted that different combinations of mono-planar and multi-planar structures with homogeneous and heterogeneous atomic distributions could be obtained as special cases of this generic nano-structure by considering the appropriate value of out-of-plane angle and necessary arrangement of atoms. H,I) Deformation induced energy components associated with the in-plane (x--y plane) and out-of-plane (x-z plane) mechanisms. The directions x and y are indicated in the figure, while the direction z is perpendicular to the x-y plane. For a heterogeneous atomic distribution, A and B indicate two different atoms, while A and B are same for a homogeneous atomic distribution. multi-planar 2D materials (such as stanene and MoS 2 ). The unit cell [49][50][51][52] of a 2D material lattice could either be formed using multiple atoms (such as hBN and MoS 2 ) or one single atom (such as graphene and stanene), based on which the 2D materials could further be classified to have heterogeneous and homogeneous atomic distribution, respectively. Depending on the mono or multi-planar configuration and homogeneity or heterogeneity in the atomic distribution, the nanostructures of 2D materials can be classified into four different groups as shown in Figure 1A-D. Accordingly, the nano-heterostructures, made of such 2D materials, can be classified into three broad categories: heterostructures having only mono-planar nanostructures (e.g., graphene-hBN heterostructure [17][18][19] ), heterostructures having both mono-planar and multi-planar nanostructures (e.g., heterostructures made of graphene-MoS 2 , [12,20] graphenestanene, [21] phosphorene-graphene, [22] phosphorene-hBN, [22]  . D) Typical representation of twist (rotation of one layer by and angle with respect to the other layer) in a two layer heterostructure (note that more than two layers with a single 2D material or multiple 2D materials could be developed). We essentially consider 1-2 as the reference global axis system of the entire heterostructure, while x-y is the local axis system for each of the layers. Here the direction -3 (or z) is along thickness of the heterostructure. and multi-layer graphene-hBN-TMDC heterostructure [23] ), and heterostructure having only multi-planar nanostructures (e.g., stanene-MoS 2 heterostructure, [24] MoS 2 -WS 2 heterostructure [25] ). For twisted nano-heterostructures, it is interesting to further classify the heterostructure from the viewpoint of multi-layered structures made of a single 2D materials and multiple 2D materials. It is noteworthy that the mechanical properties can be modulated significantly even in case of twisted multi-layer nanostructures composed of a single 2D material depending on the relative angle of twists of each layer.
Common practices for investigating the 2D materials and their heterostructures include molecular dynamics simulation and first principle studies/ ab-initio, which are able to reproduce the experimental results closely in the expense of computationally intensive and time consuming supercomputing facilities. In case of molecular dynamics simulations involving built-up nanostructures with multiple materials, the availability of interatomic potentials could be an inevitable practical barrier for the analysis. Molecular-level mechanics based formulae for the Young's moduli have been developed recently for single and multi-layer untwisted 2D materials. [26][27][28][29][30] This approach of elastic property characterization appears to be computationally very efficient, yet accurate and physically insightful. However, the molecular mechanics based analytical models concerning 2D nano-structures presented so far are limited to conventional 2D materials and their heterostructures; efficient analytical approaches have not yet been attempted for twisted heterostructures. Considering the upcoming potential of research in the field of twisted heterostructures covering their multi-functional characteristics, there is a strong rationale to develop computationally efficient closed-form formulae for the elastic moduli of these built-up artificial multilayers. Such efficient models are particularly relevant for twisted multi-layered 2D materials, since there exists a plethora of combinations corresponding to each possible twisting angle that would have different mechanical properties. It is practically impossible to carry out expensive molecular dynamics simulations or ab initio study corresponding to each of these combinations in order to fully characterize the material properties. The closed-form analytical formulae that we present in this section can serve as a ready reference to obtain the Young's moduli without carrying out expensive and time consuming computer simulations or laboratory experiments. The analytical framework for elastic moduli of twisted nano-heterostructures would be applicable to any possible combination of the twisting angles and any number of different constituent 2D material mono-layers with multi-planar or mono-planar hexagonal nanostructures. We would first discuss the mechanical equivalence of atomic bonds and thereafter the effective elastic properties of the twisted nano-structures would be determined following a multi-stage bottom-up approach. The molecular mechanics based mechanical equivalence of atomic bonds is widely reported in scientific literature. [27,36,37] Therefore, the contribution of this paper lies in developing computationally efficient and easy-to-implement generic analytical formulae for twisted nano-heterostructures and thereby presenting new insightful results corresponding to different twisting angles for various stacking sequence of the heterostructures belonging to different classes as described in the preceding paragraphs (made of mono-planar and multiplanar 2D materials with homogeneous and heterogeneous arrangement of the atoms such as graphene, MoS 2 , hBN, and stanene).
For atomic-level deformation behavior of nano-materials, the total potential energy (E) of an idealized lattice-like structure (refer to Figure 1) is expressed as a sum of various energy terms concerning bonding (such as bending energy [E b ], stretching energy [E s ], and torsion energy [E t ]) and non-bonding (E nb , such as van der Waals attraction, coulombic energy, and core repulsions) atomic interactions [36] (refer to Figure 1H,I).
Among all the bonding and non-bonding energy components, energy contributions due to bending and stretching are most predominant in case of small deformation analyses. [27,37] For multiplanar hexagonal nano-lattices (such as MoS 2 ), the bending energy consists of in-plane deformation component (E bI ) and out-of-plane deformation component (E bO ). The out-of-plane contribution of bending energy reduces to zero for monoplanar nano-lattices (such as hBN). Thus, for a general small deformation analysis, the effective total interatomic potential energy (E) is given as where Δ , Δ , and Δl are the change in-plane angle, out-of-plane angle, and bond length, respectively. Here k r and k denote the atomic force constants concerning stretching and bending of bonds, respectively. An atomic bond can be idealized as a beam with circular cross-section, leading to a lattice-like structure for the 2D materials. To find out equivalent elastic properties of a lattice-like structure, deformation behavior of each of the constituting beam elements needs to be analyzed. Deformation analysis of a beam requires stiffness parameters of the beams. Thus, it is required to establish a mechanical equivalence between the molecular mechanics parameters (k r and k ) and structural mechanics parameters (such as axial and bending stiffness, EA and EI) so that the molecular mechanics parameters can be utilized in deformation analysis of the beams using structural mechanics. The molecular and structural mechanics parameters of a circular beam with Young's modulus E, length l, uniform cross-sectional area A, and area moment of inertia I, can be bridged as: K r = EA l and k = EI l . [27,36,37] Based on this mechanical equivalence, the closed-form analytical expressions for effective Young's moduli of multilayer twisted heterostructures can be derived following a multi-stage bottom-up idealization scheme enforcing the conditions of force equilibrium and deformation compatibility (refer to Figure 3 and the Supporting Information for detailed derivation). Two in-plane effective Young's moduli of twisted heterostructures are expressed as where (•) = (•) × t in the above equation represents rigidity of the 2D material mono-layer in longitudinal, transverse and shear directions (i.e., The subscript i is used to indicate parameters corresponding to the i th layer. Here, total thickness t = ∑ n i=1 t i ; t i represents the thickness of i th layer (i = 1, 2, 3 … n). The parameter n denotes the total number of layers and i is the twisting angle of i th layer. The expressions of E xi , E yi , xyi, and G xyi in the above expressions are given as Thus, the effective Young's moduli of a twisted hexagonal heterostructure can be obtained using the closed-form analytical formulae described here based on the interatomic bond parameters such as k r , k , bond length (l), in-plane bond angle ( ), and out-of-plane angle ( ). These parameters are well-documented in the scientific literature of chemistry and material science, dealing with molecular mechanics (refer to Experimental Section for specific values of these parameters in case of the 2D materials under consideration in the present manuscript). The analytical formulae derived in this section are valid for linear analysis, that is, small deformation of the 2D material nanostructures. It can be noted that the influence of inter-layer stiffness due to Lennard-Jones potentials (that affects out-of-plane inter-layer deformation behavior) are negligible for the in-plane Young's moduli, as considered in the present study. An advantage of the developed bottom-up multi-stage approach of considering layerwise effective material properties of the twisted layers is that it allows us to ignore the effect of lattice mismatch in obtaining the effective Young's moduli for twisted heterostructures. Since the focus of this paper is on global equivalent properties like Young's moduli instead of local analysis of stresses and strains, the proposed layer-wise bottom-up approach is found to be effective in the analysis. In the derivation, we have satisfied deformation compatibility conditions of the adjacent 2D material layers, which would give rise to local strain energy at the interfaces. [12] However, effective Young's moduli being global structure-level properties as mentioned above, the influence of such local effects become negligible. It can be discerned from the derived formulae that the Young's moduli of twisted nano-heterostructures actually depend on the number and elastic properties of the layers of different constituent 2D materials rather than their exact stacking sequences, along with the individual twisting angles. An important aspect that immediately comes out from the analytical formulae of Young's moduli (refer to Equations 3 and 4) is that E 1 ≠ E 2 in the general case (note that E 1 and E 2 may be equal in certain special cases of nanostructural configuration, which is discussed in the following subsection). This outcome is on the contrary to the common notion of Young's modulus considering the standard practice for 2D materials and their heterostructures followed in literature, where the difference between E 1 and E 2 are normally not recognized and a single value of Young's modulus is reported. Introduction of twist in nano-heterostructure, in fact extends the design space for possible tailoring of anisotropy ratio ( ) as a function of the constituent 2D materials, their stacking sequence, and layer-wise twisting angle. Considering that heterostructures (including the variant of twisted heterostructures) are normally attractive for achieving multiple desired properties simultaneously, the possibility of modulating anisotropy provides a direct measure of controlling a range of direction-dependent dynamic features of the nanostructure. Thus, the derived closedform formulae of Young's moduli would be helpful in efficient characterization of different multi-functional properties simultaneously such as direction-dependent dynamic features, deformation, and stress components.
In the analytical closed-form expressions for effective Young's moduli of twisted multi-layer nano-heterostructures, each individual layer may be of multiplanar (i.e., out-of-plane angle ≠ 0) or monoplanar (i.e., out-of-plane angle = 0) configuration, and they may have any twisting angle (refer to Equations 3 and 4). It is interesting to note that the generalized closed-form formulae for the effective Young's moduli of twisted nano-heterostructures can be reduced to the analytical formulae provided by Shokrieh and Rafiee [27] for single-layer graphene when we consider n = 1, The special case of the most generalized analytical formulae for effective Young's moduli of twisted heterostructures till date, as presented above, essentially provides an exact analytical validation. For a single mono-layer of 2D material, the expressions of elastic moduli deduced from Equations 3 and 4, obtained by replacing n = 1 and = 0 • , perfectly satisfies the reciprocal theorem (i.e., E 1 21 = E 2 12 ). [38][39][40] In case of multi-layer 2D materials composed of layers having same material and without any twist in the layers (i.e., bulk material), it is expected from Equations 3 and 4 that the values of Young's moduli would reduce compared to that of mono-layer 2D materials due to the  [20] , 1 ± 0.1 [42] G/G 1.0419 1.0419 1.06 [20] , 1.04 ± 0.1 [41] M 0.1778 0.3549 0.16 [20] , 0.27 ± 0.1 [43] M/M 0.1778 0.3549 0.27 [20] , 0.2 ± 0.1 [43] G/M 0.4893 0.6025 0.53 [20] , 0.49 ± 0.05 [25] G/M/G 0.6357 0.7189 0.68 [20] , 0.56 [12] M/G/M 0.3678 0.5059 0.45 [20] presence of inter-layer spaces (which increase the effective overall thickness t). This observation has been widely reported in the scientific literature. Noteworthy is that the analytical formulae can provide such critical insights on the effective elastic properties of single-layer and bulk materials.

Validation and Numerical Predictions for the Young's Moduli of Twisted Heterostructures
Numerical results are presented for the effective Young's moduli of twisted multi-layer heterostructures using the closed-form formulae (refer to Equation 3 and 4) derived in this article.
In the preceding section, we have provided an exact analytical validation by considering a special case of the formulae for mono-layer mono-planar structures. In addition to the exact analytical validation, we would provide more numerical validations in this section before presenting numerical results on twisted multi-layer heterostructures. As investigation on such twisted nano-heterostructures is a novel and emerging field of research (specifically their mechanical characterization), the numerical results available for validating the proposed formulae of heterostructures is scarce to find in scientific literature. First, we consider multi-layer, but untwisted heterostructures composed of graphene (G) and MoS 2 (M) with different stacking sequences.
Comparative results are provided in Table 1 obtained from the developed formulae considering = 0 • and results from available literature, wherein a good agreement can be noticed for all the configurations. Further, we consider validation of twisted multi-layered heterostructures considering bi-layer configurations of graphene and MoS 2 with multiple twisting angles. Since Young's moduli of such twisted heterostructures are not available in published scientific literature, we carried out separate molecular dynamics simulations (refer to the methods section for more detail) and compared the Young's moduli obtained using the developed analytical formulae and the molecular dynamics simulations (refer to Figure 4). A good agreement can be noticed between the results obtained by two different means (considering bi-layer nanostructures with same and different materials with monoplanar and multiplanar configurations), corroborating the validity of the derived closed-form formulae. Thus, essentially we provide three different forms of validation systematically, considering single-layer of 2D materials (using exact analytical validation), multi-layered untwisted heterostructures (using numerical results from available literature), and multi-layer twisted heterostructures (using separate molecular dynamics simulations). Having adequate confidence on the proposed analytical framework, we have used the closed-form formulae for Young's moduli of twisted heterostructures to present results for different 2D materials, their stacking sequences, and angles of twist. Though the developed closed-form formulae for effective Young's moduli of twisted heterostructures are applicable to any 2D material and for every possible stacking sequences, we present the numerical results here considering four different nanostructures covering mono-planar and multi-planar configurations with homogeneous and heterogeneous atomic distributions: graphene (G), hBN (H), stanene (S), and MoS 2 (M). While most of the heterostructures considered here have been investigated for different physical and chemical properties recently, only the untwisted form of graphene-MoS 2 heterostructure has received some attention in terms of Young's moduli using molecular dynamics simulation. [12,25] Recently, the elastic properties of other heterostructures like graphene-hBN, graphene-MoS 2 , graphene-stanene, and stanene-MoS 2 have been investigated considering untwisted forms. [35] The numerical results for Young's moduli of twisted heterostructures are presented for the first time in this article.
We start presenting the numerical results by considering bi-layer twisted structures where both the layers are of same 2D material (refer to Figure 5). The variation of Young's modulus for armchair (E 1 ) and zigzag (E 2 ) directions with twisting angle are shown separately. The general trend of variation with twisting angle is found to be similar for all the considered 2D materials like graphene, hBN, MoS 2 , and stanene, albeit with different extent of fluctuation with respect to the corresponding untwisted configurations. The maximum value of Young's moduli can be achieved when one layer is rotated by 90 • with respect to the other layer. The other two peaks appear at 0 • and 180 • . The relative peak values at 0 • , 90 • , and 180 • depend on the 2D material under consideration.
It is worthy to mention here that the difference between E 1 and E 2 is not explicitly recognized in the majority of the scientific literature and the results are normally presented as E 1 = E 2 (for example, refer to the last column of Table 1). From the numerical results presented here, a general trend can be noticed for single-layer 2D materials. The Young's moduli E 1 and E 2 are different for multiplanar single-layer nanostructural forms (such as MoS 2 ), while E 1 and E 2 are same for monoplanar nanostructures (such as hBN). The numerical results for Young's moduli of graphene and MoS 2 obtained using the proposed analytical framework, as presented in Table 1, follow this trend. Similar observations are found to be reported in few published literature. [35,53] In case of multi-layered twisted nanostructures composed of same 2D material (e.g., refer to the results presented in Figure 5), corresponding to any particular twisting angle, we find E 1 = E 2 in case of monoplanar 2D materials, while E 1 ≠ E 2 in case of multiplanar 2D materials. In general, the twisted or untwisted heterostructures having at least one layer of multiplanar (i.e., ≠ 0) 2D material is expected to show different E 1 and E 2 values (corresponding to = 0 • or any particular value). From Table 1, it is noticed in case of heterostructures formed using graphene and MoS 2 that the value of E 2 is higher than E 1 . After investigating the Young's moduli of single-layer 2D materials and bi-layer nanostructures where both the layers are composed of same 2D materials and one layer is fixed while the other layer is twisted (refer to Table 1 and Figure 5), we would systematically investigate different other nanostructural configurations in the following paragraphs. It can be noted that the contour plots presented in terms of the twisting angles (refer to Figures 6-12 and Figures S3 to S13, Supporting Information) are capable of readily providing the Young's moduli corresponding to any possible combination of the angles along with the critical insights on general trend of the variation.
The results for Young's moduli of bi-layer nanostructures composed of same 2D materials are extended to the case of twist in both the layers (refer to Figure 6 and Figure S3, Supporting Information). From the figures for E 1 and E 2 concerning any combination of twisting angles corresponding to the two layers, it may be noticed that E 1 = E 2 in case of graphene and hBN bi-layers (i.e., monoplanar configuration), while E 1 ≠ E 2 in case of MoS 2 and stanene bi-layers (i.e. multiplanar configuration). A periodic 2D variation in the space of two twisting angles corresponding to the two layers can be noticed in all the plots, albeit the extent and exact nature of fluctuation depends on the nanostructural configuration.
In case of twisted bi-layer nanostructures, so far we have investigated the configurations where both the layers are made of same material. We would now explore the Young's moduli of bi-layer nanostructures made of two different 2D materials. First, we consider graphene-based twisted bi-layer heterostructures, where one layer is graphene and the other layer could be one of the three 2D materials among hBN, stanene and MoS 2 (refer to Figure 7). It can be noticed from the figures that in case of graphene-hBN www.advancedsciencenews.com www.advtheorysimul.com  heterostructure, E 1 = E 2 for any particular combination of twisting angles corresponding to the two layers. However, this observation does not hold good (i.e., E 1 ≠ E 2 ) for the other two cases where one of the layers is made of multiplanar 2D material like stanene and MoS 2 . Interestingly, the nature of variation for E 1 or E 2 along the horizontal and vertical direction varies significantly in such cases. The peak of Young's moduli is found to be realized at different possible combinations of the twisting angles 0 • , 90 • , and 180 • corresponding to the two layers. Figure 8 presents the Young's moduli of twisted bi-layer nanostructures where each of the layers is different and none of the layers is made of graphene. Since at least one layer in this form of nanostructures has a multiplanar configuration, the trend of E 1 and E 2 are found to be considerably different with the variation in twisting angles corresponding to the two layers. Figures 7C,D and 8A,B,E,F show that an extended zone of higher Young's moduli in the design space of twisting angles corresponding to bi-layer heterostructures can be obtained based on nanostructural configuration. Noting that the proposed closed-form formulae are capable of predicting the Young's moduli for any number of layers having any individual twisting angles, now we focus on investigating the Young's moduli of three-layered twisted heterostructures in the following paragraphs. First, we have considered three-layered twisted nanostructures made of a single 2D material such as graphene, hBN, stanene, and MoS 2 (refer to Figure 9 and Figure S4, Supporting Information). Even though it is readily possible to obtain the Young's moduli of the nanostructures considering simultaneous twisting of all the three layers using the proposed closed-form formulae, for presenting 2D contour plots, we have considered one of the layers to have no rotation in the figures. The other two layers are rotated (/twisted) simultaneously to present results corresponding to every possible combination of the twisting angles. Similar to the cases of bi-layer nanostructures, we can notice from the results that E 1 and E 2 become numerically equal in case of monoplanar configurations like graphene and hBN. It is interesting to note here that there is no way to modulate the Young's moduli of untwisted multi-layered nanostructures made of a single 2D material. Introduction of twist in the individual layers brings an endless possibility of having different values of Young's moduli (and other physical properties) even if the nanostructure is composed of a single 2D material. After investigating three-layered twisted nanostructures composed of a single 2D material, we have explored the Young's moduli of three-layered nanostructures made of different 2D materials in the increasing order of structural complexity (refer to Figures 10-12, and Figures S5-S13, Supporting Information). Even though it is readily possible to obtain the Young's moduli of the nanostructures considering simultaneous twisting of all the three layers using the proposed closed-form formulae, for presenting 2D contour plots, we have considered one of the layers to have no rotation in the figures. The other two layers are rotated (/twisted) simultaneously to present results corresponding to every possible combination of the twisting angles. Note that we have considered four different 2D materials in this paper, namely graphene, stanene, hBN, and MoS 2 . Let us represent the materials corresponding to the three layers of the heterostructure as L 1 , L 2 , and L 3 for the ease of understanding the nanostructural configurations. As explained above, we have kept the layer L 1 untwisted to obtain the numerical plots. Figure 10 and Figures S5, S8, and S11, Supporting Information, present the Young's moduli of heterostructures where L 1 is one of the four 2D materials considered in this paper, while L 2 and L 3 are of the same material chosen from the remaining three materials other than L 1 . Figure 11 and Figures S6, S9, and S12, Supporting Information, investigate heterostructures where L 1 is one of the four 2D materials considered in this paper and L 2 is graphene, while L 3 is chosen from the remaining two materials other than L 1 and L 2 . Figure 12 and Figures S7, S10, and S13, Supporting Information, present the most complex heterostructures where L 1 is one of the four 2D materials considered in this paper, and  three layers. Figures 11C,D and 12A,B,E,F and Figures S6C,D,  S7A,B,E,F, S9C,D, S10A,B,E,F, S12C,D, and S13A,B,E,F, Supporting Information, indicate that a rather extended zone of higher Young's moduli in the design space of twisting angles corresponding to the three layers can be obtained based on nanostructural configuration. We can notice that the Young's moduli in two perpendicular directions (longitudinal and transverse) become equal corresponding to any particular configuration when all the three layers are composed of monoplanar structures. The peak of Young's moduli is found to be realized at different possible combinations of the twisting angles 0 • , 90 • , and 180 • corresponding to the rotating layers. In general, the numerical plots including validation and new results show that it is possible to obtain the Young's moduli of the heterostructures corresponding to different combination of twist angles of the layers, leading to an efficient mechanical characterization for assessing their viability as structural components besides other multi-functional applicability. Exact numerical values of the Young's moduli corresponding to various 2D material combinations, stacking sequences, and twisting angles for any number of layers (n) can be readily obtained using the computationally efficient closed-form analytical Equations (3) and (4).

Conclusions and Perspective
The physics-based closed-form expressions for twisted heterostructures presented in this paper can obtain the Young's moduli corresponding to any twisting angle and stacking sequence of the multi-layer 2D material nano-structures. From the analytical expressions (refer to Equations (3) and (4)), it can be noted that the Young's moduli of such twisted multi-layered nanostructures are actually functions of the layer-wise material and twisting angle of the constituting individual layers. This observation is on contrary to the case of conventional untwisted nano-heterostructures where the Young's moduli depend on the number of layers of different 2D materials only rather than their exact stacking sequences. It is interesting to note here that there is no way to modulate the Young's moduli of untwisted multilayered nanostructures made of a single 2D material. Introduction of twist in the individual layers brings an endless possibility of having different values of Young's moduli (and other physical properties individually and simultaneously), even if the nanostructure is composed of a single 2D material.
Mechanical behavior of nanostructures portrayed by effective properties such as Young's moduli is normally of critical significance to access the viability of the application of nanoheterostructures as structural components in various nano-scale systems and devices. It is explained in the paper that efficient analytical formulae of the two Young's moduli (E 1 and E 2 ) could lead to anisotropy tailoring in twisted heterostructures for multifunctional characterization of several properties simultaneously, such as direction-dependent dynamic features, deformation, and stress components. The main contribution of this paper is to develop the most generalized form of closed-form analytical expressions to obtain the Young's moduli of multi-layer twisted nano-heterostructures. As special cases, these efficient formulae are also readily applicable to traditional untwisted heterostructures ( i = 0) with any stacking sequence and single-layer of materials (n = 1) with monoplanar and multiplanar configurations. Thus, the generic formulae for Young's moduli can be applied as an efficient ready-reference covering the entire spectrum of lattice-like 2D materials and the multi-layer heterostructures (twisted and untwisted) obtained by combining different layers of these materials. Such generalization in the proposed closed-form formulae, with the coupled advantage of being efficient and easy-to-implement, essentially presents a tremendous potential research scope in front of the scientific community for novel application-specific heterostructure development in a significantly expanded design space.
To achieve adequate confidence on the proposed analytical formulae, we have provided three different forms of validation systematically, considering single-layer of 2D materials (using exact analytical validation), multi-layered untwisted heterostructures (using numerical results from available literature), and multi-layer twisted heterostructures (using separate molecular dynamics simulations). In-depth novel numerical results are presented for the effective Young's moduli of heterostructures composed of different possible combinations of graphene, hBN, stanene, and MoS 2 , covering mono-planar and multi-planar configurations with homogeneous and heterogeneous atomic distributions. Even though the numerical results are obtained in this paper considering four different constituent 2D materials, the proposed formulae can be used for other nano-heterostructures containing any number and type of different 2D materials along with their layer-wise twisting angles. The bottom-up physicsbased framework leading to closed-form formulae is capable of providing a comprehensive analytical insight on the mechanical behavior of multilayer twisted heterostructures, including their quantifiable dependency on the input parameters of the design space such as bond angle, twisting angle, out-of-plane angle, bond length, nature, and quantity of the constituting 2D material layers. Noteworthy feature of such an analytical framework is the inexpensive and computationally efficient nature compared to conducting molecular dynamics simulations or nano-scale experiments. Besides in-depth deterministic analysis of Young's moduli, as presented here, the efficient analytical approach could be an attractive alternative for performing uncertainty quantification, [54][55][56] which is increasingly getting recognized as one of the most crucial aspects of nano-scale analyses. The rapid systematic scientific development in the field of 2D materials leading to the current focus on twisted heterostructures is noteworthy. After several years of exhaustive research, the scope of new exciting developments concerning pure graphene has logically come to a matured and rather saturated stage. Thus, the research focus expanded to other graphene-like 2D materials (such as stanene, hBN, MoS 2 ) in the last few years. However, the recent trend of combining single layers of different 2D materials to form a heterostructure has expanded the design space of developing new multi-functional nano-structures well beyond the scope of any single layer 2D material. Considering the prospect of multiple tunable nanoelectrochemomechanical properties of the prospective combination of so many 2D materials in different possible stacking sequences, the research activities in the field of heterostructures is growing exponentially, similar to the growth of interest in graphene a decade years ago. The attentiveness of heterostructures is expected to dramatically expand in the foreseeable future with the recent disruptive development of having twist in the 2D material layers. This essentially brings about the possibility of introducing a new dimension of multi-functional nano-scale metamaterial [57][58][59][60] development with effectively unlimited design space. Although there is a significant amount of scientific literature on the electronic, thermal, optical, and chemical properties, understanding of the mechanical properties of heterostructures are still fairly limited. The elastic properties of twisted heterostructures have not been investigated yet. Since the interest in such twisted multi-layer 2D material heterostructures is picking up rapidly, there is a strong rationale to develop simple closed-form formulae that would be able to readily obtain the Young's moduli of the heterostructure corresponding to any combination of number of layers, different 2D materials therein, stacking sequence, and the twisting angle of each layer. www.advancedsciencenews.com www.advtheorysimul.com The present article contributes significantly in this exciting research endeavor.
In summary, we have presented physics-based computationally efficient analytical formulae for predicting effective Young's moduli of twisted multi-layer 2D material based nanoheterostructures. The derived closed-form expressions are extensively validated using available scientific literature and by performing separate molecular dynamics simulations. Subsequently, a broad range of insightful new results are systematically presented for twisted heterostructures composed of different possible combinations of four different 2D materials in the order of increasing structural complexity, covering mono-planar and multi-planar configurations with homogeneous and heterogeneous atomic distributions. As the presented closed-form analytical formulae are generic in character and applicable to a wide range of 2D materials along with their heterostructures with any stacking sequence and layer-wise twisting angle, this paper can effectively serve as an efficient, yet critically insightful reference for characterizing Young's moduli in the development of various future multi-functional nano-materials.

Methods
In this section, first a brief description of the analytical framework is provided for deriving the Young's moduli of twisted nano-heterostructures that leads to the closed-form expressions (Equations (3) and (4)). Thereafter, a brief description of the molecular dynamics simulation (used for validating the analytical framework) is presented.
Analytical Approach for Effective Young's Moduli of Twisted Nano-Heterostructures: Here, a concise account of the basic principle is provided for developing the analytical approach leading to closed-form expressions for Young's moduli of twisted nano-heterostructures. A detailed derivation is provided in Supporting Information. A multi-stage bottom-up framework was followed to develop the expressions for Young's moduli. The lattice-like individual layers of 2D materials could be modeled as equivalent plate-like lamina structures at the initial stage. The multi-layer twisted heterostructure could be idealized (for computational modelling) as a layered plate-like laminated composite system with respective equivalent later-wise elastic properties and geometric dimensions (refer to Figure 3A,B). Each of the equivalent plate-like layers were assumed to be perfectly bonded with adjacent layers to satisfy deformation compatibility condition. Note that two axis systems were considered in the derivations: a local x-y-z system for each layer and a global 1-2-3 system for the entire heterostructure. The x-direction is perpendicular to the vertical atomic bonds in a hexagonal unit of a particular layer, while the y-direction is perpendicular to x-direction. The angle between x-y axis and 1-2 axis for a particular layer represent the twist, as shown in figure Figure 3C. The axis along directions z and 3 coincide with each other along thickness of the heterostructure. The effective elastic moduli (along the local axis system x-y) of each individual layer were determined first based on a molecularlevel mechanics based approach utilizing the mechanical equivalence of bond properties. The Young's moduli of each layer in the global axis system 1-2 were obtained thereafter. Finally, the equivalent Young's moduli of the entire multi-layer twisted heterostructure were determined using deformation compatibility and force equilibrium conditions. The molecular mechanics (k r and k ) and geometric parameters (l i , i , i , and t i ; here the subscript i is used to denote i th layer) of the bonds are well-documented in chemistry and materials science literature. Here, these parameters are provided for the four different materials used in obtaining the numerical results. AMBER force field can be used for graphene to obtain the molecular mechanics parameters k r and k as k r = 938 kcal mol −1 nm −2 = 6.52 × 10 −7 Nnm −1 ; k = 126 kcal mol −1 rad −2 = 8.76 × 10 −10 Nnm rad −2 . [61] The geometric parameters for graphene are: i = 0; i = 30 • ; l i = 0.142nm; t i = 0.34nm. [26] For hBN, DREIDING force model [63] can be used to obtain the molecular mechanics parameters as: k r = 4.865 × 10 −7 Nnm −1 ; k = 6.952 × 10 −10 Nnm rad −2 . [62] The geometric parameters for hBN are given as: i = 0; i = 30 • ; l i = 0.145nm; t i = 0.098nm. [28] The molecular mechanics parameters in case of stanene can be obtaind from literature as: k r = 0.85 × 10 −7 Nnm −1 ; k = 1.121 × 10 −9 Nnm rad −2 . [64,65] The geometric parameters for stanene are: i = 17.5 • ; i = 35.5 • ; l i = 0.283nm; t i = 0.172nm. [64][65][66][67] The molecular mechanics parameters in case of MoS 2 can be obtained from literature as: k r = 1.646 × 10 −7 Nnm −1 ; k = 1.677 × 10 −9 Nnm rad −2 . The geometric parameters for MoS 2 are: i = 48.15 • ; i = 48.54 • ; l i = 0.242nm; t i = 0.6033nm. [68][69][70][71] It can be noted that these molecular mechanics and geometric parameters for other 2D materials can be obtained from literature for obtaining the Young's moduli of various twisted heterostructures composed of a wide variety of materials.
Molecular Dynamics Simulation for Analyzing Twisted Nano-Heterostructures: For the purpose of validating the proposed analytical formulae, separate molecular dynamics simulations were performed to obtain the Young's moduli of twisted graphene-MoS 2 heterostructures. Here, a concise description of these molecular dynamics simulations are provided, which were carried out considering different twisting angles. A similar method was followed as reported in the published literature [12,72,73] for calculating Young's moduli of twisted graphene-MoS 2 individual bilayers and heterostructures using molecular dynamics simulation. The second-generation Brenner interatomic potential was used for modeling carbon-carbon and molybdenum-sulfur interactions, [74,75] while the heterostructures were stabilized following a standard method as per literature. [12] The adjacent graphene and MoS 2 layers were coupled by van der Waals interactions using the Lennard-Jones potential. The cut-off values in the simulations were determined by stabilizing and minimizing the heterostructures. [12]

Supporting Information
Supporting Information is available from the Wiley Online Library or from the author.