Structure Refinement Based on the Maximum-Entropy
Method from Powder Diffraction Data

Fujio Izumi and Ruben A. Dilanian
Advanced Materials Laboratory, National Institute for Materials Science,
1-1 Namiki, Tsukuba, Ibaraki 305-0044, Japan

ABSTRACTA multipurpose pattern-fitting program RIETAN2000, PRIMA for the maximum-entropy method(MEM), and VENUS for three-dimensional (3D) visualization of crystal structures and electron/nuclear
densities were substantially combined to grow into
an integrated system for MEM-based pattern fitting (MPF). In MPF, ‘observed’ structure factors,
Fo (Rietveld), estimated after Rietveld analysis are atfirst analyzed by MEM to yield 3D electron/nuclear
densities that are more or less biased in favor of astructural model in the Rietveld analysis. Then, we
evaluate structure factors, Fc (MEM), by the Fouriertransform of the densities and fit the whole powder
pattern calculated from the Fc (MEM) data to theobserved one to refine parameters irrelevant to thestructure. Fo (w.p.f.) data resulting from the whole-pattern fitting (w.p.f.) are analyzed again by MEM.
MEM analysis and whole-pattern fitting are alternately repeated until R factors in the latter no longer
decrease. These iterative processes make it possible
to minimize the bias imposed by the structural model
on final electron/nuclear densities. With this sophisticated methodology, crystal structures are represented not by structure parameters such as fractional
coordinates and atomic displacement parameters but
by 3D electron/nuclear densities. MPF is, therefore, capable of modeling static and dynamic disorder, chemical bonds, nonlocalized electrons, and anharmonic thermal motion more adequately than theconventional Rietveld method. Applications of MPF
to structure refinements of several compounds from(synchrotron) X-ray and neutron powder diffraction
data will be demonstrated, showing the excellent performance and effectiveness of MPF.

1. INTRODUCTION
Sakata et al. have been applying the maximumentropy method (MEM) [1] to determination ofelectron and nuclear (strictly speaking, coherentscattering length bc [2]) densities from powder diffraction data [3,4]. Diffraction intensities scattered in thereciprocal space are converted into electron/nuclear
densities filling the real space; that is, the densities
are visualized in three dimensions. MEM can estimate structure factors of reflections in the high-Q region excluded in the analysis of the powder diffraction
data, which is favorable for improving the S/N ratio
in the three-dimensional (3D) densities.
When dealing with compounds exhibiting heavily
overlapping reflections, integrated intensities are often estimated on the basis of the result of Rietveld
analysis [5]. In this procedure, the observed intensity at each point is apportioned in the ratio of profiles calculated from final parameters and summed up
for each reflection. The resultant ‘observed’ structure factors are then analyzed by MEM to give 3Delectron/nuclear densities in the unit cell. Unfortunately, ‘observed’ structure factors evaluated in thisway are doubly biased towards the structural model
because both phases and calculated integrated intensities used for the intensity partitioning are takenfrom the model [6].
We have been claiming that the above undesirable
effect, which is overlooked or disregarded by Sakata etal., should be overcome if only partially [7–10]. Even
if the structural model in the Rietveld analysis from
X-ray powder diffraction data is essentially correct,
chemical bonding is hardly represented with structure
parameters such as fractional coordinates and atomic
displacement parameters. This is the main reason
for the necessity of D synthesis in studies of electron
densities by single-crystal X-ray diffraction. Refinement of isotropic atomic displacement parameters, U ,in Rietveld analysis (anisotropic atomic displacement
parameters, Uij , cannot practically be refined from Xray powder diffraction data) and that of a scale factor
on the basis of the inadequate structural model will
lower the accuracy of the ‘observed’ structure factors.
In this article, we will introduce our original technology named MEM-based pattern fitting (MPF)[7–10]. We have integrated a multipurpose pattern-fitting program RIETAN-2000 [7–11], a MEM analysis program PRIMA [12] (formerly MEED developedby Kumazawa et al. [13]), and a 3D visualization program VENUS [14] to construct a system to refine
crystal structures by MPF from powder diffraction
data. With this MPF system, the bias imposed on
final electron/nuclear densities can efficiently be reduced by running RIETAN and PRIMA alternately.
Several examples of MPF structure refinementsfrom (synchrotron) X-ray and neutron powder
diffraction data will be demonstrated to verify its
high performance and effectiveness.

2. BACKGROUND
In single-crystal X-ray analysis, we often carry out
difference (D) synthesis, i.e., the inverse Fouriertransform of differences, between observed structure
factors, Fok , and structure factors, Fck , calculated
from structure parameters:∆Fk = Fok − Fck .

(1)

D synthesis is suitable for detecting subtle differences
between real and modeled structures.
Two serious problems hinder the satisfactory application of D synthesis to powder diffraction data. One
is the appearance of ‘ripples’ due to the termination
effect even in D synthesis where Fourier coefficientswith higher orders are regarded as negligible. Thiseffect makes it difficult to extract physically meaningful residual distribution from D-synthesis maps.
The other is difficulty in obtaining accurate observed
structure factors, |Fok |, owing to the collapse of the
3D reciprocal space onto the one-dimensional diffraction pattern.
In MEM analysis, no model function is fit to theobserved pattern unlike Rietveld analysis using leastsquares methods [5]; that is, we estimate density
distribution whose information entropy is maximized
within errors in observed diffraction data. Structural
information contained in diffraction data is accordingly extracted by MEM and reflected on the 3D
densities provided that the diffraction data have beenmeasured in an appropriate way. Electron densities
are calculated from Fok ’s in Fourier synthesis. On theother hand, in MEM, we solve a kind of an inverse

problem where electron densities are determined prior
to calculation of structure factors, FMk . The termination effect is, consequently, far less marked in MEM
analysis than in Fourier synthesis. Further, the MEM
can estimate structure factors of reflections in thehigh-Q region excluded in the analysis of intensity
data, as described earlier. MEM analysis, therefore,
gives less noisy 3D densities than Fourier synthesis
from a limited number of observed reflections.
The distribution of bc is also determinable by the
MEM analysis of neutron diffraction data. The term
“nuclear density” may be substituted for “electron
density” in subsequent general descriptions common
to X-ray and neutron diffraction.

3. SUPER-FAST MEM ANALYSISPROGRAM, PRIMA
Till lately, we had been utilizing MEED [13] for MEM
analysis. MEED is, however, obsolete at present because of its sluggish speed. We have recently developed our own program PRIMA [12], the stateof the art in MEM analysis of structure factors obtained experimentally by X-ray and neutron diffraction. PRIMA can determine 3D electron-density distribution from X-ray diffraction data and 3D nucleardensity distribution from neutron diffraction data. In
the case of neutron diffraction, it can deal with compounds containing elements whose bc ’s are negative,e.g., H, Li, Ti, and Mn [2]. PRIMA was written in
Fortran 90, which enables us to create storage forallocatable arrays dynamically. It is thoroughly optimized for personal computers equipped with Intel
Pentium 4 processors and operated with Microsoft
Windows.
3.1 Maximum-entropy equation
The general principle of MEM analysis is to findthe maximum of the information entropy, S, under
several constraints. There are a few variations of
the MEM formalism [15–17]. Here, we will follow
Collins’s formalism [16] based on Jaynes’s expressionof the information entropy [18].
In MEM analysis of X-ray/neutron diffraction
data, electron/nuclear densities are represented withthose in pixels (parallelepipeds) whose numbers along

a, b, and c axes are Na , Nb , and Nc , respectively. Let
N (= Na Nb Nc ) be the total number of pixels in theunit cell, ρk the normalized density at the position
rk in the 3D gridded space, ρ∗k the density at rk , andτk the normalized density at rk derived from priorinformation. Then, S is formulated as
N
X

ρk ln

ρk
τk

ρ∗
ρk = PN k

.

S=−



k=1

with

∗
k=1 ρk



(2)

(3)

S is maximized by imposing the following three constraints:
ρk > 0,
(4)

