Ab initio calculations of the surface structure and properties of the (001) surface of galena (PbS) have been performed using two complementary approaches. Both Gaussian and numerical basis functions have been utilized, as embodied within the programs CRYSTAL98 and SIESTA, with their associated implementations, in order to determine to what extent computational factors are influencing the results. Relaxations are presented based on calculations employing the HF, B3LYP, LDA, and GGA Hamiltonians. Significant relaxations (up to about 0.2 Angstrom) are found to occur at the surface which decay with increasing depth into the bulk. The magnitude of the relaxation oscillates between layers, leading to pairing, with Pb ions being displaced outwards (by approximately 0.05 Angstrom) relative to the S ions. We find that the LDA, GGA, and B3LYP Hamiltonians give results in reasonable agreement with each other, but that the structures computed with HF theory are vastly different. The present calculations also predict a surface core level shift of 0.3 eV for the sulphur 2p level at the surface itself, in excellent agreement with a recent experimental determination. Copyright (C) 2003 Elsevier Science Ltd.