Skip to content

Latest commit

 

History

History
401 lines (329 loc) · 12.5 KB

examples.md

File metadata and controls

401 lines (329 loc) · 12.5 KB

Examples

Pages = ["examples.md"]
Depth = 2

Creating a Cell

To create a Cell, we first need to create a Lattice. Then we can add atoms and their positions (in crystal coordinates):

using Spglib
lattice = [
    -3.0179389205999998 -3.0179389205999998 0.0000000000000000
    -5.2272235447000002 5.2272235447000002 0.0000000000000000
    0.0000000000000000 0.0000000000000000 -9.7736219469000005
]
positions = [[2 / 3, 1 / 3, 1 / 4], [1 / 3, 2 / 3, 3 / 4]]
atoms = [1, 1]
cell = Cell(lattice, positions, atoms)

Crystallographic choice and rigid rotation

The following example of a python script gives a crystal structure of Br whose space group type is Cmce. The basis vectors \begin{bmatrix} \mathbf{a} & \mathbf{b} & \mathbf{c} \end{bmatrix} are fixed by the symmetry crystal in the standardization. The C-centering determines the c-axis, and m and c operations in Cmce fix which directions a- and b-axes should be with respect to each other axis. This is the first one choice appearing in the list of Hall symbols among 6 different choices for this space group type.

using Spglib
lattice = Lattice([[7.17851431, 0, 0],  # a
                   [0, 3.99943947, 0],  # b
                   [0, 0, 8.57154746]])  # c
positions = [[0.0, 0.84688439, 0.1203133],
             [0.0, 0.65311561, 0.6203133],
             [0.0, 0.34688439, 0.3796867],
             [0.0, 0.15311561, 0.8796867],
             [0.5, 0.34688439, 0.1203133],
             [0.5, 0.15311561, 0.6203133],
             [0.5, 0.84688439, 0.3796867],
             [0.5, 0.65311561, 0.8796867]];
atoms = fill(35, length(positions));
cell = Cell(lattice, positions, atoms)
dataset = get_dataset(cell, 1e-5)

we get

print("Space group type: ", dataset.international_symbol)
print("Space group number: ", dataset.spacegroup_number)
print("Transformation matrix: ")
dataset.transformation_matrix
print("Origin shift: ", dataset.origin_shift)

No rotation was introduced in the idealization. Next, we swap the a- and c-axes.

lattice = Lattice([[8.57154746, 0, 0],  # a
                   [0, 3.99943947, 0],  # b
                   [0, 0, 7.17851431]])  # c
positions = [[0.1203133, 0.84688439, 0.0],
             [0.6203133, 0.65311561, 0.0],
             [0.3796867, 0.34688439, 0.0],
             [0.8796867, 0.15311561, 0.0],
             [0.1203133, 0.34688439, 0.5],
             [0.6203133, 0.15311561, 0.5],
             [0.3796867, 0.84688439, 0.5],
             [0.8796867, 0.65311561, 0.5]];
atoms = fill(35, length(positions));
cell = Cell(lattice, positions, atoms)
dataset = get_dataset(cell, 1e-5)

By this, we get

print("Space group type: ", dataset.international_symbol)
print("Space group number: ", dataset.spacegroup_number)
print("Transformation matrix: ")
dataset.transformation_matrix
print("Origin shift: ", dataset.origin_shift)

We get a non-identity transformation matrix, which wants to transform back to the original (above) crystal structure by swapping a- and c-axes. The transformation back of the basis vectors is achieved by

$$\begin{bmatrix} \mathbf{a} & \mathbf{b} & \mathbf{c} \end{bmatrix} = \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix} \mathbf{P},$$

Next, we try to rotate rigidly the crystal structure by 45^\circ around c-axis in Cartesian coordinates from the first one:

lattice = Lattice([[5.0759761474456697, 5.0759761474456697, 0],  # a
                   [-2.8280307701821314, 2.8280307701821314, 0],  # b
                   [0, 0, 8.57154746]])  # c
positions = [[0.0, 0.84688439, 0.1203133],
             [0.0, 0.65311561, 0.6203133],
             [0.0, 0.34688439, 0.3796867],
             [0.0, 0.15311561, 0.8796867],
             [0.5, 0.34688439, 0.1203133],
             [0.5, 0.15311561, 0.6203133],
             [0.5, 0.84688439, 0.3796867],
             [0.5, 0.65311561, 0.8796867]];