C=

M
1 X
|Fc (hj ) − Fo (hj )|2
= 1,
M j=1
σj2

D=

N
X

ρk = 1,

(5)

(6)

k=1

where M is the total number of reflections withknown phases, Fc (hj ) is the calculated structure factor for reflection hj , Fo (hj ) is the observed structure
factor, and σj is the estimated standard deviation(e.s.d.) of |Fo (hj )|. Equation (5) is referred to as theF -constraint (χ2 statistics type), and Eq. (6) as thenormalization constraint. When dealing with X-ray
diffraction data of a compound where T is the total
number of electrons in the unit cell, Fc (hj ) is givenasFc (hj ) = TNXρk exp (2πihj rk ) .

(7)

k=1

Let λ and µ be the Lagrange multipliers for C andD, respectively. Using the method of undetermined
multipliers, we have
Q = S − λ(C − 1) − µ(D − 1).

(8)Q can be maximized by setting the partial derivativesof Q with respect to ρk , λ, and µ at zero:
∂Q
= 0,
∂ρk

(9)

∂Q
= 0,
∂λ

(10)

∂Q
= 0.
∂µ

(11)

Then, the maximum-entropy (ME) equation is derived [19]:∂Cτkρk = exp −λ ·
Z
∂ρk
with

N
X

!

(12)

to an increase in the total number of iterations. However, the 0th-order single-pixel approximation works
well in nearly every case.
In an alternative method to solve Eq. (12), we usepart of the second term in Eq. (14), i.e., only theterm for k = l in the summation, which is referred
to as nonlinear single-pixel approximation [20]. Q is
optimized by the Newton-Raphson method:

!

∂C
Z=
τk exp −λ ·
.
∂ρk
k=1

(13)

3.2 0th order and nonlinear single-pixel approximations
Equation (12) cannot be solved analytically. Then,
Kumazawa et al. [19] introduced a 0th-order singlepixel approximation, where we expand the derivativeof C with respect to ρk in a Taylor series
N
X

∂ 2C
∂C
∂C
(ρl − τk ) ·=+∂ρk∂τk l=1∂τk ∂τl(14)

and neglect the second term in the right side ofEq. (14). This linear approximation is equivalent tothe simple replacement of ∂C/∂ρk in Eq. (12) with∂C/∂τk :
!
∂C
τk
.
(15)
ρk = exp −λ ·Z∂τkThis equation means that the ME density distribution is reached through iterations n = 0, 1, 2, · · ·:
(n)
ρk
(n+1)
ρk
=
exp

Z

−λ ·

∂C
(n)

∂ρk

!

.

(16)

We start MEM iterations, adopting the prior informationT(0)ρk = ,(17)Vviz., a uniform density in the unit cell with the volumeof V .
Because the second-order derivatives in the second
term of Eq. (14) are independent of ρ and λ, the second term is regarded as proportional to only ∆ρ =ρl −τk . That is, the 0th-order single-pixel approximation as well as the iterative procedure represented as
(n+1)
(n)
Eq. (16) is correct provided that ∆ρ = ρk
− ρk isnegligibly small. The only way to solve this problem
is to choose a sufficiently small λ value, which leads

Q(ρ + ∆ρ) = Q(ρ) + ∆ρ ·∂Q (∆ρ)2 ∂ 2 Q+·. (18)∂ρ2∂ρ2

The partial derivative of Q(ρ + ∆ρ) with respect to ρis zero provided that∂ 2Q
∆ρ = −
∂ρ2

!−1

·

∂Q
∂ρ

(19)

and
ρ(n+1) = ρ(n) + ∆ρ.

(20)

The Newton-Raphson method permits prompt convergence to the solution, decreasing the total numberof iterations considerably.
3.3 How to deal with elements of negative
coherent-scattering lengths
MEM is simply applicable to the analysis of Fo datafor compounds containing elements with negative bc
values [21–23]. In such cases, S is represented as thesum of those for positive and negative densities:
S = S + + S −.

(21)

This is a good approximation to S because of thenegligible overlap of nuclear densities for two atoms.
Thus, Eq. (7) is rewritten as
Fc (hj ) =

N 
X



−
ρ+
k + ρk exp (2πihj rk )

(22)

+
ρ+
k = T ρk ,

(23)

−
ρ−
k = T ρk ;

(24)

k=1

with

T + and T − are the total numbers of positive and negative bc ’s, respectively; ρk is the nuclear density normalized with Eq. (6). We can reach the ME nucleardensity distribution by solving the following equation
iteratively:τk±∂C
exp −λ · ± .
±
Z
∂ρk
!

ρ±
k =

(25)

3.4 Calculations of structure factors
Let P be the number of symmetry operations, and Rl
the symmetry-operation matrix, then the real part,Fc (real), and the imaginary part, Fc (imag.), of Fc (hj )in Eq. (7) for X-ray diffraction are computed from [13]Fc (real) = TNXρk Wjk (real),NXρk Wjk (imag.),(27)

cos (2πhj Rl rk ) ,

(28)

k=1

where
Wjk (real) =

P
X
l=1

Wjk (imag.) =PXsin (2πhj Rl rk ) .

(31)

λ0 = B0 + B1 G + B2 G2 ,

(32)

where


,

(33)

and A0 , A1 , A2 , A3 , B0 , B1 , and B2 are constants determined empirically. With Eqs. (30)−(32), PRIMA
is mostly capable of estimating appropriate λ0 .
PRIMA has an another feature of adjusting λ automatically during MEM iterations. In this technique,λ in the (n + 1)’th iteration (n = 1, 2, · · ·), λ(n+1) , isautomatically adjusted with a coefficient, t:

(29)

l=1

Fc (imag.) is necessary only in noncentrosymmetricspace groups.
In PRIMA, Wjk (real) and Wjk (imag.) associatedwith phases are first calculated for the asymmetricunit and stored in two-dimensional arrays allocateddynamically. The continuous use of storage for the arrays allocated dynamically accelerates the execution
speed of PRIMA dramatically.
In a relatively large system containing many reflections and pixels in the asymmetric unit, thecalculations of Wjk (real) and Wjk (imag.) are ratedetermining steps in MEM analysis. We restructured
this part fundamentally, fully utilizing space-groupsymmetry, diffraction properties, and fast evaluationof sine and cosine functions without sacrificing accuracy virtually, which considerably accelerates thecalculations of the two arrays. Details in this revolutionary breakthrough in MEM analysis will be reported elsewhere [12].
3.5 Initial value of λ and its adjustment
Learned information is needed to set the initial valueof λ, λ0 , appropriately. With too small λ0 , the convergence will be slow and may not reach C = 1. On theother hand, the convergence will never be attained
when λ0 is too large. PRIMA has a feature of automatic estimation of λ0 by using one of the following
three equations:
λ0 = G ln 10,

−1

M
1
1 X
G=
M j=1 σj2

(26)

k=1

Fc (imag.) = T

A1,(G − A2 )2 − A3

λ0 = A 0 +(30)

λ(n+1) = λ(n) + tλ0(34)λ(1) = λ0 .

(35)

with
If t = 0 (no adjustment, default value), λ is not varied. When the MEM iterations tend toward divergence, λ(n+1) is decreased to λ(n) /2. Four different
options are implemented in PRIMA for the adjustment of λ:1. λ0 is input from a file, and λ remains constant.
2. λ0 is input from a file, and λ is adjusted.
3. λ0 is automatically determined, and λ remains
constant.
4. λ0 is automatically determined, and λ is adjusted.
For example, we have empirically learned thatthe second option is the best selection in the nonlinear approximation method. This finding is explained in terms of the fact that, in the nonlinear
single-pixel approximation introducing the NewtonRaphson method, λ0 must be very small, which maylead to a failure in the automatic adjustment of λ.
3.6 R factors
Two R factors, which are measures of reliability of
MEM analysis, are defined as
PM

RF =

j=1 |Fc (hj ) − Fo (hj )|PM
j=1 |Fo (hj )|

