We performed spin-polarized periodic DFT calculations with the all-electron, frozen-core, projector augmented-wave (PAW) 36 method using the Vienna Ab initio Simulation Package (VASP) version 6.3.1. We used the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional together with Grimme's D3 dispersion correction and the Becke-Johnson damping function. We applied a kinetic-energy cutoff value of 500 eV for the plane-wave basis set. We utilized a 4´4 supercell containing four atomic layers (64 atoms) to simulate surfaces when evaluating stability of SAAs and determine reaction pathways of direct methane conversion. We constructed the periodic cell with a vacuum spacing of at least 15 Å along the direction normal to the SAA surface, resulting in a unit cell length of 30 Å along the c-axis to prevent interactions between the periodic images. The atoms in the two topmost layers of the 4´4 four-layer slab, as well as all adsorbed species, were relaxed, while the atoms in the two bottommost layers were fixed at their bulk atomic positions to simulate the semi-infinite Cu bulk crystal structure. Γ-pointcentered Monkhorst-Pack 45 k-point grids of 15´15´15 and 5´5´1 were applied to sample the Brillouin zone for the bulk unit cell and the 4´4 periodic slabs, respectively. Benchmark calculations indicated that the selection of kinetic-energy cutoff (500 eV vs 1000 eV) and k-point grids (15´15´15 vs 24x24x24 for bulk and 5´5´1 vs 6x6x1 for slab) converged the total energies to 1 meV/atom. To aid the self-consistent field convergence, we used Methfessel-Paxton smearing for Brillouin zone integration with a smearing width of 0.1 eV. In geometry optimizations, forces were converged to 0.03 eV/Å. We applied dipole-field energy and potential corrections to cancel the artificial field interaction between the slabs along the direction normal to the SAA surface.