research papers
Ab initio crystal structures and relative phase stabilities for the aleksite series, Pb_{n}Bi_{4}Te_{4}S_{n+2}
^{a}School of Chemical Engineering, The University of Adelaide, North Terrace, Adelaide, South Australia 5005, Australia, and ^{b}BHP Olympic Dam, 10 Franklin Street, Adelaide, S.A. 5000, Australia
^{*}Correspondence email: jie.yao@adelaide.edu.au
Density functional theory methods are applied to crystal structures and stabilities of phases from the aleksite homologous series, Pb_{n}Bi_{4}Te_{4}S_{n+2} (n = homologue number). The seven phases investigated correspond to n = 0 (tetradymite), 2 (aleksite21R and 42R), 4 (saddlebackite9H and 18H), 6 (unnamed Pb_{6}Bi_{4}Te_{4}S_{8}), 8 (unnamed Pb_{8}Bi_{4}Te_{4}S_{10}), 10 (hitachiite) and 12 (unnamed Pb_{12}Bi_{4}Te_{4}S_{14}). These seven phases correspond to nine single or doublemodule structures, each comprising an odd number of atom layers, 5, 7, (5.9), 9, (7.11), 11, 13, 15 and 17, expressed by the formula: S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}), where M = Pb, Bi and X = Te, S, p ≥ 2, and S and L = number of short and long modules, respectively. Relaxed structures show a and c values within 1.5% of experimental data; a and the interlayer distance d_{sub} decrease with increasing PbS content. Variable Pb—S bond lengths contrast with constant Pb—S bond lengths in galena. All phases are nfold superstructures of a rhombohedral with c/3 = d_{sub}*. Electron diffraction patterns show two brightest reflections at the centre of d_{sub}*, described by the modulation vector q_{F} = (i/N) · d_{sub}*, i = S + L. A second modulation vector, q = γ · c_{sub}*, shows a decrease in γ, from 1.8 to 1.588, across the n = 0 to n = 12 interval. The linear relationship between γ and d_{sub} allows the prediction of any theoretical phases beyond the studied compositional range. The upper PbSrich limit of the series is postulated as n = 398 (Pb_{398}Bi_{4}Te_{4}S_{400}), a phase with d_{sub} (1.726 Å) identical to that of trigonal PbS within experimental error. The aleksite series is a prime example of mixed layer compounds built with accretional homology principles.
1. Introduction
Several named minerals and a number of unnamed Bi–Pb–tellurosulfide phases were initially grouped together within a homologous series with the common formula Pb_{n}Bi_{4}Te_{4}S_{n+2}, where n is homologue number (Cook et al., 2007a). Later, the series was termed the aleksite series after the first named mineral, aleksite (Pb_{2}Bi_{4}Te_{4}S_{4}), with the generalized formula revised to Pb_{(n−1)}Bi_{2}X_{n+2} (n = homologue number, X = chalcogen) (Moëlo et al., 2008). Cook et al. (2007a) postulated the existence of a hierarchical series of Pb–Bi–tellurosulfides that can be expanded from the archetypal fiveatom tetradymite unit to larger seven, nine, 11atom units, whereas Moëlo et al. (2008) considered the tetradymite (Bi_{2}Te_{2}S) archetype as a link to layered sulfosalts. A second large group of minerals and unnamed compounds, Bi_{x}X_{y} (X = chalcogen), is also derived from the same archetype, constituting the tetradymite homologous series (Cook et al., 2007b).
Ciobanu et al. (2009) affirmed that the two homologous series derived from the tetradymite archetype share structural building principles in agreement with the formulae: (i) chalcogenrich [S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}]; p ≥ 2)] for the aleksite series, and (ii) bismuthrich [S′(Bi_{2k}X_{3})·L′[Bi_{2(k+1)}X_{3}]) for the tetradymite series. Investigation of compounds from the tetradymite series in the compositional range Bi_{2}X_{3}–Bi_{8}X_{3} using (TEM) (Ciobanu et al., 2009) confirmed that these are all nfold (N = layers in the stacking sequence) superstructures of a rhombohedral with c/3 = d_{0} ∼ 2 Å. Electron diffraction patterns show two brightest reflections in the centre of d_{0} and are described by two modulations vectors: q = (q ∼ homoatomic interval) and q_{F} = ; q_{F} = (i/N)d_{0}* = , i = S′ + L′.
The same basis for crystal structural modularity attributable to other mixed layer compounds (Amelinckx et al., 1989; Frangis et al., 1990) should extend to the aleksite series since their building modules follow the same accretional principle with the tetradymite compounds. This was demonstrated in a highangle annular darkfield scanning TEM (HAADF STEM) study of unnamed PbBi_{4}Te_{4}S_{3} (Cook et al., 2019). Based on the observed, Cook et al. (2019) showed that for a single homologue, n = 1 in this case, the structure could consist of combinations of multiple S and L_{m} modules, at constant p = 2.
Ab initio calculation of structures from the tetradymite series within the compositional interval Bi_{2}Te_{3}–Bi_{8}Te_{3} has confirmed crystal structural modularity using the accretional formalism as above (Yao et al., in the press). Moreover, the same study formulated a model combining the modulation parameter γ and d_{sub} to predict the upper (Birich) end of the tetradymite series.
Study of layered compounds in the system PbTe–Bi_{2}Te_{3} led to the definition of another homologous series based on units of fixed width: nPbTe·mBi_{2}Te_{3} (Shelimova et al., 2004) using a similar approach to the definition of nBi_{2}·mBi_{2}Te_{3} for compounds in the tetradymite series (Shelimova et al., 2000). Following the same ideas, Kuribayashi et al. (2019) discovered and named the third member of the aleksite series, hitachiite (Pb_{10}Bi_{4}Te_{4}S_{12}), and introduced the formula Bi_{2}Te_{2}S·nPbS to express homology in the series, an approach distinct from the accretional model described above.
Ab initio calculations of phases across an extended compositional range in a modular series provide an excellent tool for the exploration of modularity, crystal structures, phase stabilities, and the limits of the series. Using density functional theory (DFT) and structure simulations we study seven homologues from the aleksite series covering the compositional range Bi_{2}Te_{2}S–Pb_{12}Bi_{4}Te_{4}S_{14}. Our objectives are: (i) to describe their structures, bonding, structural–chemical modulation and phase stabilities, (ii) build a model for predicting the upper (Birich) limit of the series and (iii) discuss similarities and differences between the aleksite and tetradymite series.
2. data and selection of input files
Table 1 lists the seven phases under investigation (four minerals and three unnamed phases) and published data relating to their crystal structures. These are also shown on a diagram of (Pb/Pb+Bi) versus Te/(Te+S) (Fig. 1). They represent seven discrete homologues with evennumbered values of n (0, 2, 4, 6, 8, 10 and 12) using the formula Pb_{n}Bi_{4}Te_{4}S_{n+2} given by Cook et al. (2007a, 2019). The corresponding structures are given by the formula S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}), where p = 2, S = the 5atom layer (for simplicity termed `5layer' hereafter), and six different L modules (7, 9, 11, 13, 15 and 17layers). In addition, homologues n = 2 and n = 4, corresponding to the minerals aleksite (7layer) and saddlebackite (9layer), respectively, are each represented by simple, doublemodule [aleksite42R with the stacking sequence (5.9), and saddlebackite18H with (7.11) stacking sequence]. The aleksite54H polytype (Spiridonov, 1995), with much longer stacking sequence (77.11.77.15), was not included due to the much longer computation time required.

All phases are trigonal, but the N = N_{1} × 3 for R phases and N = N_{1} × 1 in H. The d_{sub} value is calculated from experimental data using c/N_{total.} In the unnamed phases, the c parameter is calculated assuming an interlayer distance of d_{0} ∼ 2 Å. We found that a remains constant at around 4.23 Å whereas c shows large variation depending on the and the stacking sequence of individual Nevertheless, their interlayer distances (d_{sub} = c/N_{total}) are directly comparable with one another and decrease systematically with increasing Pb and S content.
changes from (R) to (H) whenever the total number of atoms in the explicit formula is divisible by 3. The total number of layers in each structure is3. Methods
3.1. Ab initio calculations
To understand the connection between ab initio total energy calculations and structure relaxations based on density functional theory (DFT) (Hohenberg & Kohn, 1964; Kohn & Sham, 1965). We used the VASP simulation package (Kresse & Furthmüller, 1996) based on the projector augmented wave (PAW) method (Blöchl, 1994). The exchange and are treated with the generalized gradient approximation (GGA) within the Perdew, Burke and Ernzerhof (PBE) scheme (Perdew et al., 2008). The gammacentered dense k points were used to sample the (Table 2) and plane waves are expanded at cutoff energy 600 eV. van der Waals interactions (Te–Te) were included using the method of Grimme et al. (2010). Structures were relaxed with energy convergence of less than 10^{−6} eV for each ionic step and forces on each atom are within 0.02 eV Å^{−1}.
and chemistry in aleksite mixed layer compounds, we performed

In order to obtain the equation of state (EOS) for each structure, relaxations were carried out at different volumes with lattice vectors scaled from 95 to 101%. The relationship between volume and total energy was fitted using the Murnaghan (1944) equation of state:
where K_{0} and are the bulk modulus and its pressure derivative, V_{0} is the equilibrium volume and E_{0} is the reference energy. For each structure, the relaxed unitcell parameters are obtained by calculating structure relaxations at equilibrium volume.
Upon completing the structure relaxations for each phase, we calculate the formation energy (ΔE_{f}) to evaluate the relative phase stability. Applying a similar approach to that used by Woodcox et al. (2019), we establish a simple relation between ΔE_{f}, the energy for each phase (E_{phase}) and the energy of single atoms (E_{Bi}, E_{Te}, E_{Pb} and E_{S}) in equation (2):
where a, b, c and d represent the number of atoms of Bi, Te, Pb, and S, respectively, within each structure. When ΔE_{f} ≤ 0, the phase is considered potentially stable. An alternative approach to establishing the relative stability of a phase considers the energy difference to endmembers (Park et al., 2021), i.e. tetradymite (Bi_{2}Te_{2}S) and galena (PbS), using equation (3):
where , and are the total energies of each mixed phase and N_{atom} = total number of atoms.
4. Results
4.1. relaxation
We adopted the experimental a and c unitcell parameters in Table 1 as input for the total energy calculations. The same procedure used by Yao et al. (in the press) was applied to obtain the initial atomic coordinates for the N number of atoms in each The z coordinates are at equal intervals of 1/N along c, and the corresponding x, y coordinates are at , and 0 values repeating for a group of three atoms. To obtain the relaxed structures, we firstly constrain the equilibrium volume for each phase by fitting the total energy volume curves (Fig. 2) using the Murnaghan equation of state with EOS parameters tabulated in Table 3. The calculated V_{0} values agree with published data (Table 3) within 3.6% for all available structures.
The final structure parameters are obtained from the DFT calculations at the V_{0} values for all phases. The relaxed a and c unitcell parameters are within 1.5% difference with the published data (Tables 1 and 4). Comparison with experimental data shows a slight overestimation in the a parameter [Fig. 3(a)] and a good fit for d_{sub} values [Fig. 3(b)]. Both a and d_{sub} parameters show a smooth decreasing trend with increasing PbS across the compositional interval investigated. Notably, the doublemodule of aleksite and saddlebackite yield values for a and d_{sub} that are very similar to those of their respective singlemodule (Fig. 3).

4.2. models
The and 5). We note that the building modules are centred onto a slab of S–Pb–S…Pb–S flanked on each side by Bi–Te atoms. The increment of the central slab can be expressed as: Pb_{k}S_{k+1} (k = 0–6) for the six homologues discussed here. The modules are always separated by Te–Te layers (van der Waals gaps). A trigonal PbS structure (PbS_{R}) obtained by transformation from cubic galena (Noda et al., 1987) is included for comparison (Fig. 4). This shows the atomic arrangement in PbS_{R} is very similar to the central slab in aleksite structures when viewed on the zone axes.
models obtained using the relaxed unitcell parameters are plotted on the to illustrate the incremental increase in width of each structure with addition of Pb and S atoms (Figs. 4The models show that the tetradymite unit is no longer preserved as such within the singlemodule structures, although these are required to form all homologues with n < 2 such as the unnamed n = 1 phase with composition PbBi_{4}Te_{4}S_{3} (Cook et al., 2019). Additionally, 5atom tetradymite modules are constituent building blocks in the doublemodule aleksite polytype, aleksite42R (n = 2) considered in this contribution, which contains a (5.9) sequence (Fig. 4). Fiveatom tetradymite modules may potentially exist in other configurations with n > 2 (e.g. 5.13 saddlebackite), although these will not be considered here. Alongside singlemodule 9H saddlebackite (n = 4), we do, however, consider the 18H doublemodule polytype structure of saddlebackite, which features a (7.11) sequence. We note that the simplest representing the n = 2k+2 (k = integer) building modules are composed of two singlemodule units with n = 2k and n = 2k + 4.
Slight e.g. the shift between the 5 and 9atom layer stacks in the aleksite42R polytype. Variations in atom arrangements within the singlemodule structures are, in contrast, negligible.
is observed within structures containing two different modules,4.3. Bond analysis
The bond types and their variation in length across the studied phases are shown within the and 7). The three types of bonds in tetradymite (Bi—Te, Bi—S and Te—Te) are complemented by Pb—S bonds in all other compounds from the aleksite series. The Bi—Te bond decreases in length from tetradymite (3.047 Å) to the 17atom layer structure (3.036 Å) whereas the length of the Bi—S bond remains relatively constant at ∼3.019 Å. Likewise, Te—Te bond lengths increase from 3.882 Å in tetradymite to a maximum of 3.919 Å in the 11atom layer. The Te—Te bonds in the two aleksite are constant and close to those in tetradymite (3.884 Å and 3.886 Å for the 5.9 and 7.11atom layer sequences, respectively). Although the average Pb—S bond length is nearly constant ∼2.990 Å, there is a small variation within the middle Pb_{k}S_{k+1} slab, e.g. from 2.986 to 2.995 Å in the 17atom layer. The Pb—S bond lengths within the aleksite and saddlebackite doublemodule are nearly identical to those in the corresponding singlemodule polytypes.
cell for all structures (Figs. 6In Figs. 6 and 7, bond lengths are projected onto the c axis to calculate the contribution towards the d_{sub} value in each structure. Te—Te bond projections have the highest values on the c axis whereas Bi–Te and Bi–S projections are only slightly larger than Pb—S bond projections. In all singlemodule structures, there is one Te—Te bond, two Bi—Te bonds and two Bi—S bonds, whereas the number of Pb—S bonds increases from 0 in tetradymite to 12 in the 17atom layer, with an incremental step of 2. We thus divide the bond types into two groups: variable number (Pb—S) and fixed number (Bi—Te, Bi—S and Te—Te).
We have calculated the bond contribution to the d_{sub} parameter from cumulative projection values and their abundance across the compositional range studied (Fig. 8). This plot shows two opposing trends, an increase in Pb—S contribution and decrease in contribution from other bonds from tetradymite to the 17atom layer, the two lines intersecting at the 11atom layer. The increase in the contribution to d_{sub} from the Pb—S bonds is however more moderate than the contribution decreases from the other bonds, resulting in a modest decrease of d_{sub} with increased PbS concentration [Fig. 3(b)].
5.4. STEM simulation and electron diffraction
In Figs. 9 and 10 we show the relaxed structures in STEM simulations and electron diffraction (ED) patterns on . The signal intensity (I) in (HAADF) STEM imaging is proportional to Z^{2} of an element along an atomic column (Pennycook & Jesson, 1990; Ishizuka, 2002). We note that S (Z = 16) is not displayed on the simulations when the neighbouring atom columns are Bi (Z = 83) and Pb (Z = 82). A better visualization of the number of atoms in each structure can be assembled using Se (Z = 34) instead of S (inset, top right in Fig. 9 and overlays on each STEM simulation in Fig. 10). Such simulations agree very well with the atomic arrangement models displayed in Figs. 4 and 5 and are concordant with HAADF STEM images of phases from the aleksite series (Cook et al., 2019).
 Figure 9 Unitcell parameters 
The interval (cropped from the ED patterns) is essential for constraining structural modulation in terms of the increase in module width and module combinations. Within each interval, there are N_{1} − 1 number of reflections equally distributed. The asymmetric unitcell length ( d_{N1}) is correlated with layer stacks for each structure. This can be also indicated as the smallest interval, d_{N1}^{*} between two neighbour reflections along .
Typical of all phases in the series is the fact that the ED patterns show the two brightest reflections at the centre of . This interval, underpinned by the modulation vector q_{F} = , where γ_{F} = i/N_{1}·d_{sub} ∼ 1/N_{1}; i = S + L. γ_{F} values are within the range 0.2–0.059 for the analysed structures and this shows a monotonic decrease with increase in PbS concentration. In cases where there are multiple although the number of divisions is doubled relative to their single unit structure, e.g. 14 and 7 divisions for the 42R and 21R aleksite the q_{F} vector remains unchanged. Nonetheless, the q_{F} interval is split into two by a satellite reflection of lesser intensity (Fig. 9). On the other hand, the between chalcogen (S, Te, Se) and Pb and Bi atoms is underpinned by a second vector: q = (Lind & Lidin, 2003). The q modulation is depicted up to thirdorder reflections along c* (ED patterns in Figs. 9 and 10). Values of γ (1.8–1.588 for the 5 to 17atom layer module range) are calculated as 3[(N_{1} + 1)/2]/N_{1} for single modules. In the doublemodule γ = 3[(N_{1} + 2)/2]/N_{1} giving the same values of γ as the corresponding singlemodule structures (Fig. 9).
This formalism is in agreement with the crystal S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}), X = chalcogen, where S and L are the number of shorter and longer modules (Cook et al., 2019) but not the formula nPbTe·mBi_{2}Te_{3} of Shelimova et al. (2004). For example, aleksite21R would have n = m = 1, requiring two distinct modules instead of only one. The 42R polytype will have 4 modules (n = m = 2) instead of the `5' and `9' modules considered here. Such a strong correlation between electron diffraction patterns and chemical modules in a homologous series is typical for mixed layer compounds (Amelinckx et al., 1989).
5. Discussion
5.1. Phase stability and energy mixing
Formation energies [equation (2)] for the studied phases in the interval Bi_{2}Te_{2}S–Pb_{12}Bi_{4}Te_{4}S_{14} (n = 12) are given in Table 5. Calculation of the formation energy and energy of mixing requires the DFT reference energies (E_{0}) of all elements (Bi, Pb, Te and S) and endmembers (Bi_{2}Te_{2}S and PbS). The reference energies for Pb, S, Bi_{2}Te_{2}S (Table 3), and PbS are calculated from equation of state fitting [equation (1)] in this study, those for the elements Bi and Te are adopted from Yao et al. (in the press). Their parameters are summarized in Table 6. All reference energies are calculated based on the GGA functional.


The calculated ΔE_{f} values are negative for all nine phases and decrease as the PbS component increases, implying they are relative stable to the endmembers. The larger doublemodule of both aleksite (42R) and saddlebackite (18H) show the same formation energy as their corresponding singlemodule units (21R and 9H, respectively), implying they are equally stable.
Phase stability can also be evaluated from the energy of mixing (E_{mixing}), which is calculated using values of the two endmember phases, tetradymite and galena [equation (3)]; Table 5). This defines a convex hull between tetradymite and galena (PbS) with aleksite at the lowest energy point (Fig. 11). The other five studied homologues plot along or slightly below the branch between aleksite and galena. Such a distribution indicates that all studied phases can be relative stable compared with the endmembers and thus do not readily decompose into tetradymite and galena endmembers. However, whether the studied phases are thermodynamically stable may require further phonon calculations to investigate the thermal effects and contributions (e.g. Belmonte et al., 2014).
Instead of adopting the formula nPbS·mBi_{2}Te_{2}S as a working model, the energy of mixing can also be defined using the accretional model:
where S = 5atom layer, L_{1–3} represent longer 7, 9 and 11modules; M = Bi, Pb, and X = Te, S. The energy of mixing for aleksite42R (5.9) and saddlebackite18H (7.11) are found at 0 and 0.2 meV per atom, respectively. This shows ideal mixing when using the accretional model and indicates that the derived and, indeed, other multiplemodule structures in the series can be stable relative to their singlemodule components. Further calculations may, however, be required to fully validate these findings.
5.2. The γ–d_{sub} relationship: a model for the extension of the aleksite series
Preliminary work shows that homologues of the aleksite series with still greater PbS content (n = 18 and n = 30, representing 23 and 35atom layers, respectively) are present in assemblages buffered by galena (Cook et al., 2021 and unpublished data). Theoretical phases from the PbSrich end of the series, such as 403, 205, 71 and 51atom layers (corresponding to homologues with n = 398, 200, 66 and 46), can also be considered based on their chemistry, which is close to, but distinct from, PbS (Fig. 1).
Our model describes a quasilinear relationship between γ and d_{sub} (Fig. 12), which allows the prediction of d_{sub} for any phase across the 17 to 403atom layer structure range (γ = 1.588–1.504), with d_{sub} values over this interval lying in the range 1.806 to 1.726 Å. The theoretical 403atom layer phase, Pb_{398}Bi_{4}Te_{4}S_{400}, with Pb/(Pb + Bi) = 0.99 shows identical d_{sub} values as our DFTmodelled predictions for PbS_{T}, which is also within 0.8% difference of that for trigonal PbS transformed from the experimental cubic structure (Noda et al., 1987). As a result, our model is suitable to approximate d_{sub} values for aleksite series homologues across the entire compositional range from tetradymite to the PbS_{T} endmember.
5.3. Modularity and comparison with the tetradymite series
Noting the possibility of multiple et al., 2019), we introduce a modified formula:
for many, if not all, homologues in the aleksite series (Cookwhere S represents the number of 5atom layers, and L_{1}, L_{2},…L_{m} are the numbers of longer, 7, 9,…2m + 5 modules; m > 0, integer; and S, L ≥ 0. This formula is useful for expressing the range of for each homologue within the series. Therefore, applying formula (6) to Pb_{n}Bi_{4}Te_{4}S_{n+2} from Cook et al. (2007a, 2019), we can calculate the homologue number (n) by relating the total number of cations and chalcogens within the component modules:
n + 4 = 2S + 3L_{1} +…(m + 2)L_{m}, leading to n = 2S + 3L_{1} +…(m + 2)L_{m} − 4 for the number of cations, and n + 6 = 3S + 4L_{1} +…(m + 3)L_{m}, leading to n = 3S + 4L_{1} +…(m + 3)L_{m} − 6 for the number of chalcogens.
The theory of mixed layer compounds stipulates that structures built by modules which are distinct in size and chemical composition are related to one another by characteristics of electron diffraction patterns thus underpinning the modularity within a homologous series (Amelinckx et al., 1989). Both the aleksite and tetradymite series are formed by modular structures derived from the same 5atom archetype but with distinct compositional ranges, i.e. extending towards PbS (aleksite series) and Bi endmembers (tetradymite series).
The individual building modules in each series are composed of an uneven number of atoms, 7, 9, 11…2 k+ 1, but with different topology between cations (Bi, Pb) and chalcogens (Te, S, Se), i.e., symmetrical in the aleksite series and asymmetrical in the tetradymite series. Despite this, the electron diffractions of relaxed structures from the aleksite series (Figs. 9 and 10) show identical modulation vectors as corresponding phases in the tetradymite series with the same number of layers and/or building modules (Ciobanu et al., 2009; Yao et al., in the press). Such characteristics provide a strong link between the two series and prove their affiliation to a single class of mixed layer compounds built by the same accretional homology principles. The alternative homology proposed for the two series involving units of the same size, 2 and 5atom layers (Shelimova et al., 2000, 2004; Kuribayashi et al., 2019) is not supported by the crystal structures, even though it may be conceptually useful to depict chemical variation within each of the two series.
Bond analysis shows marked differences between the two series, whereby the longer Te—Te bonds are present in all homologues of the aleksite series and may be responsible for the extensive et al., in the press). Construction of incremental symmetrical modules by addition of Pb–S in the aleksite series and by asymmetrical modules involving Bi–Bi pairs in the tetradymite series leads to linear versus nonlinear features in the respective γ–d_{sub} relationships. As a result, for the same γ interval (1.8–1.5) the range of d_{sub} is greater for the aleksite series compared with the tetradymite series, i.e. ∼2 to 1.726 Å, and ∼2 to 1.973 Å, respectively (Fig. 12).
In contrast, the Te—Te bonds are only present in Terich members of the tetradymite series (Yao6. Conclusions and implications
The crystal structures and stabilities of phases from the aleksite homologous series, Pb_{n}Bi_{4}Te_{4}S_{n+2}, where n = homologue number (Cook et al., 2019), were calculated using DFT methods. The study addressed four named minerals (tetradymite, aleksite, saddlebackite and hitachiite) and three compounds yet to be described in natural specimens (Pb_{6}Bi_{4}Te_{4}S_{8}, Pb_{8}Bi_{4}Te_{4}S_{10} and Pb_{12}Bi_{4}Te_{4}S_{14}). The seven phases represent homologues where n = 0, 2, 4, 6, 8, 10 and 12. Each homologue corresponds to a singlemodule type with an uneven number of atoms (5, 7, 9, 11, 13, 15 and 17, respectively), expressed by the formula: S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}), where M = Pb, Bi, and X = Te, S), p ≥ 2, S = fiveatom layer, and L1–6 = 7, 9, 11, 13, 15 and 17atom layers. The n = 2 and n = 4 homologues are also represented by twolayer (aleksite42R and saddlebackite18H), which have structures comprising two differently sized modules, (5.9) and (7.11), respectively. Other multilayer are predicted to exist for phases across the series.
The relaxed structures show the unitcell parameters a and c within 1.5% of available experimental data. Both a and the interlayer distance d_{sub} show decrease with increasing PbS component in the relaxed structures. models and STEM simulations show that the six single modules (for structures with n > 0) are centred onto a Pb_{k}S_{k+1} slab (k = 1–6), with S–Pb–S…Pb–S arrangement flanked by Bi–Te atoms. We show variable Pb—S bond lengths in the aleksite homologues, representing an important structural difference compared to the constant Pb—S bond lengths in galena.
Electron diffraction patterns show N_{1} intervals of equal length along demonstrating that all phases are nfold superstructures of a rhombohedral with c/3 = . The modulation vector q = shows a decrease in γ, from 1.8 to 1.588, with increasing PbS component across the compositional range studied (n = 0 to 12). The ED patterns have two brightest reflections at the centre of , which are described by the modulation vector q_{F} = (γ_{F} = 0.2–0.059). The number of divisions within this central interval corresponds to the number of modules, i.e. 1 for single, and 2 for double modules. This result proves that the homologous structures can be described by the formula S(M_{p}X_{p+1})·L(M_{p+1}X_{p+2}), and not the formula nPbS·mBi_{2}Te_{2}S, involving 2 and 5atom building units (Shelimova et al., 2000).
The DFT method is also used to obtain the formation energies and energy of mixing for all seven compositions. The seven singlemodule structures and the two doublemodule
show negative formation energies, implying they can be relative stable to their endmembers.We established a linear γ and d_{sub} model which allows the calculation of d_{sub} for any phase beyond the compositional range studied, e.g. phases with n values of 46, 66, 200 and 398. The model predicted a d_{sub} value of 1.726 Å for the phase Pb_{398}Bi_{4}Te_{4}S_{400} (n = 398). This can be considered as the upper end of the series, as this is the same value obtained for PbS_{T} in DFT calculations, and lies within 0.8% of experimental data.
The aleksite and tetradymite series represent excellent examples of mixed layer compounds built by accretional homology principles derived from a shared 5atom layer archetype. This study illustrates how DFT calculations can not only support predictive models for crystal and chemical modularity, but also represent a tool to expand and ultimately constrain the limits of modular series. Potential applications exist to model other, chemically different, mixed layer structures.
Supporting information
Bi_{2}Te_{2}S  c = 30.02341 Å 
M_{r} = 705.23  V = 477.50 Å^{3} 
Trigonal, R3m  Z = 3 
a = 4.28540 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.10902  
Te1  0.00000  0.00000  0.71648  
S1  0.00000  0.00000  0.50000 
Pb_{2}Bi_{4}Te_{4}S_{4}  c = 40.41287 Å 
M_{r} = 1888.98  V = 638.16 Å^{3} 
Trigonal, R3m  Z = 3 
a = 4.27010 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.58139  
Pb1  0.00000  0.00000  0.00000  
Te1  0.00000  0.00000  0.87066  
S1  0.00000  0.00000  0.29134 
Pb_{2}Bi_{4}Te_{4}S_{4}  c = 80.73010 Å 
M_{r} = 1888.98  V = 1276.06 Å^{3} 
Trigonal, R3m  Z = 6 
a = 4.27221 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  −0.06385  
Bi2  0.00000  0.00000  0.14517  
Pb1  0.00000  0.00000  0.64551  
Te1  0.00000  0.00000  0.58074  
Te2  0.00000  0.00000  0.78975  
S1  0.00000  0.00000  0.29120  
S2  0.00000  0.00000  0.00000  
S3  0.00000  0.00000  0.50000 
Pb_{4}Bi_{4}Te_{4}S_{6}  c = 16.94310 Å 
M_{r} = 2367.51  V = 266.02 Å^{3} 
Trigonal, P3m1  Z = 1 
a = 4.25788 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.19477  
Pb1  0.33333  0.66667  0.60123  
Te1  0.33333  0.66667  −0.08940  
S1  0.33333  0.66667  0.29817  
S2  0.00000  0.00000  0.50000 
Pb_{4}Bi_{4}Te_{4}S_{6}  c = 33.81797 Å 
M_{r} = 2367.51  V = 531.62 Å^{3} 
Trigonal, P3m1  Z = 2 
a = 4.26052 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.33333  0.66667  0.20348  
Bi2  0.33333  0.66667  0.39789  
Pb1  0.33333  0.66667  0.89865  
Pb2  0.00000  0.00000  0.00000  
Pb3  0.00000  0.00000  0.50000  
Te1  0.33333  0.66667  0.74380  
Te2  0.00000  0.00000  0.34516  
S1  0.00000  0.00000  0.84818  
S2  0.33333  0.66667  0.05060  
S3  0.33333  0.66667  0.55025 
Pb_{6}Bi_{4}Te_{4}S_{8}  c = 61.24237 Å 
M_{r} = 2846.04  V = 958.12 Å^{3} 
Trigonal, R3m  Z = 3 
a = 4.25029 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.05416  
Pb1  0.00000  0.00000  0.77739  
Pb2  0.00000  0.00000  0.50000  
Te1  0.00000  0.00000  0.69161  
S1  0.00000  0.00000  0.41621  
S2  0.00000  0.00000  0.13877 
Pb_{8}Bi_{4}Te_{4}S_{10}  c = 71.45943 Å 
M_{r} = 3324.57  V = 1115.07 Å^{3} 
Trigonal, R3m  Z = 3 
a = 4.24478 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.54632  
Pb1  0.00000  0.00000  −0.07188  
Pb2  0.00000  0.00000  0.30940  
Te1  0.00000  0.00000  0.85456  
S1  0.00000  0.00000  0.23765  
S2  0.00000  0.00000  0.61879  
S3  0.00000  0.00000  0.00000 
Pb_{10}Bi_{4}Te_{4}S_{12}  c = 27.31123 Å 
M_{r} = 3803.10  V = 425.78 Å^{3} 
Trigonal, P3m1  Z = 1 
a = 4.24282 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.12157  
Pb1  0.33333  0.66667  0.75138  
Pb2  0.33333  0.66667  0.37407  
Pb3  0.00000  0.00000  0.50000  
Te1  0.33333  0.66667  −0.05593  
S1  0.33333  0.66667  0.18609  
S2  0.00000  0.00000  0.31150  
S3  0.33333  0.66667  0.56288 
Pb_{12}Bi_{4}Te_{4}S_{14}  c = 92.27622 Å 
M_{r} = 4281.63  V = 1437.84 Å^{3} 
Trigonal, R3m  Z = 3 
a = 4.24175 Å  T = 0 K 
x  y  z  B_{iso}*/B_{eq}  
Bi1  0.00000  0.00000  0.03594  
Pb1  0.00000  0.00000  0.74023  
Pb2  0.00000  0.00000  0.44421  
Pb3  0.00000  0.00000  0.14810  
Te1  0.00000  0.00000  0.68316  
S1  0.00000  0.00000  0.38839  
S2  0.00000  0.00000  0.09224  
S3  0.00000  0.00000  0.79617  
S4  0.00000  0.00000  0.50000 
Acknowledgements
We acknowledge access to the Phoenix highperformance computer at the University of Adelaide. Fabien Voisin and Mark Innes are thanked for their assistance with VASP installation and HPC configuration. We thank two anonymous reviewers for the insightful comments and editorial handling by Michal Dušek.
Funding information
This is a contribution to the Australian Research Council Linkage Project LP200100156 (Critical Metals from Complex Copper Ores) cosupported by BHP Olympic Dam.
References
Adenis, C., Langer, V. & Lindqvist, O. (1989). Acta Cryst. C45, 941–942. CrossRef ICSD CAS Web of Science IUCr Journals Google Scholar
Amelinckx, S., Van Tendeloo, G., Van Dyck, D. & Van Landuyt, J. (1989). Phase Transit. 16, 3–40. CrossRef Web of Science Google Scholar
Belmonte, D., Ottonello, G. & Zuccolini, M. V. (2014). Am. Mineral. 99, 1449–1461. Web of Science CrossRef Google Scholar
Blöchl, P. E. (1994). Phys. Rev. B, 50, 17953–17979. CrossRef Web of Science Google Scholar
Ciobanu, C. L., Pring, A., Cook, N. J., Self, P., Jefferson, D., Dima, G. I. & Melnikov, V. (2009). Am. Mineral. 94, 517–534. Web of Science CrossRef CAS Google Scholar
Clarke, R. M. (1997). Aust. J. Mineral. 3, 119–124. Google Scholar
Cook, N. J., Ciobanu, C. L., Liu, W., Slattery, A., Wade, B. P., Mills, S. J. & Stanley, C. J. (2019). Minerals, 9, 628. Web of Science CrossRef Google Scholar
Cook, N. J., Ciobanu, C. L., Slattery, A., Wade, B. P., Ehrig, K. (2021). Proc. 3rd Eur. Mineral. Conf. 29 August–2 September 2021, Krakow, Poland. Abstracts, p. 100. Google Scholar
Cook, N. J., Ciobanu, C. L., Stanley, C. J., Paar, W. H. & Sundblad, K. (2007a). Can. Mineral. 45, 417–435. Web of Science CrossRef CAS Google Scholar
Cook, N. J., Ciobanu, C. L., Wagner, T. & Stanley, C. J. (2007b). Can. Mineral. 45, 665–708. Web of Science CrossRef CAS Google Scholar
Frangis, N., Kuypers, S., Manolikas, C., Van Tendeloo, G., Van Landuyt, J. & Amelinckx, S. (1990). J. Solid State Chem. 84, 314–334. CrossRef CAS Web of Science Google Scholar
Grimme, S., Antony, J., Ehrlich, S. & Krieg, H. (2010). J. Chem. Phys. 132, 154104. Web of Science CrossRef PubMed Google Scholar
Hohenberg, P. & Kohn, W. (1964). Phys. Rev. 136, B864–B871. CrossRef Web of Science Google Scholar
Ishizuka, K. (2002). Ultramicroscopy, 90, 71–83. Web of Science CrossRef CAS Google Scholar
Kresse, G. & Furthmüller, J. (1996). Phys. Rev. 54, B11169. CrossRef Web of Science Google Scholar
Kohn, W. & Sham, L. J. (1965). Phys. Rev. 140, A1133–A1138. CrossRef Web of Science Google Scholar
Kuribayashi, T., Nagase, T., Nozaki, T., Ishibashi, J., Shimada, K., Shimizu, M. & Momma, K. (2019). Mineral. Mag. 83, 733–739. Web of Science CrossRef ICSD CAS Google Scholar
Lind, H. & Lidin, S. (2003). Solid State Sci. 5, 47–57. Web of Science CrossRef CAS Google Scholar
Lipovetskiy, A. G., Borodayev, Y. S. & Zav'yaiov, Y. N. (1979). Int. Geol. Rev. 21, 1223–1228. CrossRef Google Scholar
Littlewood, P. B. (1980). J. Phys. C. Solid State Phys. 13, 4855–4873. CrossRef CAS Web of Science Google Scholar
Liu, H. & Chang, L. L. Y. (1994). Am. Mineral. 79, 1159–1166. CAS Google Scholar
Moëlo, Y., Makovicky, E., Mozgova, N. N., Jambor, J. L., Cook, N., Pring, A., Paar, W., Nickel, E. H., Graeser, S., KarupMøller, S., BalicŽunic, T., Mumme, W. G., Vurro, F. & Topa, D. (2008). Eur. J. Mineral. 20, 7–62. Google Scholar
Murnaghan, F. D. (1944). Proc. Natl Acad. Sci. USA, 30, 244–247. CrossRef PubMed CAS Google Scholar
Noda, Y., Masumoto, K., Ohba, S., Saito, Y., Toriumi, K., Iwata, Y. & Shibuya, I. (1987). Acta Cryst. C43, 1443–1445. CrossRef ICSD CAS Web of Science IUCr Journals Google Scholar
Padaki, V. C., Lakshmikumar, S. T., Subramanyam, S. V. & Gopal, E. S. R. (1981). Pramana J. Phys. 17, 25–32. CrossRef CAS Web of Science Google Scholar
Park, S., Ryu, B. & Park, S. (2021). Appl. Sci. 11, 11341. Web of Science CrossRef Google Scholar
Pauling, L. (1975). Am. Mineral. 60, 994–997. CAS Google Scholar
Pennycook, S. J. & Jesson, D. E. (1990). Phys. Rev. Lett. 64, 938–941. CrossRef PubMed CAS Web of Science Google Scholar
Perdew, J. P., Ruzsinszky, A., Csonka, G. I., Vydrov, O. A., Scuseria, G. E., Constantin, L. A., Zhou, X. & Burke, K. (2008). Phys. Rev. Lett. 100, 136406. Web of Science CrossRef PubMed Google Scholar
Rettig, S. J. & Trotter, J. (1987). Acta Cryst. C43, 2260–2262. CSD CrossRef CAS Web of Science IUCr Journals Google Scholar
Schiferl, D. & Barrett, C. S. (1969). J. Appl. Cryst. 2, 30–36. CrossRef ICSD CAS IUCr Journals Web of Science Google Scholar
Shelimova, L. E., Karpinskii, O. G., Kosyakov, V. I., Shestakov, V. A., Zemskov, V. S. & Kuznetsov, F. A. (2000). J. Struct. Chem. 41, 81–87. Web of Science CrossRef CAS Google Scholar
Shelimova, L. E., Karpinskii, O. G., Svechnikova, T. E., Avilov, E. S., Kretova, M. A. & Zemskov, V. S. (2004). Inorg. Mater. 40, 1264–1270. Web of Science CrossRef CAS Google Scholar
Spiridonov, E. M. (1995). Zap. Vses. Mineral. Ova, 124, 24–39. CAS Google Scholar
Woodcox, M., Young, J. & Smeu, M. (2019). Phys. Rev. B, 100, 104105. Web of Science CrossRef Google Scholar
Wyckoff, R. W. (1963). Crystal Structures, Vol. 1, p. 312. New York: Interscience. Google Scholar
Yao, J., Ciobanu, C. L., Cook, N. J., Ehrig, K., Dima, G., SteinleNeumann, G. (2023). Am. Mineral. https://doi.org/10.2138/am20239018. Google Scholar
This is an openaccess article distributed under the terms of the Creative Commons Attribution (CCBY) Licence, which permits unrestricted use, distribution, and reproduction in any medium, provided the original authors and source are cited.