(36)

and
P

RwF = 


2 1/2
−2
M
j=1 σj |Fc (hj ) − Fo (hj )| 
.
PM
2
−2
j=1 σj |Fo (hj )|(37)

3.7 Appropriate pixel numbers
Under favor of space-group symmetry [24], densities
in the asymmetric unit are calculated with Eq. (16).
Na , Nb , and Nc have to be set so as to satisfy the
space-group symmetry. For example, fractional coordinates (x, y, and z) in the asymmetric unit for spacegroup Fm 3̄m (No. 225) has the following ranges [24]:10≤x≤ ,
2(38)

10≤y≤ ,(39)
4
1
0≤z≤ .
(40)
4
Then, Na must be even while Nb and Nc must be
multiples of 4; for example, Na = 82, Nb = 128, and
Nc = 128. Fractional coordinates in the asymmetricunit for space group Pm (No. 6, unique axis b) [24]are confined to
0 ≤ x ≤ 1,(41)
10≤y≤ ,
20 ≤ z ≤ 1.

(42)
(43)

In this case, Nb should be even whereas Na and Nc
may be any integers large enough for MEM analysis. It is further preferable that symmetry elements
such as rotation axes and inversion centers are located on intersections of grid lines. For further guidelines concerning the setting of pixel numbers, refer toRefs. [25, 26].

4. DERIVATION OF INTEGRATEDINTENSITIES FROM POWDERDIFFRACTION DATAIn the MEM analysis of relatively simple structuresfrom powder diffraction data, |Fok | is usually obtained by pattern decomposition such as individual
profile fitting [27], the Pawley method [28], and the
Le Bail method [29]. However, most values of |Fok | forcompounds exhibiting heavily overlapped reflections
can hardly be estimated by pattern decomposition
without introducing any structural models.
After Rietveld analysis, the ‘observed’ integrated
intensity, Iok , for reflection k is approximately computed by apportioning the net diffraction intensity,yi − yb (2θi ), at step i in the ratio of contributions ofoverlapping reflections (K = 1, 2, 3, · · ·) [5]:
Iok =

Xh

i

Yik
,
K YiK

yi − yb (2θi ) P

i

(44)

where yi is the observed intensity; yb (2θi ) is the background function; Yik and YiK are the contributions ofreflections k and K to the net diffraction intensity atPP
step i, respectively. The summations i and K are
carried out over all the data points contributing to
the profile of reflection k and over all the reflections
contributing to yi − yb (2θi ), respectively. In the casesof isolated reflections, Eq. (44) reduces to
Iok =

Xh

i

yi − yb (2θi ) .

(45)

i

Consequently, fairly accurate Iok is obtainable in amodel-free fashion as long as the background can beadequately approximated by yb (2θi ). Equations (44)and (45) are also used to derive Iok after each cycle
in Le Bail refinement [29].
Yik is evaluated from scale factor s, multiplicity mk ,preferred-orientation factor Pk , Lorentz-polarization
factor Lk , structure factor calculated from structure
parameters Fck , and profile function f (2θi − 2θk ):Yik = smk Pk Lk |Fck |2 f (2θi − 2θk ).

(46)

Substituting K for k in this equation yields YiK . Fckis calculated from structure parameters (fractional
coordinates, occupancies, and isotropic/anisotropic
atomic displacement parameters) refined in Rietveld
analysis. The observed scattering amplitude, |Fok |, isestimated with
Iok
|Fok | =smk Pk Lk


1/2

.

(47)

Takata et al. [3, 4] applied the above procedure(MEM/Rietveld method) to the estimation of |Fok |’sanalyzed by MEM. However, electron densities derived from the |Fok |’s are subject to the structuralmodel in Rietveld analysis because Eq. (46) contains
Fck . This undesirable effect enlarges with increasing
degree of overlap of reflections, lowering the accuracy of final densities necessarily. MEM can extract
structural details from the |Fok | data estimated in
the above way because they involve contributions neglected in the structural model. Furthermore, |Fok |for an isolated reflection can simply be determined
without any structural model, as described earlier.
Visualization of diffraction data by the combinationof Rietveld and MEM analyses, therefore, serves themodification of imperfect structural models [30–33].
Nevertheless, the dependence of electron-density distribution on the structural model casts a gloom overthe combined use of the Rietveld and maximumentropy methods.

5. PROCEDURE OF MEM-BASEDPATTERN FITTING
To minimize the bias imposed on electron densities
by the structural model, we devised the following iterative procedure of MPF [8, 9]:
1. Intensity data are analyzed by the Rietveld
method.
2. |Fok |’s are approximately evaluated with Eqs.
(44)−(47).
3. The |Fok |’s with phases, ψ, derived from the Rietveld analysis, are analyzed by MEM to yield
3D electron densities, ρ(x, y, z). Then, FMk iscomputed by the Fourier transform of ρ(x, y, z)
in the unit cell:
FMk = V

ZZZ

ρ(x, y, z) exp[2πi(hx + ky + lz)]dx dy dz.
(48)

4. The resulting ρ(x, y, z) data are visualized inthree dimensions.
5. After close examination of the electron-density
image, return to step 1 to modify the structuralmodel if necessary.
6. The model function is fit to the whole observed
diffraction pattern by refining parameters irrelevant to the crystal structure: scale factor, pro-file parameters, peak-shift parameters, and background parameters. Each structure factor in themodel function is not calculated from any structure parameters but fixed at FMk .
7. Terminate unless decreases in R factors in step
6 are significant compared with those in the previous whole-pattern fitting.
8. Return to step 2.
Steps 2−8 repeated till the convergence is attained
is referred to as REMEDY cycles, serving to minimize
the bias imposed by the structural model in Rietveld
analysis. The influence of the structural model in
the Rietveld analysis on ρ(x, y, z)’s diminishes by repeating REMEDY cycles. In other words, intensity
partitioning for overlapping reflections becomes moreaccurate with increasing number of iterations owing
to derivation of additional structural information by
MEM. Such an iterative method is somewhat similar
to the Le Bail method [29] but differs from it in the
point that structure factors in whole-pattern fitting
are fixed at values of FMk resulting from the previous
MEM analysis.
Figure 1 illustrates the above sequence schematically. ‘Rietveld’, ‘MEM’, and ‘w.p.f.’, each in a pair ofparentheses, denote analyses whereby structure factors are derived. Fo (Rietveld) and Fo (w.p.f.) in thisfigure correspond to Fok , and Fc (MEM) is equivalent
to FMk . Contributions of anomalous dispersion to
structure factors must be excluded in MEM analysis
in the same way as in Fourier synthesis whereas theyshould be included in whole-pattern fitting. In otherwords, structure parameters are used only to add
anomalous dispersion corrections to structure factors
in MPF.
Structural details are changed in step 3 while thegoodness-of-fit is improved by fixing each structure
factor at FMk in step 6. The effect of the structural model in the Rietveld analysis on FMk diminishes with increasing number of REMEDY cycles. We
can therefore approach the final structure reflecting
the observed intensity data more closely. Step 6 is
not Rietveld analysis but MPF to obtain values of
Fok containing additional information about structural details.

Powder diffraction data

Initial structural model

Rietveld analysis
RIETAN-2000

Revisedstructural modelFo (Rietveld)

Visualization ofdensity distribution

MEM analysis
PRIMA

Electron/nuclear
densities

VENUS

6. SOFTWARE FOR MPF

Fc (MEM)i = 1, 2, .....
Whole-pattern fitting
RIETAN-2000

Fo (w.p.f.)but also by 3D electron/nuclear densities. The establishment of MPF leads to the availability of an ultrahigh-resolution ‘X-ray/neutron microscope’ that allows us to visualize powder diffraction data as 3D
density images. X-Ray and neutron diffraction should
be utilized complementarily for (a) the more adequate
expression of chemical bonds and (b) the analysis ofanharmonic thermal vibration and disordered structures, respectively.

