{'fraction_non_alphanumeric': 0.03835417597141581, 'fraction_numerical': 0.015464493077266637, 'mean_word_length': 4.44468085106383, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The structure of thin boron nanowires predicted using evolutionary computations

This work describes the implementation of a genetic algorithm-based strategy combined with first-principles computations for identifying the structure of the most stable boron 1D structures. We focus our attention on the structure of ultrathin 1D boron structures given the lack of previous experimental and theoretical work on this topic. Our methodology yields reasonable structural candidates for further optimizations at the DFT level with tighter convergence criteria. The simulations involved 1D structures with up to 8 atoms per unit cell. We have identified four main groups of structures: flat nanowires (monatomic-height stripes) with triangular or triangular and"square"motifs, stripes with larger holes, nanowires with an open tubular shape, and regular nanowires. The diameter-dependent structural changes are discussed.


Introduction
Genetic algorithms (GAs) serve as a framework for solving optimization problems. One such problem is crystal structure prediction (CSP) based solely on stoichiometry. Brute force CSP is prohibitively expensive [1] and often lacks proper solution [1,2].
Boron is one of the least understood chemical elements in the periodic table [3]. It was discovered in 1808 [4] and first obtained in pure form in 1909 [5]. Ultrapure boron single crystals for the semiconductor industry are obtained from the decomposition of diborane assisted with Czochralski or zone melting purification processes [6]. Boron has two naturally occurring isotopes, 10 B and quantitatively dominant 11 B (both stable) [7], and over a dozen radioisotopes up to 21 B [8], which in this particular example decays through a two-neutron emission channel. Boron in the Universe and on the Earth comes from the spallation process [9]. Its abundance in the continental crust is at the level of 0.001%, which ranks this element at position number 37, after lead and before thorium [10]. Boron is important to life on Earth-it co-creates plant cell walls [11]. This element was discovered on Mars in Gale Crater in way suggesting that water was there in the past [12]. On Earth, on the west edge of the Mojave Desert (USA), in a census-designated place with not accidental name of Boron (California, USA) is located the greatest borax (Na 2 B 4 O 7 · 10H 2 O) mine on the planet. Borax can be easily converted into boric acid and then into boron trioxide, and finally into elemental boron.
Similarly to carbon, boron can form low-dimensional structures. In the two-dimensional (2D) case, theoretically, the most stable neutral structure is the α-sheet [17], however in electron-rich environments, e.g., during production on metallic substrates, the negatively charged δ-sheet would be preferred [18]. This structure was obtained experimentally on a silver substrate by molecular beam epitaxy (MBE) [19]. The first 2D boron structures were obtained in 2015 chemical vapor deposition (CVD) [20] and MBE [21]. 2D boron production is now standard procedure [22,23,24].
Boron can form materials with other chemical elements (alloys). One such example is boron nitride nanotube (BNNT) [40,41]. Paper made of BNNTs, contrary to ordinary paper or made of carbon nanotunes, is flame resistant [42]. BNNTs can also spontaneously reconstruct after bending [43].
Scientific and engineering breakthroughs of the next two decades, especially the increase in computing power, enable-at least partially-the realization of CSP. At the beginning of this century, Oganov et al. [1] realized a computer program called USPEX, which merges the GA approach with density functional theory (DFT) calculations. This method, despite being computationally quite expensive, is fruitful-boron allotrope γ-B 28 [44,45] and highpressure sodium [46] are two examples where numerical simulations were experimentally confirmed.
The USPEX program is not the first case of an evolutionary computation (EC) software package applied to the structure prediction task, a similar approach was previously used for the protein folding problem [47]. On the other hand, it seems that Oganov's work was the first attempt which encouraged the condensed matter physics community for this method-sometime later Zurek et al. [48] presented the XtalOpt program that is using similar GA/DFT combination, while scientists working on the CALYPSO project [49,50] implemented the particle swarm optimization (PSO) algorithm to determine the crystal structure of materials.
Boron in its 1D form remains the least studied, both theoretically and experimentally [3,51]. In this work, we find the structure of boron nanowires with diameters smaller than 1 nm. We are interested to study the structural changes of these ultrathin nanowires as a function of their diameter. To achieve this goal, the most stable forms of all-boron nanowires are found using a developed for this purpose software package called Quilë [52].

Theoretical approach
An evolutionary algorithm is a non-deterministic optimization method that in each run can generate different solutions, even if the initial structures are the same. The task to be solved in this work is to find-with the use of GA and DFT-the optimal crystal structure of monoatomic nanowires. Two groups of 1D structures will be considered. In the first group, the boron structures are confined to planarity. This kind of structures are commonly known as atomic stripes. In the second group of structures, the out of plane position of atoms is allowed and as a result of the structure optimization we can end up having buckled stripes, nanowires, and even nanotubes. The first step in our nanowire CSP's problem is to define the unit cell in terms of predicate formalism.

Unit cell definition in terms of predicate Q
The unit cell of our 1D structure is defined in terms of predicate Q. Without loss of generality one can make the following assumptions: • the 1D structures are periodic in the z direction, • the length of the unit cell in the z direction is equal to h, • atomic positions in the unit cell are described with a sequence of vectors ( where N a is the number of atoms in the unit cell, • the first atom in the unit cell is at the origin of the coordinate system, i.e., r 0 = 0, • flat structures are described using Cartesian coordinates (0, y, z), while nanowires are described using cylindrical coordinates (ρ, φ, z), where we choose φ 1 = 0. It can be easily seen that for flat structures N a ≥ 1, while for nanowires N a ≥ 3. The Q predicate should state two things (Q ⇔ Q 0 ∧ Q 1 ). Firstly, no two atoms should be too close to each other (Q 0 predicate). Secondly, for each atomic pair, the path of interatomic bonds connecting this pair should exist (Q 1 predicate). In order to fulfill periodic boundary conditions both predicates (i.e. Q 0 and Q 1 ) should be formulated based on two adjacent unit cells. These conditions demand the use of two parameters, labeled as d min and d max , describing the allowed range of interatomic distances, beyond which the bond cannot be formed because of too small or too big interatomic distance, respectively. The d min and d max parameters should be chosen for each material separately and a priori knowledge of the chemistry of the material might be helpful to constrain the permitted range of distances.
The Q 0 predicate can be formulated as follows: To formulate the Q 1 predicate, elements of graph theory will be used. Let V be the set of vertices consisting of atomic positions taken from two adjacent unit cells, i.e., V = { r i , r i + hê z | i ∈ ι N a }. Let E be the set of edges (undirected), which vertices fulfill the condition of maximum distance d max : By taking a pair of both objects one obtains graph H = (V, E). The Q 1 predicate might be now formulated as Q 1 ≡ ω(H) = 1, where ω(H) is the number of connected components of graph H, i.e., as a statement that graph H is connected. Graph connectivity can be checked with depth-or breadth-first search, or by analyzing the adjacency matrix. It should be noted that the aforementioned Q 1 predicate should be treated as a first approximation of the solution, i.e., every nanowire representation fulfilling Q 1 predicate is proper (nanowire is connected), but there are some proper nanowire representations, for which Q 1 predicate does not hold ("false negative"). This approximation of Q 1 is computationally cheaper than its exact version with checks employing translations along the z axis.

Encoding of the unit cell traits. Fitness function
The information about the phenotype, i.e., atomic positions and unit cell height, should be somehow encoded in the genotype. For this purpose a floating-point representation in the form of g = (x i ) c−1 i=0 ∈ [0, 1] c ⊂ R c is used. The genotypes for flat and non-flat cases are shown in Fig. 1, where the scheme of encoding between the phenotype and genotype is sketched.
One can notice that the genotype length, c, depends on the nanowire type: 3 The same is for the unit cell height, h, and atomic positions ( r i ) N a −1 i=0 . For flat nanowires h = d max (x 0 + N a −1 j=1 x 2 j ) and the atomic positions for i ≥ 1 in Cartesian coordinates are: For non-flat nanowires h = d max (x 0 + N a −1 j=1 x 3 j−1 ) and the atomic positions for i ≥ 1 in cylindrical coordinates are: The optimization problem is to find such unit cell height and atomic positions which minimize the total energy. This energy is obtained from self-consistent DFT calculations using the Quantum ESPRESSO (QE) software package [53,54].

The case of boron nanowires
The self-consistent DFT calculations are performed using the projector augmented wave (PAW) technique and the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional. A Monkhorst-Pack Brillouin zone sampling with a 1 × 1 × 8 k-point mesh is used with a half-step offset in the z direction. Broyden electron density mixing with mixing parameter β = 0.7 and a Methfessel-Paxton method with a smearing width of 0.01 Ry for the electron smearing are used. The kinetic energy cut-off of the plane wave basis is set to 60 Ry. The PAW pseudopotential is taken from the QE library [55]. The projection of the distance between extreme atoms of two adjacent unit cells on x and y axes (buffer thickness) is equal to 10 Å, which should be enough to minimize the effects of nanowire replicas placed periodically in the xy-plane.
For the purpose of random crystal generation [56] at the first GA iteration (first generation of genotypes), a mechanism of drawing genes' values from a uniform distribution is used. The genotypes violating the Q predicate are rejected and the drawing procedure is repeated until Q is fulfilled. Single arithmetic recombination with a randomreset mutation with parameter p = 1/c applied stochastically is used with recombination and mutation probabilities equal to p r = 1 and p m = 0.5, respectively. The generation size µ is set to 100, while the parents' multiset size is set to 2k = 64. The stochastic universal sampling (SUS) mechanism with linear ranking selection (lin-RS, s = 2) is used for parents' selection and selection to the next generation. The genetic search is finished after reaching the fitness function plateau, or-more precisely-if after 10 last generations the fitness function maximum has not improved more than ∆E with respect to the whole process. The ∆E value depends on the number of atoms in the unit cell: ∆E = ∆E b · N a , where ∆E b = 1 meV/atom is a DFT precision of binding energy per atom calculations resulting from the assumed approximations [57]. The bond lengths are constrained between d min = 1.54 Å and d max = 2.10 Å inclusively-these values arise from statistical analysis of known boron structures [58].
A complete implementation of boron nanowires CSP is available as an example program in the directory examples/evenstar/ of the Quilë repository. Numerical simulations for flat and non-flat nanowires for N a up to 8 are performed with the above-mentioned recipe. In order to avoid potential non-full effectiveness of the genetic search, calculations for each N a case, for both flat and non-flat nanowires, are performed 5 times and calculation logs are available online [59] (please see evolution/ directory).
The structures obtained from our GA search are further optimized using QE. Only one or, in case of degeneracy into the same energy, several best structures are further proceeded. For each of these structures, we have optimized the atomic positions and the lattice constants with an increase to a 1 × 1 × 16 k-point grid and an increase to a 15 Å buffer thickness. The total energy (E tot, isolated B atom ) for an isolated boron atom is obtained with calculations including spin polarization-boron has an unpaired electron on valence shell (results of calculations are available online [59] in the single_atom/ directory). The results of simulations are visualized with the use of the XCrySDen program [60]. Animated visualization of nanowires is presented online [59] in the results/ directory, where all structures are presented, both before and after relaxation.  The structure of flat and non-flat nanowires obtained by evolutionary computations followed by DFT relaxation. In both cases, the following data is shown (from top): number of atoms in the unit cell N a , evolutionary found structure, the consequently DFT-relaxed structure, the structure label (symbols δ i Na or κ i Na for i ∈ N), and the binding energy per atom, E b (in units of eV/atom), for the DFT-relaxed nanowires. Structures with labels in bold have the highest E b from all nanowires with a given number of atoms in the unit cell. For the non-flat nanowires, beside the side view (alongside the z axis), the top view is shown preserving the scale. Note that the scale between different structures is not preserved.

Results and discussion
The collection of the evolutionary found and DFT-relaxed 1D structures is shown in Fig. 2. They are divided into two groups: nanowires obtained through GA-evolution of flat δ i N a (gr. διαμέρισμα -flat) and non-flat κ i N a (gr. κυρτός -convex) structures, where i ∈ N is an index labeling structures (δ * N a or κ * N a ) having the same number of atoms in the unit cell. The DFT-relaxed structures are presented in Fig. 2 below the corresponding evolutionary found counterparts. The most energetically favorable nanowires within each group (δ * N a and κ * N a ) are highlighted in the figure using bold labels. The binding energy, E b , is calculated using the formula: Analyzing the results presented in Fig. 2, we can identify two groups of flat structures: stripes with triangular or triangular and "square" motifs, and stripes with larger holes. Looking at non-flat structures, we can find, among others, nanowires with an open tubular shape (e.g. κ 0 5 , κ 1 6 , or κ 1 8 ) and regular nanowires (e.g. κ 0 4 , κ 2 5 , or κ 2 6 ). The δ 0 2 and κ 1 4 stripes represent the same structure. Doubling the unit cell size (in the z direction) does not result in any structural changes what confirms the robustness of the BDC structure against Peierls distortions. The δ 0 3 and δ 0 4 stripes are also found in the κ * 3 (and also δ * 6 ) and κ * 4 (and also δ * 8 ) series of structures, respectively, and are shown only once in Fig. 2. The widest considered 1D structure is δ 3 7 (9 Å in diameter) and the thinnest is the boron chain δ 0 1 . From small-diameter structures, the most stable are fully planar stripes with boron triangular motifs. Among all flat structures, the largest value of E b belongs to the δ 0 3 nanowire (5.817 eV/atom) while among non-flat structures-to the κ 1 8 nanowire (5.840 eV/atom). Boron can also form nanotubes and the open shape of κ 1 8 suggests it as well. The value of E b for κ 1 8 is close to that of the buckled triangular boron sheet (see Fig. 3). 6 The values of E b for 1D structures are usually lower than those for 2D structures of the same material. It is worth noting that E b for some of the obtained boron nanowires is larger than that for 2D boron triangular flat sheet (s 6 -ft) [61]. This is shown in Fig. 3 where we plot E b for each structure from Fig. 2 as a function of N a . Figure 3 also shows three insets presenting the average number of generations across 5 computational series for each N a value, the respective average computational time of GA-evolution on a single CPU [62] machine, and normalized energy difference ∆E b /E b, final as a function of the generation number during GA-evolution resulting in the most energetically favorable structures, δ 0 3 and κ 1 8 . For both flat and non-flat structures, the number of generations needed to obtain the results increases with the number of atoms in the unit cell. On the other hand, at the logarithmic scale, the average computational time increases linearly with N a for flat structures and rises rapidly for non-flat structures.

Summary
In summary, we study the structure of ultrathin boron 1D periodic structures using a DFT-based GA approach. Based on our simulations, we have identified four main groups of structures: stripes with triangular or triangular and "square" motifs, stripes with larger holes, nanowires with an open tubular shape, and regular nanowires. The most stable thin structures are fully planar stripes with boron triangular motifs. However, thicker stripes tend to form tubular-like structures that can be viewed as precursors of nanotubes or nanowires. Finally, several of the studied 1D structures are more stable than the boron hexagonal (triangular) sheet which gives hope that such thin structures may be experimentally obtained. Although the whole methodology developed in this work is tested on boron 1D structures, we believe that it has general predictive capabilities and can be used to study e.g. nanowires made of silicon or other monoatomic or polyatomic materials.