atoms = fill(35, length(positions));
cell = Cell(lattice, positions, atoms)
dataset = get_dataset(cell, 1e-5)

and

print("Space group type: ", dataset.international_symbol)
print("Space group number: ", dataset.spacegroup_number)
print("Transformation matrix: ")
dataset.transformation_matrix
print("Origin shift: ", dataset.origin_shift)

The transformation matrix is kept unchanged even though the crystal structure is rotated in Cartesian coordinates. The origin shift is different but it changes only the order of atoms, so effectively it does nothing.

Transformation to a primitive cell

There are infinite number of choices of primitive cell. The transformation from a primitive cell basis vectors to the other primitive cell basis vectors is always done by an integer matrix because any lattice points can be generated by the linear combination of the three primitive basis vectors.

When we have a non-primitive cell basis vectors as given in the above example:

using Spglib
lattice = Lattice([[7.17851431, 0, 0],  # a
                   [0, 3.99943947, 0],  # b
                   [0, 0, 8.57154746]])  # c

This has the C-centring, so it must be transformed to a primitive cell. A possible transformation is shown at Transformation to the primitive cell, which is \mathbf{P}_\text{C}. With the following script:

Pc = [
        1//2 1//2 0
        -1//2 1//2 0
        0 0 1
    ]
primitive_lattice = lattice * Pc

we get the primitive cell basis vectors primitive_lattice.

find_primitive gives a primitive cell that is obtained by transforming standardized and idealized crystal structure to the primitive cell using the transformation matrix. Therefore by this script we get:

positions = [
    [0.0, 0.84688439, 0.1203133],
    [0.0, 0.65311561, 0.6203133],
    [0.0, 0.34688439, 0.3796867],
    [0.0, 0.15311561, 0.8796867],
    [0.5, 0.34688439, 0.1203133],
    [0.5, 0.15311561, 0.6203133],
    [0.5, 0.84688439, 0.3796867],
    [0.5, 0.65311561, 0.8796867],
];
atoms = fill(8, length(positions));
cell = Cell(lattice, positions, atoms)
find_primitive(cell).lattice

This is same as what we manually obtained above. Even when the basis vectors are rigidly rotated as:

new_lattice = Lattice([[5.0759761474456697, 5.0759761474456697, 0],
                       [-2.8280307701821314, 2.8280307701821314, 0],
                       [0, 0, 8.57154746]])

the relationship of a, b, c axes is unchanged. Therefore the same transformation matrix to the primitive cell can be used. Then we get:

new_lattice * Pc

However applying find_primitive rigidly rotates automatically and so the following script doesn't give this basis vectors:

lattice = Lattice([
    [5.0759761474456697, 5.0759761474456697, 0],
    [-2.8280307701821314, 2.8280307701821314, 0],
    [0, 0, 8.57154746],
])
positions = [
    [0.0, 0.84688439, 0.1203133],
    [0.0, 0.65311561, 0.6203133],
    [0.0, 0.34688439, 0.3796867],
    [0.0, 0.15311561, 0.8796867],
    [0.5, 0.34688439, 0.1203133],
    [0.5, 0.15311561, 0.6203133],
    [0.5, 0.84688439, 0.3796867],
    [0.5, 0.65311561, 0.8796867],
];
atoms = fill(8, length(positions));
cell = Cell(lattice, positions, atoms)

but gives those with respect to the idealized ones:

find_primitive(cell).lattice

To obtain the rotated primitive cell basis vectors, we can use standardize_cell as shown below:

standardize_cell(cell, to_primitive=true, no_idealize=true).lattice

which is equivalent to that we get manually. However, using standardize_cell, distortion is not removed for the distorted crystal structure.

Computing rigid rotation introduced by idealization

This example is from here.

In this package, rigid rotation is purposely introduced in the idealization step though this is unlikely as a crystallographic operation.

using StaticArrays, Spglib
lattice = Lattice([
    [5.0759761474456697, 5.0759761474456697, 0],
    [-2.8280307701821314, 2.8280307701821314, 0],
    [0, 0, 8.57154746],
]);
positions = [
    [0.0, 0.84688439, 0.1203133],
    [0.0, 0.65311561, 0.6203133],
    [0.0, 0.34688439, 0.3796867],
    [0.0, 0.15311561, 0.8796867],
    [0.5, 0.34688439, 0.1203133],
    [0.5, 0.15311561, 0.6203133],
    [0.5, 0.84688439, 0.3796867],
    [0.5, 0.65311561, 0.8796867],
];
atoms = fill(35, length(positions));
cell = Cell(lattice, positions, atoms)
dataset = get_dataset(cell, 1e-5)
dataset.international_symbol
dataset.spacegroup_number
dataset.transformation_matrix