REMEDYcyclesFig. 1 Structure refinement by MEM/Rietveld
analyses followed by iterative MEM-based
pattern fitting.
Takata et al. [3, 4] have been utilizing not all theprocesses shown in Fig. 1 but only the part of iterative MEM/Rietveld analyses (inside the upper box
consisting of broken lines) for modifying structuralmodels. As described in Sect. 4, the MEM/Rietveld
method is a tepid approach to the determination ofelectron/nuclear densities because no structural information extracted by MEM from powder diffraction data is reused by executing w.p.f. after the MEM
analysis. The MEM/Rietveld method should be used
only to modify structural models imperfect with respect to sites (defects) unfindable and overlooked in
the models.
MPF is capable of modeling disordered atomic configurations [32,34], chemical bonds, nonlocalized electrons, and anharmonic thermal motion [36] more
adequately than the conventional Rietveld method.
These purposes are attainable more satisfactorily by
MPF than with programs for conventional Rietveld
analysis because crystal structures are represented
not by structure parameters (fractional coordinates,
occupancies, and atomic displacement parameters)

We have developed a state-of-the-art system for structure refinement by MPF. This free software comprises
a multi-purpose pattern fitting program, RIETAN2000, for angle-dispersive X-ray and neutron diffraction [9,11], MEM analysis program and PRIMA [12],and VENUS [14] for 3D visualization of crystal structures and electron/nuclear densities.1
The goodness-of-fit must reach a high level to minimize the influence of the structural model on densities resulting from MPF. We have been utilizing ouroriginal technique of partial profile relaxation [11] inboth Rietveld analysis and whole-pattern fitting. At
present, MPF would not function well without thisfeature, which was certainly effective in MPF structure refinements described later.
MEED [13] is now obsolete; our original program,
PRIMA, can carry out MEM analyses incomparably
faster than MEED thanks to the full optimization of
its code, as described in Sect. 3. RIETAN-2000 provides PRIMA with structural information, and viceversa, through two files storing Fok data (*.mem) and
FMk data (*.fba).
VENUS bears the bell in 3D visualization, rendering, and manipulation of crystal structures andelectron/nuclear densities determined not only by Xray/neutron diffraction but by electronic-structure
calculations using Gaussian 98 [37], SCAT [38], and
WIEN2k [39]. VENUS fully utilizes the OpenGL
Application Programming Interface (API) with thegraphical user interface constructed with GLUT and
GLUI. Video cards equipped with GPU’s supporting
hardware acceleration for the OpenGL are desired to
1

RIETAN-2000, PRIMA, and VENUS can be downloaded
from http://homepage.mac.com/fujioizumi/index.html.

rotate, expand, and translate objects fast in three dimensions.
The VENUS system consists of two independent programs: VICS (VIsualization of
Crystal Structures) and VEND (Visualizationof Electron/Nuclear Densities). That is, VENDand VICS deal with electron/nuclear densities andcrystal structures, respectively VEND also supports
the 3D visualization of electrostatic potentials
and wave functions that have both positive andnegative values, which allows us to obtain additional
information about electronic states.

7. DATA USED IN MPF AND ITSKNOW-HOW
7.1 Three reflection types
The algorithm for the MEM analysis of Fok ’s is,
per se, the same regardless of diffraction procedures.
Some special techniques are, however, required to analyze Fok ’s derived from powder data. The MEM
deals with |Fok |’s of reflections whose phases havebeen estimated by some means. A reflection whose|Fok | is computed from Eqs. (44)−(47) at the end of
Rietveld analysis is referred to as a type 1 reflection.
Fok ’s of overlapped reflections evaluated according to
Rietveld’s procedure [5] are influenced by the refinedstructure in a varying degree. The sum of |Fok |’s isthen assigned to a group of overlapped reflections belonging to type 2 to enhance the ability of estimation
by MEM. In the case of type 3 reflections, their peak
positions are so near to 2θmax that parts of their profiles extend beyond 2θmax . Nevertheless, we should
also calculate the FMk of a type 3 reflection because
excluding them necessarily brings about a bad profile
fit near 2θmax . FMk for the type 3 reflection is thus estimated from electron densities determined from only
Fok for reflections of types 1 and 2. Type 2 reflections are decomposed to afford their individual Fok ’s.
Fok ’s, sums of Fok ’s, and only hkl indices are output
for reflections with types 1, 2, and 3, respectively.
7.2 Structural data input/output during MPF
Prior to MEM analysis, we need to estimate Fok ’s andthe e.s.d.’s of |Fok |’s: σ(|Fok |) at any rate. It should beemphasized that the Fok and σ(|Fok |) data of all the
low -Q (large-d) reflections are indispensable to MEM
analysis. The lack of these data is fatal, necessarily
leading to a failure in the MEM analysis.
Let σ(Iok ) and σ(s) be the e.s.d.’s of Iok and s,
respectively, then the estimated standard deviation(e.s.d.), σ(|Fok |), of |Fok | for a type 1 reflection isestimated from"

|Fok |  σ(Iok )σ(|Fok |) =2  Iok#2#2 1/2σ(s) 
+
.

s


"

(49)

Note that σ(|Fok |) is not actually observed but calculated from |Fok | and Iok , suffering from errors dueto the bias toward a structural model in Rietveld
analysis in the same manner as with |Fok |. Further,Eq. (49) is based on only counting statistics, neglecting the overlap of reflections. We must, hence, keep
in mind that σ(|Fok |) values estimated from powder
diffraction data are not very accurate.
The sum of |Fok |’s, Go , for a group of type 2 reflections is given by
P

Go =

k Iok

s

!1/2

P

k mk Pk Lk

(50)

P

with the summation k carried out over all the overlapped reflections. The e.s.d. of G0 is estimated by
" P
#2 "
#2 1/2
Goσ( k Iok )σ(s) 
σ(Go ) =
+
P

2 
s
k Iok

(51)

analogous to Eq. (49).
At the end of the whole-pattern fitting, Iok is evaluated in a similar manner to Eq. (44):
Iok =

Xh
i

i Y0
yi − yb (2θi ) P ik 0 .
K YiK

(52)

Yik0 is the contribution of reflection k to the net
diffraction intensity and given by
Yik0 = smk Pk Lk |FMk |2 f (2θi − 2θk ).

(53)

In this equation, Fck in Eq. (46) is replaced with FMk .
On substitution of K for k in this equation, we obtain Yik0 . |Fok | and σ(Iok ) are respectively computedwith Eqs. (47) and (14) in the same manner as theircalculations after the Rietveld analysis [1]. Because

overlapped reflections are decomposed during the first
MEM analysis without any structural model, they are
regarded as type 1 reflections in subsequent processes.
When analyzing neutron diffraction data or X-ray
diffraction data for centrosymmetric structures, |Fok |for a profile-relaxed reflection can directly be refinedin whole-pattern fitting with RIETAN-2000. This
unique feature is expected to accelerate the convergence of structure refinement, provided that calculated profiles are in good agreement with observed
ones for the relaxed reflections.
7.3 Distribution of (Gauss)
The use of the χ2 statistics in the constraint function(5) is based on the assumption that the experimental
errors, σ(|Fok |), are random with a Gaussian distribution:|FMk | − |Fok |(Gauss) =,(54)
σ(|Fok |)where (Gauss) is a sample of the random variablewith normalized Gaussian distribution [40]. No satisfactory result can be obtained when the distributionof (Gauss)’s is far from Gaussian. PRIMA has a feature to output a text file storing (Gauss) values witheither of two formats: a raw data file, *.raw, and an
Igor text file, *.err. We can easily examine the distribution of (Gauss) using this file.
7.4 Strategy for stable convergence
We should pay attention to electron-density maps
obtained from FMk data in addition to the convergence test in each iteration. There is particularly
a fair chance for structural estimation to be misdirected at the first iteration. Accordingly, we need torepeat steps 2−8 while comparing density maps resulting from MEM analysis with those obtained from
Fck ’s. This comparison is needed to judge whether
the electron-density distribution is affected by the
structural model or estimated properly by MEM. For
reflections with unsatisfactory fits between observedand calculated intensities, e.s.d.’s obtained with Eqs.
(49) and (51) should appropriately be increased to enhance the degree of freedom in estimation by MEM;
such large residuals usually arise from imperfect representation of the real structure on condition that theintensity data are collected appropriately. If this partwere partitioned among overlapped reflections in proportion to their values of YiK , MEM would estimate
the structure in a way biased in favor of the structuralmodel.
The efficiency of REMEDY cycles depends on thevalidity of the starting structural model. If the structure refined in step 1 differs considerably from thereal one, MEM analysis may fail in structure refinement owing to the inaccurate partition of observed
Bragg intensities among overlapped reflections.
R factors in the whole-pattern fitting are used as
the measure of the convergence in REMEDY cycles.
R factors in model-free MEM analysis are comparable regardless of the progress of the REMEDY cycles
because the convergence criterion is a constraint function being equal to unity. Accordingly, R factors for
MEM, Eqs. (36) and (37), should not be used to judge
the convergence of the iterations.

