Crystal structure refinement using Bloch-wave method for precession electron diffraction

ARTICLE IN PRESS

Ultramicroscopy 107 (2007) 474–482 www.elsevier.com/locate/ultramic

Crystal structure refinement using Bloch-wave method for precession electron diffraction A.P. Dudkaa,, A.S. Avilova,, S. Nicolopoulosb,c a

Institute of Crystallography of Russian Academy of Sciences, Leninsky prosprect 59, Moscow 119333, Russia b Universidad Politecnica de Valencia, Anenida de los Naranjos s/n 46022 Valencia, Spain c Nanomegas sprl, Blvd Edmond Machtens 79, B-1080 Brussels, Belgium Received 11 January 2006; accepted 28 March 2006

Abstract Procedure for crystal structure refinement using precession electron diffraction data and Bloch-wave method for accounting multibeam scattering is described. Refinement procedure takes into account features of precession geometry. Refining model consists of structural and reduced parameters determining dynamic diffraction. Difference between measured and calculated dynamic intensities of reflections is minimized with application of a nonlinear least squares method. As test example we used Si single nanocrystals. The influence of the reduced parameters on the quality of the obtained model is discussed. Refinement procedure is a part of ASTRA software. r 2006 Elsevier B.V. All rights reserved. Keywords: Precession electron diffraction; Bloch wave refinement; Multibeam scattering approximation

1. Introduction High-energy electron diffraction (hundred of keV) is the actual powerful tool for crystal structural determination. Due to strong interaction with scattering matter, this technique allows us to study very small objects: micro- and nanocrystals, thin films and surface layers. It is a rapidly developing technique with important applications for studying nanomaterials and a very useful tool for nanotechnology development. Present electron diffraction techniques can be divided into two groups. Techniques using electron diffraction cameras with wide plane-parallel electron beam belong to the first group. Actually those techniques became a basis of pioneer work on electron crystallography, originally done by Russian researchers [1]. Those techniques have received the name ‘‘electron diffraction structure analysis’’ (EDSA) and usually refer to all relevant research of atomic

Corresponding authors. Tel.: +7 495 1351020; fax: +7 495 1351011.

E-mail addresses: [email protected] (A.S. Avilov), [email protected] (S. Nicolopoulos). 0304-3991/$ - see front matter r 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.ultramic.2006.03.009

structure by methods of electron diffraction and electron microscopy. In this important early experiments performed by Russian scientists the samples that were studied were usually thin polycrystalline films having a different degree of crystalline orientations. The advantage of EDSA method consists that in the same experiment is possible to receive almost full three-dimensional set of reflections, sufficient for carrying out high-resolution structure determination, as orientation of crystals within the film is random and beam size of the diffraction cameras is considerable (size of order of several hundred microns). The appearance and common use today of transmission electron microscopes with improved electron optics characteristics allowing the use of small electron beam down to nanometer size, has opened the way for the study of the atomic structure of individual nanocrystals. Among leading diffraction techniques it is necessary to mention micro- and nanodiffraction techniques and convergent beam electron diffraction (CBED) for crystal symmetry and structure determination. In 1994 Vincent and Midgley [2] have suggested—as particular way to obtain diffraction patterns—to rotate the

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

