Re: [Spglib-users] symmetrizing positions/cell vectors without altering cell?

[Atsushi T. <atz...@gm...>]

Hi,

Thanks for making your script simple.
I see the documentation is confusing. dataset['std_lattice'] includes
rotation, but the transformation matrix doesn't rotate.

I answer to your questions 1 and 2. 1 yes. 2, the mathematical
expression must contain rotation matrix.

I wrote a script to do this job.

import spglib
import numpy as np

# Original basis vectors in row vectors
orig_basis = np.array([[0,1.41421356237309504880,0],
                       [.70710678118654752440,.70710678118654752440,1],
                       [-.70710678118654752440,.70710678118654752440,1]]).T
positions = [[0,0,0]]
atomic_numbers = [32]

cell = (orig_basis.T, positions, atomic_numbers)
dataset = spglib.get_symmetry_dataset(cell, symprec=0.1)
space_group_type = dataset['international']
centring = space_group_type[0]
t_mat = dataset['transformation_matrix']
F_centre_t_mat = np.array([[0, 1, 1],
                           [1, 0, 1],
                           [1, 1, 0]]) * 0.5 # P_F in spglib doc
rotated_std_basis = dataset['std_lattice'].T # in row vectors
std_basis = np.dot(orig_basis, np.linalg.inv(t_mat)) # in row_vectors
R = np.dot(rotated_std_basis, np.linalg.inv(std_basis)) # rotation matrix

print("space group type %s: %s-centring" % (space_group_type, centring))
print("orig basis (row vectors)\n", orig_basis)
print("transformation  matrix\n", t_mat)
print("standardized basis (row vectors)\n", std_basis)
print("standardized basis with rotation in Cartesian coordinates "
      "(row vectors)\n", rotated_std_basis)
print("Rotation matrix (R)\n", R)
print("det(R) = %f" % np.linalg.det(R))

# Tranformation back
# Rotation matrix has to be a rigit rotation matrix.
# For example, I used Euler angles to polish the above rotation matrix.
# I didn't examine carefully if my code is general or not.
# I installed http://matthew-brett.github.io/transforms3d/
from transforms3d.euler import mat2euler, euler2mat
R_new = euler2mat(*mat2euler(R))

recovered_std_basis = np.dot(np.linalg.inv(R_new), rotated_std_basis)
# Transformation matrix may be distorted.
# So an additional treatment is applied.
recovered_std_prim_basis = np.dot(recovered_std_basis, F_centre_t_mat)
t_mat_to_orig_basis_from_prim_basis = np.dot(
    np.linalg.inv(recovered_std_prim_basis), orig_basis)
# The matrix elements of this matrix have to be integers, so let's cast to int
t_mat_int = np.rint(t_mat_to_orig_basis_from_prim_basis).astype(int)
assert abs(abs(np.linalg.det(t_mat_int)) - 1) < 1e-10

final_basis = np.dot(recovered_std_prim_basis, t_mat_int)
print("final basis (row vectors)\n", final_basis)


The output is as follows:

space group type Fm-3m: F-centring
orig basis (row vectors)
 [[ 0.          0.70710678 -0.70710678]
 [ 1.41421356  0.70710678  0.70710678]
 [ 0.          1.          1.        ]]
transformation  matrix
 [[  0.00000000e+00  -5.00000000e-01  -5.00000000e-01]
 [  5.00000000e-01   5.00000000e-01  -7.51373070e-17]
 [ -5.00000000e-01   1.96261557e-17  -5.00000000e-01]]
standardized basis (row vectors)
 [[  0.00000000e+00   1.41421356e+00   1.41421356e+00]
 [ -2.22044605e-16   1.41421356e+00  -1.41421356e+00]
 [ -2.00000000e+00   0.00000000e+00   0.00000000e+00]]
standardized basis with rotation in Cartesian coordinates (row vectors)
 [[ 2.  0.  0.]
 [ 0.  2.  0.]
 [ 0.  0.  2.]]