8. EXAMPLES OF APPLYING MPF
TO VARIOUS MATERIALS
8.1 Positional disorder of K+ ions interlayered
in Kx Ti2−x/3 Lix/3 O4
A layered compound, Kx Ti2−x/3 Lix/3 O4 (space group:Cmcm) contains (a) lepidocrocite-like Ti2−x/3 Lix/3 O4
layers consisting of edge-sharing (Ti,Li)O6 octahedra
and (b) charge-balancing K+ ions in the interlayer domain [41]. Rietveld analysis from neutron diffraction
data showed K+ ions to be highly disordered between
two layers. With a split-atom model, the positional
disorder of K+ ions could not adequately be expressedin view of an extraordinarily large isotropic atomic
displacement parameter, B, of ca. 10 Å2 for K.
We determined electron densities in Kx Ti2−x/3 Lix/3 O4 (x = 0.8) to visualize how K+ ions are distributed in between Ti2−x/3 Lix/3 O4 layers [42]. Its Xray powder diffraction data were measured at room
temperature with Cu Kα radiation on a MAC Science
MXP 3TZ diffractometer with a vertical θ:θ configuration and a pair of long Soller slits with an aperture angle of 1◦ . The same diffractometer was usedin X-ray powder diffraction experiments that will be
described in 8.1−8.3.
The diffraction data were analyzed by the Rietveld

differences in the oxidation states and effective ionicradii, r, between Ti4+ (r = 0.605 Å) and Li+ (r =0.76 Å) ions sharing the same octahedral site. Thus,
our first challenge to the determination of spatial distribution for highly disordered species was quite successful.
MPF is expected to be effective for intercalationand clathrate compounds where guests have often
disordered configurations. In these compounds, splitatom models may fail in providing reasonable atomic
displacement parameters or sufficiently low R factors.
We will further apply MPF to such systems extensively.

b

K+
O2

8.2 Atomic arrangement of K clusters in a Kdoped zeolite K-LTA(Ti,Li)O1a0cFig. 2 Electron-density distribution determinedby MPF for Kx Ti2−x/3 Lix/3 O4 (x = 0.8).
Isosurfaces in three unit cells were drawn
for an equidensity level of 1.2/Å3 .
method on the basis of space group Cmcm [41]. We
assumed the random occupation of the same octahedral site (4c) by Ti4+ and Li+ ions. Final R factors [43] were Rwp = 11.51 % (S = Rwp /Re = 1.08)
and RB = 5.75 %. Subsequent three REMEDY cycles lowered Rwp to 10.59 % (S = 1.00) and RB to
1.91 %. The dramatic decrease in RB is ascribable tothe better representation of the structural details by
calculating structure factors from electron densities
than from structure parameters in Rietveld analysis.
Figure 2 shows isosurfaces of electron densities(number densities of electrons) in Kx Ti2−x/3 Lix/3 O4(x = 0.8). Neither negative densities nor distinct
ripples appeared at all, which is a great advantage of
MEM over Fourier synthesis. A winding stream of K+
ions along the [001] direction and (Ti,Li)−O bonds
were clearly visualized by virtue of MPF. The positional disorder of K+ ions is partly ascribable to large

The synthetic zeolite, Linde Type A (LTA), consistsof α- and β-cages arrayed alternately like atoms in
the CsCl-type compound with inner diameters of 11Å and 7 Å, respectively. Alkali metal ions are usually incorporated into the cages to neutralize negative charges of the LTA framework. K-type LTA(K-LTA) doped with excess K (Kx @K-LTA) is of interest in connection with its ferromagnetism in theabsence of magnetic atoms below 7 K [44]. Because
the excess K atoms release 4s electrons to change
into K+ ions, Kx @K-LTA is regarded as a stronglycorrelated electron system. To understand the mechanism of the anomalous magnetic property, nuclear
and electron densities in Kx @K-LTA were determinedby MPF from neutron and X-ray diffraction data, respectively [9, 45]. The neutron diffraction data of
K5.8 @K-LTA were measured at 40 K on HRPD at
JRR-3M (JAERI) with λ = 1.824 Å, and the X-ray
diffraction data of K5.2 @K-LTA at room temperature.
Rietveld refinements revealed that space groupFm 3̄c in K-LTA changes to F 23 with lower symmetryin Kx @K-LTA [46]. However, the disorder of K atoms
cannot sufficiently be represented in these structurerefinements. After one cycle of MEM-based wholepattern fitting with the neutron diffraction data, Rwp
slightly decreased from 5.03 % (S = 1.48) to 4.92
% (S = 1.45), and RB from 2.05 % to 1.91 %, respectively. i.e., from 8.25 % (S = 1.08) to 7.30 %(S = 0.96) for Rwp On the other hand, the goodness-

c

K

K

(a)

c

cupy the center of an α cage. K+ ions near the centerof the six-membered ring form two kinds of tetrahedral assemblies occupied statistically. Adjacent αcages contain nonequivalent K clusters with different
magnetic moments. The difference in the number of
K atoms between the two types of the α-cages was ca.
2 in a sample with the highest degree of spontaneous
magnetization, which will give rise to the ferromagnetism with antiferromagnetic components. Adjacentα-cages contain nonequivalent K clusters with different magnetic moments. The difference in the numberof K atoms between the two types of the α-cages was
ca. 2 in a sample with the highest degree of spontaneous magnetization, which will give rise to the ferromagnetism with antiferromagnetic components.
8.3 Nuclear/electron-density distribution in
the superconductor HgBa2 CuO4+δ

(b)
Fig. 3 (a) Nuclear- and (b) electron-density distribution on the (110) plane in K5.8 @K-LTA
and K5.2 @K-LTA, respectively. Contours
were plotted for (a) 0.3 × 2.1n fm/Å3 (n= 1−9) and (b) 0.3−10/Å3 with a step of
0.4/Å3 .
of-fit was improved more significantly with X-ray
diffraction data, and from 2.70 % to 1.30 % for RB
after two REMEDY cycles. Such differences are believed to reflect (a) much wider spatial distributionof electrons in comparison with atomic nuclei and (b)the highly covalent nature of Si−O and Al−O bondsin the framework.
Nuclear and electron densities on the (110) plane
are plotted in Fig. 3, where each quarter corresponds
to an α-cage. K+ ions are situated in the neighborhood of 4, 6, and 8-membered rings. Additional
K+ ions are mainly included in the four- and eightmembered rings of the framework and partially oc-