We can see the transformation matrix from the given lattice to the standardized lattice is the identity matrix, i.e., the given lattice is already a standardized lattice.

std_lattice_before_idealization = convert(Matrix{Float64}, lattice) * inv(dataset.transformation_matrix)

This is based on formula in Transformation matrix \mathbf{P} and origin shift \mathbf{p} and Passive/forward/alias transformation:

$$\begin{align} \begin{bmatrix} \mathbf{a} & \mathbf{b} & \mathbf{c} \end{bmatrix} &= \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix} \mathbf{P},\\\ \therefore \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix} &= \begin{bmatrix} \mathbf{a} & \mathbf{b} & \mathbf{c} \end{bmatrix} \mathbf{P}^{-1}. \end{align}$$

Here, \mathbf{P} is the dataset.transformation_matrix.

Note that in contrast to the Python code:

std_lattice_before_idealization = np.dot(
    np.transpose(lattice),
    np.linalg.inv(dataset['transformation_matrix'])).T

where there are multiple transpose operations, we do not have to do that in our Julia code since we choose a column-major order of stacking lattice vectors as described in Basis vectors, and we return transformation matrix in column-major order, too.

Now, we obtain the standardized basis vectors after idealization \begin{bmatrix} \bar{\mathbf{a}}_\text{s} & \bar{\mathbf{b}}_\text{s} & \bar{\mathbf{c}}_\text{s} \end{bmatrix}:

std_lattice_after_idealization = dataset.std_lattice

This is different from the standardized basis vectors before idealization \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix}. Unless this crystal structure is distorted from the crystal structure that has the ideal symmetry, this means that the crystal was rotated rigidly in the idealization step by

$$\begin{bmatrix} \bar{\mathbf{a}}_\text{s} & \bar{\mathbf{b}}_\text{s} & \bar{\mathbf{c}}_\text{s} \end{bmatrix} = \mathbf{R} \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix},$$

as stated in Rotation introduced by idealization. where \mathbf{R} is the rotation matrix. This is computed by

$$\mathbf{R} = \begin{bmatrix} \bar{\mathbf{a}}_\text{s} & \bar{\mathbf{b}}_\text{s} & \bar{\mathbf{c}}_\text{s} \end{bmatrix} \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix}^{-1}$$

In Julia code, this is

𝐀 = convert(Matrix{Float64}, std_lattice_after_idealization)
𝐁 = convert(Matrix{Float64}, std_lattice_before_idealization)
𝐑 = 𝐀 * inv(𝐁)

Note also the transpose is not applied here in contrast to the Python code:

R = np.dot(dataset['std_lattice'].T, np.linalg.inv(std_lattice_before_idealization.T))

This equals to

$$\begin{bmatrix} \cos \theta & -\sin \theta & 0\\\ \sin \theta & \cos \theta & 0\\\ 0 & 0 & 1 \end{bmatrix}$$

where \theta = -\pi/4 and \det(\mathbf{R}) = 1 when no distortion:

θ = -π/4
[
    cos(θ) -sin(θ) 0
    sin(θ) cos(θ) 0
    0 0 1
]

Compared to dataset.std_rotation_matrix:

dataset.std_rotation_matrix

we have approximately the same result.

In summary of the two steps,

$$\begin{align} \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix} = \begin{bmatrix} \mathbf{a} & \mathbf{b} & \mathbf{c} \end{bmatrix} \mathbf{P}^{-1} &= \mathbf{R}^{-1} \begin{bmatrix} \bar{\mathbf{a}}_\text{s} & \bar{\mathbf{b}}_\text{s} & \bar{\mathbf{c}}_\text{s} \end{bmatrix},\\\ \begin{bmatrix} \bar{\mathbf{a}}_\text{s} & \bar{\mathbf{b}}_\text{s} & \bar{\mathbf{c}}_\text{s} \end{bmatrix} \mathbf{P} &= \mathbf{R} \begin{bmatrix} \mathbf{a}_\text{s} & \mathbf{b}_\text{s} & \mathbf{c}_\text{s} \end{bmatrix}. \end{align}$$