electron beam along a conical surface having a common axis with the microscope optical axis. Such diffraction geometry allows to obtain reflections simultaneously from several reciprocal lattice planes (emerging from ZOLZ and even higher zone axis) in one diffraction pattern. With beam precession three-dimensional sets of reflections can be collected, and their number is increased with increasing precession angle. Through precession, as beam rotates through the optical axis, integrated diffraction intensities are obtained due to complete scanning of spots in reciprocal space by the Ewald sphere. This fact has some important consequences, among the most interesting— opens the opportunity for exact structural amplitudes calculation. Another important result related with beam precession is the weakening of many-beam effect interaction in the resulting diffraction pattern. As it is well known, the strong interaction of electrons with a crystal frequently results in a deviation from the kinematical law of diffraction. In the precession case, at diffraction from one microcrystal ˚ and more) multihaving important thickness (500–1000 A beam effects still remain, but those effects are averaged due to integration on all directions within Bragg width. It is also important to note, that when electron beam enters the crystal under a significant (precession) angle in relation with the optical axis, this results to an excitation of high-order reflections where their ‘‘extinction distance’’ is much more important in comparison with the low angle reflections. As a result, intensity of most reflections come closer to kinematical values and can be described by twobeam approximation. For the first time Vincent and Bird in 1984 [3] have paid attention to this fact. In particular for big crystal thickness and neglecting absorption effects, two-beam angle integrated diffraction intensity value does not depend on crystal thickness (as shown by Blackmann [4]). Other authors more recently [5] studied this question in detail and have shown (by using as example structure of La4Cu3MoO12 oxide) that intensity in a precession mode is close to two-beam approximation. According to that study as two beam approximation provides a simple model for simulate precession intensities, this method is more attractive for reliable structure determination in comparison with usual selected area diffraction—SAED technique. Therefore is possible to conclude as first approximation, that for routine structural determinations good results can be obtained from measured precession intensities either using kinematical theory or applying corrections based on two-beam approximation. It is very important to note that electron diffraction has a number of advantages over e.g. X-Ray diffraction, like high sensitivity at low scattering angles to valence electron distribution and higher accuracy in localization of light atoms in presence of heavy ones. Electrons are scattered on the electrostatic potential distribution, so the latter can be reconstructed by precise electron diffraction experiments as it was in Ref. [6]. But in order to refine such fine details in a

475

crystal structure, we need previously to extract exact kinematical intensity values (especially in case of thick crystals) taking in account multibeam dynamical diffraction effects. However, such multibeam calculation of intensities can be carried out only if an approximate model of the structure is known. As precession diffraction beam interacts with a single microcrystallite, is possible to apply many-beam theory for ˚ intensity calculations. In case of thin crystals (10–200 A) this has already been done with application of the multislice method in the work of Marks and Own [7], Stadelman and Oleynikov [8,9]. Using multislice approximation for intensity calculations has the advantage of speed calculations; however, in thick crystals case ˚ Bloch wave method is the most suitable (300–41000 A) for precise intensity calculations, as for instance was shown for CBED intensity profile refinements [10]. Available software that include many-beam intensity calculations in CBED can be divided in two different groups: The first group of existing software programs, refinement of kinematic amplitudes and phases of several electron diffraction reflections is performed in combination with precise structural amplitudes extracted from X-ray diffraction data [10]. In the approach used in Refs. [10,11] refinement of amplitudes and phases of all reflections present on CBED-pattern was successful only by considering sufficient number of experimental points in the CBED disks. The ‘‘experimental’’ kinematic structural factors thus received can be considered as the measurements corrected from many-beam diffraction effect, i.e. the reduced measurements. For the correct refinement it is necessary to have precise estimation of measurement errors on the reduced parameters. As has been shown in Ref. [12] such measurement errors should be evaluated statistically after repeating several measurements of those parameters; those error values should not be estimated through refinement, as their magnitude could be 10–100 times less in relation with statistical errors. Is important that experimental data have to be corrected from all type of errors (including possible systematic errors), otherwise all reduced parameters will be distorted after refinement. Nevertheless, if all experimental conditions (with possible source of errors) are well under control, is possible then that the final reduced parameters would be of high accuracy, as was shown in Ref. [12]. Other type of programs that are carrying out structural parameters refinement directly from dynamic experimental data [13,14] belong to second group. This approach is more promising as it allows to refine structure parameters directly from experimental data. Software described in Ref. [15] uses multislice method for structural parameters refinement. One important parameter that needs to be calculated during many-beam scattering refinement is the experimental crystal thickness. Therefore, the program for manybeam calculations should calculate it together with a number of other parameters determining conditions of

ARTICLE IN PRESS 476

A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

diffraction (geometry, etc.). Elements of dynamical matrices (DM) which are necessary for Bloch wave method also depend on the absorption in a sample, mainly from value of thermal diffuse scattering (TDS). The description of TDS effect on electron diffraction in crystal is given in Refs. [16,17]. The effect of isotropic thermal atomic vibrations can also be taken into account, for instance in ATOM software [18]. In Refs. [19–21] the effect of anisotropic thermal atomic vibrations is also considered. Thermal motion refinement in the TDS presence only requires modification of DM off diagonal elements. Therefore we may conclude that all refined parameters can be divided into three groups: (1) structural ones (coordinates, thermal, multipole), (2) amplitudes and phases of reflections, and (3) reduced parameters (thickness of a sample, geometrical and other experimental parameters). Refinement quality depends on the complexity of the refined model. One of the most important software features is the use of minimization method. The character of dependence of dynamical intensities on reduced parameters is neither linear nor monotonic (for example thickness change causes rocking of the minimizing functional and unit cell orientation change results in zigzag shape dependence) so the latter can present several sets of local minima. The problem can be solved by use of a method of global minimization [22], but computing time can considerably increase in comparison with use of traditional least squares (LSQ) method minimization. Other method (a local simplex method) capable to ignore presence of local minima has been used in Ref. [23] but also calculations with that method are slow. In Ref. [24] it was proposed to use the steepest-descents method (step-by-step descent) for minimization but results of application are not studied in extensively. To keep high speed of convergence of LSQ method and to overcome its weak ability to ignore local minima, it was proposed to use a quasi-Newton method in Refs. [25,26] and nonlinear LSQ in program MBFIT [13]. The most reasonable strategy seems to divide refining parameters into two groups, depending on the degree of smoothness of functional corresponding to them and to carry out refinement by different methods. Structural parameters are quickly and reliably refined using linear LSQ method. Reduced parameters can be refined by a simplex—method. In program EXTAL [27] such division of parameters was carried out, where the simplex—method is used only for refinement. Is important to stress, that all such methods described above were already used for CBED calculations. However, to our knowledge no such refinement procedure has ever been carried out using precession electron diffraction data. The following contribution is devoted to development of algorithms and refinement procedure for several crystal structure parameters (including chemical bonding and electrostatic potential) using precession electron diffraction data [28] and making use of the new dedicated software ASTRA [29]. In this work we shall study the refinement

specificity caused by use of precession technique. We will study the influence of different geometrical factors, like the use of parallel beam and the influence of the complete set of reduced parameters. The software we used (ASTRA) can refine full structural model by minimizing the differences between experimental and calculated dynamic integrated intensities for all measured reflections by using nonlinear LSQ. In general this program is different from programs described previously for dynamical CBED refinements, being quite similar with usual single crystal X-rays and neutron diffraction refinement programs. 2. Algorithm of calculations In precession electron diffraction the measured dynamic intensity of reflections is equal to the total integrated intensity resulting from continuous rotation of the primary beam. When simulating precession, calculations must proceed with many incident beam tilts, which will cover the full precession circle ð360 Þ. For our calculations the full circle is divided into a number of fixed azimuthal angular steps, where each step correspond to a fixed position of primary beam (frames). For example in our work for silicon (see below) the number of frames was set from 300 to 700. The algorithm of the refinement consists of eight basic steps (Fig. 1): (1) Assuming the structural model is approximately known kinematical structural factors U g are calculated; (2) given a chosen selected frame geometrical characteristics of reflections ‘‘excitation errors’’ are calculated; (3) further the dynamic matrix (DM) is formed using ‘‘excitation errors’’, calculated kinematic structural factors and some other (reduced) parameters; (4) using standard subroutines (see for example Ref. [30]) eigenvalues and eigenvectors of DM are calculated; (5) using current value for crystal thickness, eigenvectors and eigenvalues of DM dynamic intensities are calculated for all nonzero dynamic reflection intensities; (6) a summation of calculated dynamic intensities coming from several incident tilts (from different frames) is performed and their comparison with the measured experimental values is made; (7) during refinement, further differences ‘‘model–experiment’’ are minimized and new values of structural and reduced parameters are determined; (8) steps 1–7 are repeated iteratively in order to have optimal goodness of fit between experimental and theoretical values. The calculation program of manybeam intensities is a part of the software package ASTRA [29] which can specify structural models of crystals for any symmetry and structure. In this refinement not only positional parameters can be calculated, but is also possible to take into account thermal harmonic and anharmonic parameters up to sixth order

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

477

χ

z Rotation on ϕ

zone n

y

x

Primary beam Fig. 2. Description of coordinates for precession diffraction experiment.

Fig. 1. Flow chart for R minimization using Bloch-wave approximation.

related to atomic vibrations; again other atomic parameters relative to multipole models (up to eighth order) of Stewart–Hansen–Coppens [31], merohedric micro-twinning and extinction effect can also be refined. For calculation of the atomic factor of electron scattering the Mott-formula was used. 3. Geometrical considerations of the precession mode The results of geometrical calculations are used for diagonal elements calculation in the DM. In ASTRA, these calculations are based on the concept of a matrix of parameters and orientation UB [32]. Matrix UB is calculated starting from the orientation of a zone axis. Let us define a laboratory system of coordinates as follows: z-axis along the precession axis, x-axis to the right and y-axis forming the right-hand triple and pointing perpendicularly into the drawing plane from the observer (Fig. 2). The origin of coordinates is fixed on the top surface of the sample. Without restriction of a generality, we shall consider that the z-axis is directed vertically upwards. These coordinates coincide with ones in Ref. [10]. We shall define w0 as the beam precession angle and we shall designate as j the beam rotation angle (angle between x-axis and the beam projection onto plane xOy). Each value of the angle j has its own associated frame (set of excited reflections) and its corresponding DM. The positive direction of rotation for angle j is assigned when turning

counterclockwise from x- to y-axis and looking from a positive direction of the z-axis. The positive direction of rotation for angle w0 is from the positive direction of the z-axis onto the plane xOy. The usual set-up for calculations is fixed with a crystal zone axis ðhklÞzone selected close to the z-axis. The zone axis index is determined during the calculation and is set as input in the ASTRA program. When zone axis and elementary cell parameters are known, it is possible to calculate a matrix of parameters and orientations U 0 B, where B is the parameter matrix of the reciprocal cell, and U 0 is the preliminary orientation matrix of the reciprocal cell. Let the crystal have direct cell parameters a; b; c, a, b, g and reciprocal parameters a ; b ; c , a , b , g . Then: 0  1 a b cos g c cos b  B  c sin b cos a C B ¼ @ 0 b sin g A. 0 0 1=c The coordinates of a zone axis with origin (0 0 0) are equal to: ðxyzÞTzone ¼ B  ðhklÞTzone . Top index T specifies that the vector column is used. Let jz and wz be the Euler angles when the cell is away from the initial position, and the zone axis will coincide with the z-axis, then: 1 0 10 1 0 0 cos jz  sin jz 0 B CB cos jz 0 C U 0 ¼ XF ¼ @ 0 cos wz  sin wz A@ sin jz A. 0 sin wz cos wz 0 0 1 Quantities jz and wz are determined by ratios: tgjz ¼ xzone =yzone ;

tgwz ¼ ½ðx2zone þ y2zone Þ1=2 =zzone .

Rotation around the third axis can be determined from the two-dimensional diffraction pattern. However, with big enough number of frames it is not needed to be calculated.

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

478

The orientation error in one step j can be specified as a reduced parameter, j0 . If a zone axis does not coincide exactly with the z-axis of laboratory coordinates (precession axis), the errors of an experimental setup of crystal (angles ax , ay , az ) can be specified as the parameters of the corresponding Euler matrixes M x , M y , M z . Therefore, the specified matrix UB will be equal: UB ¼ M z M y M x  U 0 B ¼ f ðax ; ay ; az ; ðhklÞzone ; a; b; c; a; b; gÞ.

The orthogonal coordinates of the measured reflection ðhklÞ expressed in the laboratory coordinates with origin (0 0 0) of the reciprocal cell are defined by the relation g ¼ ðxyzÞT ¼ UB  ðhklÞT . The origin of the primary beam vector, K, moves around a circle and, for each rotation angle jpri (for each frame), we have K ¼ ðxyzÞTpri ¼ X Fð0; 0; 1ÞT ¼ ðjKj sin w0 cos jpri ; jKj sin w0 sin jpri ; jKj cos w0 ÞT . The sign ‘‘’’ has appeared because the vector K is approximately antiparallel to the z-axis. Angle w0 differs from the precession angle w0 because of refraction on vacuum—crystal boundary. If the normal to a crystal plane coincides with the z-axis, then the w0 angle is determined by the ratio (analogue of Snell’s formula): sin w0 ¼ sin w0 = ½1 þ l2 U 0 1=2 , where U 0 is the average inner potential ˚ 2 , and l is the relativistic wavelength of a crystal in A dependent on the value of the accelerating voltage of the microscope E 0 . Furthermore, with respect to the origin of our laboratory coordinates, the coordinates for reflection g are K þ g ¼ ðxyzÞT þ ðxyzÞTpri . Hence, the expression necessary for calculation of a diagonal element of a DM in line number g is equal to: 2KS g ¼ K 2  ðK þ gÞ2 ¼ f ðE 0 ; U 0 ; ax ; ay ; az ; j0 ; w0 ; ðhklÞzone ; a; b; c; a; b; g; ðhklÞ; jÞ, where Sg is the excitation error of a reflection. In the functional dependence f ð. . .Þ seven refining parameters, nine input parameters, the reflection indices and the angle j for an individual frame are listed. Evaluation of DM eigenvalues and eigenvectors is determined by boundary conditions. In a simple case, the normal to a surface has coordinates (0,0,1). In some cases, the normal to a surface of a crystal can deviate from direction of corresponding node line of reciprocal cell. It is possible to consider errors of orientation of a sample (angles dx , dy , dz ) as refining parameters and to process them in the same way as has been shown above for orientation error evaluation of an elementary cell.

The specified coordinates of a normal will be equal: N ¼ ðxyzÞTN ¼ M z ðdz Þ  M y ðdy Þ  M x ðdx Þ  ð0; 0; 1ÞT , where M z , M y , M x is the Euler matrixes. For the calculation of dynamic intensities (with known eigenvectors and eigenvalues of DM for each frame) is necessary to know the value K n ¼ ðK  NÞ ¼ f ðE 0 ; U 0 ; dx ; dy ; dz ; j0 ; w0 ; jÞ. Crystal plate thickness is the most important refining reduced parameter. If we assume that the plate is not exactly plane-parallel, the problem of DM diagonalization becomes considerably complicated because of the complexity of boundary conditions. However, in this situation it is possible to apply the phenomenological approach similar to one used in work [33] for description of anisotropic extinction. We allow a crystal to have the wedge-shaped form with the oblique bottom surface, where the orientation of a wedge is unknown. Then after rotation of a primary beam for each frame we shall obtain a different path length of beam into the crystal, simulating symmetry of an ellipse. Thickness variation is depending on coordinates ðx; yÞ of primary beam and can be simulated by square-law formula having three refining parameters. Then the thickness of a plate for some private direction jpri will be determined by scalar production of two threecomponential vectors ðtx ; ty ; txy Þnðx2pri ; y2pri ; 2xpri ypri Þ: t ¼ tx  x2pri þ ty  y2pri þ 2  txy  xpri ypri , Coordinates of a primary beam are determined above. The last refining parameter R in the given version of the program is the average factor of absorption m0 of a crystal: R ¼ ðtx ; ty ; txy ; m0 ; E 0 ; U 0 ; ax ; ay ; az ; dx ; dy ; dz ; j0 ; w0 Þ. Therefore, the general number of refining parameters is equal 12 for a plane-parallel plate and 14 for wedge-shaped one. This way of calculations for geometrical characteristics determination (for precession experiment) allows to use matrix operations in a modern Fortran, and it results fast enough for calculations. 4. DM formation One important question that determines the quality of adjustment of calculated dynamic to experimental intensities, is the way of reflection selection into a DM. Previous calculations with CBED data have studied this question in a number of works [23,27]. Time of calculations also depends on the size and amount of the generated DM. For the precession technique the number of DM is equal to the number of frames. In the precession experiment it should be taken in account into the refinement process all experimental reflection intensities and neglecting ‘‘forbidden’’ reflections (for which measured intensity are distinct from zero), otherwise it is impossible to receive good agreement between experimental and theoretical values. For DM

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

formation it is necessary to select reflections with small excitation errors 2KS g . Other criteria related with the quantity U g =2KS g used for CBED line profile dynamical calculations are not applicable for this purpose. It is known [10], that application of Bethe potentials allows to reduce the dimension of DM without decrease calculations accuracy. In a method of Bethe potentials DM will consist only of strongly excited reflections having big values U g =2KS g , where all elements are corrected from influence of other poorly excited reflections. It is interesting to note that precession diffraction does not lead to DM dimension decrease. In fact, forbidden reflections do appear as poorly excited where U g =2KS g 0, and they will not take part in formation of DM. On the contrary, attempts to include them in calculations lead to incorrect DM diagonalization. During the refinement procedure, is possible to include in DM reflections which participated in dynamic redistribution of scattered energy, but have not been measured experimentally for various reasons. This can raise accuracy of our refinement, though not significantly. Doing the contrary, (e.g. by decreasing the DM dimension by considering in refinement only a limited number of total observed reflections) worsens considerably the goodness of fit. Therefore, in our refinements we did not reduced the DM dimension in comparison with the number of the observed reflections in order not to loose accuracy of results. It is interesting to observe that is possible to vary quite widely the number of different angular frames into the DM and refinement quality can still remain within reasonable limits. Minimal number of frames which was used in the present work with Si was set equal to 360, corresponding to an angular step of j ¼ 1 . A small decrease in R-factors is only observed when the number of frames is doubled (angular step of 0:5 Þ. Further reducing angular step does not have any influence on refinement. In general off-diagonal elements of DM should take into account correction related with TDS (see. previous chapter). In our software ASTRA this correction is not carried out but should be made in the future. 5. Peculiarity of refinement by nonlinear LSQ In ASTRA software, during refinement of structure parameters taking in account many-beam electron scattering, the following functional is minimized: FðModel; k; t; m0 ; RÞ ¼

N obs X g¼1

wg

I dyn g obs k  cg



!2 I dyn g calc ðModel; t; m0 ; RÞ

! min ,

where N obs is the number of the measured not averaged reflections, including all forbidden reflections in kinematical approximation; I obs and w – the measured (dynamic) intensity and its statistical weight; k – the scaling factor; c

479

the possible corrections relative to geometrical and other experimental conditions; I calc – the calculated dynamic intensity; Model – the parameters of structural model; and R – the stand for reduced parameters, among them thickness t and average absorption factor m0 . The corresponding R factors should be calculated during refinement: R2w ðIÞ  R2w ðF 2 Þ " , #1=2 N obs N obs X X dyn dyn 2 dyn 2 ¼ wg ðI g obs  I g calc Þ wg ðI g obs Þ g¼1

2

R1ðIÞ  R1ðF Þ ¼

g¼1

N obs X g¼1

, jI dyn g obs



I dyn g calc j

N obs X

I dyn g obs .

g¼1

Same expressions forR-factors on modules are obtained by replacing intensities by structure factors e.g. I by F. In ASTRA the algorithm of adapted nonlinear minimization is used on the basis of Lavrent’ev–Levenberg–Marquardt [34,35] regularization. This algorithm is proved to give stable refinements using real experimental data where exists high correlation between parameters [36, p. 600]. As it is known, LSQ refinement is effective enough with refinements in a vicinity of local (or global) minimum. In our refinement, at least two groups of parameters affect dramatically smoothness of minimized functional; dependence of dynamic intensity on thickness of a plate crystal and zone axis orientation has a complex character. If refinement procedure begins in an any point and if correlations between parameters are strong enough, it is possible for refined structure model parameters to shift excessively towards wrong direction, and as result refinement may stop or may not reach a global minimum. In nonlinear LSQ (used in ASTRA), those problems are essentially solved: by a choice of a direction of search and by correction of distorting influence of correlations on a degree of conditionality of matrixes. Firstly, the second derivatives of calculated structure factor are taken into account during decomposition of minimized function. Secondly, Hessian matrix is calculated by modified iterative Davidon–Fletcher–Powell method which does not demand any matrix inverting. Thirdly, to avoid unpleasant consequences in case of possible Hessian singularity, method of Lavrent’ev–Levenberg–Marquardt is applied. Moreover, to avoid slow methods of global minimization, scanning of parameters is applied in ASTRA. Therefore it is possible to build one, two and threedimensional sections of functional having a set of local minima. In those program versions matrices of first and second derivatives of minimizing function are calculated by the finite differences method. Those algorithms may slow down slightly the calculation speed, since DM diagonalization is time consumable.

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

480

6. An experimental test for silicon crystal As preliminary crystal structure model we used high accuracy X-ray experimental data from Ref. [37]. Experimental set of X-ray data consisted of 21 reflections (including forbidden reflections 222 and 422) symmetry averaged with standard errors. Refinement parameters for the structure model (S. G. Fd-3 m) included one thermal parameter B, one parameter of spherical deformation of valence shells k [31] and the scale factor. Atomic scattering curves have been taken from [38] considering average valence shell ð3s2 þ 3p2 Þ=4 for Silicon atom. Structural model from X-ray data is shown in Table 1. At rejection of two forbidden reflections considerably decrease R-factors, but this practically does not affect final results. Electron precession experiment has been carried out in Free University of Brussels (ULB) on a Philips TEM model CM20 working at 200 kV and equipped with special precession unit Spinning Star (Nanomegas Company) [39]. The size of the specimen area illuminated by the electron beam in the precession mode is estimated to be approximately 200 nm. Precession electron diffraction from Si [1 1 0] ZA has been registered by photoplate (Fig. 3). Precession angle was set at 2:9 . Integrated intensity for the precession experiment have been estimated by using ImgView-software from EDSA software package [40]. Electron experimental data contained 152 not averaged reflections (including 33 forbidden reflections). The number of independent reflections in a set was equal to 52 (with 12 forbidden among them). R-factors of the equivalent reflections were R1ave ðIÞ ¼ 9:03%; R2wave ðIÞ ¼ 8:47%. It should be noted that by averaging equivalent Friedel reflections that there is still a quite important intensity difference between Fridel’s pairs: R1ave ðIÞ ¼ 8:8%. A refinement of structure model in kinematical approximation has been carried out (Table 2). Results reveal low relative accuracy (high R-factors) and show quite clearly that in kinematic approximation is impossible to refine any acceptable values of thermal parameterB and parameter k describing a spherical deformation of the valence shells (compare variants 1 and 2, Table 2). Rejection of forbidden reflection from the data set considerably lowers R-factors values (quite expected result as intensities of forbidden reflections have no satisfactory description in kinematical approximation).

Table 1 Refinement results for Si based on precise X-ray diffraction data N obs

B; A2

k

R1 ðIÞ; %

R2w ðIÞ; %

R1 ðF Þ; %

R2w ðF Þ; %

21 19

0.46 0.457

1.07 1.07

1.687 1.680

1.251 0.471

1.027 0.824

1.183 0.234

N obs —stands for the number of observed experimental reflections (including forbidden reflections—21 and excluding them—19).

Fig. 3. Precession electron diffraction pattern for Si crystal, plane (1 1 0).

A refinement in dynamical approximation has led to significant decrease of R-factors, both considering full data set, and on reflections that do not include the forbidden ones (Table 2). DM has been formed considering: ð2KSg Þmin ¼ ˚ 1 and (U g =2KSg Þ 22 A max ¼ 0:0001, where DM included all measured reflections; angular step j was set equal to 1 (360 steps/frames and DM-s); beam divergence of primary beam is considered as Dw ¼ 0 (a beam is considered as ideally collimated); in fact, in general case a primary beam can have a small angle divergence. Considering spread of the primary beam into the dynamical refinement is done by using additional DM; however, the inclusion of this factor does lead to any significant improvement of the refinement. Prior to dynamical refinement, structural model has been fixed according to the X-ray data, and DM were constructed on the basis of a total 152 measured reflections. The reduced parameters also have been fixed (see above), and the thickness has been refined (Fig. 4); exact zone axis orientation has been refined. The addition of 32 reflections (close to Ewald sphere, but not seen on a photo-plate) did not improve significantly the results. At further refinement step reduced parameters have been refined considering full-matrix refinement and transition to a wedge-shaped plate. At the last stage full-matrix refinement of all parameters (including structural ones) was made. Structure parameters were refined step by step: first were refined thermal parameter and then kappa, anharmonic parameters and extinction parameter [41]. Thus, refinement for dynamic scattering has been performed in a big number of different steps, where three of them are shown in Table 2. Initial R values refinement were: R1ðIÞ=R2w ðIÞ26=26% for a full set of 152 reflections including the forbidden ones.

ARTICLE IN PRESS A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

481

Table 2 Refinement results in both kinematical and dynamical approximation from precession electron diffraction data B; A2

k

N obs

R1 ðIÞ; %

R2w ðIÞ; %

R1 ðF Þ; %

R2w ðF Þ; %

Kinematical approximation 1* 1 2 3

0.46 1.95

1.07 5.21

152 152

65.130 32.453

74.599 37.560

43.920 31.769

49.781 41.506

Dynamical approximation 11 3*

0.46

1.07

4

13

0.57

1.22

5

17

0.56

1.21

152 119 152 119 152 119

18.220 17.813 17.097 16.408 14.164 12.472

19.761 20.024 17.281 17.169 12.475 11.582

12.349 11.279 12.138 10.890 11.264 9.697

12.681 11.271 12.397 10.620 11.005 9.140

No.

N par

N par —number of refined parameters. For steps 3–5 refinement is made with full set of data (upper line, N obs ¼ 152 reflections). (Lower line, N obs ¼ 119 reflections) R-factors were calculated without accounting for forbidden reflections * B and k—parameters were fixed in values, obtained from precise X-ray experiment.

Fig. 4. Crystal thickness refinement for Silicon.

During the third step of refinement (see Table 2) structural parameters were fixed; dimensions of DM was increased up to 184. During this step thickness, zone axis orientation, exact precession angle and other parameters were refined having as factors R1ðIÞ=R2w ðIÞ18=20%. Further refinement (fourth step in Table 2) relative to structural and other parameters (including anisotropic thickness of a plate) was made, improving reliability factors R1ðIÞ=R2w ðIÞ ¼ 17:0=17:3%. During the final refinement stage, a full-matrix refinement has been performed including all parameters described before like structural ones, anharmonic thermal vibrations of thirdand fourth-order and secondary extinction parameter as described in [41], R1ðIÞ=R2w ðIÞ ¼ 14:2=12:5%. After final refinement stage the following values were obtained: ˚ precession thickness of wedge-shaped plate 841–1081 A,  angle 2:87 , zone axis orientation estimated error 0:39 , 0:01 and 0:082 , thermal parameter (Debye–Waller factor) ˚ , k-parameter—1:21. for Si—0:56 A 7. Discussion Significant decrease in reliability factors for Si refinement after transition from kinematical to dynamic approximation proves, that dynamical effects in Si are quite important.

Another conclusion of our work is that after use of a Blochwave method, structural parameters refined values are close to the values obtained from X-Ray data. For the adequate description of many-beam effects it is necessary to carry out refinement of set of reduced (mainly, geometrical) parameters. Among reduced parameters that are refined, sample thickness and zone axis orientation have a great influence on R-factors. Dependence of minimizing function (R-factors) on sample thickness presents several local minima, although it is smoother in comparison with similar dependence of zone axis orientation. Obtaining correctly refined orientation of a zone axis with the help of LSQ algorithm seems difficult. In our case we have obtained zone axis orientation by scanning corresponding parameters. The influence on the refinement from transition to a plate having anisotropic thickness is very significant. However, the corresponding analytical dependence of dynamical intensities is such, that within the framework of one experiment is impossible to separate completely thickness variation effect and effect of a deviation from a normal to a plate from the precession axis. Application of dynamic calculations ‘‘brings together’’ symmetrical ‘‘corrected’’ intensity of equivalent reflections. The R factor after averaging such intensities can be estimated in range 6.5–7% against 9.03% in kinematical approximation. However, a significant part of intensity difference between Friedel’s pair reflections still remains not compensated also for dynamical approximation case. This effect is still under investigation and at present we have not any complete satisfactory explanation. The dimension of DM in precession technique is determined mainly by the number of the measured reflections. Addition into the data set of other reflections close to Ewald sphere (to an amount of 10–20% in relation with the overall reflections) leads to R-factors decrease, although not significantly. Taking into account beam divergence also raises relative accuracy of our results, although not significally and at the expense of computing time.

ARTICLE IN PRESS 482

A.P. Dudka et al. / Ultramicroscopy 107 (2007) 474–482

In those cases where a precession angle is not large and a diffraction pattern retains approximately its symmetry (Ewald sphere is precessing in the limits of Bragg width of each reflection from DM), calculation time using Bloch wave method can be reduced considerably due to crystal symmetry consideration as was shown in Ref. [42]. Is interesting to observe that in our work, just as in Ref. [12], extinction parameter has been used. At present is impossible to have a strict theoretical account of all experimental features that may influence our experiment, therefore for compensation of differences between experimental and calculated dynamical intensities refinement of phenomenological extinction parameter has been used. Differences from X-ray refinement: Refined bonding effects after introduction to k-model are more appreciable than in X-ray case, despite the quite high value of R-factors. Is important to comment that in X-ray refinement of anharmonic movement parameter of Si atom is impossible to perform, although with electron diffraction experiment this quantity is rather significant and can be refined successfully as it has been shown in our calculations. 8. Conclusion In our work is the first time that is described advanced dynamical refinement by Bloch wave calculations using precession electron diffraction intensities. Bloch wave method is applied for accounting many beam scattering; the final refined structural parameters have values that are close to the values refined from precise X-ray diffraction data. Using multibeam scattering approximation in the refinement, reduces considerably reliability factors in comparison with using only kinematic approximation. It should be noted, that the used experimental data have not allowed us to show the full opportunity of the advanced technique. However, the proof that classical LSQ structural refinement program can be used for such calculations with precession electron diffraction data, will certainly stimulate further development in EDSA with increased accuracy of experimental data. Acknowledgements Authors thank M. Nikolsky and A. Kuligin for obtention and treatment of experimental data by precession technique, J.C.H. Spence and K. Tsuda for important and useful comments to our manuscript. References [1] B.K. Vainshtein, B. Zvyagin, A.S. Avilov, in: Electron Diffraction Techniques, vol. 1, Oxford University Press, Oxford, 1992. pp. 216–312 (Chapter 6).

[2] R. Vincent, P.A. Midgley, Ultramicroscopy 53 (1994) 271. [3] R. Vincent, D.M. Bird, Philos. Mag. Lett. 53 (1984) L35. [4] J.C.H. Spence, High Resolution Electron Microscopy, Oxford University Press, Oxford, 2003. [5] C.S. Own, A.K. Subramanian, L.D. Marks, Microsc. Microanal. 10 (2004) 96. [6] A.S. Avilov, Z. Kristallogr. 218 (2003) 247. [7] C.S. Own, System design and verification of the precision electron diffraction technique, Ph.D. Dissertation, Northwestern University, 2005. [8] P. Oleynikov, Exploiting reciprocal space, Ph.D. Dissertation, Stockholm University, 2006. [9] P. Stadelman, JEMS—ENS CIME EPFL, Lausanne, 2005. [10] J.C.H. Spence, J.M. Zuo, Electron Microdiffraction, Plenum Press, New York, 1992. [11] K. Tsuda, M. Tanaka, Acta Crystallogr. A 51 (1995) 7. [12] B. Jiang, J.M. Zuo, Q. Chen, J.C.H. Spence, Acta Crystallogr. A 58 (2002) 4. [13] K. Tsuda, M. Tanaka, Acta Crystallogr. A 55 (1999) 939. [14] Y. Ogata, K. Tsuda, Y. Akishige, M. Tanaka, Acta Crystallogr. A 60 (2004) 552. [15] J. Jansen, D. Tang, H.W. Zandbergen, H. Schenk, Acta Crystallogr. A 54 (1998) 91. [16] H. Yoshioka, J. Phys. Soc. Japan 12 (1957) 618. [17] C.R. Hall, F.R.S. Hirsh, R. Soc. London A 286 (1964) 158. [18] D.M. Bird, Q.A. King, Acta Crystallogr. A 46 (1990) 202. [19] A.L. Weickenmeier, H. Kohl, Acta Crystallogr. A 47 (1991) 590. [20] G.R. Anstis, Acta Crystallogr. A 52 (1996) 3. [21] L.M. Peng, Acta Crystallogr. A 54 (1998) 481. [22] C. Deininger, G. Necker, J. Mayer, Ultramicroscopy 54 (1994) 15. [23] J.M. Zuo, J.C.H. Spence, Ultramicroscopy 35 (1991) 185. [24] K. Martiensen, R. Hoier, L.N. Bakken, in: L.D. Peachy, D.B. Williams (Eds.), Proceedings of the XIIth ICEM, vol. 2, San Francisco Press, USA, 1990, pp. 492–493. [25] B.M. Bird, M. Saunders, Acta Crystallogr. A 48 (1992) 555. [26] M. Saunders, B.M. Bird, N.J. Zaluzec, Ultramicroscopy 60 (1995) 311. [27] J.M. Zuo, Mater. Trans. Jpn. Inst. Met. 39 (1998) 938. [28] A.P. Dudka, A.S. Avilov, M.S. Nikolsky, Fifth National Conference on Applying X-rays, Synchrotron Irradiations, Neutrons and Electrons for the Investigation of Nanomaterials and Nanosystems, RSNE-NANO—2005, Moscow, 2005, Abstracts, p. 242. [29] A.P. Dudka, Crystallogr. Rep. 47 (2002) 152. [30] IMSL Fortran 90 MP Library, Visual Numerics, Inc. 1999. [31] N.K. Hansen, P. Coppens, Acta Crystallogr. A 34 (1978) 909. [32] W.R. Busing, H.A. Levy, Acta Crystallogr. A 22 (1967) 457. [33] P.J. Becker, P. Coppens, Acta Crystallogr. A 31 (1975) 417. [34] J.E. Dennis, D.M. Gay, R.E. Welsch, ACM Trans. Math. Softw. 7 (1981) 348. [35] J.E. Dennis, D.M. Gay, R.E. Welsch, ACM Trans. Math. Softw. 7 (1981) 369. [36] International Tables of Crystallography, Kluwer Academic Publishers, Dordrecht, 1992, p. 600. [37] Y.W. Yang, P. Coppens, Solid State Commun. 15 (1974) 1555. [38] Z. Su, P. Coppens, Acta Crystallogr. A 54 (1998) 646. [39] T. Weirich, J. Portillo, G. Cox, H. Hibst, S. Nicolopoulos, Ultramicroscopy 106 (2006) 164. [40] A.K. Kuligin, A.S. Avilov, XIV Russian symposium on scanning electron microscopy and analytical methods of study of solids. SEM2005, Chernogolovka, 2005, Abstracts, p. 272. [41] P.J. Becker, P. Coppens, Acta Crystallogr. A 30 (1974) 129. [42] E. Kogiso, J. Phys. Soc. Japan 42 (1977) 223.