Hex, Bugs and More Physics | Emre S. Tasci

A very simple code I wrote while studying the Python language that fetches the information from a VASP POSCAR file and prepares (and feeds) the input for Harold Stokes’ findsym (you should actually download the executable included in the ISOTROPY package instead of the referred web version) program to find the unit cell symmetry.

#!/usr/bin/python
# -*- coding: utf-8 -*-

"""
POSCAR2findsym.py
/* Emre S. Tasci <e.tasci@tudelft.nl>                    *
   A simple python code that parses the POSCAR file,
   prepares an input file for Harold Stokes' findsym
   (http://stokes.byu.edu/isotropy.html)
   executable.

   Then, checks the current system to see if findsym is
   accessible: if it is, then feeds the prepared file 
   and directly outputs the result; if findsym is not 
   executable then outputs the generated input file.

   If you are interested only in the results section of 
   the findsymm output, you should append "| grep "\-\-\-\-\-\-" -A1000"
   to the execution of this code.

   Accepts alternative POSCAR filename for argument
   (default = POSCAR)

   If a second argument is proposed then acts as if
   findsym is not accessible (ie, prints out the input
   file instead of feeding it to findsymm)

                                             10/03/09 */
"""

import sys
import re
from numpy import *
from os import popen
import commands

force_output = False

if len(sys.argv)<2:
    filename="POSCAR"
elif len(sys.argv)<3:
    filename=sys.argv[1]
else:
    filename=sys.argv[1]
    force_output = True

f = open(filename,'r')

title = f.readline().strip()
tolerance = 0.00001
latt_vec_type = 1 # We will be specifying in vectors

f.readline() # Read Dummy
lat_vec = array([])
for i in range(0,3):
    lat_vec = append(lat_vec,f.readline().strip())

# Read atom species
species_count = f.readline()
species_count = re.split('[\ ]+', species_count.strip())
species_count = map(int,species_count)
number_of_species = len(species_count)
total_atoms = sum(species_count)

initial_guess = "P" # For the centering

species_designation = array([])
for i in range(1,number_of_species+1):
    for j in range(0,species_count[i-1]):
        species_designation = append(species_designation,i)

species_designation = map(str,map(int,species_designation))
# print species_designation
sep = " "
# print sep.join(species_designation)

# Read Coordinates
f.readline() # Read Dummy
pos_vec = array([])
for i in range(0,total_atoms):
    pos_vec = append(pos_vec,f.readline().strip())

# print pos_vec

findsym_input = array([])
findsym_input = append(findsym_input, title)
findsym_input = append(findsym_input, tolerance)
findsym_input = append(findsym_input, latt_vec_type)
findsym_input = append(findsym_input, lat_vec)
findsym_input = append(findsym_input, 2)
findsym_input = append(findsym_input, initial_guess)
findsym_input = append(findsym_input, total_atoms)
findsym_input = append(findsym_input, species_designation)
findsym_input = append(findsym_input, pos_vec)
# print findsym_input
sep = "\n"
findsym_input_txt = sep.join(findsym_input)
# print findsym_input_txt

# Check if findsym is accessible:
status,output =  commands.getstatusoutput("echo \"\n\n\"|findsym ")
if(output.find("command not found")!=-1 or force_output):
    # if it is not there, then just output the input
    print findsym_input_txt
elif(status==6144):
    # feed it to findsym
    pipe = popen("echo \""+findsym_input_txt+"\" | findsym").readlines()
    for line in pipe:
        print line.strip()
quit()

Example Usage:

sururi@husniya POSCAR2findsym $ cat POSCAR
Si(S)-Au(Pd) -- Pd3S
   1.00000000000000
     4.7454403619558345    0.0098182468538853    0.0000000000000000
    -1.4895902658503473    4.5055989020479856    0.0000000000000000
     0.0000000000000000    0.0000000000000000    7.2191545483192190
   6   2
Direct
  0.1330548697855782  0.5102695954022698  0.2500000000000000
  0.4897304045977232  0.8669451452144267  0.7500000000000000
  0.5128136657304309  0.1302873334247993  0.2500000000000000
  0.8697126665752007  0.4871863042695738  0.7500000000000000
  0.0013210250693640  0.9986789749306360  0.0000000000000000
  0.0013210250693640  0.9986789749306360  0.5000000000000000
  0.6856021298170862  0.6825558526447357  0.2500000000000000
  0.3174442073552622  0.3143978111829124  0.7500000000000000

sururi@husniya POSCAR2findsym $ python POSCAR2findsym.py
FINDSYM, Version 3.2.3, August 2007
Written by Harold T. Stokes and Dorian M. Hatch
Brigham Young University
 
Si(S)-Au(Pd) -- Pd3S
Tolerance:    0.00001
Lattice vectors in cartesian coordinates:
4.74544   0.00982   0.00000
-1.48959   4.50560   0.00000
0.00000   0.00000   7.21915
Lattice parameters, a,b,c,alpha,beta,gamma:
4.74545   4.74545   7.21915  90.00000  90.00000 108.17579
Centering: P
Number of atoms in unit cell:
8
Type of each atom:
1  1  1  1  1  1  2  2
Position of each atom (dimensionless coordinates)
1   0.13305   0.51027   0.25000
2   0.48973   0.86695   0.75000
3   0.51281   0.13029   0.25000
4   0.86971   0.48719   0.75000
5   0.00132   0.99868   0.00000
6   0.00132   0.99868   0.50000
7   0.68560   0.68256   0.25000
8   0.31744   0.31440   0.75000
------------------------------------------
Space Group  40  C2v-16    Ama2
Origin at    0.00000   0.00000   0.50000
Vectors a,b,c:
0.00000   0.00000   1.00000
-1.00000  -1.00000   0.00000
1.00000  -1.00000   0.00000
Values of a,b,c,alpha,beta,gamma:
7.21915   5.56683   7.68685  90.00000  90.00000  90.00000
Atomic positions in terms of a,b,c:
Wyckoff position b, y = -0.17834, z =  0.31139
1   0.75000   0.17834   0.31139
2   0.25000   0.82166   0.31139
Wyckoff position b, y =  0.32155, z =  0.19126
3   0.75000   0.17845   0.69126
4   0.25000   0.82155   0.69126
Wyckoff position a, z =  0.00132
5   0.50000   0.00000   0.00132
6   0.00000   0.00000   0.00132
Wyckoff position b, y = -0.31592, z =  0.00152
7   0.75000   0.81592   0.50152
8   0.25000   0.18408   0.50152