We prepared a high-Tc superconductor HgBa2 CuO4+δ
(Tc = 97 K; space group P 4/mmm) where the content
of CO3 2− ions substituting for [O2−Hg−O2]2− ionswas minimized to distribute excess oxygen atoms very
uniformly by O2 annealing for a long time [47]. Its
neutron and X-ray diffraction data were collected atroom temperature on HRPD (λ = 1.823 Å) and the
MXP 3TZ diffractometer, respectively.
In structure refinements from both the neutron and
X-ray powder diffraction data, MPF afforded lower
R factors than conventional Rietveld analysis [9, 48].
The decrease in RB achieved in the MPF analysisof the X-ray diffraction data was particularly pronounced, i.e., from 1.66 % in the final Rietveld analysis to 0.89 % after two REMEDY cycles. No bonding
electrons are explicitly taken into account in the Rietveld method whereby structure parameters are refined. On the other hand, MEM analysis enables us
to express chemical bonds with electron densities in ahigher degree of freedom, which is responsible for theclose approach of Fc (MEM) towards Fo (w.p.f.). In
fact, densities of bonding electrons between metal andoxygen atoms varied appreciably during the REMEDY cycles. Figure 4 shows substantial changes inelectron densities corresponding to Cu−O1 bonds on
the CuO2 conduction sheet in three refinement stages.
This fact offers unambiguous evidence for introduction of additional structural information extracted by

Cu

O1

O1

a
b

(a)

(b)

(c)

Fig. 4 Electron-density distribution on the CuO2 sheet (z = 1/2 plane) in HgBa2 CuO4+δ after (a) Rietveld
analysis, (b) first REMEDY cycle, and (c) second REMEDY cycle. Contours were plotted for 0−3/Å3with a step of 0.3/Å3 .

MEM into intensity partitioning for overlapping reflections.
Figure 5 gives (a) nuclear and (b) electron densities on the (100) plane of HgBa2 CuO4+δ . Interposition of bonding electrons between Hg and O2 andbetween Cu and O1 is clearly recognized in Fig. 5b.
The presence of a region where electron densities are
negligibly low between the Cu and O2 atoms indicatesnearly electrostatic interaction between them. Comparison of Fig. 5a with Fig. 5b reveals that thermal
vibrations of Cu and Hg predominate along the [001]direction perpendicular to the CuO2 sheet and on thez = 0 plane perpendicular to a linear [O2−Hg−O2]2−
ion, respectively. These findings are reasonable considering the highly covalent character of the Cu−O1
and Hg−O2 bonds.
Though harmonic thermal motion was assumed in
the Rietveld analysis prior to the REMEDY cycles,
the isosurface of nuclear densities for O2 considerably
deviated from an ellipsoid, as can be appreciated from
Fig. 5a. In fact, the analysis of anharmonic thermal
vibration [36] for O2 afforded significant cubic andquartic anharmonic terms.
Interstitial oxygen defects, O3, at (1/2, 1/2, 0) areresponsible for hole doping in HgBa2 CuO4+δ . In conventional Rietveld analysis, B for a slightly occupied
site like O3 is often fixed arbitrarily at a typical value,e.g., 1 Å2 , to refine only its occupancy, g, on accountof very strong correlation between B and g. We took
more sophisticated approaches to obtain more reliable values of g(O3) and B(O3). That is, g(O3) was
estimated at 0.14 by integrating nuclear densities of
O3, and B(O3) at 1.5 Å2 by evaluating mean squaredisplacements. The g(O3) value of 0.14 corresponds
to δ = 0.14. This technique is useful for determining B and g for partially occupied sites from nuclear
densities derived by MEM.
8.4 Statistic distribution of H atoms in thehigh-temperature phase of KH2 PO4
The order−disorder-type phase transition in KH2 PO4
was studied by MPF from its neutron powder diffraction data measured at room temperature and 10 K
on HRPD (λ = 1.16 Å) at JRR-3M [49]. The procedure described in 3.4 was used to deal with H atoms
with a negative bc of −3.7390 fm [2]. The ferroelectric (space group Fdd 2) and paraelectric (space groupI 4̄2d ) phases of KH2 PO4 are stable below and above
its Curie temperature, TC , of 123 K, respectively.

(b)

(a)
Hg
O2

c

Cu

O1

b
Fig. 5 (a) Nuclear- and (b) electron-density maps for the (100) plane in HgBa2 CuO4+δ . Contours wereplotted for (a) 0.1×3n fm/Å3 (n = 1−10) and (b) 0−3/Å3 with a step of 0.3/Å3 .

In Fig. 6, nuclear densities resulting from MPFof the paraelectric phase are projected onto the (001)plane. This projection makes it possible to see the nuclear densities of all the atoms included in KH2 PO4 .
H atoms are statistically situated at a pair of positions with a distance between the twin peaks of ca.
0.36 Å, which is in good agreement with H−H distances reported in the literature. On the other hands,
H atoms in the ferroelectric phase at 10 K tend to occupy a single position near to an oxygen atom. Figure
6 illustrates that the thermal vibrations of P and O
atoms in the PO4 3− ion are of the comparable degree
but of different directional character. This finding
supports the idea that the thermal vibration of thePO4 3− ion is dominated by rigid-body motion, still
exhibiting significant internal motion.
The MEM/Rietveld method failed in visualizing
the symmetric potential valleys in the paraelectricphase, which provides evidence for the predominanceof MPF over the MEM/Rietveld method.
8.5 Orientational disorder of CO3 2− ions in
Ba4 CaCu2 O6 CO3
The tetragonal structure of Ba4 CaCu2 O6 CO3 (spacegroup P 4/mmm) can be related to the perovskitetype compound AB O3 by locating larger Ba atoms
at an A site and smaller Ca, Cu, and C atoms at B
sites. Ca atoms and Cu1O4 squares lie on the z = 0
plane, Ba and O2 atoms on the z = 1/4 plane, andCu2 (1b; 0, 0, 1/2) and C (1d; 1/2, 1/2, 1/2) atoms
on the z = 1/2 plane. Of course, three O atoms forma triangular CO3 2− ion with the C atom at its center.
Rietveld analysis from neutron powder diffraction
data measured at room temperature [50] showed that
one O3 atom is split into four pieces to occupy a 4k
site on a line connecting C and Cu2 atoms. However,the static orientational disorder of CO3 2− ions is so

a whole. Thus, MPF released us from the preconception that C fully occupies the 1d site.

O

P

010

8.6 Dynamic disorder of Cu+ ions in a superionic conductor Rb4 Cu16 I7.2 Cl12.8

H
K

O

1Å
Fig. 6 Nuclear-density distribution projectedon the (001) plane in the paraelectricphase of KH2 PO4 at room temperature.
Contours were plotted for 2.5 × 1.5n
fm/Å3 (n = 0 − 20).
pronounced that positions of remaining two O4 atoms
could not be determined unambiguously.
We have recently reanalyzed neutron powder
diffraction data measured on the BT-1 diffractometer at NIST by MPF and successfully visualized what
was overlooked in the previous Rietveld analysis (Fig.
7) [51]. MPF decreased Rwp (S), Rp , RB , and RF
from 5.17 % (1.32), 4.14 %, 5.59 %, and 3.05 % in
the final Rietveld analysis to 4.26 % (1.08), 3.46 %,2.35 %, and 1.78 %, respectively.
In addition to O3 and O4 (2h; 1/2, 1/2, z; z ≈0.43) atoms distant from the 1d site by ca. 1.29 Å,
we found four positions apart from the 1d site byca. 1.85 Å along the [110] direction (Fig. 7a). Because 1.85 Å corresponds to a Cu−O bond length,
Cu3 atoms are believed to substitute for part of C
atoms at 1d. Four O5 atoms to which each Cu3 atom
at 1d is coordinated are clearly visualized by MEM.
In Fig. 7b, peaks for O5 are much higher than those
for O3, which is explained in terms of high densitiesof O3 nuclei in the space excluding the z = 1/2 plane.
The content of O3 is much higher than that of O5 as

