Structural and chemical mechanisms governing stability of inorganic Janus nanotubes

One-dimensional inorganic nanotubes hold promise for technological applications due to their distinct physical/chemical properties, but so far advancements have been hampered by difficulties in producing single-wall nanotubes with a well-defined radius. In this work we investigate, based on Density Functional Theory (DFT), the formation mechanism of 135 different inorganic nanotubes formed by the intrinsic self-rolling driving force found in asymmetric 2D Janus sheets. We show that for isovalent Janus sheets, the lattice mismatch between inner and outer atomic layers is the driving force behind the nanotube formation, while in the non-isovalent case it is governed by the difference in chemical bond strength of the inner and outer layer leading to steric effects. From our pool of candidate structures we have identified more than 100 tubes with a preferred radius below 35 {\AA}, which we hypothesize can display unique properties compared to their parent 2D monolayers. Simple descriptors have been identified to accelerate the discovery of small-radius tubes and a Bayesian regression approach has been implemented to assess the uncertainty in our predictions on the radius.


Introduction
In the last decades miniaturization of devices has been a main trend driving the electronics industry. In addition to reducing the usage of raw materials, nanomaterials often show improved properties compared to their larger counterparts. Among these nanomaterials are two-dimensional (2D) sheets, one-dimensional (1D) structures such as nanotubes and nanoribbons, and zero-dimensional (0D) nanoparticles.
Since their discovery, nanotubes have shown promise for a wide range of applications including gas separation and capture, catalysis, solid lubrication and controlled drug delivery. [1] In addition to the well-known carbon nanotubes [2] numerous inorganic nanotubes have been synthesized . [3,4] Although the first successful synthesis of single-wall MoS 2 nanotubes has been reported [5], such structures usually appear together with numerous multi-wall tubes showing a distribution of radii and wall thicknesses. [5] These multi-wall structures alleviate the built-in strain energy through van der Waals interactions in between the layers leading to an increase in stability. [1] Overall, this has made it difficult to establish an experimental synthesis pathway to produce single-wall tubes with a specific radius and controllable physio-chemical properties.
A possible solution to this problem is the approach of considering asymmetric sheets, which can naturally wrap and form nanotubes with a well defined size. Due to the asymmetry, the unsupported sheet is expected to be unstable compared to other curled shapes, such as tubes or scrolls. Pauling, already in the 1930s, mentioned that the driving force of sheets to curve is related to the lattice-mismatch between the two inner and outer atomic layers. [6] Single-wall inorganic nanotubes with well-defined diameters hold promise for technological applications, not only because of their reduced dimensionality, but also for their unique properties, often inherently different from the ones of the corresponding asymmetric sheets. An example of a small-radius, single-wall nanotube formed from an asymmetric sheet is imogolite (Al 2 SiO 3 (OH) 4 ) which was first discovered in volcanic ash soil [7] and later synthesized. [8,9] Other tubular minerals include chrysotile (Mg 3 Si 2 O 5 (OH) 4 ) and halloysite (Al 2 Si 2 O 5 (OH) 4 ) that however occur as multi-wall tubes. [10,11] Besides naturally occurring nanotubes, "misfit-layer" compounds, composed of two separate sheets, make use of the lattice mismatch between the two sheets to induce a natural driving force to form a tube. [12] One of the possible classes of materials forming asymmetric 2D monolayers are Janus sheets, like MoSSe [14] or BiTeI [15], which can be wrapped to form 1D tubes. [16,17] A recent work [16] has shown that radii well below 35 A are needed to create single-wall Janus transition metal dichalcogenide (TMD) tubes, which have significantly different (electronic) properties from the corresponding asymmetric sheet. Although facile synthesis routes for the production of single-wall inorganic nanotubes has long been actively researched, not much attention has been paid to the question of which materials would be able to make such a structure avoiding the creation of multi-wall tubes. Consequently, a high throughput study on the stability of a wide range of Janus-based nanotubes, would provide valuable information for guiding future synthesis of small-radius single-wall nanotubes.
In this work we present a comprehensive screening study in the framework of Density Functional Theory (DFT) on the stability of 135 different inorganic nanotubes generated from the rolling of asymmetric 2D Janus sheets. The calculations focus on the stability and strain energy of the chosen nanotubes. The total number of DFT relaxations performed in this work amounts to approximately 4500. We show that for pure chalcogen or halogen tubes (isovalent anions), the wrapping mechanism is mostly governed by the lattice-mismatch between the two inner and outer atomic layers, while for mixing anions (non-isovalent anions) this is dominated by the difference in valency between the X/Y elements. These findings provide a physical foundation for designing Janus nanotubes with optimal (small) radii.

