Example system: Cu(111)-Te ad-chains

The following example covers a ViPErLEED analysis and structure optimization for a copper-telluride layer on Cu(111) based on the paper by Kißlinger et al. [18]. The experimental data and the outline for the analysis steps described below were provided by the authors of the publication. This system deals with a (5 × √3)_rect superstructure and has 10 atoms per bulk layer. It serves as an example for how a rather challenging system can be treated with ViPErLEED. The analysis presented in the original publication was performed using an early development version of ViPErLEED [18].

Fig. 16 Experimental LEED pattern of Cu(111)-(5 × √3)_rect-Te at an incident electron energy of 60 eV. Diffraction spots are labeled by the ViPErLEED spot tracker [2].

All input files needed to follow along on your own machine can be found here.

Tip

If you are new to ViPErLEED, check out other tutorials first, which explain the basics in more detail.

Introduction

The system we are analyzing is a (5 × √3)_rect superstructure of copper telluride on the hexagonal Cu(111) substrate. The tellurium coverage is 0.40 monolayers, which corresponds to four tellurium atoms per surface unit cell.

From experiments, we have a set of 79 beams in our energy range (20 to 500 eV), stored in the EXPBEAMS.csv file. These beams correspond to a total energy range of around 17500 eV, of which only around 1800 eV are from integer-order beams. Fig. 16 shows a snapshot of the experimental LEED pattern.

For our example analysis, we start from a qualitatively correct structure model (POSCAR file) with correct layer stacking [18]. However, the initial atom positions in the POSCAR are taken from bulk Cu, so we should expect significant relaxation during structure optimization. Fig. 17 shows a rendering of the initial structure (produced using VESTA [12]).

We will need to run multiple delta-amplitude calculations and structure optimizations. The individual DISPLACEMENTS files are not all shown in full below, but you can download all input files here.

Fig. 17 POSCAR of Cu(111)-(5 × √3)_rect-Te rendered in VESTA, viewed along \(\mathbf{a}\) (left), \(\mathbf{b}\) (center), and \(\mathbf{c}\) (right). Cu atoms are shown in blue, Te atoms in yellow.

PARAMETERS

Similar to the previous examples, we start by setting up a PARAMETERS file as in Listing 16.

Listing 16 PARAMETERS file for Cu(111)-(5 × √3)_rect-Te.

! ####### GLOBAL PARAMETERS #######
RUN = 0

THEO_ENERGIES = 20 500 3
LMAX = 8-14

! ####### PARAMETERS FOR INTERPRETING THE POSCAR #######

BULK_LIKE_BELOW = 0.6

SITE_DEF Te = surf top(2)
SITE_DEF Cu = surf top(2)

SUPERLATTICE M = 5 0, 1 2

! ####### PARAMETERS FOR VIBROCC #######

T_DEBYE = 343
T_EXPERIMENT = 100
VIBR_AMP_SCALE = *surf 1.3

Most of the parameters in this file have been explained in previous examples, so we will skip some details here. Worth mentioning in this particular case is that we use SITE_DEF to define the two topmost copper and tellurium atoms as explicit sites. We also give the in-plane relationship between the surface and bulk via the SUPERLATTICE parameter, to match the convention used in EXPBEAMS.csv. We use the explicit matrix notation as this particular reconstruction cannot be expressed in Wood’s notation.

Tip

For non-obvious relationships between surface and bulk unit cells it is usually better to explicitly specify the SUPERLATTICE parameter. This ensures consistency with the experimental data. This is especially the case when multiple symmetry-induced domains coexist on the surface because of a difference in the plane groups (or unit-cell sizes) of the bulk and the superstructure. Since all domains are symmetry equivalent, ViPErLEED has no way to choose among them.

As usual, the IVBEAMS and PHASESHIFTS files are generated automatically during initialization (RUN = 0). The VIBROCC file will also be generated by ViPErLEED based on the T_DEBYE, T_EXPERIMENT, and VIBR_AMP_SCALE parameters provided.

Rough DISPLACEMENTS

For the delta-amplitudes calculation and the structure optimization (RUN = 1-3), we will start out with a very rough grid of geometric optimizations (DISPLACEMENTS file). Refer to Listing 17. As is generally recommended, we start with optimizations normal to the surface (\(z\)). Immediately following that comes an in-plane optimization. Due to the complex structure, one in-plane direction is not sufficient: we need to run an optimization along both \(x\) and \(y\) directions. For details on the syntax used, see the entry on the DISPLACEMENTS file.

Listing 17 DISPLACEMENTS (30 pm range) for Cu(111)-(5 × √3)_rect-Te.