Rotation matrix (R)
 [[  0.00000000e+00   0.00000000e+00  -1.00000000e+00]
 [  7.07106781e-01   7.07106781e-01  -7.85046229e-17]
 [  7.07106781e-01  -7.07106781e-01   7.85046229e-17]]
det(R) = 1.000000
final basis (row vectors)
[[  3.33066907e-16   7.07106781e-01  -7.07106781e-01]
[  1.41421356e+00   7.07106781e-01   7.07106781e-01]
[  0.00000000e+00   1.00000000e+00   1.00000000e+00]]


Togo

On Tue, Dec 19, 2017 at 10:49 PM, Noam Bernstein
<noa...@nr...> wrote:
> On Dec 18, 2017, at 7:02 PM, Atsushi Togo <atz...@gm...> wrote:
>
> Your sample script is busy, and I want to see only essential
> information. In addition I don't have ASE in my computer.
>
>
> Sorry - I tried to make sure it was printing everything that I thought was
> helpful to understand the behavior, but here’s a version with the extra
> print statements removed and without ASE:
>
> import spglib, numpy as np
> orig_cell=np.array([[0,1.41421356237309504880,0],
>                     [.70710678118654752440,.70710678118654752440,1],
>                     [-.70710678118654752440,.70710678118654752440,1]])
> positions=[[0,0,0]]
> atomic_numbers=[32]
>
> dataset = spglib.get_symmetry_dataset((orig_cell, positions,
> atomic_numbers), symprec=0.1)
> print "orig cell\n", orig_cell
> print "supposed to be orig cell\n", np.dot(dataset['std_lattice'].T,
> dataset['transformation_matrix']).T
>
>  output:
>
> orig cell
> [[ 0.          1.41421356  0.        ]
>  [ 0.70710678  0.70710678  1.        ]
>  [-0.70710678  0.70710678  1.        ]]
> supposed to be orig cell
> [[  0.00000000e+00   1.00000000e+00  -1.00000000e+00]
>  [ -1.00000000e+00   1.00000000e+00   3.92523115e-17]
>  [ -1.00000000e+00  -1.50274614e-16  -1.00000000e+00]]
>
> Obviously, the original cell and what I expected would be the same as the
> original cell are not in fact the same, so clearly I’m misunderstanding
> something.
>
>
> By the way, do you understand the difference between transformation
> and rotation? The former is a change of basis and the later is rigid
> rotation of whole crystal in Cartesian coordinates. The former is
> applied from the right and the later is applied from the left if I am
> correct.
>
>
> Thanks for the clarification.  I’m just trying to understand why the
> documentation (at
> https://atztogo.github.io/spglib/api.html#api-origin-shift-and-transformation)
> says
>
> (a b c) = (a_s b_s c_s) P
>
> where a, b, and c are the original cell vectors, a_s, b_s, and c_s are the
> standardized cell vectors returned by get_symmetry_dataset, and P is the
> transformation matrix returned by get_symmetry_dataset.  However, when I try
> to do that multiplication by transforming the cell 3x3 matrices which have
> each lattice vector as a row into column vectors (which is my understanding
> of the expression in the docs), and right multiplying by the transformation
> matrix, the equality is not satisfied.
>
> 1. Do you agree that it should be possible to get two identical 3x3 matrices
> (assuming the initial symmetry is in fact perfect), one of which is the
> original cell vectors and the other is a product of some simple functions
> (e.g. transpose and/or inverse) of dataset[‘std_lattice’] and
> dataset[‘transformation matrix’]?
> 2. Do you agree that what my code above implements is indeed the
> mathematical expression from the documentation?  If not, what is the correct
> operation?
>
>
> thanks,
> Noam



-- 
Atsushi Togo
Elements Strategy Initiative for Structural Materials, Kyoto university
E-mail: atz...@gm...