Results
The nanotubes considered in this study consist of three layers (MXY, as illustrated in Figure   1) composed of different mid-layer elements (M = {Ti, Zr, Hf, V, Nb, Ta, Cr, Mo, W, Fe, Ge, Sn, As, Sb, Bi}) decorated with inner X and outer Y chalcogen and halogen atoms (X,Y = {O, S, Se, Te, Cl, Br, I}). Here we denote all group 16 elements including oxygen as chalcogens. For the three pnictogens (As, Sb and Bi) having 3+ as one of their possible oxidation states, we mix chalcogens and halogens in the structures. For the remaining 12 elements, the inner X and outer Y elements are either chalcogens or halogens. The idea of mixing chalcogens and halogens to form 2D MXY Janus sheets was recently explored in Ref. [18] , but has so far not been pursued in the context of nanotubes. We construct the tubes by rolling up 2D layers in both the T-and H-phase crystal structures, corresponding to the crystal structures found for the experimentally synthesized MoSSe [14] and BiTeI [15] 2D sheets, along both the armchair and zigzag directions.

Nanotube strain energy
Two main quantities that are needed to characterize an asymmetric nanotube are the optimal radius, which defines the most stable nanotube size, and the strain energy, which defines the energy associated with the wrapping of a 2D sheet into a nanotube (negative strain energies indicate a spontaneous wrapping).
The strain energy is defined as the difference between the energy of the nanotube and the corresponding 2D sheet. [19] In formula: , where E tube is the energy of a nanotube with N tube atoms and E M XY is the energy of the corresponding 2D Janus sheet with N M XY atoms in the unit cell. In the infinite limit R → ∞, the strain energy is zero, since the energy per atom of a tube is equal to the energy per atom of an infinite 2D Janus sheet.
It has been shown that for symmetric tubes (carbon, for example) the nanotube strain energy follows a 1/R 2 dependence. [20,21] This relationship does not hold for asymmetric tubes in which the strain energy curve exhibits a minimum. [22,23,24,25] Instead, it can be more accurately described using the equation: [22] E strain (R) = a Extrapolating the function using the obtained DFT data and evaluating the function at the minimum strain energy E strain−min leads to the optimal tube radius R opt (see also Figure   2). We note that, although Eq. 2 fits well in the region around the optimal tube radius, large strain energies can lead to a deviation. [26] We take this into consideration during the screening using Bayesian statistics (details in the SI) which helps to identify cases where the function chosen does not capture the observed data points well across all tube radii.
As an example, Figure 2 shows the strain energy as a function of the tube radius for three different materials, comparing a symmetric MoS 2 tube with the two studied asymmetric tubes NbSSe and BiSI. The symmetric MoS 2 tube shows a 1/R 2 dependence of the strain energy over the tube radius indicating the single-wall nanotube is less stable than the infinite sheet. This is not the case for the asymmetric NbSSe and BiSI tubes, where the strain curve exhibits a minimum. The strain energy curve for BiSI shows a strain energy minimum of -31 meV/atom at the optimal radius of ∼ 10 A while the strain energy curve for NbSSe is instead very shallow (minimum at -1.1 meV/atom at 85 A) due to the minor lattice-mismatch of 0.96 between the two parent sheets NbS 2 and NbSe 2 . Such a shallow strain energy curve makes it difficult to establish an optimal radius.  Figure 3 shows the optimal tube radii for all studied materials and its associated uncertainty.

