Synthesis of BeN4 polymorphs leading to a twodimensional material with anisotropic Dirac fermions Maxim Bykov,*,1,2 Timofey Fedotenko,3 Stella Chariton,4 Dominique Laniel,3 Konstantin Glazyrin,5 Michael Hanfland,6 Jesse S. Smith,7 Vitali B. Prakapenka,4 Mohammad F. Mahmood,2 Alexander F. Goncharov,1 Alena V. Ponomareva,8 Ferenc Tasnádi,9 Alexei I. Abrikosov,10 Talha Bin Masood,10 Ingrid Hotz,10 Alexander N. Rudenko,11,12,13 Mikhail I. Katsnelson,12,13 Natalia Dubrovinskaia,3,9 Leonid Dubrovinsky,*,14 Igor A. Abrikosov*,9 1 The Earth and Planets Laboratory, Carnegie Institution for Science, Washington, DC 20015, USA College of Arts and Science, Howard University, Washington, DC 20059, USA 3 Material Physics and Technology at Extreme Conditions, Laboratory of Crystallography, University of Bayreuth, 95440 Bayreuth, Germany 4 Center for Advanced Radiation Sources, University of Chicago, 60637 Chicago, Illinois, USA 5 Photon Sciences, Deutsches Electronen Synchrotron (DESY), D-22607 Hamburg, Germany 6 European Synchrotron Radiation Facility, 38043 Grenoble Cedex 9, France 7 HPCAT, X-ray Science Division, Argonne National Laboratory, Argonne, IL 60439, USA 8 Materials Modeling and Development Laboratory, National University of Science and Technology ‘MISIS’, 119049 Moscow, Russia 9 Department of Physics, Chemistry and Biology (IFM), Linköping University, SE-58183 Linköping, Sweden 10 Department of Science and Technology (ITN), Linköping University, Norrköping, Sweden 11 Key Laboratory of Artificial Micro- and Nano-Structures of Ministry of Education and School of Physics and Technology, Wuhan University, Wuhan 430072, China 12 Radboud University, Institute for Molecules and Materials, 6525AJ Nijmegen, The Netherlands 13 Department of Theoretical Physics and Applied Mathematics, Ural Federal University, 620002 Ekaterinburg, Russia 14 Bayerisches Geoinstitut, University of Bayreuth, 95440 Bayreuth, Germany 2 maks.byk@gmail.com, Leonid *Correspondence to: Maxim Bykov leonid.dubrovinsky@uni-bayreuth.de, Igor A. Abrikosov igor.abrikosov@liu.se Dubrovinsky Abstract: High pressure chemistry is known to inspire the creation of unexpected new classes of compounds with exceptional properties. However, the application of pressure is counterintuitive for the synthesis of layered van der Waals bonded compounds, and consequently, in a search for novel twodimensional materials. Here we report the synthesis at ~90 GPa of novel beryllium polynitrides, monoclinic and triclinic BeN4. The triclinic phase, upon decompression to ambient conditions, transforms into a compound with atomic-thick BeN4 layers consisting of polyacetylene-like nitrogen chains with conjugated π-systems and Be atoms in square-planar coordination. These layers are interconnected via weak van der Waals bonds which makes their exfoliation possible. Theoretical calculations for a single BeN4 layer show the presence of Dirac points in the electronic band structure. The BeN4 layer, i.e. beryllonitrene, represents a new class of 2D materials that can be built of a metal atom and polymeric nitrogen chains, leading to a high degree of electronic structure anisotropy. MAIN TEXT Introduction The use of non-conventional methods of materials synthesis or established synthesis techniques in a non-trivial way may lead to breakthrough discoveries that revolutionize many research fields. One of the most striking examples is graphene. It was considered as an interesting two-dimensional (2D) model, though predominantly theoretical, for about 50 years, before high quality monocrystalline graphitic films were derived by Novoselov et al. by mechanical exfoliation of small mesas of highly oriented pyrolytic graphite using adhesive tape [1]. Besides multiple exciting effects of interest for condensed matter physics and intriguing technological opportunities, e.g. in electronics and optoelectronics, the possibility to study phenomena of the quantum field theory in a condensed-matter experiment was quickly recognized [2,3], suggesting that one could mimic quantum relativistic phenomena, even those which could not be observed in high-energy physics [4]. Graphene is a material with gapless Dirac cones near the Fermi level with a linear energy-momentum dispersion near the Dirac points, and exhibiting ultrahigh carrier mobility. While the Dirac cones in graphene and many other 2D compounds are isotropic, their anisotropy could result in anisotropic carrier mobility, making it possible to realize direction-dependent quantum devices. The desired anisotropy can be achieved using mechanical stress, external periodic potentials or heterostructures [5–9]. The simplest design route of novel 2D materials is based on the substitution of elements in a known 2D material with its neighboring elements in the periodic table, for example the replacement of the carbon atom in graphene by boron and nitrogen leads to h-BN. Even more structural diversity and possible anisotropy may be achieved by considering the second-order neighbors, e.g. the compound BeN2 would be isoelectronic to graphene, but intrinsically anisotropic in nature. The search of new 2D materials commonly starts with selecting of a known 3D precursor phase, most often thermodynamically stable or very long-lived metastable at ambient conditions, typically with a layered crystal structure with hexagonal motives. Such precursors are not currently available for the Be-N system, although, several compounds in the Be-N system were explored theoretically as potential 2D materials. In particular, the predicted Be3N2 and BeN2 monolayers [10,11] are expected to be dynamically stable may possess high carrier mobilities. At the same time, the electropositive character of the beryllium metal leads to its positive charge and potential excellent adsorption ability, that might be useful in catalysis [12]. Nitrides have been predicted to possess the largest thermodynamic scale of metastability in terms of the energy differences between stable and metastable structures (~190 meV/atom), which may allow the kinetic stabilization of potentially useful compounds [13]. The synthesis of nitrogen-rich compounds is often complicated by the great stability of the dinitrogen molecule N2: the nitrides tend to decompose at relatively low temperatures, that may be below the temperatures that are required for the synthesis. One way to overcome this problem is to use the high-pressure conditions, which prevent the nitride decomposition with the evolution of N2 [14– 17]. Furthermore, high-pressure is a useful tool to affect the bonding types, coordination environments and physical properties of compounds, which may lead to completely unexpected results. For example, beryllium, known to be tetrahedrally coordinated in inorganic compounds, was recently reported to incrementally change its coordination number from four to six in the crystal structure of CaBe2P2O8 under compression [18]. Thus, potentially, the high-pressure synthesis of a beryllium nitride with octahedrally coordinated Be may provide room for unexpected structural changes on decompression. In this Letter we report the exploration of this pathway through the investigation of beryllium polynitrides. The result obtained in the present study overwhelmed our expectations Results and Discussion High-pressure synthesis High-pressure high-temperature chemical reactions in the Be-N system were studied in laser-heated diamond anvil cells (LHDACs) using synchrotron powder and single-crystal X-ray diffraction (Table S1, Materials and Methods). Laser heating of the Be metal embedded in nitrogen (N2) in the pressure chamber at 42 and 60 GPa resulted in the synthesis of previously known α-Be3N2 compound with the anti-bixbyite structure type (cubic space group Ia-3, No. 206, Fig. S1) [19]. The laser-heating of Be in N2 at 98 GPa and ~2200 K leads to the formation of a novel BeN4 compound, as determined from the single-crystal XRD analysis (Table S2). The compound has the CdP4 structure type with the monoclinic space group P21/c (No. 14) and is hereafter referred to as m-BeN4 (Fig. 1). The interconnected nitrogen atoms form 2D networks constructed of condensed corrugated N10 rings with all N atoms in tetrahedral coordination, while the Be atoms possess a distorted octahedral coordination. Such a single-bonded polymeric nitrogen network has never been observed in polynitride compounds, but has been recently predicted independently by Zhang et al. [20] and Wei et al. [21]. At slightly lower pressure of 84 GPa and ~2000 K, laser-heating of Be and N2 resulted in the synthesis of a compound with the chemical composition similar to that of m-BeN4, but with a triclinic structure, further referred to as tr-BeN4. The lattice parameters of tr-BeN4 were determined to be a = 2.4543(14), b = 3.5427(14), c = 3.5110(17) Å3, α = 103.88(4), ß = 111.52(5), γ = 97.39(4)º, V = 26.77(2) Å3. The full crystallographic information is given in the Table S3. The Be atom is coordinated by six nitrogen atoms forming an elongated octahedron. The covalently bonded nitrogen atoms form infinite zigzag chains running along the [001] direction (Fig. 2a). If viewed along the [001] direction (Fig. 2e), the structure can be understood as consisting of slightly corrugated BeN4 layers, interconnected with each other by long Be-N bonds (~1.83 Å). The tr-BeN4 is isostructural with the FeN4 compound, which was previously synthesized at ~106 GPa from Fe and N2 [22,23]. The metal-nitrogen bonds in BeN4 have a slightly more ionic character than in FeN4, which leads to a more uniform electron distribution within the polymeric chains (Fig. S2). Theoretical calculations of the structural properties of tr-BeN4 in a framework of density functional theory showed an excellent agreement between the theory and experiment for the lattice parameters, both at the synthesis pressure and on decompression (Fig. S3). Moreover, theory suggests that the structure is dynamically stable at the synthesis pressure as indicated by the absence of the imaginary frequencies in the phonon dispersion calculations (Fig. S4a). Decompression of BeN4 To investigate if m-BeN4 and tr-BeN4 are recoverable at ambient conditions, the DACs were decompressed in a few pressure steps. The m-BeN4 produced single-crystal XRD diffraction patterns and retained its monoclinic P21/c structure down to ~33 GPa. The compound was not detectable at the next pressure point of ~20 GPa (Fig. S5). The tr-BeN4 produced excellent single-crystal diffraction patterns down to ~60 GPa, but at lower pressures the crystal quality started to deteriorate. Since the theoretical model perfectly matches the experimental structure at high pressures and the evolution of the strongest and well-separated (010) reflection well agrees with that suggested theoretically (Fig. S3d), the evolution of the structure below 60 GPa was tracked by powder XRD and compared with the results of calculations. The weakening of the diffraction signal and the significant overlapping of the strongest reflections of tr-BeN4 with those of ε-N2 upon pressure decrease prevented us from unambiguous refinement of all six unit cell parameters at 20 and 34 GPa. However, the refinement was successfully made at ambient pressure (Figs. S3, S6) after the DAC was fully decompressed and nitrogen released. The decompressed sample was found to contain a significant amount of α-Be3N2, which was not observed at higher pressures. Therefore, we conclude that at ambient conditions tr-BeN4 slowly decomposes into α-Be3N2 but could potentially be completely kinetically stabilized at slightly lower temperature or by keeping the sample in a nitrogen atmosphere. The structural evolution of tr-BeN4 on decompression involves breaking of two longer Be-N bonds in the BeN6 octahedra which leads to a change of the coordination of the Be atom from octahedral to square planar. Consequently, all nitrogen atoms become sp2-hybridised and form infinite planar polyacetylenelike chains (Fig. 3 and Movie S1). The process of rearrangement becomes apparent below 50 GPa, as it can be judged from the interplanar distance (Fig. S3). The evolution of the structure to the layered motif, which is dynamically stable at ambient pressure (Fig. S4b), is accompanied by a significant depletion of the electron density between the layers (Fig. 3 and Movie S1). The topological data analysis (see Materials and Methods) confirms that the charge density value corresponding to the layers’ separation decreases with decreasing pressure. Thus, the BeN4 layers are more weakly connected at ambient pressure than at high pressure, indicating the formation of a van der Waals-bonded solid, as confirmed by our theoretical analysis of the interlayer bonding in tr-BeN4 and the exfoliation energy at ambient pressure (see Materials and Methods). Electronic structure calculations (Fig. S7) show that both the high pressure and the ambient pressure modifications of tr-BeN4 are metallic. The calculated electron density maps and electron localization functions (ELF) at a pressure of P ~ 85 GPa (Fig. S8) indicate a non-uniform distribution of electron density between the nitrogen atoms of the chains, in agreement with the crystal-chemical analysis (Fig. S2). During decompression, the electron density and interatomic distances become almost the same between all the atoms of the nitrogen chains, and the chains become planar (Fig. S8). Similar planar polyacetylene-like nitrogen chains with conjugated π-systems were recently observed in a series of metalinorganic frameworks Hf4N22, WN10, ReN10, Os5N34 [24,25] and MgN4 [26]. The high-pressure synthesis is a well-established path for the discovery of novel materials. A normal highpressure synthesis route includes pressurizing and heating in the stability field of a new compound followed by its recovery to ambient pressure, hoping to preserve the material’s interesting properties [27,28]. In this respect tr-BeN4 is an exciting exception, as its layered ambient-pressure modification obtained upon decompression appears to possess even more intriguing properties than one could expect from the high pressure one. Its bulk electronic structure is modified in a remarkable fashion upon the decompression: one can clearly see the development of a pronounced pseudogap at the Fermi energy in tr-BeN4 at low pressure. A closer examination of the bulk band structure and the Fermi surface (Fig. S7, Movie S2) allows us to conclude that there are two nodal points in the Brillouin zone (BZ) of trBeN4 located at the Fermi energy. Given that both time-reversal and inversion symmetry are preserved, the two points are four-fold degenerate, indicating that they correspond to the Dirac points. The breaking of the inversion symmetry would split each of these points into Weyl points [29]. Though the electronic DOS in bulk tr-BeN4 is finite, it originates from the well-defined Fermi surface pockets (electron and hole). Going from 3D to 2D, those pockets vanish, leaving only the Dirac points in the whole BZ. This makes the layered modification of tr-BeN4 obtained upon decompression a precursor for a new anisotropic 2D Dirac material, the beryllonitrene. Our estimation of the exfoliation energy of the single monolayer for tr-BeN4 (see Materials and Methods) indicates that it is indeed possible. Theoretical study of beryllonitrene The main interest in fundamental studies and potential applications of 2D materials is related to their peculiar electronic properties. Therefore, we now turn to the electronic structure of the beryllonitrene. The optimized crystal structure of single-layer BeN4 is shown in Fig. 4(a). It belongs to the oblique crystal system (the two-dimensional space group p2). The corresponding electronic band structure and density of states (DOS) are shown in Fig. 4(b). One can see that unlike bulk tr-BeN4, the beryllonitrene is a Dirac semimetal with linear dispersion in the vicinity of the Fermi energy and two Dirac points coinciding with the Fermi energy, similar to graphene. Around the Fermi energy, electron and hole states are essentially symmetric over an energy range of roughly 2 eV. It is this property rather than the conical point itself which opens a way to deep relations to high energy physics within its symmetry between particles and antiparticles [4]. Unlike graphene, the Dirac points are not located at the corners of the BZ, but rather along the Γ-A direction (Σ point) as shown in the inset of Fig. 4(b). This can be attributed to a lower symmetry of the BeN4 lattice compared to graphene. Indeed, one can see that the dispersion is different along the Σ-Γ and Σ-A directions, demonstrating a pronounced anisotropy of the band structure even in the long-wavelength limit. From the calculated band structure one can immediately estimate the carrier velocities in single-layer BeN4 as 𝑣𝑣(𝒌𝒌) = 1 ∂E(𝐤𝐤) ℏ ∂𝐤𝐤 . Along the Σ-Γ direction, we find the Fermi velocity to be vy = 0.8×106 m/s, which is 1.3 times smaller than the Fermi velocity in graphene [2,3]. Along the Σ-A direction, the Fermi velocity is essentially smaller (vx = 0.3×106 m/s). At small charge doping the Fermi surface is elliptical with the 1/2 dispersion given by 𝐸𝐸(𝒌𝒌) = ℏ�𝑣𝑣𝑥𝑥2 𝑘𝑘𝑥𝑥2 + 𝑣𝑣𝑦𝑦2 𝑘𝑘𝑦𝑦2 � . The anisotropy of the Fermi velocities is also reflected in the Fermi surface shown in the inset of Fig. 4(b). At sufficiently high doping, the Fermi surface is not elliptical, but is represented by two oppositely oriented semi-circles. It is worth noting that in this case the Fermi surface is perfectly nested by the wave vectors connecting two parallel parts of the two Fermi surface pockets. Therefore, the intervalley scattering in doped single-layer BeN4 is expected to be strongly anisotropic. To gain more insights into the physics of beryllonitrene, we construct a minimal lattice electronic model which describes its band structure (see Materials and Methods). The corresponding “effective lattice” (electronic structure model) is shown in Fig. 4(a) superimposed onto the real-space crystal structure. It is described by a slightly distorted honeycomb net reminiscent to that of graphene. Each node of the honeycomb net is located at the center of the N-N bond. The model allows us to assess the relevance of screening effects in single-layer BeN4. Calculations show (see Materials and Methods) that dimensionless interaction constant ε ≈ 8 is notably larger compared to graphene (3.5), suggesting that Coulomb interaction effects are less important in beryllonitrene. Another essential difference with graphene is the existence of van Hove singularities close enough to the Fermi energy (the separation is about 0.5 eV, contrary to approximately 2 eV in graphene). One can assume that this van Hove singularity is available by a reasonable doping; the Fermi energy shift by 0.5 eV is relatively easily reachable in graphene by chemical doping [30]. Recently, even much higher shift, about 2 eV, was reached in graphene by Gd intercalation [31]. In 2D systems with a Fermi energy close enough to the van Hove singularity one can expect a flat band formation due to many-body effects and all kinds of correlation-induced instabilities (see Ref. [32] and refs. therein), including, e.g., the emergence of superconducting or magnetic phases [33]. To conclude, applying high-pressure synthesis followed by decompression to ambient conditions, we have synthesized new layered material, BeN4. Van der Waals bonds between its layers and presence of anisotropic Dirac cones in its electronic structure show that 2D BeN4, “beryllonitrene”, has unique properties. Indeed, the high degree of electron-hole symmetry makes the 2D BeN4 system similar, in some respect, to the world of high-energy particles with its symmetry between particles and antiparticles. However, contrary to the “true” Universe, which is isotropic, the massless Dirac fermions in the beryllonitrene are essentially anisotropic (see Materials and Methods). This opens a door into the whole new world and the question arises on modification of “conventional” quantum relativistic effects such as Zitterbewegung, chiral tunneling, relativistic collapse at supercritical charges etc. [4] for the anisotropic massless fermions. Materials and Methods Synthesis In every experiment the Be metal was placed in the sample chamber of a BX90 diamond anvil cell equipped with Boehler-Almax type diamonds (Re gasket, 120-250 µm culet diameter). The DACs were loaded with molecular nitrogen which served as a pressure-transmitting medium and as a reagent. The DACs were compressed up to the target pressures and heated with the double-sided laser-heating systems of the P02.2 (Petra III, DESY, Hamburg, Germany) and 13IDD (APS, Argonne, USA) beamlines. The summary of the heating experiments is shown in the Table S1. X-ray diffraction XRD measurements were performed at the beamlines P02.2 beamline of Petra III (DESY, Hamburg, Germany, λ = 0.2891 Å, 2×2 µm2 focusing, Perkin Elmer XRD1621 detector), the 13IDD beamline (APS, Argonne, USA, λ = 0.2952 Å, 3×3 µm2, CdTe Pilatus 1M detector ) , the 16-IDB beamline (APS, Argonne, USA, λ = 0.34453 Å, 3×3 µm2, Pilatus 1M detector) and the ID15B (ESRF, Grenoble, France, λ = 0.4117 Å, 10×10 µm2, MAR555 flat panel detector). For the single-crystal XRD measurements, the samples were rotated around a vertical ω-axis in the ±35° range. The diffraction images were collected with an angular step Δω = 0.5° and an exposure time of 5s or 10s/frame. For the analysis of the single-crystal diffraction data (indexing, data integration, frame scaling and absorption correction) we used the CrysAlisPro software package. To calibrate an instrumental model in the CrysAlisPro software, i.e., the sample-to-detector distance, detector’s origin, offsets of goniometer angles, and rotation of both X-ray beam and the detector around the instrument axis, we used a single crystal of orthoenstatite ((Mg1.93Fe0.06)(Si1.93, Al0.06)O6, Pbca space group, a = 8.8117(2), b = 5.18320(10), and c = 18.2391(3) Å). Powder diffraction measurements were performed either without sample rotation (still images) or upon continuous rotation in the ±20°ω range. The images were integrated to powder patterns with the DIOPTAS software [34]. Le-Bail fits of the diffraction patterns were performed with the Jana2006 and Topas 4.2 software packages [35]. The structures were solved with the ShelXT structure solution program [36] using intrinsic phasing and refined with SHELXL using Olex2 [37,38]. The equations of state (EOS) of the synthesised materials were obtained by fitting pressure-volume dependence data using EoSFit7-GUI software [39]. CSD 2005719-2005720 contain the supplementary crystallographic data for this paper. These data can be obtained free of charge from FIZ Karlsruhe via www.ccdc.cam.ac.uk/structures. First principles calculations of tr-BeN4 The ab-initio calculations were performed using the projector-augmented-wave (PAW) method [40] as implemented in the VASP code [41–43]. The sampling for BZ integrations was performed using the Gamma scheme with 18×18×18 k-point grids. In the calculation of the Fermi surface denser k-points sampling, 31×31×31 was employed. The energy cutoff for the plane waves included in the expansion of wave functions was set to 800 eV. For calculation of exchange–correlation energy we used the PerdewBurke-Ernzerhof GGA (PBE) [44] and optB88b-vdW [45] methods with nonlocal contributions, where the van der Waals interactions are included and the exchange functional is optimized for the correlation part [46]. The convergence criterion for the electronic subsystem was chosen to be equal to 10-4 eV for two subsequent iterations, and the ionic relaxation loop within the conjugated gradient method was stopped when forces became of the order of 10-3 eV/Å. Finite displacement method was used to calculate the phonon dispersions of tr-BeN4 at 0 GPa with Phonopy [47]. Converged phonon dispersions were achieved using a (4 × 4 × 4) supercell with 320 atoms. Electronic properties of tr-BeN4 Analysis of the electronic density of states of tr-BeN4 (Fig. S7) demonstrates that at high pressure the material is metallic and the main contribution to the DOS at the Fermi level comes from px electrons of nitrogen chains, forming π-system. Upon decompression the DOS becomes more localized, the electronic bands become narrower and small energy gaps and pseudogaps between states appear, including a pseudogap in the energy range from -6 to -5 eV due to the break of multicenter Be-N2-Be bonds. We can also observe the pronounced pseudogap at the Fermi energy due to the relaxation of the π-system. Topological data analysis of BeN4 charge density To elucidate the possibility of the synthesis of a single layer of BeN4 we investigated the nature of interlayer bonding of the layered tr-BeN4 considering the evolution of the charge density between the layers upon decompression. From Fig. 3 and Supplementary video 1 one sees that the charge density between the layers decreases significantly with decreasing pressure, and at ambient pressure the charge density becomes almost zero. To quantify the pressure evolution of the charge density, we performed a topological data analysis [48,49] of the valence charge density of tr-BeN4 at seven different pressure values, 83.6 GPa, 40.5 GPa, 25.6 GPa, 17.8 GPa, 11.8 GPa, 6.8 GPa and 0 GPa. A 2x2x2 supercell was created from the unit cell. Merge tree based topological data analysis was done for the supercells at volumes corresponding to the above specified pressures with periodic boundary conditions. Let f:D→R be the charge density field. For a given charge density iso-value ρ, the super level set is defined as the set of all points p∈D with f(p)≥ρ. A merge tree tracks the changes in the topology of the super level set as the iso-value is decreased from high to low values. Thus, it captures the birth, merge and death events of the super level set components. The analysis allows us to identify exact iso values at which qualitative changes occur in the topology of the charge density and to classify the changes. Since we are using the valence charge density for our analysis, Be atoms that donate their valence electrons to the N atoms, become effectively invisible in the valence charge density plots. Thus, the analysis focuses on the behavior of charge density in associated with the N atoms. Let us consider an example of the merge tree constructed at zero pressure. It captured following significant events in the topology of the charge density iso surfaces when charge density value is decreased from high to low (see Fig. S9). 1. The super-level set components are born at the maxima in the charge density (single point, Fig. S9a). These maxima are due to the high number of electrons that N atoms attract. If one would continue to lower the iso value, discrete iso surfaces form around N atoms (Fig. S9 b). 2. These discrete iso surfaces merge at the saddle points in the charge density that correspond to critical values at N-N bonds (according to Bader’s definition [50]) and eventually enclose N chains (Fig. S9c). 3. The charge density enclosing N chains within the same BeN4 layer continue to grow, and eventually merge into one charge density layer (Fig. S9 d,e). 4. All the charge density layers merge into a single connected component (Fig. S9f), meaning that we can no longer distinguish between different BeN4 layers by just looking at the charge density at the given iso value. Fig. S10a shows the obtained critical charge density iso values where super-level set components are born, while Fig. S10b shows the obtained critical charge density iso values where topological change in the connectivity of super level sets (points 3 and 4 above) happens. These values are plotted as a function of pressure so that one can observe the difference in topology of charge density as the pressure changes. Figure S10a shows that at synthesis pressure the distribution of charge is unequal between the N atoms with two different maxima meaning that one pair of N atoms getting a bit more electrons than the other pair, this of course can be an effect of the system being in a six coordination state and two N atoms interact with neighboring Be atoms from different layers and attract electrons from them. As the pressure decreases to ambient pressure the maxima for the different N atoms converge suggesting a uniform interaction with the electron giving Be atom and consequently better separation between BeN4 layers. Figure S10b shows that the layer separation charge density iso value to decreases as the pressure decreases to ambient pressure, that is, the BeN4 layers are more strongly connected at synthesis pressure. In contrast, we observe an increase in the chain separation values with decreasing pressure, which suggests the strengthening of interactions between N chains within a layer at ambient pressure. Combined, we observe a stronger interaction within layers and weaker interaction between layers as the pressure is decreased. In fact, at 83.6 GPa, the chain separation and layer separation occur at similar values of the chare density, suggesting a somewhat uniform distribution of electrons at synthesis pressure and a more discreet distribution at ambient pressure. Van der Waals nature of interlayer bonding of layered tr-BeN4 and exfoliation energy at ambient pressure To elucidate the possibility of the synthesis of a single layer of BeN4 we investigated the nature of interlayer bonding of the layered tr-BeN4. First, we studied the importance of Van der Waals (VdW) interactions in calculations of the structural properties of BeN4. To do so we compare the results of theoretical calculations within VdW-optB88 functional that includes the interactions and semi-local PBE calculations. The results show excellent agreement between PBE, vdW and experimental results at high pressure (Fig. S3a). However, as one sees in Fig. S3a upon decompression the agreement of the results with vdW correction with experiment becomes increasingly better in comparison with the semi-local PBE calculations. This provides additional argument that at ambient pressure we deal with a van der Waals material. Second, we calculated the exfoliation energy, the energy required to remove one atomic layer from the surface of the bulk material, as the difference in the ground-state energy between a slab of N atomic layers and a slab of (N -1) atomic layers plus an atomic layer separated from the slab: 𝐸𝐸𝑒𝑒𝑥𝑥𝑒𝑒 = −[𝐸𝐸𝑁𝑁−𝑙𝑙𝑙𝑙𝑦𝑦𝑒𝑒𝑙𝑙 𝑠𝑠𝑙𝑙𝑙𝑙𝑠𝑠 − 𝐸𝐸(𝑁𝑁−1)−𝑙𝑙𝑙𝑙𝑦𝑦𝑒𝑒𝑙𝑙 𝑠𝑠𝑙𝑙𝑙𝑙𝑠𝑠 − 𝐸𝐸𝑠𝑠𝑠𝑠𝑠𝑠𝑠𝑠𝑙𝑙𝑒𝑒 𝑙𝑙𝑙𝑙𝑦𝑦𝑒𝑒𝑙𝑙 ], where E - is the energy of supercells per surface atom. Two series of calculations with supercells N=10 and N=15 give Eexf = 79.8 meV/atom and 79.3 meV/atom, respectively. We find that the vdW results for the tr- BeN4 exfoliation energy are in the range of values corresponding to the experimental [51,52] and theoretical [53] results for the graphite exfoliation energy. First-principles calculations of single-layer BeN4 Calculations were carried out using density functional theory (DFT) within the projected augmented wave method [40] as implemented in the Vienna ab initio simulation package (VASP) [42,43]. A generalized gradient approximation [44] was employed to describe exchange-correlation effects. The kinetic energy cutoff was set to 400 eV and the BZ sampled by a (14×16) k-point mesh. The vacuum thickness of 25 Å was used to avoid spurious interactions between the supercell periodic images in the direction perpendicular to the 2D plane. An energy window of at least 50 eV was used in the polarizability calculations. The construction of the Wannier functions and tight-binding parametrization of the DFT Hamiltonian was performed using the maximal localization procedure [54] by means of the Wannier90 code [55]. Minimal electronic model of BeN4 band structure In our model calculations, the electronic lattice is described by a slightly distorted honeycomb structure reminiscent of the graphene lattice. Each lattice site is located at the center of the N-N bond Fig. 4a, and is represented by a Wannier function in the form of a superposition of pz-like orbitals, |𝑤𝑤⟩ ≈ |𝑝𝑝z1 ⟩ − |𝑝𝑝z2 ⟩. The electronic properties of such effective lattice can be described by a tight-binding Hamiltonian 𝐻𝐻 = ∑𝑠𝑠𝑖𝑖 𝑡𝑡𝑠𝑠𝑖𝑖 𝑐𝑐𝑠𝑠† 𝑐𝑐𝑖𝑖 , where 𝑡𝑡𝑠𝑠𝑖𝑖 is the hopping integral between lattice sites i and j, and 𝑐𝑐𝑠𝑠† �𝑐𝑐𝑖𝑖 � is the creation (annihilation) operator of an electron at site i (j). The leading interaction terms correspond to the nearest- and next-nearest-neighbor terms with magnitudes t1=1.35 eV and t2=−0.49 eV, respectively. The difference in the hopping parameters is attributed to the lattice anisotropy. The constructed model accurately describes the electronic structure around the Fermi energy, as can be seen from Fig.4(b). To provide further details, we assess the relevance of screening effects in single-layer BeN4. In the longwavelength limit, the dielectric function of a pristine 2D Dirac system in Random Phase Approximation 2 𝑒𝑒 (RPA) reads ε = 1 + π2 α, where α = ℏ𝑣𝑣 is an effective fine-structure constant [4]. Note that lattice Quantum Monte Carlo simulations [56] show that RPA results for the dielectric constant of massless Dirac fermions are pretty close to the exact ones. In anisotropic case we have the same expression but with 𝑣𝑣 = �𝑣𝑣𝑥𝑥 𝑣𝑣𝑦𝑦 [57]. Substituting our calculated velocities, we obtain ε ≈ 8. This may, oppositely, weaken the many-body effects for the limit of low doping but, in any case, perturbation treatment is not enough. Beyond the long-wavelength limit, the screening effects can be quantified by the strength of the effective Coulomb interaction defined in terms of the Hubbard model [4]. Practically, it can be estimated within the constrained random phase approximation (cRPA) [58]. The cRPA assumes static dielectric screening of the bare Coulomb interaction V(rij) between Wannier orbitals defined in terms of the tight-binding model by all high-energy states non included in the model: U(q)=[1−V(q)П(q)]V(q), where U(q) is the effective Coulomb interaction in q-space, and П(q) is the polarizability capturing relevant screening processes. Our estimate yields U0=5.3 eV for the on-site interaction, and U1=3.8 eV and U2=2.2 eV for the nearest- and next-nearest-neighbor interactions, respectively. The obtained interactions are up to two times smaller compared to graphene [59]. At the same time, the half-bandwidth in single-layer BeN4 (W=3.1 eV) is comparable to that in graphene (2.8 eV). The degree of the local electronic correlations can be estimated by the ratio U0/W. We arrive at the following upper limit U0/W ~ 1.7, which suggests that single-layer BeN4 is less correlated material compared to graphene. Particularly, pristine single-layer BeN4 is expected to be further away from the Mott metal-insulator transition, magnetic instabilities, and other exotic quantum effects [60]. Massless Dirac particles with Coulomb interaction For massless Dirac particles with Coulomb interaction the Fermi velocity is renormalized by many-body corrections logarithmically divergent at 𝑘𝑘𝐹𝐹 → 0, 𝑘𝑘𝐹𝐹 is the Fermi wave vector [4,61]. For the case of anisotropic massless Dirac fermions these corrections were calculated in Ref. [57]. The result is: 𝛿𝛿𝑣𝑣𝑥𝑥 = 𝑣𝑣𝑥𝑥2 𝑒𝑒 2 𝑣𝑣𝑥𝑥 �𝐾𝐾(𝑚𝑚) − 2 𝑅𝑅(𝑚𝑚)� 𝜉𝜉, 𝑣𝑣𝑦𝑦 𝜋𝜋ℏ 𝑣𝑣𝑦𝑦 𝛿𝛿𝑣𝑣𝑦𝑦 = 𝑒𝑒 2 𝜋𝜋ℏ [𝐾𝐾(𝑚𝑚) − 𝑃𝑃(𝑚𝑚)]𝜉𝜉 , 2 where 𝜉𝜉 = ln 𝑘𝑘1𝑙𝑙, a is the lattice constant, 𝑚𝑚 = 1 − 𝑣𝑣𝑣𝑣𝑥𝑥2 (we assume 𝑣𝑣𝑥𝑥 ≤ 𝑣𝑣𝑦𝑦 ), 𝐹𝐹 𝑦𝑦 𝜋𝜋/2 𝐾𝐾(𝑚𝑚) = � is the full elliptic integral of the first kind, 0 𝑑𝑑𝑑𝑑 √1 − 𝑚𝑚𝑠𝑠𝑠𝑠𝑠𝑠2 𝑑𝑑 𝜋𝜋/2 𝑃𝑃(𝑚𝑚) = ∫0 With our parameters, we have 𝛿𝛿𝑣𝑣𝑥𝑥 𝑣𝑣𝑥𝑥 𝑑𝑑𝑑𝑑𝑠𝑠𝑠𝑠𝑠𝑠2 𝑑𝑑 (1−𝑚𝑚𝑠𝑠𝑠𝑠𝑠𝑠2 𝑑𝑑)3/2 = 1.9𝜉𝜉, 𝛿𝛿𝑣𝑣𝑦𝑦 𝑣𝑣𝑦𝑦 2 𝑑𝑑𝑑𝑑𝑐𝑐𝑐𝑐𝑠𝑠 𝑑𝑑 , 𝑅𝑅(𝑚𝑚) = ∫0𝜋𝜋/2 (1−𝑚𝑚𝑠𝑠𝑠𝑠𝑠𝑠 2 𝑑𝑑)3/2 . = −3.0𝜉𝜉. The difference in sign reflects the fact that renormalization tends to reduce the anisotropy of carrier velocities, in agreement with Ref. [57]. Keeping in mind that 𝑘𝑘𝐹𝐹 𝑎𝑎 ≫ 1 and 𝜉𝜉 ≳ 1, the obtained corrections turn out to be larger than the order of the calculated velocities. On the one hand, this result means that the first-order perturbative corrections are insufficient to correctly describe the velocity renormalization. On the other hand, this suggests an important role of the many-body effects in the physics of anisotropic Dirac fermions in single-layer BeN4. Anyway, the conclusion is that the many-body effects due to the long-range Coulomb interaction are much stronger in ≈ 0.5𝜉𝜉 [4]. BeN4 than in graphene where in the same approximation 𝛿𝛿𝑣𝑣 𝑣𝑣 References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] K. S. Novoselov, A. K. Geim, S. V. Morozov, Y. Jiang, S. V. Dubonos, I. V. Grigorieva, and A. A. Firsov, Science 306, 666 (2004). K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, M. I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and A. A. Firsov, Nature 438, 197 (2005). Y. Zhang, Y.-W. Tan, H. L. Stormer, and P. Kim, Nature 438, 201 (2005). M. I. Katsnelson, The Physics of Graphene, 2nd ed. (Cambridge University Press, Cambridge, 2020). S.-M. Choi, S.-H. Jhi, and Y.-W. Son, Phys. Rev. B 81, 081407 (2010). C.-H. Park, L. Yang, Y.-W. Son, M. L. Cohen, and S. G. Louie, Nat. Phys. 4, 213 (2008). S. Katayama, A. Kobayashi, and Y. Suzumura, J. Phys. Soc. Japan 75, 054705 (2006). Y. Suzumura, J. Phys. Soc. Japan 85, 053708 (2016). Y. Zhao, X. Li, J. Liu, C. Zhang, and Q. Wang, J. Phys. Chem. Lett. 9, 1815 (2018). C. Chen, B. Huang, and J. Wu, AIP Adv. 8, 105105 (2018). C. Zhang and Q. Sun, J. Phys. Chem. Lett. 7, 2664 (2016). Y. Ding, Y. Ji, H. Dong, N. Rujisamphan, and Y. Li, Nanotechnology 30, 465202 (2019). W. Sun, S. T. Dacek, S. P. Ong, G. Hautier, A. Jain, W. D. Richards, A. C. Gamst, K. A. Persson, and G. Ceder, Sci. Adv. 2, e1600225 (2016). B. A. Steele, E. Stavrou, J. C. Crowhurst, J. M. Zaug, V. B. Prakapenka, and I. I. Oleynik, Chem. Mater. 29, 735 (2017). M. Bykov, K. R. Tasca, I. G. Batyrev, D. Smith, K. Glazyrin, S. Chariton, M. Mahmood, and A. F. Goncharov, Inorg. Chem. 59, 14819 (2020). K. Niwa, T. Terabe, K. Suzuki, Y. Shirako, and M. Hasegawa, J. Appl. Phys. 119, 065901 (2016). M. Bykov, K. V. Yusenko, E. Bykova, A. Pakhomova, W. Kraus, N. Dubrovinskaia, and L. Dubrovinsky, Eur. J. Inorg. Chem. 2019, 3667 (2019). A. Pakhomova, G. Aprilis, M. Bykov, L. Gorelova, S. S. Krivovichev, M. P. Belov, I. A. Abrikosov, and L. Dubrovinsky, Nat. Commun. 10, 2800 (2019). O. Reckeweg, C. Lind, A. Simon, and F. J. DiSalvo, Zeitschrift Für Naturforsch. B 58, 159 (2003). S. Zhang, Z. Zhao, L. Liu, and G. Yang, J. Power Sources 365, 155 (2017). [21] S. Wei, D. Li, Z. Liu, W. Wang, F. Tian, K. Bao, D. Duan, B. Liu, and T. Cui, J. Phys. Chem. C 121, 9766 (2017). [22] M. Bykov, E. Bykova, G. Aprilis, K. Glazyrin, E. Koemets, I. Chuvashova, I. Kupenko, C. McCammon, M. Mezouar, V. Prakapenka, H.-P. Liermann, F. Tasnádi, A. V. Ponomareva, I. A. Abrikosov, N. Dubrovinskaia, and L. Dubrovinsky, Nat. Commun. 9, 2756 (2018). [23] M. Bykov, S. Khandarkhaeva, T. Fedotenko, P. Sedmak, N. Dubrovinskaia, and L. Dubrovinsky, Acta Crystallogr. Sect. E Crystallogr. Commun. 74, 1392 (2018). [24] M. Bykov, S. Chariton, E. Bykova, S. Khandarkhaeva, T. Fedotenko, A. V. Ponomareva, J. Tidholm, F. Tasnádi, I. A. Abrikosov, P. Sedmak, V. Prakapenka, M. Hanfland, H.-P. Liermann, M. Mahmood, A. F. Goncharov, N. Dubrovinskaia, and L. Dubrovinsky, Angew. Chemie Int. Ed. 59, 10321 (2020). [25] M. Bykov, E. Bykova, E. Koemets, T. Fedotenko, G. Aprilis, K. Glazyrin, H.-P. P. Liermann, A. V. Ponomareva, J. Tidholm, F. Tasnádi, I. A. Abrikosov, N. Dubrovinskaia, and L. Dubrovinsky, Angew. Chemie Int. Ed. 57, 9048 (2018). [26] D. Laniel, B. Winkler, E. Koemets, T. Fedotenko, M. Bykov, E. Bykova, L. Dubrovinsky, and N. Dubrovinskaia, Nat. Commun. 10, 4515 (2019). [27] M. Bykov, S. Chariton, H. Fei, T. Fedotenko, G. Aprilis, A. V. Ponomareva, F. Tasnádi, I. A. Abrikosov, B. Merle, P. Feldner, S. Vogel, W. Schnick, V. B. Prakapenka, E. Greenberg, M. Hanfland, A. Pakhomova, H.-P. Liermann, T. Katsura, N. Dubrovinskaia, and L. Dubrovinsky, Nat. Commun. 10, 2994 (2019). [28] P. F. McMillan, Nat. Mater. 1, 19 (2002). [29] S. M. Young, S. Zaheer, J. C. Y. Teo, C. L. Kane, E. J. Mele, and A. M. Rappe, Phys. Rev. Lett. 108, 140405 (2012). [30] R. R. Nair, I.-L. Tsai, M. Sepioni, O. Lehtinen, J. Keinonen, A. V. Krasheninnikov, A. H. Castro Neto, M. I. Katsnelson, A. K. Geim, and I. V. Grigorieva, Nat. Commun. 4, 2010 (2013). [31] S. Link, S. Forti, A. Stöhr, K. Küster, M. Rösner, D. Hirschmeier, C. Chen, J. Avila, M. C. Asensio, A. A. Zakharov, T. O. Wehling, A. I. Lichtenstein, M. I. Katsnelson, and U. Starke, Phys. Rev. B 100, 121407 (2019). [32] D. Yudin, D. Hirschmeier, H. Hafermann, O. Eriksson, A. I. Lichtenstein, and M. I. Katsnelson, Phys. Rev. Lett. 112, 070403 (2014). [33] V. Y. Irkhin, A. A. Katanin, and M. I. Katsnelson, Phys. Rev. B 64, 165107 (2001). [34] C. Prescher and V. B. Prakapenka, High Press. Res. 35, 223 (2015). [35] V. Petricek, M. Dusek, and L. Palatinus, Zeitschrift Für Krist. 229, 345 (2014). [36] G. M. Sheldrick, Acta Crystallogr. Sect. A Found. Adv. 71, 3 (2015). [37] G. M. Sheldrick, Acta Crystallogr. Sect. A 64, 112 (2008). [38] O. V. Dolomanov, L. J. Bourhis, R. J. Gildea, J. A. K. Howard, and H. Puschmann, J. Appl. Crystallogr. 42, 339 (2009). [39] J. Gonzalez-Platas, M. Alvaro, F. Nestola, and R. Angel, J. Appl. Crystallogr. 49, 1377 (2016). [40] P. E. Blöchl, Phys. Rev. B 50, 17953 (1994). [41] G. Kresse and J. Furthmüller, Comput. Mater. Sci. 6, 15 (1996). [42] G. Kresse and J. Furthmüller, Phys. Rev. B 54, 11169 (1996). [43] G. Kresse and D. Joubert, Phys. Rev. B 59, 1758 (1999). [44] J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). [45] J. Klimeš, D. R. Bowler, and A. Michaelides, Phys. Rev. B 83, 195131 (2011). [46] J. Klimeš, D. R. Bowler, and A. Michaelides, J. Phys. Condens. Matter 22, 022201 (2010). [47] A. Togo and I. Tanaka, Scr. Mater. 108, 1 (2015). [48] H. Carr, J. Snoeyink, and U. Axen, Comput. Geom. 24, 75 (2003). [49] C. Heine, H. Leitte, M. Hlawitschka, F. Iuricich, L. De Floriani, G. Scheuermann, H. Hagen, and C. Garth, Comput. Graph. Forum 35, 643 (2016). [50] R. F. W. Bader, Atoms in Molecules: A Quantum Theory (Clarendon Press, Oxford, 1994). [51] W. Wang, S. Dai, X. Li, J. Yang, D. J. Srolovitz, and Q. Zheng, Nat. Commun. 6, 7853 (2015). [52] J. Wang, D. C. Sorescu, S. Jeon, A. Belianinov, S. V. Kalinin, A. P. Baddorf, and P. Maksymovych, Nat. Commun. 7, 13263 (2016). [53] Y. Han, K. C. Lai, A. Lii-Rosales, M. C. Tringides, J. W. Evans, and P. A. Thiel, Surf. Sci. 685, 48 (2019). [54] N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Rev. Mod. Phys. 84, 1419 (2012). [55] A. A. Mostofi, J. R. Yates, Y.-S. Lee, I. Souza, D. Vanderbilt, and N. Marzari, Comput. Phys. Commun. 178, 685 (2008). [56] N. Y. Astrakhantsev, V. V. Braguta, M. I. Katsnelson, A. A. Nikolaev, and M. V. Ulybyshev, Phys. Rev. B 97, 035102 (2018). [57] V. K. Dugaev and M. I. Katsnelson, Phys. Rev. B 86, 115405 (2012). [58] F. Aryasetiawan, M. Imada, A. Georges, G. Kotliar, S. Biermann, and A. I. Lichtenstein, Phys. Rev. B 70, 195104 (2004). [59] T. O. Wehling, E. Şaşıoğlu, C. Friedrich, A. I. Lichtenstein, M. I. Katsnelson, and S. Blügel, Phys. Rev. Lett. 106, 236805 (2011). [60] M. V Ulybyshev, P. V Buividovich, M. I. Katsnelson, and M. I. Polikarpov, Phys. Rev. Lett. 111, 056801 (2013). [61] J. González, F. Guinea, and M. A. H. Vozmediano, Nucl. Phys. B 424, 595 (1994). [62] S. Anzellini, A. Dewaele, F. Occelli, P. Loubeyre, and M. Mezouar, J. Appl. Phys. 115, 043511 (2014). Acknowledgments Parts of this research were carried out at the Extreme Conditions Beamline (P02.2) at DESY, a member of Helmholtz Association (HGF). Portions of this work were performed on beamline ID15 at the European Synchrotron Radiation Facility (ESRF), Grenoble, France. Portions of this work were performed at GeoSoilEnviroCARS (The University of Chicago, Sector 13) and at HPCAT (sector 16) of the Advanced Photon Source (APS), Argonne National Laboratory Funding: Research was sponsored by the Army Research Office and was accomplished under the Cooperative Agreement Number W911NF-19-2-0172. N.D. and L.D. thank the Deutsche Forschungsgemeinschaft (DFG projects DU 954-11/1 and, DU 393-9/2, and DU 393-13/1) and the Federal Ministry of Education and Research, Germany (BMBF, grant no. No. 05K19WC1) for financial support. D.L. thanks the Alexander von Humboldt Foundation for financial support. Theoretical analysis of chemical bonding was supported by the Russian Science Foundation (Project No. 18-12-00492). Calculations of the phonon dispersion relations were supported by the Ministry of Science and High Education of the Russian Federation in the framework of Increase Competitiveness Program of NUST MISIS (No. K2-2019-001) implemented by a governmental decree dated 16 March 2013, No. 211. Support from the Knut and Alice Wallenberg Foundation (Wallenberg Scholar Grant No. KAW2018.0194), the Swedish Government Strategic Research Areas in Materials Science on Functional Materials at Linköping University (Faculty Grant SFO-Mat-LiU No. 2009 00971) and SeRC, the Swedish Research Council (VR) grant No. 2019-05600 and the VINN Excellence Center Functional Nanoscale Materials (FunMat-2) Grant 2016–05156 is gratefully acknowledged. The computations were enabled by resources provided by the Swedish National Infrastructure for Computing (SNIC) partially funded by the Swedish Research Council through grant agreement no. 2016-07213. The work of MIK was supported by the JTC-FLAGERA Project GRANSPORT. GeoSoilEnviroCARS is supported by the National Science Foundation – Earth Sciences (EAR – 1634415) and Department of Energy-Geosciences (DE-FG0294ER14466). HPCAT operations are supported by DOE-NNSA’s Office of Experimental Sciences. Advanced Photon Source is U.S. Department of Energy (DOE) Office of Science User Facility operated for the DOE Office of Science by Argonne National Laboratory under Contract No. DE-AC0206CH11357 Author contributions: M. B., L. D. and N. D. initiated and designed the study. M.B., T.F., S.C., D.L., K.G., M.H., J.S.S., V.B.P., M. F. M., A.F.G, N.D. and L.D. contributed to synchrotron experiments. M.B. analyzed single-crystal and powder diffraction data. T.F. analyzed powder diffraction data for Be3N2. A. N. R. and M. I. K. developed the theoretical model of BeN4 monolayer. A. I. A., T.B.M. and I. H. performed scientific visualization and analysis of charge density. A. V. P., F. T., I. A. A. performed first principles calculations. All authors discussed the results and contributed to the final manuscript. Competing interests: Authors declare no competing interests. Data and materials availability: CSD 2005719-2005720 contain the supplementary crystallographic data for this paper. These data can be obtained free of charge from FIZ Karlsruhe via www.ccdc.cam.ac.uk/structures. Figures and Tables Fig. 1. Crystal Structure of m-BeN4 at 98 GPa. (a) View along the [100] direction; (b) view along the [010] direction; (c) part of the crystal structure showing the smallest repeating units of four nitrogen atoms, each unit is shown in its own color; (d) 10-member ring of N atoms with interatomic distances indicated in Å; (e) BeN6 coordination polyhedron; (f) Lewis formula of the smallest repeating nitrogen unit N42-. Grey and blue represent Be and N atoms, respectively. Fig. 2. Crystal structure of tr-BeN4 at 84 GPa. (a) The view along the [100] direction; (b) coordination of a Be atom by one of the polymeric nitrogen chains. (c) BeN6 octahedron with Be-N distances given in Å. (d) Coordination geometry of the nitrogen chain. (e) Scheme showing the transformation of tr-BeN4 upon decompression. (f) Single layer of BeN4 at ambient pressure. (g) Stacking sequence of polymeric nitrogen chains belonging to neighboring layers of tr-BeN4 at ambient pressure; viewing direction is perpendicular to the layer. Beryllium and nitrogen atoms are represented as grey and blue spheres, respectively. Fig. 3. Evolution of the crystal structure and the charge density of tr-BeN4 upon decompression from the synthesis pressure to ambient. At the synthesis pressure, Be and N atoms form corrugated nets and the Be atoms (green) are octahedrally coordinated by six N atoms (blue) (left figure). Upon decompression, below about 50 GPa, the structure relaxes in such a way that the corrugated nets of Be and N atoms flatten, their increasing mutual separation reduces the coordination number of Be atoms to 4, thus turning the structure motif to layered (right figure). The grey planes between the layers indicate the position of the charge density slice and visualize the depletion of the charge density between two layers of tr-BeN4 at zero pressure relative to that at synthesis pressure. Topological data analysis of the charge density (see Materials and Methods, Figs. S9-S10) confirms a stronger interaction within layers and weaker interaction between layers as the pressure is decreased. Experimental lattice parameters for trBeN4 at 1 bar: a = 3.275(2), b = 4.212(2), c = 3.704(2) Å, α = 103.43(3), β = 105.61(4), γ = 111.86(4)º, V = 42.4 Å3. Optimized atomic positions: Be (0.5 0 0), N1 (0.1592 0.66170 0.5144), N2 (0.15840 0.6630 0.14780). Fig. 4. Electronic structure of the beryllonitrene(a) Effective electronic lattice (green) superimposed onto the real-space crystal structure of single-layer BeN4 (a). Green balls centered on the N-N bond correspond to the wave functions shown on the right side. Fourfold-coordinated gray balls depict Be atoms, dirty-blue three-fold coordinated balls correspond to N atoms. Arrows denote two main interaction parameters defined in terms of a tight-binding Hamiltonian. (b) Electronic band structure and the density of states calculated using DFT (red) and the electronic lattice model (blue) shown in (a). The inset shows the Brillouin zone with the high-symmetry points, as well as the Fermi surface corresponding to the Fermi energy of -0.5 eV Supplementary Materials Fig. S1. (a) Powder diffraction pattern of as-synthesised Be3N2 at 42 GPa (λ = 0.29 Å). (b) The crystal structure of α-Be3N2. Green and blue spheres represent the positions of Be and N atoms, respectively. The structure can be understood as a defect fluorite-type structure with distorted cubic closest packing of the nitrogen atoms, while the metal atoms occupy 75% of the tetrahedral holes in an ordered manner. The [𝟒𝟒] [𝟔𝟔] coordination pattern can be expressed as ∝𝟑𝟑[𝑩𝑩𝑩𝑩𝟑𝟑 𝑵𝑵𝟐𝟐 ], where the polyhedra around Be and around the nitrogen atoms are both distorted. (c) The pressure dependence of the unit cell volume of Be3N2. Experimental points are shown by circles. The solid curve is the fit of the experimental P-V data using the 3rd order Birch-Murnaghan EOS with the following parameters: V0= 542.3(±1.3) Å3 𝑲𝑲𝟎𝟎 = 269(14) GPa, and 𝑲𝑲′𝟎𝟎 = 1.6(3). Errors in experimental data are within the symbol size. Green symbols and points show the calculated EoS of Be3N2 obtained within vdW-optB88 functional. (d) F-f plot of α-Be3N2 based on the experimental data. Fig S2. Comparison of the bonding within the polymeric nitrogen chains in a series of tetranitrides FeN4, BeN4 and MgN4. On the first glance, the connectivity within the nitrogen chains of BeN4 is similar to that in FeN4: The N1 atom has a trigonal–planar coordination, while N2 is tetrahedrally coordinated, suggesting sp2 and sp3 hybridization, respectively. Therefore, the nitrogen atoms should form catenapoly[tetraz-1-ene-1,4-diyl] anions [–N–N–N=N–]2-∞. However, the N chains in BeN4 are much more planar than in FeN4 with N2N2N1N1 torsion angle of 17.1º vs 36º in FeN4. BeN4 is in the intermediate situation between FeN4 and MgN4, [22,26] with the latter containing planar polyacetylene-like N chains with N2N2N1N1 angle of 0º. Polyacetylene-like chains with conjugated π-system and polytetrazene chains are isoelectronic with an ionic formula [N4]2-, however the polyacetylene chains are expected to be a more stable structural unit due to the stabilization by conjugation. Nevertheless, significant polar covalent character of the Fe-N bonds imposes geometric constraints, that leads to the redistribution of electron density in the nitrogen chains and to the formation of distorted chains with alternating double and single bonds. Indeed, the electronegativity of metal decreases in the row Fe – Be – Mg suggesting that the latter should form ionic compounds with nitrogen. Calculated Bader charges on Be and Fe in BeN4 and FeN4 (+1.77 and +1.22 respectively) support more ionic character of BeN4. We therefore suggest that BeN4 contains unique N-chains that represent the intermediate state between polytetrazene and polyacetylene-type chains. Fig. S3. (a-c) Calculated and experimental lattice parameters of tr-BeN4 as a function of pressure. (d) Pressure dependence of the d-spacing of the reflection (010) compared with the theoretical model. (e) Calculated d-d interplanar distance of tr-BeN4 as a function of pressure. Fig. S4. Phonon dispersion relations calculated for tr-BeN4 at volumes 26.8 Å3 (a) and 43.05 Å3 (b), which correspond to pressures of P ~ 85 GPa and P ~ 0 GPa. Fig. S5. Calculated EoS of m-BeN4 (P21/c) (a) and pressure dependence of enthalpy difference between tr-BeN4 (P-1) and m-BeN4 ΔH= H(tr-BeN4)-H(m-BeN4) (b). The transformation between tr-BeN4 and mBeN4 leads to a conversion of some double N=N bonds to more energetic single N-N bonds, that should increase the overall energy of the nitrogen network. However, this loss in energy is compensated by the gain in volume with V(tr-BeN4) = 26.47 Å3/fu and V(m-BeN4) = 25.53 Å3/fu at the same pressure of 93 GPa, which makes the PΔV term for the transformation tr-BeN4 → m-BeN4 as high as -0.55 eV. Theoretical calculations (Fig. S3b) show that m-BeN4 become more energetically favorable polymorph of BeN4 at a pressure of ~105 GPa. Fig. S6. (a) Powder diffraction pattern of decompressed sample #2 containing layered tr-BeN4 and αBe3N2. (b) Powder diffraction pattern of sample #2 at the synthesis pressure (88 GPa). Wavelength 0.2952 Å for both patterns. e) Fig. S7. Electronic band structure (a, b) and density of states (c, d) of tr-BeN4 at pressures P~85 GPa (a, c) and P~0 GPa (b, d), which correspond to volumes of 26.8 Å3 and 43.05 Å3, respectively . (e) The Fermi surface tr-BeN4 at volume 43.05 Å3 (pressure P~0 GPa) The position of the Dirac cone (D point) is calculated to be at (-0.37, -0.18, -0.12) in reciprocal coordinates. See Movie S2 for a close view of the anisotropy of the cone. Fig S8. (a, b) Calculated charge density map of tr- BeN4 at pressures P ~ 85 GPa and P ~ 0 GPa in the plane of the nitrogen chain. (c,d) Calculated electron localization function with isosurface value 0.8 (ELF) of tr-BeN4 at pressures P~85 GPa and P~0 GPa, which correspond to volumes of 26.8 Å3 and 43.05 Å3. Calculated electron localization function indicates strong covalent bonding between nitrogen atoms at high and ambient pressures. Besides, at P ~ 85 GPa we can see attractors localized near Be and N2 atoms corresponding to multicenter bonds, which allow beryllium to have six-fold coordination at such a compression. The bond between the Be and N1 atoms is two-center polar covalent. Upon decompression the multicenter Be-N2-Be bonds break. (a) (d) (b) (e) (c) (f) Fig S9. Topological data analysis of tr-BeN4 charge density at zero pressure. Merge tree is constructed by following evolution of the charge density representation as the iso value decreases from high to low values. The exact value of the charge density corresponding to topological changes of the obtained iso surfaces is captured by the merge tree analysis. (a) High iso values give rise to an iso surface that appears first around the N atoms. (b) This surface grows as the iso value is lowered. (c) The first topological change occurs as the charge density iso surfaces merge into N chains. (d) As we continue to lower the iso value the iso surfaces keep expanding strictly within a BeN4 layer. (e) The second topological change is observed as the charge density iso surfaces surrounding the N chains grow large enough to merge and cover the entire BeN4 layer. Charge density value corresponding to this topological change will be referred to as chain separation. (f) Finally, the third topological change occurs when the iso value is so low that the iso surfaces covering individual BeN4 layers merge into one. Charge density value corresponding to this topological change will be referred to as layer separation. (a) (b) Fig S10. (a) The difference in the max value between two types of N atoms as a function of pressure. One can see that at synthesis pressure there is a clear difference between the two N atoms, while the atoms become similar at ambient pressure. (b) The critical charge density values corresponding to chain and layer separation (see Fig. S9 e and f, respectively) as a function of pressure. At low pressure, the values are far apart, which means that the charge density is significantly more localized within physical BeN4 layers. On the contrary, at the synthesis pressure the chain and layer separation values of the charge density are close to each other, indicating a more uniform distribution of the charge density between BeN4 layers. Table S1. Summary of syntheses in Be-N system. Reaction Pressure, GPa* Temperature, K Products mixture Detection method DAC # 1 Be + N2 41.8(5) 2100(200) α-Be3N2 PXRD DAC # 1 Be + N2 + Be3N2 60.5(5) > 3500(300) α-Be3N2 PXRD DAC # 2 Be + N2 84(4) 2020(200) P-1 BeN4 SCXRD DAC # 3 Be + N2 98(4) 2145(200) P21/c BeN4 SCXRD DAC # 4 Be + N2 93(4) 2000(200) P-1 BeN4 + P21/c SCXRD BeN4 * In the experiments 2-4 pressure was determined based on the Equation of State of Re (Ref. [62]). For the pressure determination, the diffraction pattern from the very edge of the rhenium gasket was measured. According to Anzellini et al. [62], such pressure determination method gives pressure with ~5% accuracy. Table S2. Crystal structure, data collection and refinement details of mBeN4 at 96.8 GPa Crystal data Chemical formula BeN4 Mr 65.05 Crystal system, space group Monoclinic, P21/c Temperature (K) 293 Pressure (GPa) 96.8 a, b, c (Å) 3.283(2), 3.2765(10), 4.8185(19) ß (°) 99.67 (6) V (Å3) 51.10 (4) Z 2 Radiation type Synchrotron, λ= 0.34453 Å μ (mm-1) 0.10 Crystal size (mm) 0.002 × 0.002 × 0.002 Data collection Diffractometer APS 16IDB, Pilatus 1M detector Absorption correction Multi-scan CrysAlis PRO 1.171.40.45a (Rigaku Oxford Diffraction, 2019) Empirical absorption correction using spherical harmonics, implemented in SCALE3 ABSPACK scaling algorithm. Tmin, Tmax 0.442, 1.000 No. of measured, independent and observed [I > 2s(I)] reflections 102, 60, 53 Rint 0.024 (sin θ/λ)max (Å-1) 0.732 Refinement R[F2 > 2s(F2)], wR(F2), S 0.097, 0.229, 1.21 No. of reflections 60 No. of parameters 20 Δρmax, Δ ρmin (e Å-3) 0.71, -0.63 Crystal Structure Wyckoff Site Coordinates (x, y, z) Uiso (Å2) Be 2d 0.5, -0.5, 0 0.018(4) N1 4e 0.095(3), -0.892(2), -0.1048(17) 0.015(5) N2 4e 0.238(3), -0.1947(17), -0.2559(18) 0.010(3) Table S3. Crystal structure, data collection and refinement details of trBeN4 at 84 GPa Crystal data Chemical formula BeN4 Mr 65.05 Crystal system, space group Triclinic, 𝑃𝑃1� Temperature (K) 293 Pressure (GPa) 83.5 a, b, c (Å) 2.4062 (8), 3.5119 (9), 3.5060 (14) α, β, γ (°) 104.29 (3), 111.39 (3), 96.71 (3) V (Å3) 26.02 (2) Z 1 Radiation type Synchrotron, λ = 0.2885 Å µ 0.09 (mm-1) Crystal size (mm) 0.002 × 0.002 × 0.002 Data collection Diffractometer LH@P02.2 Perkin Elmer XRD 1921 Absorption correction Multi-scan CrysAlis PRO 1.171.40.54a (Rigaku Oxford Diffraction, 2019) Empirical absorption correction using spherical harmonics, implemented in SCALE3 ABSPACK scaling algorithm. Tmin, Tmax 0.728, 1 No. of measured, independent and observed [I > 2σ(I)] reflections 167, 107, 94 Rint 0.032 (sin θ/λ)max (Å-1) 1.041 Refinement R[F2 > 2σ(F2)], wR(F2), S 0.074, 0.188, 1.17 No. of reflections 107 No. of parameters 10 ∆ρmax, ∆ρmin (e Å-3) 0.41, -0.48 Crystal Structure Wyckoff Site Coordinates (x, y, z) Uiso (Å2) Be 1d 0.5000 0.0000 0.0000 0.0085(8) N1 2i 0.1120(10) 0.6776(5) 0.5034(13) 0.0066(5) N2 2i 0.0574(10) 0.6894(6) 0.1292(13) 0.0069(6)