The conductivity of Rb4 Cu16 I7 Cl13 (space group
P 41 32) at room temperature is 0.34 S cm−1 , whichis the highest of all the superionic conductors. Its
structure changes without any phase transitions froma low-temperature phase where Cu+ ions hardly migrate to the high ionic conductor because lattice defects are generated by thermal excitation. Figure
8 illustrates results of MPF structure refinements
for Rb4 Cu16 I7.2 Cl12.8 at 11 K and 280 K from neutron powder diffraction data measured on the HB-4
diffractometer at HFIR [34]. In the case of the datataken at 280 K, two REMEDY cycles lowered Rwp(S) from 6.06 % (1.05) to 5.98 % (1.03) and RB from6.25 % to 5.18 %.
The most significant difference between Figs. 8aand 8b is that Cu+ ions are localized at definite positions at 11 K but delocalized spirally at room temperature. Rb4 Cu16 I7 Cl13 has three Cu sites in theasymmetric unit: Cu1 and Cu2 at 24e and Cu3 at
8c. From Fig. 8b, we confirmed that [Cu1−Cu2−]∞
chains are main conduction pathways in this material, as previously suggested from Rietveld analysisof time-of-flight neutron powder diffraction data [35].
Visualization of distribution of Cu+ ions in the unit
cell with VEND is very useful to understand the positional disorder in this system three-dimensionally,
as can be appreciated from Fig. 8.
8.7 Positional and orientational disorder in asolid solution of Sr9.3 Ni1.2 (PO4 )7
To apply MPF to high-resolution powder diffraction data, we measured synchrotron X-ray powder
diffraction data on a powder diffractometer (beamline BL15XU at SPring-8) with the Debye−Scherrer
geometry using capillary tubes rotated at a speedof 6.3 rad/s. Incident beams from an undulator
were monochromatized with inclined double-crystal
monochromators of Si(111). An instrumental resolution at the top international level is attained in thisdiffractometer.

c
Cu2
C/Cu3
Cu2

O3

O5

O5
O3

O2
O4

Ba
Ca

0

O1

Cu1

a

(a)

(b)

Fig. 7 (a) Isosurfaces of nuclear densities at a 0.6 fm/Å3 level and (b) a bird’s eye view of nuclear
densities up to 8 % of the maximum (663.9 fm/Å3 ) on the z = 1/2 plane in Ba4 CaCu2 O6 CO3 .

Cu1
Cl
I

Cl
I

Cl

I
Rb

Cu2

Cu3
I

(a)

(b)

Fig. 8 Isosurfaces of nuclear densities in Rb4 Cu16 I7.2 Cl12.8 at (a) 11 K and (b) 280 K. Equidensity
level: 0.2 fm/Å3 . In (b), sections at unit-cell edges clearly display [Cu1−Cu2−]∞ pathways.

8.8 Structure refinement of LTL by MPF

P1

Fig. 9 Three-dimensional electron-density
image of a P1O4 3− ion in Sr9.3 Ni1.2 (PO4 )7 with an equidensity level of
1.6/Å3 .

A dehydrated zeolite, Linde Type L (LTL,
K9 Al9 Si27 O72 , space group P 6/mmm) [52], is a microporous material with a disordered structure including
large one-dimensional channels parallel to the c axis.
Its synchrotron X-ray powder diffraction data (λ =1.2000 Å) were collected at room temperature with
the Ge(111) analyzer. Rwp (S), Rp , RB , and RF wererespectively 4.90 % (2.47), 3.19 %, 1.60 %, and 1.59
% in the final Rietveld refinement and 4.46 % (2.24),
2.93 %, 0.69 %, and 0.84 % after two REMEDY cycles. Figure 11 shows the final profile fit to be very
excellent. Electron-density images clearly showed aframework consisting of highly covalent (Si,Al)−O
bonds. We also found that electron-density distribution around K+ ions incorporated in the framework
is not spherical but anisotropic.

9. CLOSING REMARKS
X-Ray powder diffraction data of a whitlockite-like
phosphate, Sr9.3 Ni1.2 (PO4 )7 (space group R 3̄m) [32],were taken at 100 K with a wavelength of 0.8000 Å
and a Ge(111) analyzer. Very sharp profiles were observed with FWHM’s ranging from 0.008◦ to 0.019◦ .
In Sr9.3 Ni1.2 (PO4 )7 , part of Sr2+ and PO4 3− ionsexhibit highly disordered arrangements. Then, we
found positions of these atoms by the MEM/Rietveld
method, building a split-atoms model. Subsequent
two REMEDY cycles decreased Rwp (S), Rp , RB , andRF a great deal from 5.39 % (2.30), 3.93 %, 2.80 %,and 1.88 %, in the final Rietveld refinement adopting the split-atom model to 5.18 % (2.20), 3.71 %,1.03 %, and 0.82 %, respectively. As usual, RB andRF dropped dramatically in the first cycle. Figure
9 illustrates isosurfaces of electron densities determined by MPF for P1O4 tetrahedron. This figure
revealed that such pronounced orientational disordercan hardly be represented with the split-atom model,
which is the main reason why MPF lowered the R
factors surprisingly. MPF afforded an excellent fit
between observed and calculated patterns (Fig. 10).

These results proved that the MPF method is effective even in the structure refinement from theultrahigh-resolution X-ray diffraction data.

The revolutionary MPF technology was devised toovercome limitations of Rietveld analysis. RIETAN2000 is distinguished from other programs for classical Rietveld analysis by this elegant methodology. In
MPF, electron/nuclear-density distribution changes
noticeably with accompanying decreases in RB and
RF during REMEDY cycles, which is a strong piece ofevidence for the reduction in the bias of the structuralmodel and the serious flaw of the MEM/Rietveld
method on its use to determine electron/nuclear densities [3, 4].
MEM is, per se, model-free whereas the MEM/
Rietveld method is more or less model-dependent.
The MEM/Rietveld method should never be confounded with MEM. This expedient method is, therefore, unsound on its application to determination ofelectron/nuclear-density distribution though it is certainly effective for the rough visualization of electron/nuclear densities just after Rietveld analysis.
We agree that the MEM/Rietveld method is more
effective in modifying imperfect structural models
than Fourier/D synthesis thanks to higher S/N ratios due to estimation of structure factors for unobserved reflections and non-negative densities obtained by MEM. However, it cannot give accurate 3D
densities particularly when reflections heavily overlap

300000

Intensity

250000200000150000100000
50000
0
10

15

20

25

2 θ / deg

30

35

40

45

Fig. 10 Observed (+ symbols), calculated (upper solid line), and difference (lower solid line) patterns resulting from the final MPF structure refinement from the synchrotron X-ray powder diffraction
data for Sr9.3 Ni1.2 (PO4 )7 . Tick marks denote the peak positions of possible Bragg reflections.

60000040000
50000030000

Intensity

40000020000

30000010000
0200000100000

10

15

20

25

30

35

40

45

50

55

60

10

15

20

25

30

35

40

45

50

55

60

0

5

2 θ / deg

Fig. 11 Observed, calculated, and difference patterns resulting from the final MPF structure refinement from the synchrotron X-ray powder diffraction data for LTL. The
patterns between 10◦ and 60◦ are magnified in the inset.

with each other because of low symmetry and insufficient instrumental resolution. Frankly speaking, the
MEM/Rietveld method is nothing but iterative conventional Rietveld refinements where structural models are modified, if necessary, after visualizing electron/nuclear densities by MEM.
The number of REMEDY cycles required for convergence generally increases with lowering resolution
in a diffraction pattern because observed intensities
of Bragg reflections are repartitioned with Fc (MEM)
calculated from Fo (Rietveld) or Fo (w.p.f.). Even withhigh-resolution synchrotron X-ray powder diffraction
data, at least one cycle of MPF proved to be indispensable for the determination of more accurate electron densities. In MPF, parameters other than structure parameters are refined, and partition of observed
intensities is rapidly improved because of the use ofFc (MEM). The MPF method is, hence, so much superior to the MEM/Rietveld method, where the analytical method in the Fourier synthesis of Fo (Rietveld)
data is merely replaced with MEM. From the amazing results presented above for the several materials,
we conclude MPF to be a powerful method of structure refinement where imperfect structural representations in conventional Rietveld analysis are supplemented with model-free MEM.

ACKNOWLEDGMENTS
The technology of MPF owes its existence to Shintaro
Kumazawa’s proposal to implement it in RIETAN2000. The development of both RIETAN-2000 and
VENUS was financially supported by the Project forSuperconductivity Research at NIMS.