Stability and optimal radius
For around 20 structures the uncertainty on the radius is estimated to be larger than 30 A (blue shaded triangles in upper right corner in Figure 3). We employ Bayesian regression (details in the SI) to automatically asses the uncertainty associated with fitting the obtained DFT data to equation 2. This makes it possible to spot data points during the screening study that might require a more detailed investigation, while the data points that fit the underlying fitting function show reduced uncertainties. We have identified three situations where the uncertainties are large: (1) The model is not able to fully capture the variation of strain energies across different tube radii. Only few structures, such as SbTeBr, do not fit the proposed model (similar to imogolite, Figure 2 in the SI), leading to large uncertainties in the prediction. (2) The strain energy curves are rather flat due to a less pronounced energy minimum making it difficult to estimate a precise optimal tube radius. This is, for

R opt
T-a Uncertainty Figure 3: Overview of the 135 investigated materials and their extrapolated optimal radius according to equation 2. It is ensured that at least three data points after filtering according to the criteria described in the SI exist (otherwise marked with a hatched box). The asterisk indicates that the difference in energy between the H-/T-phase is less than 10 meV/atom and that the smaller optimal tube radius is chosen given the extrapolated optimal radius from both prototypes. The lowercase letters a/z indicate which wrapping direction (armchair/zigzag) is preferred.
example, the case for NbSSe, TiSSe or VSSe. (3) Too few data points are available, as it happens for BiSI, BiSBr and BiSCl. By looking at the strain energy curve, we note that in these cases adding additional data points to the curve is not necessary, as the optimal radius has already been found using only three data points. However, this leads to poor statistics for the Bayesian regression. In general the values for the optimal tube radii are similar to what has been found in literature for a much smaller set of TMD Janus nanotubes.
[17] radius as shown in Figure 3. To compare the energy of the tube to its most stable 3D bulk structures, we also calculate the convex hull energy of the tube E CH at its optimal radius (E CH−min ). Except for Fe, all studied nanotubes show good stability against the decomposition into competing bulk structures (taken from the Materials Project database [27]) for at least one of the calculated combinations (here we define a combination stable when the energy of the candidate compound is within 0.2 eV/atom above the convex hull to account for a possible metastability [28,29,30]). This is in good agreement with published studies on the stability of 2D Janus monolayers. [18] Given the mid-layer element is in its preferred oxidation state, the resulting tube shows higher stability compared to the case of an unfavored oxidation state as expected. Ge, for example, generates more stable tubes when combined with two halogens (Ge 2+ ) than with two chalcogens (which would need a Ge 4+ , instead). We do not find stable nanotubes for 23 combinations (hatched boxes in Figure 4).
18 of these materials can be attributed to the transition metal in its unfavored 2+ oxidation state, when paired with two halogens. The remaining 5 materials contain either Fe, Sn or Ge which form in general less stable nanotubes for most of the studied combinations (i.e. more than 0.1 eV/atom above the convex hull). Because of their stability in the 3D form, oxygen containing tubes are in general more prone to decompose compared to the pure chalcogenide ones.
For almost 90 % of the materials, the energy difference between the armchair and zigzag wrapping direction is below 10 meV/atom, which indicates that there is only a weak driving force causing wrapping up around a specific direction. Although the armchair and zigzag wrapping directions only have a minor impact on the stability of the tube, it can be expected that the wrapping direction has a larger impact on the electronic properties due to the difference in bond distances in these tube configurations. We have recently demonstrated this for MoSTe, but a more comprehensive study focusing specifically on the electronic properties would be needed to establish general design rules. [13] The combination of the three metals As, Sb and Bi mixed together with a chalcogen sitting inside and a halogen element sitting outside of the structure generates stable and small-radius nanotubes with a rather small strain energy associated with its optimal radius.
For instance, AsSI has a minimum strain energy of -80 meV/atom which is ∼ 40 meV smaller than that of the experimentally observed imogolite nanotube when compared to computational reference data found in literature. [31] The MXY nanotubes share with imogolite, which is known to exist as a single-wall nanotube, the shape of the strain energy curve (see BiSI in SI Figure 2).  Janus sheet (inter-layer distance). A similar plot showing the optimal radius versus the convex hull energies can be found in the SI in Figure 3. Figure 5 (b) indicates that very low strain energies are not necessary to obtain small-radius nanotubes. Nanotubes with a diameter smaller than 15 A are predicted to be found in a range of strain energy minima from -15 meV/atom (BiSeBr) to -280 meV/atom (TaOTe, which has the most negative strain energy in our dataset). The optimal radius increases sharply when approaching the limit of no lattice-mismatch between the MX 2 and MY 2 parent sheets (i.e. lattice-mismatch a M X 2 /a M Y 2 = 1, see Figure 5 (a) -blue dashed curve). Additionally, the lattice-mismatch is correlated with the ratio of ionic radii when mixing two chalcogens (Janus-Cal-Cal class following trend OTe > OSe > OS > STe > SeTe > SSe) or two halogens (Janus-Hal-Hal class following ClI > BrI > ClBr) with each other.
The lattice-mismatch between the two parent sheets can give an estimate on the optimal tube radius in the case of isovalent anions, while it fails for mixing halogens with chalcogens  in the BiTeI layer is similar to the one found in bismuth tellurides. [32] To illustrate this, we take the BiSI structure in Figure 6 (Figure 5 (d)). Additionally, we observe for the non-isovalent structures a shortening in M-X bond lengths upon wrapping the infinite 2D Janus sheet into a nanotube leading to the effective thickness of the MX layer in a tube being even thinner than in the sheet.
In order to further clarify this point we compare the stable BiSI tube and the energetically unstable inverse wrapped counterpart BiIS in Figure 6 on the right. In the inverse wrapped case steric effects in between I-I elements even lead to a bond expansion of the Bi-Bi bond, d M −M , from 4.2 A to 4.6 A. Given both tube configurations differ by 0.14 eV/atom in energy in favor of the tube in which the chalcogen element sits inside, the importance of the (BiTe) layer for the stability is underlined. Plotting the stability difference in between the stable (chalcogen inside) and the inverse wrapped tube (chalcogen outisde) for nine selected non-isovalent tubes we find the points to follow a positive linear correlation with the ratio in Figure 4 in the SI). We conclude that the steric effects caused by a difference in layer thickness significantly influences the stability of the tubes. Steric effects can also explain the reduced minimum strain energy for BiXI < BiXBr < BiXCl structures, as well as the reduced optimal radius for BiSY < BiSeY < BiTeY structures. In these cases, the element forming stronger bonds determines the wrapping direction.

