I would like to study the stability of different surface orientations of wurtzite/zinc blende materials which exhibit asymmetric surfaces for certain lattice orientations and therefore would like to use the periodic surface model (electrode) to compute the relative stabilities of these surfaces.

I do not yet have the full licence in my disposal, so I ran some test calculations on my desktop for Au(111) surface. However, when I compute the total energy of surfaces with varying thickness, in some instances I get negative surface energies.

For example, when using surfaces with 6 and 8 layers of gold, medium basis set and SG15 pseudo, I get energies

6 -22519.98379

7 -26273.25555

8 -30026.51107

Making a linear fit gives a slope of -3751.264 eV, close to that of a bulk calculation (-3753.231 eV). However, the intercept is negative and so is the surface energy.

If I use a HGH pseudo, I get a positive surface energies that vary depending on the basis set used.

What I would like to know:

Is this a correct approach to take and are the results some kind of an artifact resulting from a thin surface model or am I making a fundamental error somewhere or how should I approach the issue? In the original paper ("First-principles Greenâ€™s-function method for surface calculations: A pseudopotential localized basis set approach") the authors discuss developing the procedure for SG15 pseudos. However, QuantumATK allows me to set from a variety of different pseudos. Is the Green's function method suitable for all the available pseudos and basis sets?

- Weckman