References[1] Sakata, M., Mori, R., Kumazawa, S., Takata,M., Toraya, H. 1990, J. Appl. Crystallogr. 23,
526.
[2] Sears, V.F. 1999, International Tables for Crystallography, Vol. C, 2nd ed., A.J.C. Wilson,
E. Prince (Eds.), Kluwer Academic Publishers,
Dordrecht, pp. 440−450.
[3] Takata, M., Umeda, B., Nishibori, E., Sakata,
M., Saito, Y., Ohno, M., Shinohara, H. 1995,
Nature (London) 377, 46.

[4] Takata, M., Nishibori, E., Sakata, M. 2001, Z.
Kristallogr. 216, 71.
[5] Rietveld, H.M. 1969, J. Appl. Crystallogr. 2, 65.
[6] McCusker, L.B., Von Dreele, R.B., Cox, D.E.,
Louër, D., Scardi, P. 1999, J. Appl. Crystallogr.
32, 36.
[7] Izumi, F., Kumazawa, S., Ikeda, T., Ida, T.
1998, Powder Diffraction, Proc. Int. School
Powder Diffr. (ISPD’98), S.P. Sen Gupta (Ed.),
Allied Publishers, New Delhi, pp. 24−36.
[8] Izumi, F. 2000, Rigaku J. 17, No. 1, 34.
[9] Izumi, F., Kumazawa, S., Ikeda, T., Hu, W.-Z.,
Yamamoto, A., Oikawa, K. 2001, Mater. Sci.
Forum 378−381, 59.
[10] Izumi, F., Ikeda, T. 2001, Commission on Powder Diffraction, IUCr, Newslett. No. 26, 7.
[11] Izumi, F., Ikeda, T. 2000, Mater. Sci. Forum
321−324, 198.
[12] Dilanian, R.A., Izumi, F. 2003, unpublishedwork.
[13] Kumazawa, S., Kubota, Y., Takata, M., Sakata,M., Ishibashi, Y. 1993, J. Appl. Crystallogr. 26,
453.
[14] Dilanian, R.A., Izumi, F. 2002, unpublishedwork.
[15] Gull, S.F., Daniel, G.J. 1978, Nature (London)
272, 686.
[16] Collins, D.M. 1982, Nature (London) 298, 49.
[17] Bricogne, G. 1988, Acta Crystallogr., Sect. A
44, 517.
[18] Jaynes, E.T. 1957, Phys. Rev. 106, 620.
[19] Kumazawa, S., Takata, M., Sakata, M. 1995,
Acta Crystallogr., Sect. A 51, 47.
[20] Wilkins, S.W., Varghese, J.N., Lehmann, M.S.
1983, Acta Crystallogr., Sect. A 39, 47.
[21] Sakata, M., Uno, T., Takata, M. 1993, J. Appl.
Crystallogr. 26, 159.
[22] Burger, K. 1998, MEND available for download from http://www.uni-tuebingen.de/uni/
pki/maxent/maxent.html.
[23] Tanaka, H., Takata, M., Nishibori, E., Kato, K.,
Iishi, T., Sakata, M. 2002, J. Appl. Crystallogr.
35, 282.
[24] Hahn, Th. (Ed.) 2002, International Tables for
Crystallography, Vol. A, 5th ed., Kluwer Academic Publishers, Dordrecht, Part 7.

[25] Eyck, L.F.T. 1973, Acta Crystallogr., Sect. A
29, 183.
[26] Rowicka, M., Kudlicki, A., Otwinowski, Z. 2002,
Acta Crystallogr., Sect. A 58, 574.
[27] Toraya, H. 1995, The Rietveld Method, R.A.
Young (Ed.), Oxford University Press, Oxford,
Chap. 14.
[28] Pawley, G.S. 1981, J. Appl. Crystallogr. 14, 357.
[29] Le Bail, A., Duroy, H., Fourquet, J.L. 1988,
Mater. Res. Bull. 23, 447.
[30] Ikeda, T., Izumi, F., Kodaira, T., Kamiyama,T. 1998, Chem. Mater. 10, 3996.
[31] Ikeda, T., Akiyama, Y., Izumi, F., Kiyozumi,
Y., Mizukami, F., Kodaira, T. 2001, Chem.
Mater. 13, 1286.
[32] Belik, A.A., Izumi, F., Ikeda, T., Morozov,
V.A., Dilanian, R.A., Torii, S., Kopnin, E.M.,
Lebedev, O.I., Van Tendeloo, G., Lazoryak, B.I.
2002, Chem. Mater. 14, 4464.
[33] Takada, K., Sakurai, H., Takayama-Muromachi,E., Izumi, F., Dilanian, R.A., Sasaki, T. 2003,
Nature (London) 422, 53.
[34] Oikawa, K., Kamiyama, T., Kannno, R., Izumi,
F., Ikeda, T., Chakoumakos, B.C. Mater. Sci.
Forum, in press.
[35] Kanno, R., Ohno, K., Kawamoto, Y., Takeda,
Y., Yamamoto, O., Kamiyama, T., Asano, H.,
Izumi, F., Kondo, S. 1993, J. Solid State Chem.
102, 79.
[36] Kumazawa, S., Takata, M., Sakata, M. 1995,
Acta Crystallogr., Sect. A 51, 651.
[37] Frisch, Æ., Frisch, M.J. 1999, Gaussian 98
User’s Reference, 2nd ed., Gaussian, Inc., Pittsburgh.
[38] Adachi, H., Tsukada, M., Satoko, C. 1978, J.
Phys. Soc. Jpn. 45, 875.
[39] Blaha, P., Schwarz, K., Madsen, G., Kvasnicka,
D., Luitz, J. 2001, WIEN2k, an Augmented
Plane Wave + Local Orbitals Program for Calculating Crystal Properties, K. Schwarz (Ed.),
Technische Universität Wien, Wien.
[40] Palatinus, L., van Smaalen, S. 2002, Acta Crystallogr., Sect. A 58, 559.
[41] Sasaki, T., Kooli, F., Iida, M., Michiue, Y.,
Takenouchi, S., Yajima, Y., Izumi, F., Watanabe, M., Chakoumakos, B.C. 1998, Chem.
Mater. 10, 4123.

[42] Izumi, F., Ikeda, T., Sasaki, T., Kumazawa, S.
2000, Mol. Cryst. Liq. Cryst. 341, 253.
[43] Young, R.A. 1995, The Rietveld Method, R.A.
Young (Ed.), Oxford University Press, Oxford,
Chap. 1.
[44] Nozue, Y., Kodaira, T., Goto, T. 1992, Phys.
Rev. Lett. 68, 3789.
[45] Ikeda, T., Kodaira, T., Izumi, F., Kumazawa,S. 2000, Mol. Cryst. Liq. Cryst. 341, 447.
[46] Ikeda, T., Kodaira, T., Izumi, F., Kamiyama,T., Ohshima, K. 2000, Chem. Phys. Lett. 318,
93.
[47] Yamamoto, A., Hu, W.-Z., Izumi, F., Tajima,
S. 2001, Physica C (Amsterdam) 351, 329.
[48] Izumi, F., Yamamoto, A., Khasanova, N.R.,
Kumazawa, S., Hu, W.-Z., Kamiyama, T. 2000,
Physica C (Amsterdam) 335, 239.
[49] Kumazawa, S., Yamamura, S., Nishibori, E.,
Takata, M., Sakata, M., Izumi, F., Ishii, Y.
1999, J. Phys. Chem. Solids 60, 1407.
[50] Kikuchi, M., Izumi, F., Kikuchi, M., Ohshima,E., Morii, Y., Shimojo, Y., Syono, Y. 1995,
Physica C (Amsterdam) 247, 183.
[51] Dilanian, R.A. et al., unpublished work.
[52] Barrer, R.M., Villiger, H. 1969, Z. Kristallogr.
128, 352.