Optimal radius descriptors
All studied tubes are more stable if the chalcogen element is placed inside. In fact, for some of the studied non-isovalent nanotubes the optimal tube radius is independent of the halogen sitting on the outside of the tube (e.g. SbSCl, SbSBr and SbSI). Thus, only considering MX 2 and MY 2 parent structures provides an inaccurate picture when studying non-isovalent nanotubes. Additionally, the parent structures of non-isovalent tubes will each have a Melement oxidation state different from the 2D Janus sheet (e.g. BiI 2 , BiS 2 and BiSI). Two models from literature use the parent structure lattice constants as descriptors for predicting the optimal tube radius of isovalent nanotubes. The first approach is based on the plate theory as described by Timoshenko [33,34], while the second approach makes use of the Poisson ratios of the parent sheets (Poisson model) [16]. Both models require the calculation of the stiffness tensors of the parent structures. Our approach, which we refer to as Innerbond model, is instead based on solely geometrical parameters, which can be readily obtained from existing 2D databases (e.g. [35]). It uses the lattice constants of the 2D MXY Janus and MX 2 parent sheets as well as the 2D Janus t M −X layer thickness, as discussed above. It is initially assumed that the optimal tube radius is determined by the the lattice mismatch between the MXY Janus sheet and the corresponding MX 2 parent sheet. The derivation of the Inner-bond model and the used formulas for the two reference models are provided in the SI.  Figure 7: (a) Optimal tube radius using the three models based on the Timoshenko plate theory [33], considering the Poisson ratios (Poisson) [16] and using lattice mismatches between Janus MXY and MY 2 parent sheets (Inner-bond model). The optimal radius based on DFT calculations R opt is also shown which has been obtained extrapolating the function given in equation 2. The asterisk indicates that due to prototype changes when straining the SbBr 2 sheet, it was not possible to get the elastic tensor required for the Timoshenko and Poisson model. (b) The radius plotted against the lattice mismatch MX 2 /MXY.
All three models capture the rolling mechanism for isovalent anions for a selected set of materials (with discrepancies up to 30 % compared to the calculated DFT radius - Figure 7 (a)), while only the Inner-bond model predicts the optimal radii of non-isovalent nanotubes to be small. Nevertheless, it fails to capture that the optimal DFT radius for some tubes, in which a chalcogen sits inside the tube, is almost independent of the chosen halogen element sitting outside of the tube (see orange bars for SbSCl, SbSBr and SbSI in Figure 7 (a)).
Here, the layer thickness t M −X is so thin, that the steric effects of the halogen atoms at small radii do not decisively impact the stability of the tube. In this picture, the wrapping mechanism is solely governed by the steric effects of the chalcogen elements. Controversely, the 2D Janus sheet lattice constants increases for SbSCl < SbSBr < SbSI due to the size difference of the halogen atoms, leading to a mismatch in predicted values for the Inner-bond model.
One feature all the studied materials share is that the lattice mismatch between the parent MX 2 and the Janus MXY sheets is always smaller than 1 (Figure 7 (b)). This is different from the idea of comparing MX 2 and MY 2 lattice constants, for which, in the case of mixing chalcogen and halogen X/Y elements, the lattice mismatch exceeds 1 (Figure 5 (a)).
Although only giving quantitative predictions, the Inner-bond model shows that the bond distances of the elements in the inner MX 2 sheet together with the MXY lattice constant are good descriptors for predicting the optimal radius within the studied compound space.
In the case of isovalent anions the lattice-mismatch of the parent structures can also be used as a descriptor.