== SEARCH z
    = GEO_DELTA
        Te_surf       z = -0.30 0.30 0.02
        Te_def        z = -0.30 0.30 0.02
        Cu_surf       z = -0.20 0.20 0.02
        Cu_def L(2-4) z = -0.10 0.10 0.01

== SEARCH xy
    = GEO_DELTA
        Te_surf xy         = -0.30 0.30 0.02
        Te_def  xy         = -0.30 0.30 0.02
        Cu_surf xy         = -0.20 0.20 0.02
        Cu_def  xy L(2-4)  = -0.10 0.10 0.01

Note

Note that we allow displacements of up to 0.3 Å for tellurium to speed up the convergence in this example. Normally, this is not recommended, because the tensor-LEED approximation will lead to substantial errors for such large displacements.

As usual, we can see a visualization of the optimization convergence in the Search-report.pdf file of the OUT directory. Fig. 18 shows the first page of the Search-report.pdf file. The upper panel shows the \(R\) factor as a function of the search progress (search generations). The lower panel reports the deviation of the structural parameters over time. Together, the two plots give an idea of the quality of the convergence of the optimization.

Fig. 18 First page of the Search-report.pdf file produced by viperleed.calc for the first, rough structure optimization for Cu(111)-(5 × √3)_rect-Te.

The initial reference calculation yields an \(R\) factor \(R_\mathrm{P} \approx 0.82\) since our starting configuration is very far from the ideal positions. Over the rough optimization, parameter values are shifted by up to 0.24 Å compared to the initial model (this is a lot!). The \(R\) factor drops to \(R_\mathrm{P} \approx 0.47\). This is still quite poor, but the progress is encouraging.

Remember to keep the best fit-structure by calling the bookkeeper utility with the --cont flag before proceeding:

$ python3 bookkeeper.py --cont #[or ./bookkeeper --cont]

The \(R\) factor resulting from an additional reference calculation, is decreased to \(R_\mathrm{P} \approx 0.33\) compared to the \(R_\mathrm{P} \approx 0.47\) value at the end of the optimization run — a big difference. This comes from the mentioned tensor-LEED error.

Fine DISPLACEMENTS and search parameters

Since in the former fit no parameter value reached the edge of the respective variation range, we should follow up by running a second, less coarse, optimization. For this stage, a choice of 10 pm range with 1 pm steps (0.1 Å range with 0.01 Å steps) should be reasonable (see Listing 18).

Listing 18 DISPLACEMENTS (10 pm range) for Cu(111)-(5 × √3)_rect-Te.

== SEARCH z
    = GEO_DELTA
        Te_surf   z        = -0.10 0.10 0.01
        Te_def    z        = -0.10 0.10 0.01
        Cu_surf   z        = -0.10 0.10 0.01
        Cu_def    z L(1-4) = -0.05 0.05 0.005

== SEARCH xy
    = GEO_DELTA
        Te_surf xy         = -0.10 0.10 0.01
        Te_def  xy         = -0.10 0.10 0.01
        Cu_surf xy         = -0.10 0.10 0.01
        Cu_def  xy L(1-4)  = -0.05 0.05 0.005

Hint

In preparing this example we found that the default search parameters lead to rather slow convergence in this and the following steps. To speed up the process, we recommend using these settings for SEARCH_START and SEARCH_CONVERGENCE (simply append the lines to PARAMETERS):

SEARCH_START = centered
SEARCH_CONVERGENCE gaussian = 0.05 0.5
SEARCH_CONVERGENCE dgen dec = 50 1.5

After the previous search steps, the parameters are close to the optimum. Therefore, the danger of getting trapped in a local \(R\)-factor minimum close to the starting position is low, and we can initialize the search at the previously determined values. The SEARCH_CONVERGENCE dgen parameter ensures that the search range shrinks more rapidly than with standard parameters.

The optimization on the 1 pm (0.01 Å) grid allows us to further reduce the \(R\) factor to about \(R_\mathrm{P} \approx 0.23\), which is again a good improvement on the previous value of \(R_\mathrm{P} \approx 0.33\).

Full-dynamic optimization

If we now visually compare experimental with calculated I(V) curves, we already notice a good qualitative agreement. However, we find that the peak widths in the experimental dataset seem to be consistently narrower than in our calculation. This is generally a sign that the imaginary part of the inner potential (\(V_{0\mathrm{i}}\)) of our calculation is off. To be more precise, \(V_{0\mathrm{i}}\) is likely too large, as larger \(V_{0\mathrm{i}}\) increases peak widths and smooths out the curves.

\(V_{0\mathrm{i}}\) strongly affects the \(R\) factor, but is hard to estimate for an unknown system. Since we did not specify a value for \(V_{0\mathrm{i}}\) in PARAMETERS, ViPErLEED took the default value of 4.5 eV (see V0_IMAG). The parameter \(V_{0\mathrm{i}}\) is not accessible in the tensor-LEED approximation, but we can use a full-dynamic optimization to find an optimal value. To do this, we set RUN = 6 and add this line to PARAMETERS:

OPTIMIZE V0i = step 0.5

Warning

Always make sure that the optimized value used for \(V_{0\mathrm{i}}\) is (i) reasonable (\(V_{0\mathrm{i}} \lesssim 7\) eV), (ii) forms an actual minimum, and (iii) gives qualitatively correct I(V) curves.

Since larger \(V_{0\mathrm{i}}\) values smooth out the I(V) curves, it is possible to deceptively decrease the \(R\) factor by arbitrarily increasing \(V_{0\mathrm{i}}\). This usually happens in cases of very high \(R\)-factor values, where \(V_{0\mathrm{i}}\) optimization is not the main concern anyhow.

Fig. 19 File FD_Optimization.pdf showing parabolic fit and minimum value for \(V_{0\mathrm{i}}\).

The optimized value for \(V_{0\mathrm{i}}\) will be output in the log file and automatically added to the PARAMETERS file for subsequent runs. Furthermore, ViPErLEED produces the files FD_Optimization.pdf (Fig. 19) and FD_Optimization_beams.pdf (Fig. 20) in the OUT directory. FD_Optimization_beams.pdf shows the calculated diffraction intensities for different values of the optimized parameter, while FD_Optimization.pdf shows the \(R\) factors corresponding to each trial value.

Fig. 20 Part of FD_Optimization_beams.pdf showing the effects of \(V_{0\mathrm{i}}\) on the (1|0) beam. Note that the \(V_{0\mathrm{i}}\) variation only leads to minor changes of the spectral appearance.

Refined structure fit

As usual, we can now perform some final structure fits over a fine-grained grid with subpicometre step. In particular, we should also optimize the vibrational amplitudes, which we have skipped so far. We recommend starting with the vibrational amplitudes here, since we have not touched them at all in the previous optimization step. Listing 19 shows the DISPLACEMENTS file for the refined structure optimization.

Listing 19 DISPLACEMENTS (0.5 pm step) for Cu(111)-(5 × √3)_rect-Te.

== SEARCH u
    = VIB_DELTA
        Te_surf   = -0.05 0.05 0.005
        Te_def    = -0.05 0.05 0.005
        Cu_surf   = -0.05 0.05 0.005

== SEARCH z
    = GEO_DELTA
        Te_surf       z = -0.05 0.05 0.005
        Te_def        z = -0.05 0.05 0.005
        Cu_surf       z = -0.05 0.05 0.005
        Cu_def L(1-4) z = -0.03 0.03 0.005

== SEARCH xy
    = GEO_DELTA
        Te_surf xy       = -0.05 0.05 0.005
        Te_def  xy       = -0.05 0.05 0.005
        Cu_surf xy       = -0.05 0.05 0.005
        Cu_def L(1-4) xy = -0.03 0.03 0.005

You may want to finish up with a last “fine tuning” of the vibrational amplitudes and (\(z\)) positions. See, for instance, DISPLACEMENTS_fine_2 in the input files, but feel free to play around with the setting yourself, to get a feeling for the available options. Altogether, this should bring us to a Pendry \(R\) factor \(R_\mathrm{P} \approx 0.19\), which is already a good agreement for such a heavily corrugated surface. For more details concerning corrugated surfaces see the discussion in Ref. 18.

Error calculation

We can also perform an error calculation for this system to gauge how sensitive our result is to minor changes of structural parameters. This can also give us an estimate of the uncertainty of the fit parameters. To run an error calculation, we need a DISPLACEMENTS file specifying the requested steps.

Listing 20 shows an example for displacements along the \(x\) direction (i.e., parallel to the \(\mathbf{a}\) unit-cell vector for the POSCAR file used here). Examples for \(y\), \(z\), and vibrational amplitudes are provided in the input files.

Listing 20 DISPLACEMENTS for error calculation in \(x\) direction.

== SEARCH error
    = GEO_DELTA
        Te_surf xy[1 0] = -0.10 0.10 0.004
        Te_def  xy[1 0] = -0.10 0.10 0.004
        Cu_surf xy[1 0] = -0.10 0.10 0.004
        Cu_def  xy[1 0] = -0.10 0.10 0.004

The results are plotted in the Errors.pdf file, shown in Fig. 21. We see that displacements for atoms in all layers have a drastic impact on the \(R\) factor. Thus, we can be fairly confident that all varied atoms are indeed present within the true surface structure.

Fig. 21 Page 1 of file Errors.pdf for (left) displacements in \(x\) direction and (right) changes of vibrational amplitudes.