pyextal.adp¶

Anisotropic displacement parameters and their symmetry constraints.

The Bloch engine has always been anisotropic. The elastic scattering factor in sfsub.f90 applies

\[T(\mathbf{h}) = \exp(-\mathbf{h}^{T} \beta \mathbf{h})\]

with a full symmetric \(\beta\) read from blochwave.cryst.atpar rows 3..8 (Fortran atpar(4:9)), stored in the order (b11, b22, b33, b12, b13, b23). What was missing was a Python path to it: pyextal.optimize.CoarseOptimize.optimizeDWF() used to write only the three diagonal components, from a single isotropic B per atom.

blochwave.cryst.dw is not an independent quantity. It is the equivalent isotropic \(B_{eq} = \frac{4}{3}\mathrm{Tr}(\beta G)\), which xtal0.f90 derives from \(\beta\) at parse time, and only the absorptive form factor uses it. Both therefore have to be written together, which is what set_beta() does.

Conventions, all taken from the Fortran rather than assumed:

Note that the isotropic mapping uses the reciprocal metric tensor \(G^{*}\) including its off-diagonal elements. xatoms in xtal0.f90 does exactly this. Using 1/G_ii for the diagonal and leaving the off-diagonals alone – what the old optimizeDWF did – agrees only for orthogonal cells.

Symmetry¶

An atom’s \(\beta\) is not free. sfsub.f90 builds the orbit of each asymmetric-unit atom by transforming the reflection rather than the atom, so the tensor an orbit-mate effectively carries is \(W \beta W^{T}\). For an operation that maps the atom onto itself, that has to reproduce the tensor it started from:

\[\beta = W \beta W^{T} \quad \forall (W, \mathbf{w}) \in \mathrm{site\ symmetry}\]

allowed_beta_basis() turns those conditions into the basis of the subspace they leave free, so a refinement can vary the site’s real degrees of freedom and nothing else. For SrTiO3 that is one parameter on Sr (1a) and Ti (1b), and two on O (3c, site symmetry 4/mmm): refining all six on the oxygen would be refining four directions the crystal’s symmetry does not have.

Attributes¶

`BETA_ROWS`
`_VOIGT`
`TWO_PI_SQ`

Classes¶

DWFModel

The free temperature-factor parameters of a refinement.

Functions¶

`beta_matrix`(vec)	Expands a 6-component beta into the symmetric 3x3 tensor.
`beta_vector`(mat)	Contracts a symmetric 3x3 tensor to the engine's 6-component order.
`isotropic_beta`(B, gmxr)	Converts an isotropic B to the equivalent beta tensor.
`uij_to_beta`(uij, gmxr)	Converts crystallographic `Uij` to the engine's beta.
`beta_to_uij`(beta, gmxr)	Inverse of `uij_to_beta()`.
`equivalent_B`(beta, gmx)	Equivalent isotropic B, matching `xtal0.f90`.
`space_group_ops`()	Reads the space group operations out of the engine.
`site_symmetry`(position[, ops, tol])	Finds the operations that map an atom onto itself.
`_transform_matrix`(W)	The 6x6 map that `beta -> W beta W^T` induces on the 6-vector.
`beta_projector`(rot)	Averages the constraint over a set of operations.
`allowed_beta_basis`(rot[, tol])	Basis of the beta components a site symmetry leaves free.
`atom_positions`()	Fractional coordinates of the asymmetric-unit atoms.
`atom_labels`()	Element label of each asymmetric-unit atom.
`get_beta`(atom)	Reads one atom's beta out of the engine.
`set_beta`(atom, beta, gmx)	Writes one atom's beta into the engine, keeping `dw` consistent.
`describe`([gmx, gmxr, ops])	Summarises what each atom's symmetry allows.
`atom_line`(label, position, beta[, occupancy])	Formats an atom line for a `.dat` written with `dw = bij`.

Module Contents¶

pyextal.adp.BETA_ROWS = (3, 4, 5, 6, 7, 8)¶

pyextal.adp._VOIGT = ((0, 0), (1, 1), (2, 2), (0, 1), (0, 2), (1, 2))¶

pyextal.adp.TWO_PI_SQ = 19.7392088¶

pyextal.adp.beta_matrix(vec)¶

Expands a 6-component beta into the symmetric 3x3 tensor.

Parameters:: vec (array-like) – (b11, b22, b33, b12, b13, b23).
Returns:: (3, 3) symmetric array.
Return type:: np.ndarray

pyextal.adp.beta_vector(mat)¶

Contracts a symmetric 3x3 tensor to the engine’s 6-component order.