Strain energy descriptors
The Janus sheets show an asymmetry between the two sides of the layer. We can make a sheet where the mirror symmetry along the mid-layer (M-layer) is restored through rearranging the X/Y elements into alternating X/Y rows (as shown in Figure 8). The graph in Figure 8 shows the energy difference between the alternating and the Janus sheet (E 2D−alt -E 2D−Janus ) as a descriptor for the strain energy minumum (E strain−min ). Few outliers ( > 50 meV/atom error in prediction) correspond to structures with a low convex hull stability, i.e. E CH > 0.2 eV/atom.
For instance, this descriptor, which does not require the calculation of the tubes, could be used to indicate if a sheet exfoliated onto a host structure might undergo spontaneous curling. The strain energy would need to be larger than the adsorption energy of the sheet on the host structure, which is for example not the case for the experimentally exfoliated BiTeI sheet [15]. Whether a large strain energy is needed to form single-wall nanotubes (in favor Different synthesis procedures can be imagined leading to the synthesis of the Janus nanotubes suggested here. Besides the exfoliation of a monolayer onto a host structure and its subsequent spontaneous wrapping, which might need lithography techniques to cut the monolayer into the required dimensions, the intercalation of ions into Janus multilayered materials, combined with ball-milling, [37] can be a viable experimental synthesis path allowing to produce larger quantities. A common issue with these methods would be preventing unwanted reactions occurring at unsaturated bonds. A possible alternative synthesis pathway could be to use Atomic Layer Deposition combined with lithography to cut the 2D layer into a desired shape before it spontaneously wraps. Regardless the synthesis procedure, the intrinsic driving force of the 2D Janus sheets to form 1D tubes (or alternatively other curled shapes) should be large enough to be observed.
A promising starting material to experimentally verify our predictions is the layered BiTeI 3D bulk crystal structure for which experimental synthesis routes have been established. [32] Other interesting candidates include tubes based on Vanadium or Titanium, which, by displaying a variety of oxidation states, are suitable for electrochemical applications. [38,39] Curled multiwall VO x structures with much larger inner radii than the tubes predicted here have already been reported. [39] Other non-isovalent anion based tubes, such as BiSI, might be feasible. However, the weaker ionic halogen bonds might drive reconstructions, in contrast to less reactive covalent bonds.

Discussion
In this work, we investigated the formation mechanism of nanotubes by studying their stability and optimal nanotube radius. 135 different Janus nanotubes have been calculated using a structure prototype approach (T-/H-phase). Each three-layered material consisted of one of 15 different cation mid-layer elements in combination with inner and outer atomic layers being occupied by either chalcogens, halogens or their mixture as anions. For isovalent anions, the wrapping mechanism could be explained by the lattice-mismatch between the two inner and outer atomic layers, while for non-isovalent anions steric effects caused through short pnictogen-chalcogen inside and longer pnictogen-halogen bonds outside of the material drive the stability. These effects are beneficial to the formation of some of the smallest identified nanotubes showing optimal radii below 10 A. We noted that in general a large minimum strain energy is not needed to find tubes with an optimal radius smaller than 35 A. Additionally, the minimum strain energy was reasonably well estimated using the energy difference between a 2D Janus and alternating sheet as a descriptor. We employed Bayesian statistics to assess the quality of our fitting in order to identify uncertainties in our predictions of the optimal radius due, for example, to a misfit between the obtained data and the underlying strain energy curve function.
Nanotubes based on BiTeI and related compositions appear to be a particular interesting starting point for experimental verification due to their synthesizability in a 3D structure and exfoliability into 2D layers, as well as their predicted stability in 1D form. Following the interest in Vanadium-and Titanium-based nanotubes for electrochemical applications, we suggest that the metastable combinations of these metals paired with either OTe, OSe, or OS can be additional interesting candidates.
The findings reported here, shed light on the mechanism behind the curling of 2D Janus sheets and define a new path for the synthesis of nanotubes with small radii, for which the lattice mismatch and the bonding character of the anions play a fundamental role.

Methods
The first step to create our library of nanotubes is to relax the 2D Janus sheets, taken from the Computational 2D Materials Database (c2db). [35] If a 2D Janus sheet is not present in the database, the 2D Janus sheet is created by averaging the lattice constants from its constituent MX 2 and MY 2 parent sheets. Subsequently, tubes are generated by repeating and wrapping the 2D sheets both along the armchair and zigzag wrapping directions, thereby obtaining tubes with various radii (similar to what is shown in Figure 1). In details, the initial number of unit cell repetitions is n = (6,8,10,12,14) for the armchair and n = (10,13,16,19,22) for the zigzag wrapping direction, which correspond to tubes with a radius smaller than 20 A. We apply a set of filters to decide whether or not the relaxed structure is accepted for further investigation. These filters include assuring that a tube retained its circular shape and that no unwanted changes into different prototypes occurred during the relaxation. The filters discard ∼ 40% of the data generated. A detailed discussion on how the data is filtered prior to visualization can be found in the SI. For consistent and reproducible calculations, we implement a workflow combining the Atomic Simulation Environment (ASE) [40] and the workflow scheduling system MyQueue [41]. Inspired by the CUSTODIAN package [27], we establish an "ASE error handler" to handle common DFT errors thus limiting the need for user intervention. A similar approach has been recently implemented to autonomously discover battery electrodes. [42] All calculations are carried out with the Vienna ab initio Simulation Package (VASP) using a plane-wave basis set with an energy cutoff of 550 eV. [43,44,45] In order to approximate the exchange-correlation effects the Perdew-Burke-Ernzerhof (PBE) form generalized gradient approximation (GGA) is used. [46] A k-point density > 4.7 / A −1 is used to sample the Brillouin zone. All forces are converged to less than 0.01 eV/ A. The structures are relaxed in a non-magnetic configuration, i.e. without applying initial magnetic moments on the elements. A minimum vacuum in between repeating images of 16 A is ensured. Dipole corrections are applied along the non-periodic direction for materials with an out-of-plane dipole moment. To assess the stability versus 3D phases, we use a convex hull analysis, where the reference structures are the ones defining the convex hull in the Materials Project database. [47] These structures are then relaxed with the matching input parameters used in this work, [27] while the reference energy of oxygen is obtained by calculating the difference in energy between water and hydrogen in the gas phase, including the zero point energy (ZPE) corrections. [48] Code availability The workflow is continuously developed and may be accessed at https://gitlab.com/ivca/Nanotubes.

Data availability
The structures presented in this work are available on DTU DATA with the identifier "doi.org/10.11583/DTU.13153355".