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

[Noam B. <noa...@nr...>]

> On Sep 25, 2017, at 9:15 AM, Atsushi Togo <atz...@gm...> wrote:
> 
> If the distortion is not that much, then the approach by
> Grosse-Kunstleve and Adams should work. If you write a python script,
> I think I can comment to it. The C-code is about here,
> https://github.com/atztogo/spglib/blob/master/src/site_symmetry.c#L244 <https://github.com/atztogo/spglib/blob/master/src/site_symmetry.c#L244>


Hmm.  I tried something like this initially, and it wasn't working at all, but I must have had a bug, since I redid it cleanly (so I could show the code), and the new version's apparently working.  Maybe because I was coding for an audience I was more careful.  In any case, while I think it would be useful to have a version of standardize_cell() that only moved atoms, without touching the cell vectors, I think my problem is solved for now (especially if I can recode it without a tight loop in python).  It should probably be redone without so many tight loops in explicit python, but here it is, in case it’s useful for anyone else:

#!/usr/bin/env python

import os, sys
import spglib, ase, ase.io
import numpy as np

at = ase.io.read(sys.argv[1])

dataset = spglib.get_symmetry_dataset(at, symprec=1.0e-5)
sys.stderr.write("precise initial symmetry group number {}, international (Hermann-Mauguin) {} Hall {} prec {}\n".format(dataset["number"],dataset["international"],dataset["hall"],1.0e-5))
dataset = spglib.get_symmetry_dataset(at, symprec=0.1)
sys.stderr.write("loose initial symmetry group number {}, international (Hermann-Mauguin) {} Hall {} prec {}\n".format(dataset["number"],dataset["international"],dataset["hall"],0.1))

p = at.get_scaled_positions()
cell = at.get_cell()

new_position = np.zeros(at.get_positions().shape)
for i in range(len(at)):
    # zero accumulators for rot and transl
    mean_r = np.zeros((3,3))
    mean_t = np.zeros((3))
    n_contrib=0
    position = p[i][:]
    # loop over opreations
    for (r, t) in zip(dataset['rotations'], dataset['translations']):
        # transformed position
        pos = np.dot(r, position) + t
        # if it's close to periodic image, contribute to accumulators
        dp = pos - position
        dp_in_cell = dp - np.round(dp)
        if np.linalg.norm(dp_in_cell) < 0.1:
            mean_r += r
            mean_t += t - np.round(dp)
            n_contrib += 1
    mean_r /= float(n_contrib)
    mean_t /= float(n_contrib)

    new_position[i,:] = np.dot(position, mean_r) + mean_t

at.set_scaled_positions(new_position)
dataset = spglib.get_symmetry_dataset(at, symprec=1.0e-5)
sys.stderr.write("manual symmetrization got symmetry group number {}, international (Hermann-Mauguin) {} Hall {} prec {}\n".format(dataset["number"],dataset["international"],dataset["hall"],1.0e-5))
ase.io.write(sys.stdout, at, format="extxyz")


thanks again,
Noam