Parameters:: mat (np.ndarray) – (3, 3) symmetric array.
Returns:: (b11, b22, b33, b12, b13, b23).
Return type:: np.ndarray

pyextal.adp.isotropic_beta(B, gmxr)¶

Converts an isotropic B to the equivalent beta tensor.

Parameters:

B (float) – Isotropic Debye-Waller factor, Angstrom^2.
gmxr (np.ndarray) – (3, 3) reciprocal metric tensor, param.gmxr.

Returns:

The 6-component beta.

Return type:

np.ndarray

Note

The off-diagonal components are non-zero whenever the cell is not orthogonal, and they matter – omitting them is not “the isotropic approximation”, it is a different tensor.

pyextal.adp.uij_to_beta(uij, gmxr)¶

Converts crystallographic Uij to the engine’s beta.

Parameters:

uij (array-like) – (U11, U22, U33, U12, U13, U23), Angstrom^2.
gmxr (np.ndarray) – (3, 3) reciprocal metric tensor.

Returns:

The 6-component beta.

Return type:

np.ndarray

pyextal.adp.beta_to_uij(beta, gmxr)¶

Inverse of uij_to_beta().

Parameters:

beta (array-like) – The 6-component beta.
gmxr (np.ndarray) – (3, 3) reciprocal metric tensor.

Returns:

(U11, U22, U33, U12, U13, U23), Angstrom^2.

Return type:

np.ndarray

pyextal.adp.equivalent_B(beta, gmx)¶

Equivalent isotropic B, matching xtal0.f90.

Parameters:

beta (array-like) – The 6-component beta.
gmx (np.ndarray) – (3, 3) direct metric tensor, param.gmx.

Returns:

B_eq = 4/3 Tr(beta @ gmx).

Return type:

float

pyextal.adp.space_group_ops()¶

Reads the space group operations out of the engine.

Returns:: (W, w) of shapes (n, 3, 3) and (n, 3), acting as x' = W x + w on fractional coordinates.
Return type:: tuple[np.ndarray, np.ndarray]

Note

Requires a parsed crystal – call after pyextal.callBloch.bloch_parse().

spgr.gp holds only the coset representatives. Two further pieces have to be folded in to get the group atmult actually generates, and leaving either out silently under-counts site symmetry:

the inversion-related half, when spgr.centro is set;
the spgr.ncsub4 centring translations in spgr.csub4 – without them Si in Fd-3m looks like a site of order 6 with two free temperature-factor components instead of order 24 with one, because the operations that return the atom to itself do so only up to a face-centring vector.

pyextal.adp.site_symmetry(position, ops=None, tol=0.0001)¶

Finds the operations that map an atom onto itself.

Parameters:

position (array-like) – (3,) fractional coordinates.
ops (tuple, optional) – (W, w) from space_group_ops(). Read from the engine when omitted.
tol (float, optional) – Tolerance on the coordinate match, in fractional units. Defaults to 1e-4.

Returns:

(m, 3, 3) array of the rotation parts, always including the identity.

Return type:

np.ndarray

pyextal.adp._transform_matrix(W)¶: The 6x6 map that beta -> W beta W^T induces on the 6-vector.

pyextal.adp.beta_projector(rot)¶

Averages the constraint over a set of operations.

Parameters:: rot (np.ndarray) – (m, 3, 3) rotation parts of a site-symmetry group.
Returns:: (6, 6) projector onto the invariant subspace.
Return type:: np.ndarray

pyextal.adp.allowed_beta_basis(rot, tol=1e-08)¶

Basis of the beta components a site symmetry leaves free.

Parameters:

rot (np.ndarray | list) – (m, 3, 3) operations of one site, or a list of such arrays to constrain several sites to a shared tensor (the intersection of their invariant subspaces).
tol (float, optional) – Singular-value cutoff. Defaults to 1e-8.

Returns:

(6, nfree) orthonormal basis. Columns span the allowed tensors; nfree is 1 for a cubic site and 6 for a site with no symmetry beyond the identity.

Return type:

np.ndarray

pyextal.adp.atom_positions()¶

Fractional coordinates of the asymmetric-unit atoms.

Returns:: (natoms, 3).
Return type:: np.ndarray

pyextal.adp.atom_labels()¶

Element label of each asymmetric-unit atom.

Returns:: natoms labels, e.g. ['SR', 'TI', 'O'].
Return type:: list[str]

pyextal.adp.get_beta(atom)¶

Reads one atom’s beta out of the engine.

Parameters:: atom (int) – Index into the asymmetric unit.
Returns:: The 6-component beta.
Return type:: np.ndarray

