- Email: [email protected]
A procedure to obtain reliable pair distribution functions of non-crystalline materials from diffraction data
NUCLEAR
INSTRUMENTS
AND M E T H O D S
143 ( 1 9 7 7 )
569-575
; ©
NORTH-HOLLAND
PUBLISHING
CO.
A P R O C E D U R E TO OBTAIN RELIABLE PAIR D I S T R I B U T I O N F U N C T I O N S OF N O N - C R Y S T A L L I N E MATERIALS F R O M D I F F R A C T I O N DATA F. YSSING HANSEN
Fysisk-Kemisk Institut, The Technical University o.I Denmark, DK 2800 Lyngby, Denmark and KIM CARNEIRO
Physics Laboratory 1, H.C. Orsted Institute, Universitetsparken 5, DK 2100 Copenhagen 0, Denmark Received 10 January 1977 A simple numerical method, which unifies the calculation of structure factors from X-ray or neutron diffraction data with the calculation of reliable pair distribution functions, is described. The objective of the method is to eliminate systematic errors in the normalizations and corrections of the intensity data, and to provide measures for elimination of truncation errors without loosing information about the structure. This is done through an iterative procedure, which is easy to program for computers. The applications to amorphous selenium and diatomic liquids are briefly reviewed.
1. Introduction
The determination of pair distribution functions from structure factor data provided either by X-ray or neutron diffraction experiments is sensitive both to proper corrections of the experimental results as well as to proper Fourier transformation of the structure factor. The two problems may at first sight seem uncorrelated, but it is very difficult to judge from the structure factor alone, whether or not the instrumental corrections, including normalization, of the diffraction data have been performed consistently. The consistency of the data treatment however is clearly revealed, when the structure factor data are Fourier transformed to give the pair distribution function. The influence of improper normalization and corrections of the diffraction data on the pair distribution function has been discussed several times (see e.g. refs. 1 and 2), and the conclusion is that all types of errors in the data treatment result in strong unphysical oscillations in the pair distributions function at small r, where the pair distribution function is known to be equal to zero due to the hard core repulsion between atoms. Therefore, this region may be used to see, if the data treatment has been performed properly. It is evident that an erroneous data treatment of the experimental results also demolishes the pair distribution function at larger r, but in this region one has no way to separate out "false" contributions from " t r u e " contributions. The Fourier transformation of the structure fac-
tor may also introduce errors in the pair distribution function. It is well known 3,4) that the Fourier transform of a function, which is not known until convergence to zero, contains artificially broadened peaks and around those, also false satellites. These satellites are mainly concentrated around the main peaks of the transformed function. In typical Xray and thermal neutron experiments, the maxim u m wavevector transfer /£m is 9--13 A - 1 . At this wavevector the structure factor of almost all non-crystalline materials still shows significant features, which indicate that the integrand of the Fourier integral has not yet converged to zero. A direct transform of the incomplete structure factors will therefore always give false features in the pair distribution function and smear out finer details. Reported pair distribution functions in the literature often show strong oscillations at small r, partly due to truncation errors in the Fourier transform and partly to improper treatment of the diffraction data. In order to remove the false satellites from the pair distribution function it has been common practice to apply convergence factors or window functions 4,s) to the measured structure factor to make the integrand in the Fourier integral converge. Such a procedure removes the false peaks, but this advantage is achieved at the expense of poorer resolution and still leaves the resulting pair distribution function quantitatively inadequate. An extension of the experimental region to higher/£m values should therefore make it possible
570
F.
YSSING
HANSEN
to avoid truncation errors in the final pair distribution function. The recent occurrence of experimental facilities in the form of either accelerators or reactor hot sources providing epithermal neutrons has, in fact, made such extensions feasible, However, the pioneering LINAC measurements on diatomic liquids 6-9) have shown that the problems in deriving reliable pair distribution functions are not solved by going to large Km'S alone. As truncation errors should be eliminated, the deficiciences of the results are likely due to improper normalization and corrections of the diffraction data and in particular, as shown by Page and Powles9), the Placzek correction plays an important role in some of these experiments. This clearly stresses the dilemma in structure studies of non-crystalline materials. On the one hand an accurate Fourier transformation requires that K"m is large, but on the other hand, since at large wavevectors the Placzek correction becomes progressively more problematic, Km should not be too large. In this paper we will show, how it is possible to choose a Km for which both requirements above are met in the sense that although none of the two problems is eliminated, both of them are kept at a manageable level. The same problems have also been discussed recently by Konnert et al. 12) and Rao et al. 2°) and their methods together with our method may be taken as a natural extension of the work of Kaplow e t al.l'19). 2. Theory The relation between intensity data I(2, K) in an X-ray diffraction experiment and the structure factor of a non-crystalline sample with random orientations may be written 4) I(2, K) = A(2) T(2, tc) [f2(~c) S(K)+ s(tc)],
(l)
AND
KIM
CARNEIRO
both 2 and K and is easy to determine, if the total scattering length is known, and if the shape of the sample is simple, like a spere, a cylinder or a platel4). For other shapes one may have to consider Monte Carlo simulations in order to evaluate T(2, K). Some examples of the T(2, K) factor are shown in fig. 1. The atomic scattering factor f ( x ) in X-ray experiments is usually taken from calculations and depends strongly on K, while the corresponding coherent scattering cross-section acoh and the incoherent scattering cross-section %~ in neutron diffraction experiments are independent of 2 and K. The contribution to the intensity TABLE 1 The multiple scattering cross-section o m calculated for the cases studied [after Bleck and Averbachl5)]. Sample Se Se N2 02
(2 (2 (2 (2
= = = =
0.852 1.705 1.06 1.06
am/aco h A) A) A) A)
0.235) 0.205) 0.206) 0.106)
I
I
1
i
I
U..R
"--.~.7
~ ( x - ~ __ .i.~" O,,,~ /
-I-•"
t
2
while in the case of a neutron diffraction experiment, eq. (1) takes the form 5,~3)
so(x:
o.8s2A)_
_ - 4 - - - -_ " -- ~0.3
I(~, ~) = A(~) T(X,K) {S(K) +f~(2,~) +
02
Br, ( k: 1.044 A )
O'inc O'm (,~)/ +-[1 +f~(,~, ,0-1 + • O'coh O'coh J
(2)
2 is the wavelength of the incident neutrons or X-rays and x - 4 n sin 0/2, where we use 0 for half the total scattering angle. Apart from the normalization constants A(2), the correction terms in
eqs. (1) and (2)
are
as follows"
The transmission coefficient T(2, K) depends on
_~ (~--5.o__%,~
Br,
-
-
N2 - 0 . t
.9,0 ( ) , : t o e k ) 1
o
I
20
I
I
40
I
6O
I
80 28 [DEGREES]
Fig. 1. The reciprocal transmission factor T(2,K) 14) for cylindrical samples with radius R and total scattering length # cm -1 as a function of the scattering angle 20. The actual
curves for the samples considered are shown.
RELIABLE PAIR DISTRIBUTION FUNCTIONS from the multiple scattering in X-ray measurements s(K) is strongly dependent on K, while the corresponding contribution in neutron scattering am(2) mainly depends on 2, and usually shows a rather modest K dependence for samples of cylindrical shape15). Some typical values are shown in table 1. A more detailed analysis of the multiple scattering may be obtained from computer simulations16). The static structure factor S(K) may be defined in terms of the dynamical structure factor S(K, co) by S(K)
=
~00
S(K, co) doJ.
(3)
In actual measurements the intensity is measured at fixed scattering angles, 20, which means that the integration in eq. (3) has not been performed correctly. This is only of practical importance for neutrons, and the correction of this effect is given by the Placzek correctionl°). It can be expressed as a series expansion in K2, with the frequency moments of S(K, co) as coefficients. Only the first and the second moments of the frequency are independent of the interactions between the atoms of the sample, implying that direct measurements of S(K, co) are necessary for the determination of the higher moments. Including only terms involving the first and the second moment of the frequency, the Placzek correction may approximately be written
fp(2, t¢)
=
(kBT
2Eo
x
571
the different magnitudes of the factors involved in the two cases are appreciable. For practical purposes this is important, since the equality of the structure factors derived from X-ray and neutron experiments serves as a good test that both results are correct. The fact that minor differences are present in the reported S(K) for many materials, points towards the care with which one has to draw definite conclusion about the details of the structure of non-crystalline matter. 3. The procedure We now proceed to the derivation of the pair distribution function g(r). The relation between this and the structure factor S(x) is given by
2 [ 9 ~c[S-(~c)-1] sin(~cr) dt¢, (5) 4 rtrpo [g (r) - 1] = -jo
7[
and
4~rrpo[ g ( r ) -
tc[S(t¢)- 1] =
1] sin(for) dr,
where P0 is the number density of the matter. Physically we know that atoms do not get closer than a certain distance rl, determined by the range of the hardcore interactions, i.e.
g(r)
=
0,
for
r < rI .
(7)
In this range we therefore obtain from eq. (5)
-
4 rrrpo = -2 f ~ /c[S(/c)-
1] sin(for) dtc.
0
fp ('~.)
tCo~)(M )
+
-=+(~ "tc2k~ 83tc~)
× 0.
+ ...
(4)
T is the temperature, kB the Boltzman constant, E0 the energy of the incident neutrons, k0 the wavevector of the neutrons, m the neutron mass and M the atomic mass of the constituent atoms of the sample. It is evident that for heavy nuclei, the significance of fp(2, x) decreases while for large K the importance o f £ ( 2 , K)increases. This is illustrated in fig. 2, where we show £ ( 2 , K) for the studied cases using the actual experimental parameters. Having discussed the corrections necessary to extract S(x), let us point out that the various terms are quite different in the cases of X-ray and neutron diffraction. In part this is reflected in the different appearences of eqs. (1) and (2), but also
(6)
S~e ~
:
o
Se
.
8
s
\\
2
1
x =toz,~,.~
-0.05
-0 l r
=1.06,~ ,
o
Fig.
1 s
x =lo63, , , ,
1. , ~o
~[~-~] ;--
2. The Placzek corrections for the cases studied.
(8)
572
F. Y S S I N G
HANSEN
Eq. (8) may be used to determine the normalization constants in either eq. (1) or (2) by a leastsquares fit of the integral in eq. (8) to the line -4rcrpo for small r. The procedure may even be used to determine other u n k n o w n parameters in eq. ( 1 ) o r (2). As an exampleS), we considered Gob and ainc of Se to be too poorly known for our purpose, and kept them as variable parameters. Indeed we found that the choice of cross-sections had an appreciable influence on the fit, and as a by-product of the c~ilculations we got an indirect determination of the scattering cross-sections, The problem is now that in ordinary X-ray and thermal neutron diffraction experiments, the integrand in eq. (8)is not known until convergence to zero, as mentioned in the introduction. Thus, a direct integration will produce oscillations in the region of small r. We therefore have to find a method, which eliminates these truncation errors, To set up such a method, let us first introduce a model pair distribution function gm(r) 4 rCpor 2 g m(r) --
i
Ni
r exp If-(r-ri)2.1
i=1 x/(2rur 2) ri. 4rtpo r2
gm(r)
--
'
r
2ai 4ztPo r2 ,
r ~> r t .
(9) (10)
According to eq. (9) it is assumed that the structure may be resolved in a series of well defined neighbour shells around a set of mean distances ri, characterized by certain a~'s, which define the range of each shell, and the number of atoms Ni per shell. The distribution functions are not exactly Gaussian due to the factor r/r~, which is introduced for mathematical convenience [cf. eq. (11)]. Significant deviations from the Gaussian form will only occur for peaks, where tr~/r~ is not small. In practice only the first few nearest neighbours are so well defined that it makes sense to use an expression for the structure like eq. (9). For greater r, the peaks become less pronounced corresponding to many terms like in eq. ( 9 ) w i t h closely spaced characteristic distances, going to be a continuous distribution giving a uniform density at large r. For our purpose we therefore assume the density to be uniform for distances beyond a certain d i s t a n c e r t as given by eq. (10). This is of course a very crude model of the matter, but the point is that we do not need an exact model for our purpose. From the theory of Fourier transformation we know that the contribution to S(K) at large K
AND
KIM C A R N E I R O
is due to short range correlations in the matter. In the model, we therefore focus our interest on the short distance correlations, as the correlations at larger distances are not going to make contributions at large K values. Satisfactory information about the contributions from these distances are directly extracted from the experimental structure factor. The summation in eq. (9)extends over the first few peaks in the pair distribution function, i.e., n is small number, rt corresponds roughly to the position of the last neighbour shell considered in the model, i.e., r,,~r t. We need not define it more exactly, nor do we have to specify in detail how the two regions are connected. Our model of the matter is sketched in fig. 3 for the case of n - 2. It is obvious that we cannot determine the parameters in eq. ( 9 ) j u s t by fitting eq. (9) to the direct Fourier transform of the experimental structure factor, as the resulting pair distribution function is contaminated by truncation errors. Consequently we eliminate these in the following way: First we calculate from eq. (9) and (10) the model structure factor sm(K) using eq. (6). With our choice of gm(r) this can be done analytically. We find sm(K)-- 1 = ~ Ni sin(~cri)exp(-1Cr2tc2) + i= 1
Kri
4rtpo [Krt c o s ( t c r , ) -
+
sin(tert)]
x
(11)
~c3
x exp(-½tr 2/(72). In eq. (11) we have introduced at to simulate the smooth transition between the two different regions in our model pair distribution function as sketched in fig. 3. Even without the exponential damping factor one finds, that the last term in eq. (11) for most non-crystalline materials has converged to zero at 8-9 A -1 . Therefore, an exact definition of the details of the transition region (degm [r)
~S m (x)
1
',\ i
/ 2
rl
r2
',__ r
j
_
~m
Fig. 3. A sketch of a typical model pair distribution function gin(r) in the case of n = 2 and the corresponding model struc-
ture
factor
SIn(K).
RELIABLE
PAIR
DISTRIBUTION
fined by rt and at) is of no interest to us, as we only intend to use eq. (11) for larger x's. The model structure factor is also shown in fig. 3. Second we use eq. ( 1 1 ) t o calculate the model pair distribution function gt(r), when the model structure factor is truncated at the experimental Xm (see fig. 4). According to eq. (5) we find
4rcrpo [ g t ( r ) - 1] = 2 re'
dt¢
I
i= 1
Ni ~
ri
×
--1
× e x p ( - ½ a 2K2) I - 4rcpor, (12) where the last term on the right-hand side of eq. (12) is the Fourier transform of the corresponding term in eq. (11), which is assumed to have converged to zero at Xm. This is not the case for the integrand of the first terms in eq. (12), so they have to be evaluated numerically. Therefore, gm(r) and gr~(r) will be different, as there are introduced truncation errors in the latter (see fig. 4). We may now compare the first peak and part of the seond peak (see fig. 5 ) o f g(r), which is calculated as the direct Fourier transform of the experimental structure factor S(~:), with the similar part of gt(r), as the same truncation errors have been introduced in both functions. Hence, the influence of truncation errors in the comparison has been eliminated as intended. It is, therefore, possible to determine g~n (r)
s~(x)
l
I
/
Vv
I 1(. m
Fig. 4. A sketch of the model pair distribution function gr~(r) calculated from the truncated m o d e l structure factor sin(K). False satellites are shown in gr~(r) due to truncation errors. g(r)
573
FUNCTIONS
the peak positions r;, the widths a~, and the number of atoms N~ in each shell by a least-squares fit of the relevant parts of g(r) and g~(r). When the parameters have been determined we may use eq. (11) to extend the experimental S(K) until convergence, and from the combined experimental and model structure factors obtain the final pair distribution function (fig. 6). This may be stated in another way. In the fit we obtain the parameters (r;, or;, N), and these are used in our final description of g(r) for the neighbour shells, which give rise to scattering intensity for wavevectors larger than the experimental cut-off Xm. Our numerical procedure for normalization and correction of diffraction data in order to obtain the structure factor and the pair distribution function of non-crystalline materials may be summarized in the following way: 1) From the intensity data we perform a calculation of the unknown constants of eq. (1) or (2) by a least-squares fit to the straight line -4zrrpo [eq. (8)]. The calculation involves only r-values less than those characteristic of the first peak in the pair distribution function. Then calculate the experimental structure factor S(K) and g(r)(fig. 5). 2) We now estimate a set of parameters r~, a;, Nj from g(r) and use these as first guesses of the final values, which are found by a least-squares fit of the relevant parts of g~(r) and g(r) using eq. (5). 3) In the calculation of the constants of eq. (1) or (2), all the influence of truncation errors may not have been eliminated. The results may therefore be refined, if the constants are redetermined in a least-squares fit calculation to 4np0r[g~(r)-1] rather than to the line -4rcp0r. The influence of the truncation errors will thus have been eliminated. (2) and (3) are then repeated until convergence in the constants and the parameters of the model. The best criteria for a consistent data treatment is, that all spurious oscillations in g(r) for small r
s(x) g(r)
1
= fi
Vv
A
V
"7
,....-
"x,m
Fig, 5. The pair distribution function g(r) calculated from the experimental structure factor S(K). g(r) shows truncation errors.
r
~.m
"IK
Fig. 6. The final structure factor S(K). Full curve is the experimental part, and dotted curve is the model structure factor, sm(x) shown here only for x > x m. g(r) is the final pair distribution function.
574
F
YSSING
HANSEN
AND
KIM
CARNEIRO
TABLE 2 Parameters of the model pair distribution f u n c t i o n s for the various s y s t e m s studied.
Se5)
N2
Br 2
02
293 K
80 K
X-ray18)
Neutron6)
X-ray 18)
Neutron 6)
Neutron8)
r1 N1 a1
2.32 1.95 0.068
2.32 1.95 0.095
1.096 1.002 0.08
1.115 0.96 0.03
1.217 1.00 0.08
1.216 0.94 0.04
2.273 0.9 0.04
~) N2 a2
3.96 6.90 0.257
3.96 7.3 0.235
have vanished. In the least square fit calculations we have used the method of Powell~7), which has proved to be very effective. It only uses the analytic expressions for the functions to be fitted, while analytic expressions of the derivatives are not needed. This makes the fitting procedure rather flexible.
width to decrease with decreasing temperature, as the atomic vibration becomes less, unless a temperature decrease also injects structural changes, which, as judged from the neutron work, seem to
4. Discussion We have described a procedure, which unifies the calculation of structure factors from intensity data with the calculation of reliable pair distribution functions. Even after the appearance of epithermal neutrons a method like this seems to be of value, as important correction terms at large wavevectors are unknown. We have applied the present method to thermal neutron diffraction data on amorphous Se 5)to X-ray data ~8) of liquid N2 and liquid 02 and to thermal neutron data on the liquids N2, 02, and Br2 6,8). The most import a n t results concerning the parametrization of the model pair distribution function are summarized in table 2. In the case of amorphous Se it is seen that we have parametrized the first two peaks. Only the first peak is well separated, but the second peak is still rather pronounced. We found that the parametrization of the second peak showed some influence on the determination of the first peak due to the mutual influence in the region between the peaks. This effect was not observed in the case of liquid N2, 02, and Br2, where only the first peak was parametrized [i.e., n = 1 in eq. (9)]. A peculiar effect was noticed in the case of Se, where the width of the first peak was increased, when the temperature was decreased. As the first peak is well separated, the parameters have a clear physical meaning. One would therefore expect the
o s
RMS
.
1.0
1=.1
112 ~ [/~]
1.0
1.5
~RMS
o.s
0.5
N~
RMS
o.s
1
,
,
,
,
2
3
4
5
_)_
2 [~] ++a+
Fig. 7. Root mean square deviation curves, RMS, for the three parameters r 1, N l, a I d e t e r m i n e d in the data t r e a t m e n t of the
X-ray data of N2.
RELIABLE PAIR DISTRIBUTION
take place. In view of our discussion above, it should however be pointed out that the conclusions drawn from the X-ray work of Kaplow et al. 19) do not agree on this point. The parameters of peaks, which are not well separated, do not have a well defined physical interpretation. They depend to a certain degree on how much of the peak has been used in the fit. In the case of Se we have thus used the second peak up to a wave vector transfer of 4 A -1. It was sufficient to describe the rising part of the peak, and the region around the maximum. That was all needed for a precise evaluation of the parameters of the first peak, which again was sufficient for an extension of the experimental structure factor until convergence. In the case of liquid N2, 02, and Br2 we were able to obtain good results for the bond lengths and for the vibrational amplitudes, in the two former cases even better results than those originally obtained from LINAC-data, mainly because these were not Placzek corrected. As the atomic vibrations in N2, 02, and Br2 are almost harmonic, the first peaks should be of Gaussian form. The deviation from that form in our model eq. (10) due to the factor r/ri is easily seen to be of no importance, as a~ir~ is small. Our experience with the method indicates that the most sensitive parameter is r;, then N;, and finally a;. This is illustrated in fig. 5, where the mean square deviations in a typical fit are shown for the three parameters. As at is determined with the greatest uncertainty, one may obtain a valuable check on the calculated parameter by a comparison with large ~c LINAC-data. Some of the parameters in eq. (1) or (2) are quite often given with great uncertainty, and it is im-
FUNCTIONS
575
portant to notice that our method also makes it possible to keep such parameters floating and to determine their values together with the normalization constant in the fit to the known parts of
g(r). References
1)
R. Kaplow, T. A. Rowe and B. L. Averbach, Phys. Rev. 168 (1968) 1068. 2) R. D. Mountain, J. Chem. Phys. 57 (1972)4346. 3) R. B. Blackmann and J. W. Tuckey, The measurement ojpower spectra (Dover, New York, 1958). 4) B. E. Warren, X-ray dijfraction (Addison-Wesley, New York) sect. 10. 5) F. Yssing Hansen, T. S. Knudsen and K. Carneiro, J. Chem. Phys. 62, no. 4 (1975) 1556. 6) j. C. Dore, G. Walford and D. I. Page, Molec. Phys. 29 (1975) 565. 7) j. H. Clarke. J. C. Dore and R. N. Sinclair, Molec. Phys. 29 (1975) 581. 8) j. H. Clarke, J. C. Dore, G. Walford and R. N. Sinclair, Molec. Phys. 31 (1976) 883. 9) D. I. Page and J. G. Powles, Molec. Phys. 29 (1975) 000. 10) G. Placzek, Phys. Rev. 86 (1953) 377. ll) j. L. Yarnell, M. J. Katz, R. G. Wenzel and S. H. Koenig, Phys. Rev. A7 (1973) 2130. 12) j. H. Konnert and J. Karle, Acta Cryst. A 29 (1973) 702. 13) A. de Graaf and B. Mozer, J. Chem. Phys. 55 (1971) 4076. 14) International tables for X-ray crystallography (Kynoch, Birmingham, England, 1972). 15) I. A. Blech and B. L. Averbach, Phys. Rev. 137 (1965) A1113. 16) j. R. D. Copley, Comp. Phys. Comm. 7 (1974) 289. 17) M. J. D. Powell, Comp. J. 7 (1965) 303. 18) H. W. Furumoto and C. H. Shaw, Phys. Fluids 7 (1964) 1026; tabulated by P. W. Schmidt and C. W. Thompson, in Simple dense fluids (ed. H. L. Frisch and Z. W. Salsburg; Academic, New York, 1968). 19) R. Kaplow, S. L. Strong and B. L. Averbach, Acta Cryst. 19 (1965) 1043. 20) K. R. Rao and B. A. Dasannacharya, Pram/ma 5, no. 6 (1975) 328.
INSTRUMENTS
AND M E T H O D S
143 ( 1 9 7 7 )
569-575
; ©
NORTH-HOLLAND
PUBLISHING
CO.
A P R O C E D U R E TO OBTAIN RELIABLE PAIR D I S T R I B U T I O N F U N C T I O N S OF N O N - C R Y S T A L L I N E MATERIALS F R O M D I F F R A C T I O N DATA F. YSSING HANSEN
Fysisk-Kemisk Institut, The Technical University o.I Denmark, DK 2800 Lyngby, Denmark and KIM CARNEIRO
Physics Laboratory 1, H.C. Orsted Institute, Universitetsparken 5, DK 2100 Copenhagen 0, Denmark Received 10 January 1977 A simple numerical method, which unifies the calculation of structure factors from X-ray or neutron diffraction data with the calculation of reliable pair distribution functions, is described. The objective of the method is to eliminate systematic errors in the normalizations and corrections of the intensity data, and to provide measures for elimination of truncation errors without loosing information about the structure. This is done through an iterative procedure, which is easy to program for computers. The applications to amorphous selenium and diatomic liquids are briefly reviewed.
1. Introduction
The determination of pair distribution functions from structure factor data provided either by X-ray or neutron diffraction experiments is sensitive both to proper corrections of the experimental results as well as to proper Fourier transformation of the structure factor. The two problems may at first sight seem uncorrelated, but it is very difficult to judge from the structure factor alone, whether or not the instrumental corrections, including normalization, of the diffraction data have been performed consistently. The consistency of the data treatment however is clearly revealed, when the structure factor data are Fourier transformed to give the pair distribution function. The influence of improper normalization and corrections of the diffraction data on the pair distribution function has been discussed several times (see e.g. refs. 1 and 2), and the conclusion is that all types of errors in the data treatment result in strong unphysical oscillations in the pair distributions function at small r, where the pair distribution function is known to be equal to zero due to the hard core repulsion between atoms. Therefore, this region may be used to see, if the data treatment has been performed properly. It is evident that an erroneous data treatment of the experimental results also demolishes the pair distribution function at larger r, but in this region one has no way to separate out "false" contributions from " t r u e " contributions. The Fourier transformation of the structure fac-
tor may also introduce errors in the pair distribution function. It is well known 3,4) that the Fourier transform of a function, which is not known until convergence to zero, contains artificially broadened peaks and around those, also false satellites. These satellites are mainly concentrated around the main peaks of the transformed function. In typical Xray and thermal neutron experiments, the maxim u m wavevector transfer /£m is 9--13 A - 1 . At this wavevector the structure factor of almost all non-crystalline materials still shows significant features, which indicate that the integrand of the Fourier integral has not yet converged to zero. A direct transform of the incomplete structure factors will therefore always give false features in the pair distribution function and smear out finer details. Reported pair distribution functions in the literature often show strong oscillations at small r, partly due to truncation errors in the Fourier transform and partly to improper treatment of the diffraction data. In order to remove the false satellites from the pair distribution function it has been common practice to apply convergence factors or window functions 4,s) to the measured structure factor to make the integrand in the Fourier integral converge. Such a procedure removes the false peaks, but this advantage is achieved at the expense of poorer resolution and still leaves the resulting pair distribution function quantitatively inadequate. An extension of the experimental region to higher/£m values should therefore make it possible
570
F.
YSSING
HANSEN
to avoid truncation errors in the final pair distribution function. The recent occurrence of experimental facilities in the form of either accelerators or reactor hot sources providing epithermal neutrons has, in fact, made such extensions feasible, However, the pioneering LINAC measurements on diatomic liquids 6-9) have shown that the problems in deriving reliable pair distribution functions are not solved by going to large Km'S alone. As truncation errors should be eliminated, the deficiciences of the results are likely due to improper normalization and corrections of the diffraction data and in particular, as shown by Page and Powles9), the Placzek correction plays an important role in some of these experiments. This clearly stresses the dilemma in structure studies of non-crystalline materials. On the one hand an accurate Fourier transformation requires that K"m is large, but on the other hand, since at large wavevectors the Placzek correction becomes progressively more problematic, Km should not be too large. In this paper we will show, how it is possible to choose a Km for which both requirements above are met in the sense that although none of the two problems is eliminated, both of them are kept at a manageable level. The same problems have also been discussed recently by Konnert et al. 12) and Rao et al. 2°) and their methods together with our method may be taken as a natural extension of the work of Kaplow e t al.l'19). 2. Theory The relation between intensity data I(2, K) in an X-ray diffraction experiment and the structure factor of a non-crystalline sample with random orientations may be written 4) I(2, K) = A(2) T(2, tc) [f2(~c) S(K)+ s(tc)],
(l)
AND
KIM
CARNEIRO
both 2 and K and is easy to determine, if the total scattering length is known, and if the shape of the sample is simple, like a spere, a cylinder or a platel4). For other shapes one may have to consider Monte Carlo simulations in order to evaluate T(2, K). Some examples of the T(2, K) factor are shown in fig. 1. The atomic scattering factor f ( x ) in X-ray experiments is usually taken from calculations and depends strongly on K, while the corresponding coherent scattering cross-section acoh and the incoherent scattering cross-section %~ in neutron diffraction experiments are independent of 2 and K. The contribution to the intensity TABLE 1 The multiple scattering cross-section o m calculated for the cases studied [after Bleck and Averbachl5)]. Sample Se Se N2 02
(2 (2 (2 (2
= = = =
0.852 1.705 1.06 1.06
am/aco h A) A) A) A)
0.235) 0.205) 0.206) 0.106)
I
I
1
i
I
U..R
"--.~.7
~ ( x - ~ __ .i.~" O,,,~ /
-I-•"
t
2
while in the case of a neutron diffraction experiment, eq. (1) takes the form 5,~3)
so(x:
o.8s2A)_
_ - 4 - - - -_ " -- ~0.3
I(~, ~) = A(~) T(X,K) {S(K) +f~(2,~) +
02
Br, ( k: 1.044 A )
O'inc O'm (,~)/ +-[1 +f~(,~, ,0-1 + • O'coh O'coh J
(2)
2 is the wavelength of the incident neutrons or X-rays and x - 4 n sin 0/2, where we use 0 for half the total scattering angle. Apart from the normalization constants A(2), the correction terms in
eqs. (1) and (2)
are
as follows"
The transmission coefficient T(2, K) depends on
_~ (~--5.o__%,~
Br,
-
-
N2 - 0 . t
.9,0 ( ) , : t o e k ) 1
o
I
20
I
I
40
I
6O
I
80 28 [DEGREES]
Fig. 1. The reciprocal transmission factor T(2,K) 14) for cylindrical samples with radius R and total scattering length # cm -1 as a function of the scattering angle 20. The actual
curves for the samples considered are shown.
RELIABLE PAIR DISTRIBUTION FUNCTIONS from the multiple scattering in X-ray measurements s(K) is strongly dependent on K, while the corresponding contribution in neutron scattering am(2) mainly depends on 2, and usually shows a rather modest K dependence for samples of cylindrical shape15). Some typical values are shown in table 1. A more detailed analysis of the multiple scattering may be obtained from computer simulations16). The static structure factor S(K) may be defined in terms of the dynamical structure factor S(K, co) by S(K)
=
~00
S(K, co) doJ.
(3)
In actual measurements the intensity is measured at fixed scattering angles, 20, which means that the integration in eq. (3) has not been performed correctly. This is only of practical importance for neutrons, and the correction of this effect is given by the Placzek correctionl°). It can be expressed as a series expansion in K2, with the frequency moments of S(K, co) as coefficients. Only the first and the second moments of the frequency are independent of the interactions between the atoms of the sample, implying that direct measurements of S(K, co) are necessary for the determination of the higher moments. Including only terms involving the first and the second moment of the frequency, the Placzek correction may approximately be written
fp(2, t¢)
=
(kBT
2Eo
x
571
the different magnitudes of the factors involved in the two cases are appreciable. For practical purposes this is important, since the equality of the structure factors derived from X-ray and neutron experiments serves as a good test that both results are correct. The fact that minor differences are present in the reported S(K) for many materials, points towards the care with which one has to draw definite conclusion about the details of the structure of non-crystalline matter. 3. The procedure We now proceed to the derivation of the pair distribution function g(r). The relation between this and the structure factor S(x) is given by
2 [ 9 ~c[S-(~c)-1] sin(~cr) dt¢, (5) 4 rtrpo [g (r) - 1] = -jo
7[
and
4~rrpo[ g ( r ) -
tc[S(t¢)- 1] =
1] sin(for) dr,
where P0 is the number density of the matter. Physically we know that atoms do not get closer than a certain distance rl, determined by the range of the hardcore interactions, i.e.
g(r)
=
0,
for
r < rI .
(7)
In this range we therefore obtain from eq. (5)
-
4 rrrpo = -2 f ~ /c[S(/c)-
1] sin(for) dtc.
0
fp ('~.)
tCo~)(M )
+
-=+(~ "tc2k~ 83tc~)
× 0.
+ ...
(4)
T is the temperature, kB the Boltzman constant, E0 the energy of the incident neutrons, k0 the wavevector of the neutrons, m the neutron mass and M the atomic mass of the constituent atoms of the sample. It is evident that for heavy nuclei, the significance of fp(2, x) decreases while for large K the importance o f £ ( 2 , K)increases. This is illustrated in fig. 2, where we show £ ( 2 , K) for the studied cases using the actual experimental parameters. Having discussed the corrections necessary to extract S(x), let us point out that the various terms are quite different in the cases of X-ray and neutron diffraction. In part this is reflected in the different appearences of eqs. (1) and (2), but also
(6)
S~e ~
:
o
Se
.
8
s
\\
2
1
x =toz,~,.~
-0.05
-0 l r
=1.06,~ ,
o
Fig.
1 s
x =lo63, , , ,
1. , ~o
~[~-~] ;--
2. The Placzek corrections for the cases studied.
(8)
572
F. Y S S I N G
HANSEN
Eq. (8) may be used to determine the normalization constants in either eq. (1) or (2) by a leastsquares fit of the integral in eq. (8) to the line -4rcrpo for small r. The procedure may even be used to determine other u n k n o w n parameters in eq. ( 1 ) o r (2). As an exampleS), we considered Gob and ainc of Se to be too poorly known for our purpose, and kept them as variable parameters. Indeed we found that the choice of cross-sections had an appreciable influence on the fit, and as a by-product of the c~ilculations we got an indirect determination of the scattering cross-sections, The problem is now that in ordinary X-ray and thermal neutron diffraction experiments, the integrand in eq. (8)is not known until convergence to zero, as mentioned in the introduction. Thus, a direct integration will produce oscillations in the region of small r. We therefore have to find a method, which eliminates these truncation errors, To set up such a method, let us first introduce a model pair distribution function gm(r) 4 rCpor 2 g m(r) --
i
Ni
r exp If-(r-ri)2.1
i=1 x/(2rur 2) ri. 4rtpo r2
gm(r)
--
'
r
2ai 4ztPo r2 ,
r ~> r t .
(9) (10)
According to eq. (9) it is assumed that the structure may be resolved in a series of well defined neighbour shells around a set of mean distances ri, characterized by certain a~'s, which define the range of each shell, and the number of atoms Ni per shell. The distribution functions are not exactly Gaussian due to the factor r/r~, which is introduced for mathematical convenience [cf. eq. (11)]. Significant deviations from the Gaussian form will only occur for peaks, where tr~/r~ is not small. In practice only the first few nearest neighbours are so well defined that it makes sense to use an expression for the structure like eq. (9). For greater r, the peaks become less pronounced corresponding to many terms like in eq. ( 9 ) w i t h closely spaced characteristic distances, going to be a continuous distribution giving a uniform density at large r. For our purpose we therefore assume the density to be uniform for distances beyond a certain d i s t a n c e r t as given by eq. (10). This is of course a very crude model of the matter, but the point is that we do not need an exact model for our purpose. From the theory of Fourier transformation we know that the contribution to S(K) at large K
AND
KIM C A R N E I R O
is due to short range correlations in the matter. In the model, we therefore focus our interest on the short distance correlations, as the correlations at larger distances are not going to make contributions at large K values. Satisfactory information about the contributions from these distances are directly extracted from the experimental structure factor. The summation in eq. (9)extends over the first few peaks in the pair distribution function, i.e., n is small number, rt corresponds roughly to the position of the last neighbour shell considered in the model, i.e., r,,~r t. We need not define it more exactly, nor do we have to specify in detail how the two regions are connected. Our model of the matter is sketched in fig. 3 for the case of n - 2. It is obvious that we cannot determine the parameters in eq. ( 9 ) j u s t by fitting eq. (9) to the direct Fourier transform of the experimental structure factor, as the resulting pair distribution function is contaminated by truncation errors. Consequently we eliminate these in the following way: First we calculate from eq. (9) and (10) the model structure factor sm(K) using eq. (6). With our choice of gm(r) this can be done analytically. We find sm(K)-- 1 = ~ Ni sin(~cri)exp(-1Cr2tc2) + i= 1
Kri
4rtpo [Krt c o s ( t c r , ) -
+
sin(tert)]
x
(11)
~c3
x exp(-½tr 2/(72). In eq. (11) we have introduced at to simulate the smooth transition between the two different regions in our model pair distribution function as sketched in fig. 3. Even without the exponential damping factor one finds, that the last term in eq. (11) for most non-crystalline materials has converged to zero at 8-9 A -1 . Therefore, an exact definition of the details of the transition region (degm [r)
~S m (x)
1
',\ i
/ 2
rl
r2
',__ r
j
_
~m
Fig. 3. A sketch of a typical model pair distribution function gin(r) in the case of n = 2 and the corresponding model struc-
ture
factor
SIn(K).
RELIABLE
PAIR
DISTRIBUTION
fined by rt and at) is of no interest to us, as we only intend to use eq. (11) for larger x's. The model structure factor is also shown in fig. 3. Second we use eq. ( 1 1 ) t o calculate the model pair distribution function gt(r), when the model structure factor is truncated at the experimental Xm (see fig. 4). According to eq. (5) we find
4rcrpo [ g t ( r ) - 1] = 2 re'
dt¢
I
i= 1
Ni ~
ri
×
--1
× e x p ( - ½ a 2K2) I - 4rcpor, (12) where the last term on the right-hand side of eq. (12) is the Fourier transform of the corresponding term in eq. (11), which is assumed to have converged to zero at Xm. This is not the case for the integrand of the first terms in eq. (12), so they have to be evaluated numerically. Therefore, gm(r) and gr~(r) will be different, as there are introduced truncation errors in the latter (see fig. 4). We may now compare the first peak and part of the seond peak (see fig. 5 ) o f g(r), which is calculated as the direct Fourier transform of the experimental structure factor S(~:), with the similar part of gt(r), as the same truncation errors have been introduced in both functions. Hence, the influence of truncation errors in the comparison has been eliminated as intended. It is, therefore, possible to determine g~n (r)
s~(x)
l
I
/
Vv
I 1(. m
Fig. 4. A sketch of the model pair distribution function gr~(r) calculated from the truncated m o d e l structure factor sin(K). False satellites are shown in gr~(r) due to truncation errors. g(r)
573
FUNCTIONS
the peak positions r;, the widths a~, and the number of atoms N~ in each shell by a least-squares fit of the relevant parts of g(r) and g~(r). When the parameters have been determined we may use eq. (11) to extend the experimental S(K) until convergence, and from the combined experimental and model structure factors obtain the final pair distribution function (fig. 6). This may be stated in another way. In the fit we obtain the parameters (r;, or;, N), and these are used in our final description of g(r) for the neighbour shells, which give rise to scattering intensity for wavevectors larger than the experimental cut-off Xm. Our numerical procedure for normalization and correction of diffraction data in order to obtain the structure factor and the pair distribution function of non-crystalline materials may be summarized in the following way: 1) From the intensity data we perform a calculation of the unknown constants of eq. (1) or (2) by a least-squares fit to the straight line -4zrrpo [eq. (8)]. The calculation involves only r-values less than those characteristic of the first peak in the pair distribution function. Then calculate the experimental structure factor S(K) and g(r)(fig. 5). 2) We now estimate a set of parameters r~, a;, Nj from g(r) and use these as first guesses of the final values, which are found by a least-squares fit of the relevant parts of g~(r) and g(r) using eq. (5). 3) In the calculation of the constants of eq. (1) or (2), all the influence of truncation errors may not have been eliminated. The results may therefore be refined, if the constants are redetermined in a least-squares fit calculation to 4np0r[g~(r)-1] rather than to the line -4rcp0r. The influence of the truncation errors will thus have been eliminated. (2) and (3) are then repeated until convergence in the constants and the parameters of the model. The best criteria for a consistent data treatment is, that all spurious oscillations in g(r) for small r
s(x) g(r)
1
= fi
Vv
A
V
"7
,....-
"x,m
Fig, 5. The pair distribution function g(r) calculated from the experimental structure factor S(K). g(r) shows truncation errors.
r
~.m
"IK
Fig. 6. The final structure factor S(K). Full curve is the experimental part, and dotted curve is the model structure factor, sm(x) shown here only for x > x m. g(r) is the final pair distribution function.
574
F
YSSING
HANSEN
AND
KIM
CARNEIRO
TABLE 2 Parameters of the model pair distribution f u n c t i o n s for the various s y s t e m s studied.
Se5)
N2
Br 2
02
293 K
80 K
X-ray18)
Neutron6)
X-ray 18)
Neutron 6)
Neutron8)
r1 N1 a1
2.32 1.95 0.068
2.32 1.95 0.095
1.096 1.002 0.08
1.115 0.96 0.03
1.217 1.00 0.08
1.216 0.94 0.04
2.273 0.9 0.04
~) N2 a2
3.96 6.90 0.257
3.96 7.3 0.235
have vanished. In the least square fit calculations we have used the method of Powell~7), which has proved to be very effective. It only uses the analytic expressions for the functions to be fitted, while analytic expressions of the derivatives are not needed. This makes the fitting procedure rather flexible.
width to decrease with decreasing temperature, as the atomic vibration becomes less, unless a temperature decrease also injects structural changes, which, as judged from the neutron work, seem to
4. Discussion We have described a procedure, which unifies the calculation of structure factors from intensity data with the calculation of reliable pair distribution functions. Even after the appearance of epithermal neutrons a method like this seems to be of value, as important correction terms at large wavevectors are unknown. We have applied the present method to thermal neutron diffraction data on amorphous Se 5)to X-ray data ~8) of liquid N2 and liquid 02 and to thermal neutron data on the liquids N2, 02, and Br2 6,8). The most import a n t results concerning the parametrization of the model pair distribution function are summarized in table 2. In the case of amorphous Se it is seen that we have parametrized the first two peaks. Only the first peak is well separated, but the second peak is still rather pronounced. We found that the parametrization of the second peak showed some influence on the determination of the first peak due to the mutual influence in the region between the peaks. This effect was not observed in the case of liquid N2, 02, and Br2, where only the first peak was parametrized [i.e., n = 1 in eq. (9)]. A peculiar effect was noticed in the case of Se, where the width of the first peak was increased, when the temperature was decreased. As the first peak is well separated, the parameters have a clear physical meaning. One would therefore expect the
o s
RMS
.
1.0
1=.1
112 ~ [/~]
1.0
1.5
~RMS
o.s
0.5
N~
RMS
o.s
1
,
,
,
,
2
3
4
5
_)_
2 [~] ++a+
Fig. 7. Root mean square deviation curves, RMS, for the three parameters r 1, N l, a I d e t e r m i n e d in the data t r e a t m e n t of the
X-ray data of N2.
RELIABLE PAIR DISTRIBUTION
take place. In view of our discussion above, it should however be pointed out that the conclusions drawn from the X-ray work of Kaplow et al. 19) do not agree on this point. The parameters of peaks, which are not well separated, do not have a well defined physical interpretation. They depend to a certain degree on how much of the peak has been used in the fit. In the case of Se we have thus used the second peak up to a wave vector transfer of 4 A -1. It was sufficient to describe the rising part of the peak, and the region around the maximum. That was all needed for a precise evaluation of the parameters of the first peak, which again was sufficient for an extension of the experimental structure factor until convergence. In the case of liquid N2, 02, and Br2 we were able to obtain good results for the bond lengths and for the vibrational amplitudes, in the two former cases even better results than those originally obtained from LINAC-data, mainly because these were not Placzek corrected. As the atomic vibrations in N2, 02, and Br2 are almost harmonic, the first peaks should be of Gaussian form. The deviation from that form in our model eq. (10) due to the factor r/ri is easily seen to be of no importance, as a~ir~ is small. Our experience with the method indicates that the most sensitive parameter is r;, then N;, and finally a;. This is illustrated in fig. 5, where the mean square deviations in a typical fit are shown for the three parameters. As at is determined with the greatest uncertainty, one may obtain a valuable check on the calculated parameter by a comparison with large ~c LINAC-data. Some of the parameters in eq. (1) or (2) are quite often given with great uncertainty, and it is im-
FUNCTIONS
575
portant to notice that our method also makes it possible to keep such parameters floating and to determine their values together with the normalization constant in the fit to the known parts of
g(r). References
1)
R. Kaplow, T. A. Rowe and B. L. Averbach, Phys. Rev. 168 (1968) 1068. 2) R. D. Mountain, J. Chem. Phys. 57 (1972)4346. 3) R. B. Blackmann and J. W. Tuckey, The measurement ojpower spectra (Dover, New York, 1958). 4) B. E. Warren, X-ray dijfraction (Addison-Wesley, New York) sect. 10. 5) F. Yssing Hansen, T. S. Knudsen and K. Carneiro, J. Chem. Phys. 62, no. 4 (1975) 1556. 6) j. C. Dore, G. Walford and D. I. Page, Molec. Phys. 29 (1975) 565. 7) j. H. Clarke. J. C. Dore and R. N. Sinclair, Molec. Phys. 29 (1975) 581. 8) j. H. Clarke, J. C. Dore, G. Walford and R. N. Sinclair, Molec. Phys. 31 (1976) 883. 9) D. I. Page and J. G. Powles, Molec. Phys. 29 (1975) 000. 10) G. Placzek, Phys. Rev. 86 (1953) 377. ll) j. L. Yarnell, M. J. Katz, R. G. Wenzel and S. H. Koenig, Phys. Rev. A7 (1973) 2130. 12) j. H. Konnert and J. Karle, Acta Cryst. A 29 (1973) 702. 13) A. de Graaf and B. Mozer, J. Chem. Phys. 55 (1971) 4076. 14) International tables for X-ray crystallography (Kynoch, Birmingham, England, 1972). 15) I. A. Blech and B. L. Averbach, Phys. Rev. 137 (1965) A1113. 16) j. R. D. Copley, Comp. Phys. Comm. 7 (1974) 289. 17) M. J. D. Powell, Comp. J. 7 (1965) 303. 18) H. W. Furumoto and C. H. Shaw, Phys. Fluids 7 (1964) 1026; tabulated by P. W. Schmidt and C. W. Thompson, in Simple dense fluids (ed. H. L. Frisch and Z. W. Salsburg; Academic, New York, 1968). 19) R. Kaplow, S. L. Strong and B. L. Averbach, Acta Cryst. 19 (1965) 1043. 20) K. R. Rao and B. A. Dasannacharya, Pram/ma 5, no. 6 (1975) 328.