## 2 Materials and Methods  
All calculations implemented the Vienna Ab initio Simulation Package [19; 20; 21; 22] and the Perdew-Burke-Ernzerhof formulation to the generalized-gradient approximation for the exchange-correlation functional [23; 24]. By using a \(3\times 3\times 3\) fcc lattice supercell (108 inequivalent atoms), we may gauge the short- and long-range effects of the stabilization solute, Ga, and the mono-vacancy defect. In addition, this supercell has numerous nearest neighboring (nn) shells complete, i.e., the first (1 nn), second (2 nn), third (3 nn), fourth (4 nn), and sixth (6 nn) are filled. The seventh (7 nn), ninth (9 nn), and 10 th (10 nn) are incomplete, and there is no atoms in the eighth (8 nn) shell due to the finite size of the supercell and periodic boundary conditions. The Ga concentrations investigated included 0.926, 1.852, and 2.778 at. % Ga, and in cases where there are multiple Ga atoms, the Ga atoms where located at 4 nn distance in a substitutional Pu site from each other since this was determined to be the most energetically favorable site and as observed in EXAFS data (Figure 1) [3; 25; 26; 27; 28]. Density functional theory has determined that for metallic Pu an anti-ferromagnetic (AFM) arrangement is the ground state [29; 30], even though it has been experimentally determined that there is no static long range ordering [31]. By using an AFM configuration, this allows for a continuous description of the electronic structure throughout the many phases of Pu [29]. The use of magnetism within DFT allows an approximation for the strongly correlated 5f electrons and provides comparable structural properties to experimental data [25; 30], which allows a first-principles approach to understanding the complex energetic potential landscape for defects as these calculations require the use of supercells. We would also like to note that a non-magnetic configuration significantly contracts the lattice [25], and test calculations with the Ga atom yielded an instability of site substitutional Ga,unless the volume is contracted. As for theoretical methods beyond DFT, such as dynamical mean field theory, there is an inherent limitation of the theory that will not allow the use of supercells that are necessary for representation of the symmetry breaking due to a defect or alloying solute [32; 33]. Theoretical work extending to the Hubbard \(U\) model (DFT+U) has recently suggested that a small U value of \(\sim\)1 eV may be used for Pu [34] instead of the higher value of \(\sim\)4 eV [35], and it has also been shown that \(\delta\)-Pu is highly sensitive to the U value in lattice parameter, bulk modulus, and density of states [36]. [PERSON] et al. [30] also calculated that the larger U value (\(\sim\)4 eV) will give the unalloyed \(\delta\)-phase an unphysical lower energy than the \(\alpha\)-phase. Thus, in general, Pu is highly sensitive to the chosen value of U, and since it is expected that local distortions may take on the form as the low symmetry phases of Pu [37; 38], which exhibit a higher degree of delocalization than the \(\delta\)-phase [30; 33], the addition of a U value could highly restrict this natural progression of possible theoretical mechanisms of destabilization. Additional theoretical works have taken into account a disordered moment as a stabilization mechanism for \(\delta\)-Pu and have explored this along with other stabilization alloying impurities in addition to Ga [39; 40; 41; 42]. For the purpose of this work, the Pu atoms were configured as a spin-polarized AFM arrangement with no spin-orbit coupling and orbital polarization added. It has been determined that the addition of relativistic effects (spin-orbit coupling and orbital polarization) compensates for the spin-polarization effects [43]. In addition, energetic trends predicted by DFT will not significantly change as it has been calculated that the energy difference between no-spin-orbit-coupling and spin-orbit-coupling for formation energies of a vacancy in 3.125 at.