## 2 Computational Methods  
### Electronic Structure Method  
The results from the electronic-structure calculations presented here come from applying density-functional theory that took shape in the important works by [PERSON], [PERSON], and [PERSON] [50; 51]. Density-functional theory is in principle correct, but it depends on practical simplifications. The most important one is the treatment of the electrons exchange and correlation interactions. For these interactions many different approximations have been proposed. Here we are applying the generalized gradient approximation (GGA) that is better for the actinides than the previously used local density approximation [52] and even some more recent formulations [53]. In addition to spin-orbit interaction we extend DFT to include effects of orbital-orbital interactions as we will discuss below. In terms of the technical details of the calculations, we are applying two approaches. For all calculations, except those related to the self-consistent phonon method (see below), we utilize an all-electron full-potential linear muffin-tin orbital (FPLMTO) method that is well established [54] and has been described in detail [55]. In this method, no approximations are made for the electron core states that lie deeper in energy than the valence states, unlike the pseudopotential method. The core approximation made in the latter technique is computationally efficient but can sometimes cause inaccuracies that are avoided in all-electron calculations. The present implementation does not make any assumptions beyond that of the electron exchange and correlation functional. Basis functions, electron densities, and potentials are calculated without any geometrical approximation and these are expanded in spherical harmonics (with a cut-off \(l_{\text{max}}=8\)) inside non-overlapping (muffin-tin) spheres surrounding each atom and in Fourier series in the region between these muffin-tin spheres. There is a choice how to define the length of the muffin-tin sphere radius (r\({}_{\text{MT}}\)) and here it is chosen so that r\({}_{\text{MT}}\)/r\({}_{\text{WS}}\)\(\sim\) 0.8, where r\({}_{\text{WS}}\) is the Wigner-Seitz (atomic-sphere) radius. The radial part of the basis functions inside the muffin-tin spheres are calculated from a wave equation for the \(l=0\) component of the potential that include all relativistic corrections including spin-orbit coupling for d and f states but not for the p states, following the arguments in Reference [56]. Likewise, orbital-orbital interaction (orbital polarization, OP) is only considered for d and f states. The formulation of the orbital polarization is obtained from [PERSON], [PERSON], and [PERSON] [57] and this interaction has been shown to be important for actinides, particularly plutonium [58]. For the two-atom per cell UN structure we applied 500 k points in the irreducible part of the Brillouin zone and each energy eigenvalue was broadened with a Fermi-Dirac distribution (T = 300 K). Uranium mononitride is an anti-ferromagnet but within the GGA a ferromagnetic ordering has lower energy and we chose to adopt ferromagnetism for most calculations. However, the energetics and bonding are very similar for both and close to experiment [21], see Table 1, and the actual choice (ferromagnetic or anti-ferromagnetic) does likely not matter much for the properties we focused on here. The reported atomic volume, bulk modulus and its pressure derivative in Table 1 are obtained from a Vinet fit [59] to the calculated DFT pressures. For the self-consistent phonon method (see below) we applied a 54-atom UN cell in the calculations of the forces that are required for this approach. As this task is computationally rather challenging, we utilized a more efficient electronic-structure code with the additional simplifications of ignoring magnetic ordering and spin-orbit interaction. Hence, we employed the Vienna ab initio simulation package VASP [60; 61; 62].