pyextal.adp.set_beta(atom, beta, gmx)¶

Writes one atom’s beta into the engine, keeping dw consistent.

Parameters:

atom (int) – Index into the asymmetric unit.
beta (array-like) – The 6-component beta.
gmx (np.ndarray) – (3, 3) direct metric tensor, for the B_eq the absorptive form factor needs.

Note

blochwave.resetugh() still has to be called before the change reaches the structure matrix.

pyextal.adp.describe(gmx=None, gmxr=None, ops=None)¶

Summarises what each atom’s symmetry allows.

Parameters:

gmx (np.ndarray, optional) – Direct metric tensor, for B_eq.
gmxr (np.ndarray, optional) – Reciprocal metric tensor, unused but accepted so callers can pass both.
ops (tuple, optional) – (W, w) from space_group_ops().

Returns:

One entry per atom with its label, position, site-symmetry order, number of free beta components and current values.

Return type:

list[dict]

pyextal.adp.atom_line(label, position, beta, occupancy=1.0)¶

Formats an atom line for a .dat written with dw = bij.

Parameters:

label (str) – Element label.
position (array-like) – (3,) fractional coordinates.
beta (array-like) – The 6-component beta.
occupancy (float, optional) – Site occupancy. Defaults to 1.0.

Returns:

One atom line.

Return type:

str

Note

A file using these lines must declare dw = bij on its crystal line, and then every atom in it needs nine numbers – the parser reads the whole block one way.

class pyextal.adp.DWFModel(groups, bases, x0, bounds, names, gmx, gmxr, aniso)¶

The free temperature-factor parameters of a refinement.

Sits between an optimizer’s flat parameter vector and the per-atom tensors the engine reads. Two things it takes care of that a plain “one B per atom” vector cannot:

Symmetry. In anisotropic mode each site is refined only along the directions allowed_beta_basis() says are free, so a cubic site keeps one parameter and a 4/mmm site gets two – never six.
Units. Every coefficient is rescaled to Angstrom^2 of equivalent B, so a vector of coefficients sums to \(B_{eq}\) for that site, the bounds mean the same thing they did for the isotropic refinement, and an optimizer’s step length is comparable across parameters. Directions that carry no trace (pure shear) are instead scaled to the magnitude of an isotropic tensor with B = 1, and bounded symmetrically about zero.

Variables:

groups (list[list[int]]) – Atom indices sharing each set of parameters.
bases (list[np.ndarray]) – (6, nk) rescaled basis per group.
x0 (np.ndarray) – Starting coefficients, read from the current crystal.
bounds (list[tuple]) – Per-coefficient bounds.
names (list[str]) – Human-readable name of each coefficient.

groups¶

bases¶

x0¶

bounds¶

names¶

gmx¶

gmxr¶

aniso¶

slices = []¶

property n¶

Number of free parameters.

Type:: int

classmethod build(param, aniso=False, group='atom', boundary=(0.0, 2.0))¶

Constructs the model for the crystal currently in the engine.

Parameters:

param (SimParams) – Supplies gmx and gmxr.
aniso (bool, optional) – Refine the symmetry-allowed anisotropic components rather than one isotropic B. Defaults to False.
group (str | list, optional) – 'atom' for one set of parameters per asymmetric-unit atom, 'type' to tie atoms of the same element together, or an explicit list of index lists. Defaults to 'atom'.
boundary (tuple, optional) – Bounds on a B-like coefficient. Defaults to (0.0, 2.0).

Returns:

The model, starting from the crystal’s current tensors.

Return type:

DWFModel

betas(x)¶

The beta tensor each group takes for a parameter vector.

Parameters:: x (array-like) – Coefficients.
Returns:: One 6-component beta per group.
Return type:: list[np.ndarray]

penalty(x)¶

Rejects unphysical tensors without running a simulation.

A beta that is not positive definite means an imaginary mean-square displacement. Catching it here costs microseconds; letting it reach the engine costs a full diagonalization and returns a meaningless number.

Parameters:: x (array-like) – Coefficients.
Returns:: A penalty above any attainable objective value, or None when every tensor is positive definite.
Return type:: float | None

apply(x)¶

Writes a parameter vector into the engine.

Parameters:: x (array-like) – Coefficients.

Note

blochwave.resetugh() must follow before the change reaches the structure matrix.

report(x)¶

A one-line-per-group summary of a parameter vector.

Parameters:: x (array-like) – Coefficients.
Returns:: Formatted table of coefficients, B_eq and beta.
Return type:: str